bpsk.m

Mtlab toolbox containing useful m files to generate LPI signals

%******************************************************************************
% bpsk.m
% Use:  Code to generated BPSK signals called by lpit.m
%
% Inputs:       Amplitude of the carrier signal
%		        Carrier signal frequency - f (Hz)
%               Sampling frequency - fs (HZ)
%               Desired SNR in dB
%               Number of bits per Barker code
%               Number of code periods
%               Number of cycles per Barker phase code (cpp)
%               Number of periods to view on graphs (CW signal)
%
% Output:        In-phase (I) and Quadarture (Q) components of the signal
%               Plots
%******************************************************************************


clc;clear all;
disp('************************************');
disp('**************** BPSK *****************');
disp('***************************************');

%	DEFAULT VARIABLES 
A  = 1;						        %Amplitud of the carrier signal
f    = 1e3;					        %Carrier signal frequency - f (Hz)
fs = 7e3;                           %Sampling frequency - fs (Hz)
SNRdb  = 0;						    %Desired SNR in dB
barker = 7;                         %Number of bits per Barker code
np   = 175;						    %Number of code periods
cpp = 1;						    %Number of cycles per Barker bit
NPV  = 55;					        %Number of periods to view on graphs (CW signal)
  
% NEW INPUT 
newvar = 1;
while newvar == 1;
    disp(' ')
    disp('WHICH PARAMETER DO YOU WANT TO SET ?  ')
    disp(' ')
    fprintf('1. Amplitude of the carrier signal - (A)= %g.\n', A)
    fprintf('2. Carrier frequency - f (Hz) = %g.\n', f)
    fprintf('3. Sampling frequency - fs (Hz)= %g.\n', fs)
    fprintf('4. Signal to noise ratio - SNRdb (dB)= %g.\n', SNRdb)
    fprintf('5. Number of bits per Barker code - barker (13/11/7)= %g.\n', barker)
    fprintf('6. Number of code periods- np= %g.\n', np)
    fprintf('7. Number of cycles per Barker phase code - cpp = %g.\n',cpp)
    fprintf('8. Number of code periods to view on graphs= %g.\n', NPV)
    fprintf('9. No changes\n')
    disp(' ')
    option= input('Select a option: ');
   
    switch option
    case 1   
        A=input('New amplitude of carrier signal (Volts) = ');       
    case 2
        f=input('New carrier signal frequency (Hz) = ');      
    case 3
        fs=input('New sampling frequency (Hz) = ');   
    case 4
        SNRdb=input('New signal to Noise Ratio (dB) = ');        
    case 5
        barker=input('New number of bits for Barker Code = ');      
    case 6
        np=input('New number of periods of carrier signal to generate = '); 
    case 7        
        cpp=input('New number of cycles per Barker phase = ');      
    case 8
        NPV=input('New number of periods to view on graphs = ');
    case 9
        newvar = 0;
    end
     clc;
end

%  Variables below are calculated based upon the above user defined variables
ti = 1/(f*100);				            %determine suitable time increment to plot carrier signal
t  = 0:ti:np/f;				               %set up time vector
xt = A*sin(2*pi*f.*t);	            %representation of carrier signal(continuous)

SAR=floor(fs/f);
n  = 0:1:SAR*np;				   %sample vector n, covers np periods of carrier frequency

xnT  = A*cos(2*pi.*n*f/fs);     %create vector of sampled function, for I
xnT2 = A*sin(2*pi.*n*f/fs);   %create shifted version of function, for Q


%************
% This section creates the modulating signal and modulates the carrier
%************
fm = f/cpp;					       %set frequency of modulating signal
pw = floor(fs/fm);						%determine # of samples (of n) per Barker bit

%create sequence of the (7, 11 or 13)-bit Barker code with pw samples of each bit,
%string several sequences end to end to match length of sampled signal

if barker==13
    brk = [ones(1,pw*5),-(ones(1,pw*2)),ones(1,pw*2),-ones(1,pw),ones(1,pw),-ones(1,pw),ones(1,pw)];%13 bit
elseif barker==11
    brk = [ones(1,pw*3),-(ones(1,pw*3)),ones(1,pw),-ones(1,pw),-ones(1,pw),ones(1,pw),-ones(1,pw)];% 11 bit
else
    brk = [ones(1,pw*3),-(ones(1,pw*2)),ones(1,pw),-ones(1,pw)];% 7 bits
end

brkseq = [];					            %initialize brkseq to hold a sequence of brk
ns = floor(length(n)/(pw*barker)); %integer number of complete 13-bit Barker sequences

for step = 1:1:ns;				            %create multiple n-bit Barker sequence vector
   brkseq = [brkseq,brk];
end

%************
% Here is where the modulation is done and the I and Q channels created
%************

I   = xnT(1:barker*ns*pw).*brkseq;	      %modulate sampled carrier with Barker sequence, I
Q   = xnT2(1:barker*ns*pw).*brkseq;	    %modulate shifted function with Barker sequence, Q
rl  = length(I);						                %get length of the new vector

%******************
% In this section we create White Gaussian Noise (WGN) and add it to both the 
% I and the Q channels.  In the original program the noise was developed using 
% the normrnd MATLAB function. To standarize noise for all the programs, we
%use the wgn Matlab function.
%******************

[a,b]=size(xnT);
SNR=10^(SNRdb/10);
power=10*log10(A^2./(2*SNR));%calculate SNR in dB for WGN function
noise=wgn(a,b,power); %calculate noise at specified SNR
IwN =I+noise(1:length(I));             %add noise to I 
QwN=Q+noise(1:length(Q));           %add noise to Q 
 xnTwN = xnT + noise;					%add WGN to sampled signal, for reference only

%Generation of plots 
disp(' ')
plt = input('Do you want to generate plots of the signal (Y/y or N/n) ?','s');
disp(' ')
if (plt == 'Y') | (plt =='y')
    disp(' ')
       %*************************
        % Plot a CW representation of the signal and a sampled version of the signal
        %*************************
        figure(1);
        %plot carrier frequency (CW) and sampled function
        subplot(2,1,1);				%CW representation
        
        if ( np > NPV );				%just show some of the data to keep screen uncluttered
            plot(t(1:1:100*NPV),xt(1:1:100*NPV));grid;
            axis([0 0.012 -A*1.5 A*1.5]);
            title(['Continuous Wave Signal, f = ',num2str(f),', #periods= ',num2str(np),', ',num2str(NPV),' shown']);
        else
            plot(t,xt);	grid;			%plot carrier frequency as continuous function

            title(['Continuous Wave Signal, f = ',num2str(f),', #periods= ',num2str(np)]);
        end

        xlabel('Time (sec)'),ylabel('Amplitude');


        subplot(2,1,2);				%show sampled version of carrier signal

        if ( np > NPV )				%keep display uncluttered
            stem(n(1:1:SAR*NPV),xnT(1:1:SAR*NPV));grid;%plot sampled carrier frequency (portion of)
            axis([1 SAR*NPV -A*1.5 A*1.5]);
            title(['Sampled signal, fs = ',num2str(fs),', #periods= ',num2str(np),', ',num2str(NPV),' shown']);
        else
            stem (n,xnT);grid;				%plot sampled carrier frequency
            title(['Sampled signal, fs = ',num2str(fs),', #periods= ',num2str(np)]);
        end

        xlabel('n (sample #)'),ylabel('x[n]');
        title(['Sampled signal, fs = ',num2str(fs),', #periods= ',num2str(np)]);  
        
        %***********************
        %Plot Modulating signal and PSD
        %***********************
        figure(2);							
        subplot(2,1,1);					%show modulating signal, distinguishing between bits
        set (gca,'YLim',[1.2*min(brk) 1.2*max(brk)]);%force y axis limits to better see signal
        xlabel('n (sample #)'),ylabel('Amplitude');
        if ( np > NPV )					%limit displayed amount
            flag = 1;						%flag tracks bit changes to change display
   
            for step=1:pw:SAR*NPV			%step through each Barker bit
                if ( flag == 1 )			%change plotting for each bit
                    hold on,plot(n(step:1:step+pw-1),brkseq(step:1:step+pw-1));grid;
                else
                    hold on,plot(n(step:1:step+pw-1),brkseq(step:1:step+pw-1),'x');grid
                end
      
            flag = -flag;				%swap flag
            end
            step=ceil(step/pw);			%determine # of bits displayed
            title([num2str(barker),'-bit Barker Sequence',' # periods shown = ' num2str(step/barker)]);
        else
            hold on, plot(n(1:rl),brkseq(1:rl),'x');grid;%plot the modulating signal
            title([num2str(barker),'-bit Barker Sequence']);
        end
        subplot(2,1,2);	%show Power Spectral Density of the modulating signal
        psd(brkseq,[],fs);
        title(['Power Spectral Density of the', num2str(barker),'-bit Barker Sequence', ' # samples = ',num2str(length(brkseq))]);

        %************************ 
        %Plot the sampled signal and the Barker sequence in subplot 1, then show the 
        % results after modulation in subplot 2
        %************************
        figure(3);					
        subplot(2,1,1);			

        if ( np > NPV )					%show only as much as the user asked for
            stem(n(1:1:SAR*NPV),xnT(1:1:SAR*NPV));
            axis([0 SAR*NPV -A*1.5 A*1.5]);
            hold on, plot(n(1:1:SAR*NPV),brkseq(1:1:SAR*NPV),'r--');%modulating signal in red
        else
            stem(n(1:rl),xnT(1:rl));						%plot the sampled carrier signal
            axis([0 SAR*NPV -A*1.5 A*1.5]);
            hold on, plot(n(1:rl),brkseq(1:rl),'r--');%plot the modulating signal in red
        end

        xlabel('n (sample #)'),ylabel('Amplitude');
        title([num2str(barker),'-bit Barker Sequence Overlayed on Sampled Signal',', cpp= ',num2str(cpp)]);

        subplot(2,1,2);					%show results after modulation (I signal)
        
        if ( np > NPV )					
            stem(n(1:1:SAR*NPV),I(1:1:SAR*NPV));grid
            axis([0 SAR*NPV -A*1.5 A*1.5]);
        else
            stem(n(1:rl),I);grid;				%plot the modulated signal
            axis([0 SAR*NPV -A*1.5 A*1.5]);
        end
        xlabel('n (sample #)'),ylabel('Amplitude');
        title(['Sampled Signal Modulated by',num2str(barker),'-bit Barker Sequence',' cpp= ',num2str(cpp)]);

        %*********************
        %PSD of the carrier frequency before and afeter modulation
        %*********************
        figure(4);							
        subplot(2,1,1);					%plot PSD of unmodulated carrier
        psd(xnT,[],fs);
        xlabel('Frequency (Hz)');
        title(['PSD of Original Sampled Function, fs = ',num2str(fs),', # samples= ',num2str(rl)]);


        subplot(2,1,2);					%plot PSD of modulated carrier (I)
        psd(I,[],fs);
        xlabel('Frequency (Hz)');
        title(['PSD of Modulated Function, fs = ',num2str(fs),', # samples= ',num2str(rl),', cpp= ',num2str(cpp)]);
        
        %**********************
        %Show statistics and portion of the noise vector
        %**********************
        figure(5);											
        subplot(2,1,1);									%show noise signal vs time (sample)

        if (length(noise) > 100 );
        plot(noise(1:100));grid;
        title(['Plot of Noise Signal, channel I, # samples= ',num2str(rl),', only 100 shown']);
        else
        plot(noise);grid;							%look at complete noise vector
        title(['Plot of Noise Signal, channel I, # samples= ',num2str(rl)]);
        end

        xlabel('Sample #'),ylabel('Magnitude');


        subplot(2,1,2);						%show histogram, mu, and sigma of noise in I
        hist(noise,30);						%show that the noise has normal distribution
        xlabel('Magnitude'),ylabel('Frequency');
        title(['Histogram of Noise Signal, channel I']);

        %**********************
        %Compare I signal without noise to I signal with noise
        %**********************
        % Figure 6 shows how the noise affects the signal.  Depending upon parameter
        % settings, the user may wish to use the plot commands instead of the stem
        % commands that are used.
        figure(6);								
        subplot(2,1,1);						%show I signal without noise

        if ( np > NPV )
            stem(n(1:1:SAR*NPV),I(1:1:SAR*NPV));grid;
            axis([0 SAR*NPV -A*1.5 A*1.5]);
            %plot(n(1:1:SAR*NPV),I(1:1:SAR*NPV),'x');   
        else
            stem(n(1:rl),I);grid;			%plot the modulated signal
            axis([0 SAR*NPV -A*1.5 A*1.5]);
        end

        xlabel('n (sample #)'),ylabel('Amplitude');
        title(['Sampled Signal Modulated by', num2str(barker),'-bit Barker Sequence',', cpp= ',num2str(cpp)]);


        subplot(2,1,2);						%show I signal after noise has been added

        if ( np > NPV )
            stem(n(1:1:SAR*NPV),IwN(1:1:SAR*NPV));grid;
            %plot(n(1:1:SAR*NPV),IwN(1:1:SAR*NPV),'x');
        else
            stem(n(1:rl),IwN); grid;					%plot the modulated signal with noise added
        end

        xlabel('n (sample #)'),ylabel('Amplitude');
        title(['Sampled Signal Modulated by', num2str(barker),'-bit Barker Sequence, WGN added',', SNR= ',num2str(SNRdb ,3),'dB']);

        figure(7);
        subplot(2,1,1);					%PSD of unmodulated signal with noise added
        psd(xnTwN,[],fs);
        title(['PSD of Unmodulated Signal with Noise, fs = ',num2str(fs),', # samples= ',num2str(rl),', cpp= ',num2str(cpp),', SNR= ',num2str(SNRdb,2),'dB']);

        subplot(2,1,2);					%PSD of modulated signal with noise added
        psd(IwN,[],fs);
        title(['PSD of Modulated Function with Noise, fs = ',num2str(fs),', # samples= ',num2str(rl),', cpp= ',num2str(cpp)]);
        
        
else
    disp('Signal not plotted')
    fprintf('\n\n')
end



% Save as I and Q with the Noise included
format short e;
IwN=IwN';
QwN=QwN';
disp(' ')

plt2 = input('Do you want to save the new signal (Y/y or N/n) ?','s');
disp(' ')
if (plt2 == 'Y') | (plt2 =='y')

    I2=I; Q2=Q;
    I=IwN; Q=QwN;
    ff=floor(f/1e3);
    ffs=floor(fs/1e3);
    save(['B_', num2str(ff), '_', num2str(ffs), '_', num2str(barker), '_', num2str(cpp),'_',num2str(SNRdb)],'I','Q');
    I=I2';
    Q=Q2';
    save(['B_' num2str(ff) '_' num2str(ffs) '_' num2str(barker) '_' num2str(cpp) '_s'],'I','Q');
    disp(' ');
    disp(['Signal and noise save as :  B_', num2str(ff), '_', num2str(ffs), '_', num2str(barker), '_', num2str(cpp),'_',num2str(SNRdb)]);
    disp(['Signal only save as :         B_', num2str(ff), '_', num2str(ffs), '_', num2str(barker), '_', num2str(cpp),'_s']);
    disp(['Directory:                       ' num2str(cd)]); 
else
    disp(' ')
    disp('Signal not saved')
    fprintf('\n\n')
end