pertone.m

ADSL经典仿真程序

   s = get(gco,'String');
   vv = eval(s,num2str(v(1)));
   % Dmax cannot be smaller than 1 and Dmin
   if vv < 1 | vv < Dmin, 
      vv = v;
   end
   set(gco,'Userdata',vv,'String',num2str(vv))
   return
   
   % set N   
elseif strcmp(action,'setN'),
   hndlList=get(gcf,'Userdata');
   v = get(gco,'Userdata');
   s = get(gco,'String');
   vv = eval(s,num2str(v(1)));
   % N has to be a power of 2
   if (floor(log(vv)/log(2))-(log(vv)/log(2))) ~= 0
      vv = v; 
   end
   set(gco,'Userdata',vv,'String',num2str(vv))
   return
   
 elseif strcmp(action,'setHDelaySearch'),
   hndlList=get(gcf,'Userdata');
   hd=hndlList(22);
   od=hndlList(23);
  set(hd,'value',get(hd,'Max'));
  set(hd,'userdata','1');
  set(od,'value',get(od,'Min'));
  set(od,'userdata','0');
         
   return
   
 elseif strcmp(action,'setODelaySearch'),
   hndlList=get(gcf,'Userdata');
   hd=hndlList(22);
   od=hndlList(23);
  set(od,'value',get(od,'Max'));
  set(od,'userdata','1');
  set(hd,'value',get(hd,'Min'));
  set(hd,'userdata','0');
   return
   
   
   % design TEQ   
elseif strcmp(action,'design'),
   set(gcf,'Pointer','watch');
   % get handles
   axHndl=gca;
   fhndlList=get(gcf,'Userdata');
   freqzHndl = fhndlList(1);
   TGHndl = fhndlList(2);
   NwHndl = fhndlList(3);
   NHndl = fhndlList(4);
   CGHndl = fhndlList(5);
   MarginHndl = fhndlList(6);
   DminHndl = fhndlList(7);
   DmaxHndl = fhndlList(8);
   PwrHndl = fhndlList(9);
   AWGNHndl = fhndlList(10);
   CSAHndl = fhndlList(11);
   viewHndl = fhndlList(12);
   helpHndl = fhndlList(13);
   closeHndl = fhndlList(14);
   methodHndl = fhndlList(15);
   methRateHndl = fhndlList(16);
   methSNRHndl = fhndlList(17);
   methSSNRHndl = fhndlList(18);
   methMSEHndl = fhndlList(19);
   methDelayHndl = fhndlList(20);
   methMaxHndl = fhndlList(21);
   HDHndl=fhndlList(22);
   ODHndl=fhndlList(23);
   
   colors = get(gca,'colororder'); 
  
   % init variables
   TG = get(TGHndl,'userdata');
   Nw = get(NwHndl,'userdata');
   N = get(NHndl,'userdata');
   codingGain = get(CGHndl,'userdata');
   margin = get(MarginHndl,'userdata');
   Dmin = get(DminHndl,'userdata');
   Dmax = get(DmaxHndl,'userdata');
   totalInputPower = get(PwrHndl,'userdata');
   AWGNpower = get(AWGNHndl,'userdata');
   loopNum = get(CSAHndl,'userdata');
   method = get(methodHndl,'value');
   HD=get(HDHndl,'value');
   OD=get(ODHndl,'value');
   if (HD==1)&(OD==0)
       delaySearch=1; %delaySearch=1 => heuristic delay search
   elseif (HD==0)&(OD==1)
       delaySearch=0; %delaySearch=0 => optimal delay search
   else
       delaySearch=1;
   end
   barflag = 1;
   chflag=1;
   %number of iteration per TEQ tap
   numIter = 1000;

   
   %generate signal,noise, and channel data
   [recNoisySig,receivedSignal,noise,channel,inputSignal,gamma,fs,NEXTnoise, AWGN, FX] = ...
      siggen_pertone(N,AWGNpower,loopNum,totalInputPower,codingGain,margin,barflag,chflag);

  %save tmp1.mat inputSignal, receivedSignal, noise;

   %estimate the power spectrums
   [inputSpecAll, noiseSpecAll, channelGainAll]=...
      specestim(inputSignal,noise,channel,N,barflag);
%     figure
%     plot(inputSpecAll);
   %calculate SNRs and usable channels
   [subMFBall,subMFB,usedSubs,noiseSpec,channelGain,inputSpec] = ...
   calcsnrs(inputSpecAll,channelGainAll,noiseSpecAll,gamma);

   %chose method and design pertone
   [subSNR, Ps, Pn, WM, D, title_str, py, pX, rateVector] = selmeth_pte(method,FX,...
   receivedSignal,noise,channel, Nw, N, barflag, AWGN, delaySearch, usedSubs, Dmin, Dmax, margin, codingGain);
   

%    Xhat = pertoneeq(WM.',py);
% 
%    [subSNR, Ps, Pn] = snr_pte(pX, Xhat, usedSubs);
   
   
   bmap(usedSubs) = ba_cal(subSNR(usedSubs), 9.8, margin, codingGain, 0, 100, 0);
   
   RDMTfinal = fs/(N+32)*sum(bmap);
   
   bmfb(usedSubs) = ba_cal(subMFBall(usedSubs), 9.8, margin, codingGain, 0, 100, 0);
   RDMTmfb = fs/(N+32)*sum(bmfb); %include the cp
   % get performance results
   %[SSNR, SNR, subSNR,geoSNRfinal,geoSNRmfb,bDMTfinal,bDMTmfb,...
    %     RDMTfinal,RDMTmfb,hw,Fh,Fw,colorNoiseaft,Fhw] = ...
     % perform(W,B,channel,D,Nb,N,inputSignal,noise,...
      %channelGainAll,inputSpecAll,noiseSpecAll,margin,codingGain,fs,subMFBall,usedSubs);
   
   % put results into table
   str = sprintf('%1.3fe6',RDMTfinal/1e6);
   set(methRateHndl,'string',str);
   str = '-';
   set(methSNRHndl,'string',str);
   str = '-';
   set(methSSNRHndl,'string',str);
   str = sprintf('%d',D);
   set(methDelayHndl,'string',str);
   str = '-';
   set(methMSEHndl,'string',str);
   str = sprintf('%1.3fe6',RDMTmfb/1e6);
   set(methMaxHndl,'string',str);
   
   
   % save results to disk in order to use it in other functions 
   save pertonetmpresults subMFBall subSNR N usedSubs title_str noiseSpec noiseSpecAll WM D py pX delaySearch rateVector Dmin Dmax fs;
   
   %TEQequalizerTaps = W;
   %TEQtargetImpulseResp = B;
   %TEQoptimalDelay = D;
   %TEQoriginalChannel = channel;
   %TEQequalizedChannel = hw;
   %TEQmatchedFilterBoundPerChannel = subMFB;
   %TEQSNRperChannel = subSNR;
   %TEQoriginalChannelFreqResp = Fh;
   %TEQequalizerFreqResp = Fw;
   %TEQequalizedChannelFreqResp = Fhw;
   %TEQchanNoisePowSpecAfterEqual = colorNoiseaft;
   %TEQchanNoisePowSpecBeforeEqual = noiseSpec;
    TEQDmin = Dmin;
    TEQDmax = Dmax;
    TEQTG = TG;
    TEQNw = Nw;
    TEQN = N;
    TEQPT = WM;
    TEQPTSubsnr = subSNR;
    TEQoriginalChannel = channel;
    TEQoptimalDelay = D;
    TEQCodingGain=codingGain;
    TEQMargin=margin;
    TEQInputPower=totalInputPower;
    TEQAWGNPower=AWGNpower;
    TEQLoopnum=loopNum;
   
    save pertoneresults TEQDmin TEQDmax TEQTG TEQNw TEQN TEQPT TEQPTSubsnr TEQoriginalChannel ...
    TEQoptimalDelay TEQCodingGain TEQMargin TEQInputPower TEQAWGNPower TEQLoopnum;

   %save teqresults TEQequalizerTaps TEQtargetImpulseResp  ...
	   %TEQoptimalDelay TEQoriginalChannel TEQequalizedChannel ...
	   %TEQmatchedFilterBoundPerChannel TEQSNRperChannel ...
	   %TEQoriginalChannelFreqResp TEQequalizerFreqResp ...
	   %TEQequalizedChannelFreqResp TEQchanNoisePowSpecAfterEqual...
	   %TEQchanNoisePowSpecBeforeEqual TEQDmin TEQDmax TEQNb TEQNw TEQN;
	   
   % plot result
   pertone('plotGraph');
   set(gcf,'Pointer','arrow');
   return
   
   % plot graph
elseif strcmp(action,'plotGraph'),
  legend off
  set(gcf,'Pointer','watch');
   % get handles
   fhndlList = get(gcf,'Userdata');
   viewHndl = fhndlList(12);
   freqzHndl = fhndlList(1);
   % set axis
   axHndl = gca;
   axes(freqzHndl);
   % which graph is going to be plotted
   graphNum = get(viewHndl,'value'); 
   %1= TIR/SIR, 2= TEQ impulse res., 3= TEQ freq. res., 4= SNR/MFB
   %5=original and short channel impulse res, 6 = noise power spec., 
   %7= delay plot, 8= equalized channel freq
   
   % load results for current TEQ
   load pertonetmpresults;
   freqAxis = 1:N/2+1;
   
   % plot chosen grap
   if graphNum == 1
      legend off
      plot(freqAxis(usedSubs),10*log10(subSNR(usedSubs)+eps))
      hold on;
      plot(freqAxis(usedSubs),10*log10(subMFBall(usedSubs)+eps),'r')
      hold off
      grid on
      axisSet = axis;
      axisSet(1) = min(freqAxis(usedSubs))-2;
      axisSet(2) = max(freqAxis(usedSubs))+2;
      axis(axisSet);
      xlabel('frequency')
      ylabel('magnitude (dB)')
      title([title_str,' achieved and upper-bound SNR'])
      legend('SNR','MFB',0);
      
  elseif graphNum == 2
      legend off
      plot(freqAxis(usedSubs),20*log10(noiseSpecAll(usedSubs)+eps),'r')
      grid on	
      axisSet = axis;
      axisSet(1) = min(freqAxis(usedSubs))-2;
      axisSet(2) = max(freqAxis(usedSubs))+2;
      axis(axisSet);
      xlabel('frequency')
      ylabel('magnitude (dB)')
      title([title_str,' Noise Power Spectrum'])
      legend('before TEQ','after TEQ',0);  
      
    elseif graphNum == 3 
      if delaySearch==0
        legend off
        plot([Dmin:Dmax],rateVector*fs/(N+32));
        xlabel('delay');
        ylabel('bit rate');
        title([title_str,' bit rate vs. delay']);
     else
        legend off
        plot([D],rateVector*fs/(N+32));
        xlabel('delay');
        ylabel('bit rate');
        title([title_str,' bit rate vs. delay']); 
      end
      grid on
  end
   
   set(gcf,'Pointer','arrow');
   return
   
     %save pertone to a file 
elseif strcmp(action,'savepertone');
   
    [filename,path]=uiputfile('*');
    if filename==0
        return;    
    end
    
    load pertoneresults;
    output=fopen([path,filename],'w');
    fprintf(output,'%s','optimal delay =');
    fprintf(output,' %d ;\n',TEQoptimalDelay);
    fprintf(output,'\n');
    [rowNum,colNum]=size(TEQPT);
    for j=1:colNum
        fprintf(output,'subcarrier %d FEQ =[\n',j);
        for i=1:rowNum
            fprintf(output,'%f \n',TEQPT(i,j));
        end
        fprintf(output,'%s\n','];');
        fprintf(output,'\n');
    end


    fclose(output);
    
    
   %save all parameters
elseif strcmp(action,'saveall')
    [filename,path]=uiputfile('*');
    if filename==0;
        return;
    end
    
    load pertoneresults;
    output=fopen([path,filename],'w');
    fprintf(output,'all the parameters are also saved into %s.mat\n',filename);
    fprintf(output,'%s','tone group =');
    fprintf(output,' %d\n',TEQTG);
    fprintf(output,'%s','FEQ length =');
    fprintf(output,' %d;\n',TEQNw);
    fprintf(output,'%s','FFT size =');
    fprintf(output,' %d;\n',TEQN);
    fprintf(output,'%s','coding gain =');
    fprintf(output,' %f;\n',TEQCodingGain);
    fprintf(output,'%s','margin =');
    fprintf(output,' %d;\n',TEQMargin);
    fprintf(output,'%s','Dmin =');
    fprintf(output,' %d;\n',TEQDmin);
    fprintf(output,'%s','Dmax =');
    fprintf(output,' %d;\n',TEQDmax);
    fprintf(output,'%s','input power =');
    fprintf(output,' %f (dBm);\n',TEQInputPower);
    fprintf(output,'%s','AWGN power =');
    fprintf(output,' %f (dBm/Hz);\n',TEQAWGNPower);
    fprintf(output,'%s','loop number =');
    fprintf(output,' %d;\n',TEQLoopnum);
    fclose(output);
    str=sprintf(['save ', filename,' TEQTG TEQNw TEQN TEQCodingGain TEQMargin TEQDmin',...
            ' TEQDmax TEQInputPower TEQAWGNPower TEQLoopnum;']);
    eval(str);
    
    
   % get help
elseif strcmp(action,'info'),
   ttlStr = get(gcf,'name');
   myFig = gcf;
   
   topic1 =  ['Per Tone design demo'];
   helptop1 = [...
     'The demonstration facilitates the design of a per tone equalizer      '
     'structure.  The per tone equalizer consists of a sliding FFT and      '
     'a linear combiner of each vector of sliding FFT values for each       '
     'subchannel (tone) with a vector of equalizer coefficients tuned       '
     'for that tone.  In essence, the time domain equalizer of Nw taps      '
     'in the conventional equalizer has been moved into the single-tap      '
     'frequency domain equalizer to form the Nw-tap linear combiner.  A     '
     'reduced complexity per tone equalizer structure exists which          '
     'replaces the sliding FFT with a single FFT and tapped delay line      '
     'of time domain difference terms.  For either structure, Nw N/2        '
     'complex coefficients need to be trained, where Nw is the equalizer    '
     'length (combiner taps per tone) and N/2 is the number of              '
     'subchannels (tones).  More information concerning the per tone        '
     'equalizer can be found in the following paper:                        '
     '                                                                      '
     'K. Van Acker, G. Leus, M. Moonen, O. van de Wiel, and T. Pollet,      '
     '"Per Tone Equalization for DMT-based systems", IEEE Trans. on         '
     'Communications, vol. 49, no. 1, pp. 109-119, Jan. 2001.               '
     '                                                                      '
     'Usage: Set parameters to desired values and hit the Calculate button. '
     '                                                                      '
     'The following design methods are supported:                           '
	 '  Least squares matching of frequency domain input-output.            '
     '  MMSE pertone method                                                 '
     '    G. Ysebaert, M. Moonen and T. Pollet, "Combined RLS-LMS           '
     '    initialization for per tone equalizers in DMT-receivers"          ' 
     '    Acoustics, Speech, and Signal Processing, 2002 IEEE International '
     '    Conference on , Volume: 3 , 2002                                  '
	 ];
   
   topic2 =  ['Design and Simulation Parameters'];
   helptop2 = [...
	 'On the top of he control window on the right is a pulldown menu       '
     'from which a design method can be chosen.                             '
	 'Below this pulldown menu are the following editable text windows      '
     'which are used to set the desired parameters,                         '
     '  Tones per group. Set the tones for each group.                      '
     '   Per tone equalizer length. Defines the number of taps of           ' 
	 '   the equalizer for each tone.                                       '
     '  Fast Fourier transform (FFT) size. Sets the FFT size used in DMT    '
     '   modulation. It is twice the number of subchannels.                 '
     '  Coding gain (dB). Defines a coding in dB used during capacity       '
     '   calculations.                                                      '
     '  Margin (dB). Sets the desired system margin in dB. This is also     '
     '   used in capacity calculations.                                     '
     '  Dmin and Dmax. The optimal delay is search in the interval between  '
     '   Dmin and Dmax.                                                     '
	 '  Input power (dBm). Defines the input signal power in dBm.           '
	 '  AWGN pow (dBm/Hz). Sets the amount of additive white Gaussian noise '
     '   in dBm/Hz. AWGN is added to the near-end crosstalk noise.          '
     '  CSA loop \# (1-8). Selects the desired channel to run the simulation'
     '   on. Currently 8 standard CSA loops are supported.                  '];
   
     topic3 =  ['Graphics'];
   helptop3 = [...
	 'Below the editable text windows is another pull-down menu which       '
     'is used to select the desired graph to be displayed.                  '
     'The following graphics can be selected:                               '
     '  SNR and MFB. The SNR and matched filter bound to the SNR is         '
     '   displayed as a function of frequency (subchannels).                '
     '  Noise Power Spectrum. Shows the power spectrum of the noise which   '
     '   consists of NEXT noise plus AWGN.                                  '
     ];
 
   topic4 =  ['Performance Measures'];
   helptop4 = [...
'The two remaining buttons in the control frame are                    '
'   Info. Displays this help.                                          '
'   Calculate. Starts the calculation and performance evaluation of the'
'    TEQ.                                                              '
'The following performance measures are calculated and listed in the   '
'table:                                                                '
'   Rate. Gives the achievable bit rate with the given channel and Per '
'    Tone Equalizers settings.                                         '
'   Max Rate. Shows the absolute maximum achievable bit rate given the '
'    channel and equalizer settings. It is basically the capacity      '
'   calculated from the MFB.                                           '];
		 
   str =  { topic1 helptop1; topic2 helptop2; topic3 helptop3; topic4 helptop4 };
   
   helpwin(str,'Topic 1','Per Tone design demo')                   
   return  
   
else
   disp(sprintf( 'Internal ERROR: ''%s'' not recognized',action))
end