+ Post New Thread
Results 1 to 20 of 31

17th March 2012, 23:02 #1
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
How to graph BER with perfect results?
Hello
i have a problem with matlab. suppose after completing a simulation with this results :
snr=0:5:35
for qpsk i have this vector: ber_qpsk=[1e3 1e4 1e5 0 0 0 0]
for 16qam i have this vector: ber_16qam=[1e2 1e3 1e4 1e5 0 0 0]
for 64qam i have this vector: ber_64qam=[0.5e2 1e2 1e3 1e4 1e5 0 0]
now, i want to graph them all in one graph. tries berfit and i got '[]'.

22nd March 2012, 22:55 #2
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
anyone can help me with that ?

27th March 2012, 00:05 #3
 Join Date
 Aug 2007
 Location
 Germany
 Posts
 64
 Helped
 9 / 9
 Points
 2,104
 Level
 10

29th March 2012, 16:26 #4
 Join Date
 Aug 2011
 Location
 New Jersey, USA
 Posts
 19
 Helped
 6 / 6
 Points
 801
 Level
 6
Re: How to graph BER with perfect results?
berfit is doing more than you want to do. To plot those three vectors in one graph just use:
figure(1);
subplot(311), plot(snr, ber_qpsk);
subplot(312), plot(snr, ber_16qam);
subplot(313), plot(snr, ber_64qam);
Then, you can use the bertool command to create vectors of the theoretical BER curves for psk and qam, and plot them next to those.
 Ryan
1 members found this post helpful.

30th March 2012, 11:56 #5
Re: How to graph BER with perfect results?
Hi,
Try this
snr=0:5:30;
ber_qpsk=[1e3 1e4 1e5 0 0 0 0];
ber_16qam=[1e2 1e3 1e4 1e5 0 0 0];
ber_64qam=[0.5e2 1e2 1e3 1e4 1e5 0 0];
figure(1);
semilogy(snr,ber_qpsk,'ko','LineWidth',1);
hold on
semilogy(snr,ber_16qam,'kx','LineWidth',1);
hold on
semilogy(snr,ber_16qam,'kd','LineWidth',1);
grid on
legend('qpsk','16qam','64qam',3);
xlabel('SNR (dB)');
ylabel('BER');
title(['BER versus SNR for different receivers'])
axis([min(snr) max(snr) 10^(6) 1]);
Best regards
2 members found this post helpful.

30th March 2012, 11:56

31st March 2012, 17:26 #6
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
hi , i'll try this solutions thx.
but in addition i asked my lecturer and he told that i need to increase the number of bits. so for my expected value for my simulation i want to be around ber = 10^6 so this demanding to transmit at least 10^7 bits. i can do it, but this can take me a day i've checked it. i know that my sim' isn't so efficient , cause i have many FOR's.
so here comes your part :)
i can reduce the for loops in some cases like :
instead of
for 1:N_sym %number of symbols usually 1000 up to 4500
x_Tx_bits=randsrc(1,2256,[0 1]);
CodedBits=encode(henc,x_Tx_Bits); %using ldpc encoder with 2256*4512 matrix
end
i want create matrix of bits:
x_Tx_Bits=randsrc(N_sym,2256,[0 1]);
is there any way i can encode this giant matrix at once ?

31st March 2012, 17:26

31st March 2012, 21:31 #7
Re: How to graph BER with perfect results?
Hi,
Yes to simulate ber=10^6, you must send 10^7 bits but the incertitude is big. If you sen 10^8 you will have ber=10^6 with 1% incertitude.
You can encode a matrix. This is an example :
N_sym=4512;
n = 6; k = 4; % Set codeword length and message length.
x_Tx_Bits=randsrc(N_sym,4,[0 1]); % Message is a binary matrix.
code = encode(x_Tx_Bits,n,k,'cyclic');
But, you don't give enough parameters of LDPC coder. henc ?
1 members found this post helpful.

31st March 2012, 21:31

31st March 2012, 22:50 #8
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
I am posting my raw code. meaning i am doing all my tests on this code. If i am getting expected results I am changing my simulation files.
The code including awgn/LTI/LTV with all the blocks that needed from ofdm scheme(not including Interleaving , * exist but not straightforward )
The LDPC matrix taken from the net and it working very well for AWGN+LTI channel I'm getting with QPSK the BER=0 (i know it's wrong write it) issue just with SNR~1.5dB.
My goal is to shorten the LTV stage with OVERLAP&ADD method, but first i want to stabilize my results.
For now the LTI+awgn working well just need to fix BER.
In the LTV(section) channel i have hard to figure the math and to normalize my parameters good enough.
many thanks to u all, it is very appreciated . keep helping me :)
Code:function NumErrors=OurLDPCEncandDec_without_MMSE(itertion,Rate,Modulation,SNRdB) % close all; clc ; clear % Rate=2;Modulation=4;SNRdB=3 % Rate: 2 for 1/2, 4 for 3/4, 6 for 5/6, 8 for 7/8 NumErrors=zeros(1,itertion); for iter=1:itertion %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Parameters Section %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% load 802_16e_4512_d5_H H mat; H1=H(:,1:2256); H2=H(:,2257:4512); H=[H2 H1]; if Rate==2 ParityLength=2256; elseif Rate==4 ParityLength=2256/3; elseif Rate==6 ParityLength=464; elseif Rate==8 ParityLength=320; elseif Rate==10 ParityLength=240; end; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Create Signal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% xt_long_with_CP=[]; xt_long_bits=[]; N_sym=100; henc = fec.ldpcenc(H); for k=1:N_sym InfoBits=randsrc(1,2256,[1 0]); % Generate 2256 information bits xt_long_bits=[xt_long_bits InfoBits]; % Encode and Modulate CodedBits= encode(henc,InfoBits); % Rate Matching ParityLocations=2256+randperm(2256); ParityLocations=ParityLocations(1:ParityLength); TotalLocations(:,k)=[1:2256,ParityLocations]; % Should be an output of the encoder CodedBits=CodedBits(TotalLocations(:,k)); Symbols=ModulateQAM(CodedBits,Modulation); FFTsize=2^(ceil(log2(length(Symbols)))); ActiveSCs=length(Symbols); PilotLocations=1:4:FFTsize; %pilot loc /for mod=4,16 spacing = 3. /for mod=64 spacing =4 SC_loc=setdiff(1:FFTsize,PilotLocations); SC_loc=SC_loc(1:ActiveSCs); %info location %%%%pilots%%%% x_f=zeros(FFTsize,1); x_f(SC_loc)=Symbols; Xp=randsrc(1,length(PilotLocations),[1 1]); %Pilot generation Pilots_PRBS(:,k)=Xp.'; % PRBS patterns x_f(PilotLocations)=Xp.'; %insertion of the pilots xt=ifft(x_f); xt_withCP=[xt(end1023:end);xt]; xt_long_with_CP=[xt_long_with_CP;xt_withCP]; end %%%%%%%%%%%%%%%%%%%%%%%%Channel%%%%%%%%%%%%%%%%%%%%%%%%%% Fs=10e6; %Sampling rate f_m=10; %Doppler spread N_samples=length(xt_long_with_CP); T=N_samples/Fs; %Sample time SNR=10^(SNRdB/10); % Taking channel parameters from ITU table TapPowerdB=[0 9.7 19.2 22.8]; % in dB Delay=[1 110 190 410]*1e9; % in Sec % ChannelLength=length(TapPower); a_samples=[]; %%%%Tap power in dB to linear%%%% TapPower=10.^(TapPowerdB/10); % TapPower=TapPower/norm(TapPower); % Normalize to unit channel power F_low=f_m*10; % usaully 10*f_m the low sampling % freq of the channel coefficients N_low=ceil(T*F_low); %%%%Gausssian random process pass throw filter(flat/Ushaped)%%%% h=firpm(64,[0 2*f_m/F_low 2*f_m/F_low*1.2 1],[1 1 0 0 ]); h=h/norm(h); % channel_length=length(h); for k=1:length(TapPower) a_low_temp=conv((randn(1,N_low+64)+1i*randn(1,N_low+64))/sqrt(2),h);%was ...,h)*TapPower(k); a_low_temp=a_low_temp(65:N_low+64); % Generate first tap sampled at low freq a_low(k,:)=a_low_temp; channel_length=6; %%%%Upsample To Fs%%%% a_samples(k,:)=resample(a_low(k,:),Fs,F_low); Sigvar=var(a_samples(k,:)); a_samples(k,:)=a_samples(k,:)*sqrt((TapPower(k)/Sigvar)); end %%%%Initilaize Delay Vector%%%% DelayInSamples=ceil(Delay*Fs); DelayInSamples=DelayInSamplesDelayInSamples(1); y=zeros(1,length(xt_long_with_CP)); %%%%Convolution%%%% for k=max(DelayInSamples)+1:length(xt_long_with_CP) for j=1:length(TapPower) a(j)=xt_long_with_CP(kDelayInSamples(j))*a_samples(j,k); end y(k)=sum(a); end % Add Noise SNR_corrected=SNR*(1366+2256)/FFTsize; %need to generic fit for qpsk 2256 infobits + 1366 pilots sigma_noise=norm(xt_long_with_CP)/sqrt(length(y)*SNR_corrected); y=y+sigma_noise*(randn(size(y))+1i*randn(size(y)))/sqrt(2); yBeforeCP_RemoveMatrix=reshape(y.',FFTsize+1024,N_sym); yCP_RemoveMatrix=yBeforeCP_RemoveMatrix(1025:end,:); Y_matrix_PostFFT=fft(yCP_RemoveMatrix); for kk=1:N_sym, % H_est = MMSE_CE(Y_matrix_PostFFT(:,kk),Pilots_PRBS(:,kk),PilotLocations,FFTsize,4,h,SNR); % h_est = ifft(H_est); h_DFT = h_est(1:channel_length); % H_DFT = fft(h_DFT,FFTsize); Y_Pilots=Y_matrix_PostFFT(PilotLocations,kk); EstTimeDomainChannel=ifft(Y_Pilots.*Pilots_PRBS(:,kk));%Taking only pilots EstTimeDomainChannel=EstTimeDomainChannel(1:(channel_length)); % Kill contibutions outside the CP duration H_est=fft(EstTimeDomainChannel,FFTsize); % est. channel demodulatedSCs=Y_matrix_PostFFT(:,kk)./H_est; Y_equalized_withoutPilots(:,kk)=demodulatedSCs(SC_loc); ppSNR(:,kk)=abs(H_est(SC_loc)).^2/(sigma_noise^2*FFTsize); end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Decode Signal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% DecodedBits_long=[]; for kk=1:N_sym, LLRs=SISO_LLR(Y_equalized_withoutPilots(:,kk),ppSNR(:,kk),Modulation); LLR_rate_matched=zeros(1,2*2256); LLR_rate_matched(TotalLocations(:,kk))=LLRs; hdec=fec.ldpcdec(H); DecodedBits=decode(hdec, LLR_rate_matched); DecodedBits_long=[DecodedBits_long DecodedBits]; end NumErrors(1,iter)=sum(DecodedBits_long ~=xt_long_bits); % Count errors end % function QAM_Symbols=ModulateQAM(CodedBits,Modulation) % QPSK mapping MapQPSK= ([1+i 1i 1+i 1i])/sqrt(2); % 16QAM mapping Map_16QAM= ((2/sqrt(10))* ... [ 0.5+0.5*i 0.5+1.5*i 0.50.5*i 0.51.5*i ... 1.5+0.5*i 1.5+1.5*i 1.50.5*i 1.51.5*i .... 0.5+0.5*i 0.5+1.5*i 0.50.5*i 0.51.5*i .... 1.5+0.5*i 1.5+1.5*i 1.50.5*i 1.51.5*i]); % 64QAM mapping Map_64QAM= ((2/sqrt(42))* ... [ 1.5+1.5*i 1.5+0.5*i 1.5+2.5*i 1.5+3.5*i ... 1.51.5*i 1.50.5*i 1.52.5*i 1.53.5*i ... 0.5+1.5*i 0.5+0.5*i 0.5+2.5*i 0.5+3.5*i ... 0.51.5*i 0.50.5*i 0.52.5*i 0.53.5*i ... 2.5+1.5*i 2.5+0.5*i 2.5+2.5*i 2.5+3.5*i ... 2.51.5*i 2.50.5*i 2.52.5*i 2.53.5*i ... 3.5+1.5*i 3.5+0.5*i 3.5+2.5*i 3.5+3.5*i ... 3.51.5*i 3.50.5*i 3.52.5*i 3.53.5*i ... 1.5+1.5*i 1.5+0.5*i 1.5+2.5*i 1.5+3.5*i ... 1.51.5*i 1.50.5*i 1.52.5*i 1.53.5*i ... 0.5+1.5*i 0.5+0.5*i 0.5+2.5*i 0.5+3.5*i ... 0.51.5*i 0.50.5*i 0.52.5*i 0.53.5*i ... 2.5+1.5*i 2.5+0.5*i 2.5+2.5*i 2.5+3.5*i ... 2.51.5*i 2.50.5*i 2.52.5*i 2.53.5*i ... 3.5+1.5*i 3.5+0.5*i 3.5+2.5*i 3.5+3.5*i ... 3.51.5*i 3.50.5*i 3.52.5*i 3.53.5*i]); if Modulation==4 Map=MapQPSK; elseif Modulation==16 Map=Map_16QAM; elseif Modulation==64 Map=Map_64QAM; end; BitsPerCarrier=log2(Modulation); CodedBitMat=reshape(CodedBits.',BitsPerCarrier,length(CodedBits)/BitsPerCarrier); SymbolIndices=bi2de(fliplr(CodedBitMat.'))+1; QAM_Symbols=Map(SymbolIndices);
 Post added at 22:50  Previous post was at 22:24 
To simulate 10^8 bits i will need to seat for a week... i think... this is way i asking how can i get rid from most of the FOR loops.
With the code i sent last post i will have to generate (1e8/2256=44000 symbols=> for i=1:44000)(with QPSK code rate=0.5)Last edited by WUID; 31st March 2012 at 22:41.

1st April 2012, 02:53 #9
Re: How to graph BER with perfect results?
Hi,
You can began with optimizing your code and tested for 10^7 bits. After that simulate with 10^8 to get more precision for 10^6 ber.
It is normal that BER simulation takes much time (more than one day) but it is possible to win time by changing code.
I have some question before i can propose some modifications
 Why you are using the iterative loop ? could be for iterative decoding of LDPC code ?!
 There a simple modification that you can do to save time like
. % Parameters Section must be before the iterative loop. it is useless to repeat loading 802_16e_4512_d5_H every time.
. hdec=fec.ldpcdec(H) instruction is repeated twice. Once at coding and one at decoding. Only one before the loop will be ok.
 encode(henc,InfoBits) doesn't support Matrix input. I will see what you can change after i analyse more the code.
Regards
Chaker
1 members found this post helpful.

1st April 2012, 03:17 #10
 Join Date
 Aug 2011
 Location
 New Jersey, USA
 Posts
 19
 Helped
 6 / 6
 Points
 801
 Level
 6
Re: How to graph BER with perfect results?
I must have been tired when I answered. This does not plot the three vectors in one graph. Mohamed's way is right.
Matlab is optimized for vector multiplication, so use the '.' operator as much as you can to avoid nested for loops. This is a godsend in communications simulation.
For instance if you want to implement a matched filter for your receiver, where received_sig and template_sig are the two signal vectors, you can multiply an integrateanddump in one line:
corr = sum(received_sig.*template_sig);
This would do the same as:
corr = 0;
for ii=1:LENGTH
corr = corr + received_sig(ii)*template_sig(ii);
end
In my experience this can really speed things up. I would like to quantify that statement, but I'm not on a computer with Matlab right now.
 Ryan
1 members found this post helpful.

1st April 2012, 05:49 #11
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
Hi Mohamed
 Why you are using the iterative loop ? could be for iterative decoding of LDPC code ?!
 encode(henc,InfoBits) doesn't support Matrix input.
as for % Parameters Section your are very right, but this is the deal with my simulation(this code that is modified to function). I want to create Data matrix in one line:
InfoBits=randsrc(1000,2256,[0 1]);
and be able to encode it at once.
the aim right now with eliminating the FOR is in the "%Create Signal" part.

1st April 2012, 12:59 #12
Re: How to graph BER with perfect results?
Hi,
I have made some change in your code by putting all fixed parameters before iteration loop and 1:N_sym loop. I was forced to add some parameters. Try the code if it give same results. I can't do this because i don't have SISO_LLR function.
function NumErrors=OurLDPCEncandDec_without_MMSE(itertion,R ate,Modulation,SNRdB)
% close all; clc ;
%clear all ;
% Rate=2;
% Modulation=4;
% SNRdB=3;
% itertion=5;
% Rate: 2 for 1/2, 4 for 3/4, 6 for 5/6, 8 for 7/8
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% Parameters Section
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
Symbols_length=2256;
CP_lenght=1024;
N_sym=100;
if Rate==2
ParityLength=2256;
elseif Rate==4
ParityLength=2256/3;
elseif Rate==6
ParityLength=464;
elseif Rate==8
ParityLength=320;
elseif Rate==10
ParityLength=240;
end;
load 802_16e_4512_d5_H H mat;
H1=H(:,1:2256);
H2=H(:,2257:4512);
H=[H2 H1];
henc = fec.ldpcenc(H);
hdec = fec.ldpcdec(H);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% IFFT/FFT Parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
FFTsize=2^(ceil(log2(Symbols_length)));
ActiveSCs=Symbols_length;
PilotLocations=1:4:FFTsize; %pilot loc /for mod=4,16 spacing = 3. /for mod=64 spacing =4
SC_loc=setdiff(1:FFTsize,PilotLocations); SC_loc=SC_loc(1:ActiveSCs); %info location
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
%Channel
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
Fs=10e6; %Sampling rate
f_m=10; %Doppler spread
N_samples=N_sym*(FFTsize+CP_lenght);
T=N_samples/Fs; %Sample time
% Taking channel parameters from ITU table
TapPowerdB=[0 9.7 19.2 22.8]; % in dB
Delay=[1 110 190 410]*1e9; % in Sec
% ChannelLength=length(TapPower);
a_samples=[];
%%%%Tap power in dB to linear%%%%
TapPower=10.^(TapPowerdB/10);
% TapPower=TapPower/norm(TapPower); % Normalize to unit channel power
F_low=f_m*10; % usaully 10*f_m the low sampling
% freq of the channel coefficients
N_low=ceil(T*F_low);
%%%%Gausssian random process pass throw filter(flat/Ushaped)%%%%
h=firpm(64,[0 2*f_m/F_low 2*f_m/F_low*1.2 1],[1 1 0 0 ]);
h=h/norm(h);
% channel_length=length(h);
for k=1:length(TapPower)
a_low_temp=conv((randn(1,N_low+64)+1i*randn(1,N_lo w+64))/sqrt(2),h);%was ...,h)*TapPower(k);
a_low_temp=a_low_temp(65:N_low+64); % Generate first tap sampled at low freq
a_low(k,:)=a_low_temp;
channel_length=6;
%%%%Upsample To Fs%%%%
a_samples(k,:)=resample(a_low(k,:),Fs,F_low);
Sigvar=var(a_samples(k,:));
a_samples(k,:)=a_samples(k,:)*sqrt((TapPower(k)/Sigvar));
end
%%%%Initilaize Delay Vector%%%%
DelayInSamples=ceil(Delay*Fs);
DelayInSamples=DelayInSamplesDelayInSamples(1);
%%% Noise parameters
SNR=10^(SNRdB/10);
SNR_corrected=SNR*(1366+2256)/FFTsize; %need to generic fit for qpsk 2256 infobits + 1366 pilots
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% Began iteration
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
NumErrors=zeros(1,itertion);
for iter=1:itertion
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% Transmitter
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
xt_long_with_CP=[];
xt_long_bits=[];
for k=1:N_sym
InfoBits=randsrc(1,2256,[1 0]); % Generate 2256 information bits
xt_long_bits=[xt_long_bits InfoBits];
% Encode and Modulate
CodedBits= encode(henc,InfoBits);
% Rate Matching
ParityLocations=2256+randperm(2256);
ParityLocations=ParityLocations(1:ParityLength);
TotalLocations(:,k)=[1:2256,ParityLocations]; % Should be an output of the encoder
CodedBits=CodedBits(TotalLocations(:,k));
Symbols=ModulateQAM(CodedBits,Modulation);
%%%%pilots%%%%
x_f=zeros(FFTsize,1); x_f(SC_loc)=Symbols;
Xp=randsrc(1,length(PilotLocations),[1 1]); %Pilot generation
Pilots_PRBS(:,k)=Xp.'; % PRBS patterns
x_f(PilotLocations)=Xp.'; %insertion of the pilots
xt=ifft(x_f);
xt_withCP=[xt(end1023:end);xt];
xt_long_with_CP=[xt_long_with_CP;xt_withCP];
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% Channel
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
y=zeros(1,length(xt_long_with_CP));
%%%%Convolution%%%%
for k=max(DelayInSamples)+1:length(xt_long_with_CP)
for j=1:length(TapPower)
a(j)=xt_long_with_CP(kDelayInSamples(j))*a_samples(j,k);
end
y(k)=sum(a);
end
%%%Noise%%%
sigma_noise=norm(xt_long_with_CP)/sqrt(length(y)*SNR_corrected);
y=y+sigma_noise*(randn(size(y))+1i*randn(size(y)))/sqrt(2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% Receiver
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
yBeforeCP_RemoveMatrix=reshape(y.',FFTsize+1024,N_ sym);
yCP_RemoveMatrix=yBeforeCP_RemoveMatrix(1025:end,: );
Y_matrix_PostFFT=fft(yCP_RemoveMatrix);
for kk=1:N_sym,
% H_est = MMSE_CE(Y_matrix_PostFFT(:,kk),Pilots_PRBS(:,kk),P ilotLocations,FFTsize,4,h,SNR);
% h_est = ifft(H_est); h_DFT = h_est(1:channel_length);
% H_DFT = fft(h_DFT,FFTsize);
Y_Pilots=Y_matrix_PostFFT(PilotLocations,kk);
EstTimeDomainChannel=ifft(Y_Pilots.*Pilots_PRBS(:, kk));%Taking only pilots
EstTimeDomainChannel=EstTimeDomainChannel(1:(chann el_length)); % Kill contibutions outside the CP duration
H_est=fft(EstTimeDomainChannel,FFTsize); % est. channel
demodulatedSCs=Y_matrix_PostFFT(:,kk)./H_est;
Y_equalized_withoutPilots(:,kk)=demodulatedSCs(SC_ loc);
ppSNR(:,kk)=abs(H_est(SC_loc)).^2/(sigma_noise^2*FFTsize);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
% Decode Signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%
DecodedBits_long=[];
for kk=1:N_sym,
LLRs=SISO_LLR(Y_equalized_withoutPilots(:,kk),ppSN R(:,kk),Modulation);
LLR_rate_matched=zeros(1,2*2256);
LLR_rate_matched(TotalLocations(:,kk))=LLRs;
DecodedBits=decode(hdec, LLR_rate_matched);
DecodedBits_long=[DecodedBits_long DecodedBits];
end
NumErrors(1,iter)=sum(DecodedBits_long ~=xt_long_bits); % Count errors
end
Finely, for me the time taken by the code is related to FFTsize=4096 and Channel effects (LTV+noise) addition witch are applied to all signal xt_long_with_CP.
For IFFT/FFT size we can't change.
For Channel effect, why you don't apply effects to one symbol. It is long enough 2256.
For changing signal to matrix form. encode and decode function will not support matrix input. It is mentioned in matlab help (help fec.ldpcenc/encode)
Regards
Chaker
1 members found this post helpful.

1st April 2012, 16:05 #13
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
For Channel effect, why you don't apply effects to one symbol. It is long enough 2256.
Look a variable T, if the symbols/samples(xt_long_with_CP) is to small T is small >N_low is small > the channel is very short and it will not affect the Tx signal.
As for the iteration, I did it just because i dont get constant results with the LTV channel. At the end of all the iteration i am getting Numerror vector length 100
with different results. and it's not what i am expecting...
***Thanks for the code, thought like that and implement what you did.

1st April 2012, 17:45 #14
Re: How to graph BER with perfect results?
Hi,
There is so many constraint in your system.
I get final suggestion.
I propose to use frame notion where 1 frame = 10 * N_sym or more. By this, you can fix y length to 5120*10 for exemple.
Then you add while loop where simulation stop after FrameNb send (FrameNb=100frames) or stop after MaxErrorBits detected (MaxErrorBits=100 error bits).
while ((i_frame<FrameNb)&&(NumErrors<MaxErrorBits))
Good luck WUID
1 members found this post helpful.

1st April 2012, 19:40 #15
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
Muhamed , i liked your suggestion.
Code:N_sym=100; Frame=10*N_sym; FrameNb=100; i_frame=1; ErrorBits=0; MaxErrorBits=100; while ((i_frame<FrameNb)&&(ErrorBits<MaxErrorBits))... for k=1:Frame
anyone else (or you :) ) on another subject...
On the first channel section where i'm generating the Gaussian random process. I'm having hard times to understand and implement the math. I mean to ask, about my figures and normalization of the filter, taps variance... is it make scene to any1 who is reading this post ?

1st April 2012, 20:30 #16
Re: How to graph BER with perfect results?
I try to code the solution but i can't simulate.
Code:function NumErrors=OurLDPCEncandDec_without_MMSE(itertion,Rate,Modulation,SNRdB) % % close all; clc ; % clear all ; % Rate=2; % Modulation=4; % SNRdB=3; % itertion=5; % Rate: 2 for 1/2, 4 for 3/4, 6 for 5/6, 8 for 7/8 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Parameters Section %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Symbols_length=2256; CP_lenght=1024; N_sym=100; Frame=10*N_sym; FrameNb=100; MaxErrorBits=1000; if Rate==2 ParityLength=2256; elseif Rate==4 ParityLength=2256/3; elseif Rate==6 ParityLength=464; elseif Rate==8 ParityLength=320; elseif Rate==10 ParityLength=240; end; load 802_16e_4512_d5_H H mat; H1=H(:,1:2256); H2=H(:,2257:4512); H=[H2 H1]; henc = fec.ldpcenc(H); hdec = fec.ldpcdec(H); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % IFFT/FFT Parameters %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% FFTsize=2^(ceil(log2(Symbols_length))); ActiveSCs=Symbols_length; PilotLocations=1:4:FFTsize; %pilot loc /for mod=4,16 spacing = 3. /for mod=64 spacing =4 SC_loc=setdiff(1:FFTsize,PilotLocations); SC_loc=SC_loc(1:ActiveSCs); %info location %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %Channel %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Fs=10e6; %Sampling rate f_m=10; %Doppler spread N_samples=N_sym*(FFTsize+CP_lenght); T=N_samples/Fs; %Sample time % Taking channel parameters from ITU table TapPowerdB=[0 9.7 19.2 22.8]; % in dB Delay=[1 110 190 410]*1e9; % in Sec % ChannelLength=length(TapPower); a_samples=[]; %%%%Tap power in dB to linear%%%% TapPower=10.^(TapPowerdB/10); % TapPower=TapPower/norm(TapPower); % Normalize to unit channel power F_low=f_m*10; % usaully 10*f_m the low sampling % freq of the channel coefficients N_low=ceil(T*F_low); %%%%Gausssian random process pass throw filter(flat/Ushaped)%%%% h=firpm(64,[0 2*f_m/F_low 2*f_m/F_low*1.2 1],[1 1 0 0 ]); h=h/norm(h); % channel_length=length(h); for k=1:length(TapPower) a_low_temp=conv((randn(1,N_low+64)+1i*randn(1,N_low+64))/sqrt(2),h);%was ...,h)*TapPower(k); a_low_temp=a_low_temp(65:N_low+64); % Generate first tap sampled at low freq a_low(k,:)=a_low_temp; channel_length=6; %%%%Upsample To Fs%%%% a_samples(k,:)=resample(a_low(k,:),Fs,F_low); Sigvar=var(a_samples(k,:)); a_samples(k,:)=a_samples(k,:)*sqrt((TapPower(k)/Sigvar)); end %%%%Initilaize Delay Vector%%%% DelayInSamples=ceil(Delay*Fs); DelayInSamples=DelayInSamplesDelayInSamples(1); %%% Noise parameters SNR=10^(SNRdB/10); SNR_corrected=SNR*(1366+2256)/FFTsize; %need to generic fit for qpsk 2256 infobits + 1366 pilots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Began iteration %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% NumErrors=zeros(1,itertion); for iter=1:itertion %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Transmitter %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% i_frame=1; ErrorBits_all_frame=0; while ((i_frame<FrameNb)&&(ErrorBits_all_frame<MaxErrorBits)) xt_long_with_CP=[]; xt_long_bits=[]; for k=1:N_sym InfoBits=randsrc(1,2256,[1 0]); % Generate 2256 information bits xt_long_bits=[xt_long_bits InfoBits]; % Encode and Modulate CodedBits= encode(henc,InfoBits); % Rate Matching ParityLocations=2256+randperm(2256); ParityLocations=ParityLocations(1:ParityLength); TotalLocations(:,k)=[1:2256,ParityLocations]; % Should be an output of the encoder CodedBits=CodedBits(TotalLocations(:,k)); Symbols=ModulateQAM(CodedBits,Modulation); %%%%pilots%%%% x_f=zeros(FFTsize,1); x_f(SC_loc)=Symbols; Xp=randsrc(1,length(PilotLocations),[1 1]); %Pilot generation Pilots_PRBS(:,k)=Xp.'; % PRBS patterns x_f(PilotLocations)=Xp.'; %insertion of the pilots xt=ifft(x_f); xt_withCP=[xt(end1023:end);xt]; xt_long_with_CP=[xt_long_with_CP;xt_withCP]; end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Channel %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% y=zeros(1,length(xt_long_with_CP)); %%%%Convolution%%%% for k=max(DelayInSamples)+1:length(xt_long_with_CP) for j=1:length(TapPower) a(j)=xt_long_with_CP(kDelayInSamples(j))*a_samples(j,k); end y(k)=sum(a); end %%%Noise%%% sigma_noise=norm(xt_long_with_CP)/sqrt(length(y)*SNR_corrected); y=y+sigma_noise*(randn(size(y))+1i*randn(size(y)))/sqrt(2); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Receiver %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% yBeforeCP_RemoveMatrix=reshape(y.',FFTsize+1024,N_sym); yCP_RemoveMatrix=yBeforeCP_RemoveMatrix(1025:end,:); Y_matrix_PostFFT=fft(yCP_RemoveMatrix); for kk=1:N_sym, % H_est = MMSE_CE(Y_matrix_PostFFT(:,kk),Pilots_PRBS(:,kk),PilotLocations,FFTsize,4,h,SNR); % h_est = ifft(H_est); h_DFT = h_est(1:channel_length); % H_DFT = fft(h_DFT,FFTsize); Y_Pilots=Y_matrix_PostFFT(PilotLocations,kk); EstTimeDomainChannel=ifft(Y_Pilots.*Pilots_PRBS(:,kk));%Taking only pilots EstTimeDomainChannel=EstTimeDomainChannel(1:(channel_length)); % Kill contibutions outside the CP duration H_est=fft(EstTimeDomainChannel,FFTsize); % est. channel demodulatedSCs=Y_matrix_PostFFT(:,kk)./H_est; Y_equalized_withoutPilots(:,kk)=demodulatedSCs(SC_loc); ppSNR(:,kk)=abs(H_est(SC_loc)).^2/(sigma_noise^2*FFTsize); end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Decode Signal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% DecodedBits_long=[]; for kk=1:N_sym, LLRs=SISO_LLR(Y_equalized_withoutPilots(:,kk),ppSNR(:,kk),Modulation); LLR_rate_matched=zeros(1,2*2256); LLR_rate_matched(TotalLocations(:,kk))=LLRs; DecodedBits=decode(hdec, LLR_rate_matched); DecodedBits_long=[DecodedBits_long DecodedBits]; end ErrorBits_one_frame=sum(DecodedBits_long ~=xt_long_bits); % Count errors ErrorBits_all_frame=ErrorBits_all_frame+ErrorBits_one_frame; i_frame=i_frame+1; end % End while loop (all frame send or 100 errors bits detected) NumErrors(1,iter)=ErrorBits_all_frame; all_transmitted_bits=(i_frame1)*N_samples; ber_values(1,iter)=ErrorBits_all_frame/all_transmitted_bits; end % End iter loop
that is why i add
Code:NumErrors(1,iter)=ErrorBits_all_frame; all_transmitted_bits=(i_frame1)*N_samples; ber_values(1,iter)=ErrorBits_all_frame/all_transmitted_bits;
For your question, concerning the AWGN noise and OFDM yes but other part no :/
best regardsLast edited by ChakerQ; 1st April 2012 at 20:37.
1 members found this post helpful.

1st April 2012, 21:47 #17
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
thx Mohamed!!!
for now I've dropped down the iteration part for now. as forCode:(i_frame1)*N_samples
Frame=10*N_sym
Totalbits=2256*Frame*FrameNb
can u tell me about the "other part"? what do you think ? or is it completely wrong ?
***attaching my simulation code with your suggestion and my modifications including SISO_LLR function :)

1st April 2012, 22:46 #18
Re: How to graph BER with perfect results?
Hi WUID,
I have done some simulation.
I get ber=4 10^2 for SNRdB=0 and ber=5 10^5 for SNRdB=1.
According to ieee article (see files attached), this is not true.
I think there is a problem with awgn noise generation.
Began with verifying SNR value.
2 members found this post helpful.

1st April 2012, 23:19 #19
 Join Date
 Mar 2009
 Posts
 45
 Helped
 0 / 0
 Points
 1,072
 Level
 7
Re: How to graph BER with perfect results?
ya i know it looks that the encoder doing some magic and in addition you can see that the SNR is being corrected in SNR_corrected to the number of fft bins i am using e.g
Code:SNR_corrected=SNRdB*(1366 pilots + 2256 data bins)/FFTsize

1st April 2012, 23:19

1st April 2012, 23:35 #20
Re: How to graph BER with perfect results?
There is no magic in digital communication lol
I think the ber you get it is not for the snr declared.
To have corresponding values, you must correct the SNR expression.
I used this code to simulate AWGN
Code:EbNo = 0:2:15; EsNo= EbNo+10*log10(Used_subcarrier/len_fft)+ 10*log10(len_fft/Symbol_interval) + 10*log10(k); snr=EsNo  10*log10(len_fft/Symbol_interval); chan_awgn = awgn(ser_data,snr(ii),'measured');
http://www.dsplog.com/2008/08/26/ofd...nnelberbpsk/
+ Post New Thread
Please login