+ Post New Thread
Page 1 of 2 1 2 LastLast
Results 1 to 20 of 31
  1. #1
    Member level 2
    Points: 1,072, Level: 7

    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=[1e-3 1e-4 1e-5 0 0 0 0]
    for 16-qam i have this vector: ber_16qam=[1e-2 1e-3 1e-4 1e-5 0 0 0]
    for 64-qam i have this vector: ber_64qam=[0.5e-2 1e-2 1e-3 1e-4 1e-5 0 0]

    now, i want to graph them all in one graph. tries berfit and i got '[]'.

  2. #2
    Member level 2
    Points: 1,072, Level: 7

    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 ?



  3. #3
    Member level 3
    Points: 2,104, Level: 10
    Achievements:
    7 years registered

    Join Date
    Aug 2007
    Location
    Germany
    Posts
    64
    Helped
    9 / 9
    Points
    2,104
    Level
    10

    Re: How to graph BER with perfect results?

    Quote Originally Posted by WUID View Post
    anyone can help me with that ?
    It seems that there is some problem in simulation as you are getting 0 BER just after 1e-5. Just plot them on one graph using plot but it wont be a water fall curve.



  4. #4
    Junior Member level 1
    Points: 801, Level: 6

    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.

  5. #5
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    Re: How to graph BER with perfect results?

    Hi,

    Try this


    snr=0:5:30;
    ber_qpsk=[1e-3 1e-4 1e-5 0 0 0 0];
    ber_16qam=[1e-2 1e-3 1e-4 1e-5 0 0 0];
    ber_64qam=[0.5e-2 1e-2 1e-3 1e-4 1e-5 0 0];
    figure(1);
    semilogy(snr,ber_qpsk,'k-o','LineWidth',1);
    hold on
    semilogy(snr,ber_16qam,'k-x','LineWidth',1);
    hold on
    semilogy(snr,ber_16qam,'k-d','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.

    •   Alt30th March 2012, 11:56

      advertising

        
       

  6. #6
    Member level 2
    Points: 1,072, Level: 7

    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 ?



    •   Alt31st March 2012, 17:26

      advertising

        
       

  7. #7
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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.

    •   Alt31st March 2012, 21:31

      advertising

        
       

  8. #8
    Member level 2
    Points: 1,072, Level: 7

    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(end-1023: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]*1e-9; % 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/U-shaped)%%%%
           
           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=DelayInSamples-DelayInSamples(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(k-DelayInSamples(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 1-i -1+i -1-i])/sqrt(2);
    
    % 16QAM mapping
    Map_16QAM= ((2/sqrt(10))* ...
        [ 0.5+0.5*i 0.5+1.5*i 0.5-0.5*i 0.5-1.5*i ...
        1.5+0.5*i 1.5+1.5*i 1.5-0.5*i 1.5-1.5*i ....
        -0.5+0.5*i -0.5+1.5*i -0.5-0.5*i -0.5-1.5*i ....
        -1.5+0.5*i -1.5+1.5*i -1.5-0.5*i -1.5-1.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.5-1.5*i 1.5-0.5*i 1.5-2.5*i 1.5-3.5*i   ...
                0.5+1.5*i 0.5+0.5*i 0.5+2.5*i  0.5+3.5*i  ...
                0.5-1.5*i 0.5-0.5*i 0.5-2.5*i 0.5-3.5*i   ...
                2.5+1.5*i 2.5+0.5*i 2.5+2.5*i  2.5+3.5*i  ...
                2.5-1.5*i 2.5-0.5*i 2.5-2.5*i 2.5-3.5*i   ...
                3.5+1.5*i 3.5+0.5*i 3.5+2.5*i  3.5+3.5*i  ...
                3.5-1.5*i 3.5-0.5*i 3.5-2.5*i 3.5-3.5*i   ...
                -1.5+1.5*i -1.5+0.5*i -1.5+2.5*i  -1.5+3.5*i  ...
                -1.5-1.5*i -1.5-0.5*i -1.5-2.5*i -1.5-3.5*i   ...
                -0.5+1.5*i -0.5+0.5*i -0.5+2.5*i  -0.5+3.5*i  ...
                -0.5-1.5*i -0.5-0.5*i -0.5-2.5*i -0.5-3.5*i   ...
                -2.5+1.5*i -2.5+0.5*i -2.5+2.5*i  -2.5+3.5*i  ...
                -2.5-1.5*i -2.5-0.5*i -2.5-2.5*i -2.5-3.5*i   ...
                -3.5+1.5*i -3.5+0.5*i -3.5+2.5*i  -3.5+3.5*i  ...
                -3.5-1.5*i -3.5-0.5*i -3.5-2.5*i -3.5-3.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 ----------

    Quote Originally Posted by Mohamed Chaker View Post
    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 ?
    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.



  9. #9
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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.

  10. #10
    Junior Member level 1
    Points: 801, Level: 6

    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?

    Quote Originally Posted by rhearty1 View Post
    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
    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 integrate-and-dump 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.

  11. #11
    Member level 2
    Points: 1,072, Level: 7

    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 ?!
    I did the iterative loop because i got non consistent results with the LTV channel. e.g for 100 iter i getting numerror vector that most of it contain "zero" errors and part of it is with different number of errors. I believe is due to my misunderstanding the math completely and some of the vectors in the LTV section aren't normalized well.

    - encode(henc,InfoBits) doesn't support Matrix input.
    the encoder/decoder function receives 'henc' sparse matrix.

    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.



  12. #12
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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]*1e-9; % 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/U-shaped)%%%%

    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=DelayInSamples-DelayInSamples(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(end-1023: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(k-DelayInSamples(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
    I still don't understand the utility of iterative loop iter=1:itertion. Based on my understanding to the code, you calculate the BER 'itertion' time for the same configuration. Only random signal, noise and interleaving are modified witch supposed not affecting the BER.

    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.

  13. #13
    Member level 2
    Points: 1,072, Level: 7

    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.



  14. #14
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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))
    This loop will eliminate the send off 1000 symbol every simulation.

    Good luck WUID


    1 members found this post helpful.

  15. #15
    Member level 2
    Points: 1,072, Level: 7

    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
    Hope that's what your meant...

    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 ?



  16. #16
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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]*1e-9; % 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/U-shaped)%%%%
           
    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=DelayInSamples-DelayInSamples(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(end-1023: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(k-DelayInSamples(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_frame-1)*N_samples; 
    ber_values(1,iter)=ErrorBits_all_frame/all_transmitted_bits;
    end % End iter loop
    Note that the number of transmitted bits is variable.
    that is why i add
    Code:
    NumErrors(1,iter)=ErrorBits_all_frame;
    all_transmitted_bits=(i_frame-1)*N_samples; 
    ber_values(1,iter)=ErrorBits_all_frame/all_transmitted_bits;
    Try to simulate and give me feedback :)

    For your question, concerning the AWGN noise and OFDM yes but other part no :/

    best regards
    Last edited by ChakerQ; 1st April 2012 at 20:37.


    1 members found this post helpful.

  17. #17
    Member level 2
    Points: 1,072, Level: 7

    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 for
    Code:
    (i_frame-1)*N_samples
    it's wrong. for all type of modulation I am transmitting 2256 InfoBits. so
    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 :)



  18. #18
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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.

  19. #19
    Member level 2
    Points: 1,072, Level: 7

    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
    -> so the SNRdB value decrease



    •   Alt1st April 2012, 23:19

      advertising

        
       

  20. #20
    Member level 1
    Points: 660, Level: 5
    ChakerQ's Avatar
    Join Date
    Mar 2012
    Location
    Tunisia
    Posts
    40
    Helped
    19 / 19
    Points
    660
    Level
    5

    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');
    I follow this tutorial
    http://www.dsplog.com/2008/08/26/ofd...nnel-ber-bpsk/



--[[ ]]--