How can I add cp and zp to code of ofdm ber and get it out?

2 ビュー (過去 30 日間)
mohamad gholami
mohamad gholami 2018 年 12 月 22 日
Hello...
I have a code about ofdm that is just ber.
I want to add cp and zp using qam16 modulation in awgn channel.
But i can not...
Please help me..
Or if you can send me the code....
I want to go to these outputs for which I send them a picture.
I will send you the code.
Thanks
% OFDM_basic.m
clc
clear all
tic
NgType=1; % NgType=1/2 for cyclic prefix/zero padding
if NgType==1, nt='CP'; elseif NgType==2, nt='ZP'; end
Ch=0; % Ch=0/1 for AWGN/multipath channel
if Ch==0, chType='AWGN'; Target_neb=100; else chType='CH'; Target_neb=500; end
figure(Ch+1), clf
PowerdB=[0 -8 -17 -21 -25]; % Channel tap power profile 'dB'
Delay=[0 3 5 6 8]; % Channel delay 'sample'
Power=10.^(PowerdB/10); % Channel tap power profile 'linear scale'
Ntap=length(PowerdB); % Chanel tap number
Lch=Delay(end)+1; %Channel length
Nbps=4; M=2^Nbps; % Modulation order=2/4/6 for QPSK/16QAM/64QAM
Nfft=64; % FFT size
Ng=3; %Nfft/4; % Ng=0: Guard interval length
Ng=Nfft/4;
Nsym=Nfft+Ng; % Symbol duration
Nvc=Nfft/4; % Nvc=0: no virtual carrier
Nused=Nfft-Nvc;
EbN0=[0:5:20]; % EbN0
N_iter=1e5; % Number of iterations for each EbN0
Nframe=3; % Number of symbols per frame
sigPow=0; % Signal power initialization
file_name=['OFDM_BER_' chType '_' nt '_' 'GL' num2str(Ng) '.dat'];
fid=fopen(file_name, 'w+');
norms=[1 sqrt(2) 0 sqrt(10) 0 sqrt(42)]; % BPSK 4-QAM 16-QAM
for i=0:length(EbN0)
randn('state',0); rand('state',0); %Ber2=ber(); % BER initialization
Neb=0; Ntb=0; % Initialize the number of error/total bits
for m=1:N_iter
% Tx______________________________________________________________
X = randi(M,1,Nused*Nframe) - 1;% bit: integer vector
Xmod= qammod(X,M,0,'gray')/norms(Nbps);
if NgType~=2, x_GI=zeros(1,Nframe*Nsym);
elseif NgType==2, x_GI= zeros(1,Nframe*Nsym+Ng);
% Extend an OFDM symbol by Ng zeros
end
kk1=[1:Nused/2]; kk2=[Nused/2+1:Nused]; kk3=1:Nfft; kk4=1:Nsym;
for k=1:Nframe
if Nvc~=0, X_shift= [0 Xmod(kk2) zeros(1,Nvc-1) Xmod(kk1)];
else X_shift= [Xmod(kk2) Xmod(kk1)];
end
x= ifft(X_shift);
x_GI(kk4)= guard_interval(Ng,Nfft,NgType,x);
kk1=kk1+Nused; kk2= kk2+Nused; kk3=kk3+Nfft; kk4=kk4+Nsym;
end
if Ch==0, y= x_GI; % No channel
else % Multipath fading channel
channel=(randn(1,Ntap)+j*randn(1,Ntap)).*sqrt(Power/2);
h=zeros(1,Lch); h(Delay+1)=channel; % cir: channel impulse response
y = conv(x_GI,h);
end
if i==0 % Only to measure the signal power for adding AWGN noise
y1=y(1:Nframe*Nsym); sigPow = sigPow + y1*y1';
continue;
end
% Add AWGN noise________________________________________________
snr = EbN0(i)+10*log10(Nbps*(Nused/Nfft)); % SNR vs. Eb/N0
noise_mag = sqrt((10.^(-snr/10))*sigPow/2);
y_GI = y + noise_mag*(randn(size(y))+j*randn(size(y)));
% Rx_____________________________________________________________
kk1=(NgType==2)*Ng+[1:Nsym]; kk2=1:Nfft;
kk3=1:Nused; kk4=Nused/2+Nvc+1:Nfft; kk5=(Nvc~=0)+[1:Nused/2];
if Ch==1
H= fft([h zeros(1,Nfft-Lch)]); % Channel frequency response
H_shift(kk3)= [H(kk4) H(kk5)];
end
for k=1:Nframe
Y(kk2)= fft(remove_GI(Ng,Nsym,NgType,y_GI(kk1)));
Y_shift=[Y(kk4) Y(kk5)];
if Ch==0, Xmod_r(kk3) = Y_shift;
else Xmod_r(kk3)= Y_shift./H_shift; % Equalizer - channel compensation
end
kk1=kk1+Nsym; kk2=kk2+Nfft; kk3=kk3+Nused; kk4=kk4+Nfft; kk5=kk5+Nfft;
end
X_r=qamdemod(Xmod_r*norms(Nbps),M,0,'gray');
Neb=Neb+sum(sum(de2bi(X_r,Nbps)~=de2bi(X,Nbps)));
Ntb=Ntb+Nused*Nframe*Nbps; %[Ber,Neb,Ntb]=ber(bit_Rx,bit,Nbps);
if Neb>Target_neb, break; end
end
if i==0
sigPow= sigPow/Nsym/Nframe/N_iter;
fprintf('Signal power= %11.3e\n', sigPow);
fprintf(fid,'%%Signal power= %11.3e\n%%EbN0[dB] BER\n', sigPow);
else
Ber = Neb/Ntb;
fprintf('EbN0=%3d[dB], BER=%4d/%8d =%11.3e\n', EbN0(i), Neb,Ntb,Ber)
fprintf(fid, '%d\t%11.3e\n', EbN0(i), Ber);
if Ber<1e-6, break; end
end
end
if (fid~=0), fclose(fid); end
disp('Simulation is finished');
plot_ber(file_name,Nbps);
toc.
i need acod that give me this two figure in output::
aaa.png

回答 (0 件)

カテゴリ

Help Center および File ExchangePropagation and Channel Models についてさらに検索

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by