ofdm_memorysaleh_memorylessgainpd.m

上传的是一个16Qam预失真处理的程序

% OFDM_memorySaleh_memorylessgainPD.m
%
% Simulation 16QAM-OFDM system with memoryless Gain PD
%
% Programmed by Ruby Lin

%******************** Preparatin part *******************************

sr = 256000.0;                 % Symbol rate
m1 = 4;                        % m1:Number of modulation levels
                               % (BPSK:m1=1, QPSK:m1=2, 16QAM:m1=4)
br = sr .* m1;                 % Bit rate
                               % (BPSK:m1=1, QPSK:m1=2, 16QAM:m1=4)
nd = 1024;                     % Number of symbols htat simulates in
                               % each loop

L = 4;                                     % IFFT(FFT)采样点数
Sdata = 2^7;                               % IFFT(FFT)每帧长度

SNR_dB = 1:12;                           % 仿真信噪比范围
SNR1_dB = 0:0.1:12;  

Ar = 2.0;                                % TWTA的Dong-Seog Han参数
Br = 1;
Ap = pi/3;
Bp = 1;

%******************** Data generation *******************************
     
     data1 = rand(1,nd*m1) > 0.5;    % rand:built in function

%******************** 16QAM Modulation ******************************

     [ich,qch] = qammod(data1,1,nd,m1);
%      figure(1);
%      plot(ich,qch,'*');
      
   QAM_mod = ich + i*qch;
   
   %**************************** OFDM调制 ****************************

OFDM_mod_group = zeros(L*length(QAM_mod),1);
OFDM_mod = zeros(L*length(QAM_mod),1);

for n = 1:length(QAM_mod)/Sdata         
    
    OFDM_mod_group((n-1)*L*Sdata+Sdata/2+1:(n-1)*L*Sdata+Sdata/2*3) =...
        QAM_mod((n-1)*Sdata+1:n*Sdata);
      
    OFDM_mod((n-1)*(L*Sdata)+1:n*(L*Sdata)) = ...
        ifft(OFDM_mod_group((n-1)*(L*Sdata)+1:n*(L*Sdata)),(L*Sdata));
    
end;

OFDM_mod_reshape = reshape(OFDM_mod,1,length(OFDM_mod));

%**************************** memoryless Gain PD ****************************

delta = 1/2^6;                                         % 查找表精度
t = [delta:delta:1];                                   % 表格长度
Fr = (1-sqrt(1-(t/2).^2))./((t/2).^2);                 % 存储的幅值
Fp = (-pi/3)*(Fr.^2)./(1+Fr.^2);
F = Fr.*exp(i*Fp);                                     % 增益预失真的函数


PD_in_Env = abs(OFDM_mod_reshape);
PD_in_Phase = angle(OFDM_mod_reshape);

num = fix(PD_in_Env/delta) + 1;                % 将信号对应到相应的表格
for n=1:length(num)
    FF(n) = F(num(n));
end;
PD_out = OFDM_mod_reshape.*FF;

%**************************** memory Saleh ****************************

FIR_c = [0.7692 0.1538 0.0769];

A_memory = filter(FIR_c,1,PD_out);

PA_in_Env = abs(A_memory);
PA_in_Phase = angle(A_memory);

PA_out_r = Ar*PA_in_Env./(1+Br*PA_in_Env.^2);                                  % 幅度非线性放大
PA_out_p = Ap*PA_in_Env.^2./(1+Bp*PA_in_Env.^2) + PA_in_Phase;               % 相位非线性放大

PA_out = PA_out_r.*exp(i*PA_out_p);  %/nf_ibo*nf/2;     % 去归一化和去功率回退


ich2 = real(PA_out);
qch2 = imag(PA_out);

for ebn0 = 1:length(SNR_dB)+1
     
%******************** START CALCULATION *****************************

nloop = 10;                       % Number of simulation loops
noe = 0;                           % Number of error data
nod = 0;                           % Number of transmitted data

   %for iii = 1:nloop     
%******************** Attenuation Calculation ***********************

      spow = sum(ich2.*ich2+qch2.*qch2)/nd;
          % sum:built in function
      attn = 0.5*spow*sr/br*10.^(-(ebn0-1)/10);
      attn = sqrt(attn);
          % sqrt:built in function
     
%************** Add White Gaussian Noise (AWGN) *********************

     [ich3,qch3] = comb(ich2,qch2,attn);
          
OFDM_AWGN = ich3 + j*qch3;
OFDM_AWGN_reshape = reshape(OFDM_AWGN,length(OFDM_AWGN),1);
%**************************** OFDM解调 ****************************

for m = 1:length(QAM_mod)/Sdata         
    
    OFDM_demod_group((m-1)*(L*Sdata)+1:m*(L*Sdata)) = ...
        fft(OFDM_AWGN_reshape((m-1)*(L*Sdata)+1:m*(L*Sdata)),(L*Sdata));
    
    OFDM_demod((m-1)*Sdata+1:m*Sdata) = ...
        OFDM_demod_group((m-1)*L*Sdata+Sdata/2+1:(m-1)*L*Sdata+Sdata/2*3);
         
end;

ich5 = real(OFDM_demod);
qch5 = imag(OFDM_demod);
figure(2);
plot(ich5,qch5,'*');
    
%******************** 16QAM Demodulation ****************************
     
     [demodata] = qamdemod(ich5,qch5,1,nd,m1);
      
%******************** Bit Error Rate (BER) **************************
  
     noe2 = sum(abs(data1-demodata));
     nod2 = length(data1);
     noe = noe + noe2;
     nod = nod + nod2;
     
      ber(ebn0) = noe/nod;
  
   
end  % for iii = 1:nloop

t1 = [0:0.1:12];
tt1=exp(t1*log(10)/10);  
B1 = 3/8.*erfc(sqrt(2/5.*tt1))-9/64.*erfc(sqrt(2/5.*tt1)).*erfc(sqrt(2/5.*tt1));
t11=[0:length(SNR_dB)];

figure(3);
semilogy(t11,ber,'r');

[pxx1,f1] = pwelch(OFDM_mod_reshape,512);                 
 pxxdB1 = 10 * log10(pxx1 / max(pxx1));

[pxx2,f2] = pwelch(OFDM_AWGN_reshape,256);                 
pxxdB2 = 10 * log10(pxx2 / max(pxx2));
figure(4);                                          % OFDM调制、解调信号功率谱密度
plot(f1,pxxdB1,'v-',f2,pxxdB2,'r');

% eyediagram((ich+qch*i),2);
eyediagram((ich5(1:500)+qch5(1:500)*i),2)

%******************** end of file ***********************************