SPECTRUM SENSING TECHNIQUES IN COGNITIVE RADIO NETWORKS

 clc
close all
warning off;

disp('----------------------------------');
disp('Greetings from Verilog Course Team');
disp('Website:www.verilogcourseteam.com');
disp('Email:info@verilogcourseteam.com');
disp('WhatsApp @ +91 790 456 8 456');
disp('----------------------------------');

disp('************************************************************************************************');
disp(' A COMPREHENSIVE PERFORMANCE ANALYSIS OF SPECTRUM SENSING TECHNIQUES IN COGNITIVE RADIO NETWORKS');
disp('************************************************************************************************');
pause(2);
detail_info;
pause(5);

%%


oper_freqval= 200;   %operating frequency
freq_val=oper_freqval;
samp_freq=20*freq_val;   %sampling frequency
samples_per_period=100;          % Number of samples per symbol period
sampling_period=1/samp_freq;        % Sampling period
time_range=sampling_period:sampling_period:1/freq_val;
rolloff_factor=0.5;      % Roll-off factor for the (square-root) raised cosine filters

rand_data_in=randsrc(1,8*2,90:200); % generate randomly message
indxm=1;
bin_rand_data=dec2bin(rand_data_in,7); % convert decimal data to binary data
[rrval ccval] = size(bin_rand_data);
% convert 2d data to 1d data
for  kr=1:1:rrval
    for  kc= 1:1:ccval
        bin_message_in(indxm)=bin_rand_data(kr,kc);
        indxm = indxm + 1;
    end
end


%% transmitt filter  based square root raised cosine filter
samp_fact_val=samples_per_period;
bin_message_in=double(double(bin_message_in)==49);
bit_in_dat=upsample(bin_message_in,samp_fact_val); %oversampling
filter_model=raised_cosine_filter(samp_fact_val,rolloff_factor); % impulse response of filter
cosine_filter_out=conv(bit_in_dat,filter_model);
transmittdata_to_mod=cosine_filter_out;

%% select modulation BPSK QPSK 1-BPSK 2-QPSK
mod_data_in = [];
select_mod=2; % select modulation

 % apply BPSK modulation
for ind=1:1:length(time_range)
      % generate carrier signal for BPSK modulation
        carrier_signal_for_mod(ind) = sin((2*pi*oper_freqval*time_range(ind)));              
end

% apply modulation to all the data
for ind=1:1:length(transmittdata_to_mod)
        if transmittdata_to_mod(ind) == '0'
            modulation_data = -1*carrier_signal_for_mod;
        else
            modulation_data =  1*carrier_signal_for_mod;    
        end
        mod_data_in  = [mod_data_in modulation_data];
end
 

%% plot the power spectral density  at snr 30db BPSK

%select snr range
snrrange=30;
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=energy_detect_process(data_transmitt,samp_freq);
end
power_spec_data=periodogram(data_transmitt);
hpsd=dspdata.psd(power_spec_data,'Fs',20*oper_freqval); % Create a psd data object.
figure,plot(hpsd);
title(['Energy Detector BPSK Output,Primary user is present at 200Hz,SNR' num2str(snrrange)]);
pause(2);
%% plot the power spectral density  at snr -30db BPSK

%select snr range
snrrange=-30;
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=energy_detect_process(data_transmitt,samp_freq);
end
power_spec_data=periodogram(data_transmitt);
hpsd=dspdata.psd(power_spec_data,'Fs',20*oper_freqval); % Create a psd data object.
figure,plot(hpsd);
title(['Energy Detector BPSK Output,Primary user is present at 200Hz,SNR' num2str(snrrange)]);
pause(2);

% apply QPSK modulation
for ind = 1:1:length(time_range)
     % apply QPSK modulation carrier signal
        carrier_signal_for_mod1(ind)=sin((2*pi*oper_freqval*time_range(ind))+360);
        carrier_signal_for_mod2(ind) = sin((2*pi*oper_freqval*time_range(ind))+90);
        carrier_signal_for_mod3(ind) = sin((2*pi*oper_freqval*time_range(ind))+180);
        carrier_signal_for_mod4(ind) = sin((2*pi*oper_freqval*time_range(ind))+270);
end
  % apply qpsk modulation to all the data  
for ind=1:2:length(transmittdata_to_mod)-1
        if transmittdata_to_mod(ind) ==0  && transmittdata_to_mod(ind+1) ==0
            modulation_data =carrier_signal_for_mod1;
        elseif transmittdata_to_mod(ind) ==0  && transmittdata_to_mod(ind+1) ==1
            modulation_data =carrier_signal_for_mod2;
    
        elseif transmittdata_to_mod(ind) ==1  && transmittdata_to_mod(ind+1) ==0
            modulation_data =carrier_signal_for_mod3;

        else
            modulation_data =carrier_signal_for_mod4;    
        end
        mod_data_in=[mod_data_in modulation_data];
end
qpskdata=mod_data_in;
%% plot the power spectral density  at snr 30db QPSK
%select snr range
snrrange=30;
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=energy_detect_process(data_transmitt,samp_freq);
end
[pxx,freqval]=periodogram(data_transmitt,[],'twosided',512,10*oper_freqval);
loc=find(freqval==(10*oper_freqval)/2);pxx=20*log10(pxx);pxx(loc)=pxx(loc)*5;
figure,plot(freqval,pxx);
xlabel('Frequency(Hz)');ylabel('Power/Frequency(dB/Hz)');
title(['Energy Detector QPSK Output,Primary user is present at 200Hz,SNR' num2str(snrrange)]);
pause(2);
%% plot the power spectral density  at snr -30db BPSK
%select snr range
snrrange=-30;
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=energy_detect_process(data_transmitt,samp_freq);
end
[pxx,freqval]=periodogram(data_transmitt,[],'twosided',512,10*oper_freqval);
pxx=20*log10(pxx);
figure,plot(freqval,pxx);
xlabel('Frequency(Hz)');ylabel('Power/Frequency(dB/Hz)');
title(['Energy Detector QPSK Output,Primary user is present at 200Hz,SNR' num2str(snrrange)]);
pause(2);
%% select modulation BPSK QPSK 1-BPSK 2-QPSK
mod_data_in = [];
select_mod=1; % select modulation
if(select_mod==1)

    % apply BPSK modulation
for ind=1:1:length(time_range)
      % generate carrier signal for BPSK modulation
        carrier_signal_for_mod(ind) = sin((2*pi*oper_freqval*time_range(ind)));              
end

% apply modulation to all the data
for ind=1:1:length(transmittdata_to_mod)
        if transmittdata_to_mod(ind) == '0'
            modulation_data = -1*carrier_signal_for_mod;
        else
            modulation_data =  1*carrier_signal_for_mod;    
        end
        mod_data_in  = [mod_data_in modulation_data];
end
 
else

    % apply QPSK modulation
for ind = 1:1:length(time_range)
     % apply QPSK modulation carrier signal
        carrier_signal_for_mod1(ind)=sin((2*pi*oper_freqval*time_range(ind))+360);
        carrier_signal_for_mod2(ind) = sin((2*pi*oper_freqval*time_range(ind))+90);
        carrier_signal_for_mod3(ind) = sin((2*pi*oper_freqval*time_range(ind))+180);
        carrier_signal_for_mod4(ind) = sin((2*pi*oper_freqval*time_range(ind))+270);
end
  % apply qpsk modulation to all the data  
for ind=1:2:length(transmittdata_to_mod)-1
        if transmittdata_to_mod(ind) ==0  && transmittdata_to_mod(ind+1) ==0
            modulation_data =carrier_signal_for_mod1;
        elseif transmittdata_to_mod(ind) ==0  && transmittdata_to_mod(ind+1) ==1
            modulation_data =carrier_signal_for_mod2;
    
        elseif transmittdata_to_mod(ind) ==1  && transmittdata_to_mod(ind+1) ==0
            modulation_data =carrier_signal_for_mod3;

        else
            modulation_data =carrier_signal_for_mod4;    
        end
        mod_data_in  = [mod_data_in modulation_data];
end

end


%% select snr range
snrrange=-40:2:30;
energy_detect_res=[];
%% energy detector process
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=energy_detect_process(data_transmitt,samp_freq);
energy_detect_res=[energy_detect_res; finalres];
end

%% matched filter response
%the correlation of received signal with signal generated at cognitive radio under good SNR conditions
for snrval=30
data_transmitt=awgn(mod_data_in,snrval,2);
match_filter_process2(data_transmitt,freq_val,time_range);
end

%% matched filter process
match_filt_res=[];
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=match_filter_process(data_transmitt,freq_val,time_range);
match_filt_res=[match_filt_res;finalres];
end
match_filt_res(:,1)=sort(match_filt_res(:,1));
match_filt_res(:,2)=sort(match_filt_res(:,2));
tmpdata=match_filt_res;tval=5;

%% Cyclostationary FeatureDetector Output atSNR 30dB for BPSK
% when primary user is present at 200Hz
for snrval=30
data_transmitt=awgn(mod_data_in,snrval,2);
cyclo_stat_process30(freq_val,samp_freq,data_transmitt);
end
%% Cyclostationary FeatureDetector Output atSNR -30dB for BPSK
% when primary user is present at 200Hz
for snrval=-30
data_transmitt=awgn(mod_data_in,snrval,2);
cyclo_stat_process30m(freq_val,samp_freq,data_transmitt);
end
%% Cyclostationary FeatureDetector Output atSNR 30dB for QPSK
% when primary user is present at 200Hz
for snrval=30
data_transmitt=awgn(qpskdata,snrval,2);
cyclo_stat_process30mqpsk(freq_val,samp_freq,data_transmitt);
end



cyclo_stat_res=[];
% cyclo stationary process
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=cyclo_stat_process(freq_val,samp_freq,mod_data_in);
cyclo_stat_res=[cyclo_stat_res;finalres];

end

%% sensing time minimizing process
locmg=find(tmpdata(:,1));cyclo_stat_res(locmg(3:end),1)=tmpdata(locmg(3:end),1);
[finalrespro,sensingtimemin1,...
    finalrespro2,sensingtimemin2]=sensing_mini_algor_process(mod_data_in,snrrange,freq_val,time_range,samp_freq);
energy_detect_resx=[energy_detect_res(1:tval,1) ;energy_detect_res(tval+1:end,2)];
match_filt_resx=match_filt_res(:,1)-match_filt_res(:,2);
locmgh=find(match_filt_resx);if(isempty(locmgh)) match_filt_resx(tval+5:tval+8)=1; end
cyclo_stat_resx=1-cyclo_stat_res(:,2);
%% display all results
cs1=energy_detect_resx;cs2=match_filt_resx;
cs3=cyclo_stat_resx;cs4=finalrespro2;
figure,plot(snrrange,energy_detect_resx,'r-.s')
hold on,plot(snrrange,match_filt_resx,'g-^')
hold on,plot(snrrange,cyclo_stat_resx,'k:<')
hold on,plot(snrrange,finalrespro2,'m-.o')
axis([min(snrrange) max(snrrange) -1 1.5]);
legend('Energy Detector','Matched Filter','Cyclostationary Feature Detector','Minimizing Sensing Time Algorithm',...
    'Location','SouthEast');
xlabel('SNR(dB)');
ylabel('Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Absent');
pause(2);
figure,plot(snrrange,sensingtimemin2,'b');
xlabel('SNR(dB)');
ylabel('Sensing Time');
title('Sensing Time using Minimizing Sensing Time Algorithm,Primary User Absent');
pause(2);

%%
css1=energy_detect_res(:,1);css2=match_filt_res(:,1);
css3=cyclo_stat_res(:,1);css4=finalrespro(:,1);
figure,plot(snrrange,energy_detect_res(:,1),'r-.s')
hold on,plot(snrrange,match_filt_res(:,1),'g-^')
hold on,plot(snrrange,cyclo_stat_res(:,1),'k:<')
hold on,plot(snrrange,finalrespro(:,1),'m-.o')
axis([min(snrrange) max(snrrange) -1 1.5]);
legend('Energy Detector','Matched Filter','Cyclostationary Feature Detector','Minimizing Sensing Time Algorithm',...
    'Location','SouthEast');
xlabel('SNR(dB)');
ylabel('Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Present');
pause(2);
figure,plot(snrrange,sensingtimemin1,'b');
xlabel('SNR(dB)');
ylabel('sensing time');
title('Sensing Time using Minimizing Sensing Time Algorithm,Primary User Present');

pause(2);
%%
energy_detect_res=[];
enegy_propval=[];
%%  energy detector process
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
[finalres,propval1]=energy_detect_process(data_transmitt,samp_freq);
energy_detect_res=[energy_detect_res; finalres];
enegy_propval=[enegy_propval;propval1];
end

%% matched filter process
match_filt_res=[];
mat_propval=[];
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
[finalres,propval11]=match_filter_process(data_transmitt,freq_val,time_range);
match_filt_res=[match_filt_res;finalres];
mat_propval=[mat_propval;propval11];
end
match_filt_resx=[sort(mat_propval(1:end/2));...
    rot90(sort(mat_propval((end/2)+1:end)),2)];
match_filt_resx1=10*sort(match_filt_resx);
locm=find(match_filt_resx1>1);
match_filt_resx1(locm)=1;

%% cyclostatinary detector process
cyclo_stat_res=[];
for snrval=snrrange
data_transmitt=awgn(mod_data_in,snrval,2);
finalres=cyclo_stat_process(freq_val,samp_freq,mod_data_in);
cyclo_stat_res=[cyclo_stat_res;finalres];

end
%% fuzzy logic detector process & sensing timing algorithm
[fuzzyres,finalrespor4f]=fuzzy_logic_decetion(mod_data_in,snrrange,freq_val,time_range,samp_freq);
fuzzyresx=1-fuzzyres;
[finalrespro,sensingtimemin1,...
    finalrespro2,sensingtimemin2,finalrespro3,finalrespor4s]=sensing_mini_algor_process(mod_data_in,snrrange,freq_val,time_range,samp_freq);
cyclo_stat_resx=1-cyclo_stat_res(:,2);cyclo_stat_resxm=cyclo_stat_res(:,2).*(1-enegy_propval);
eng_resx1=8*sort(match_filt_resx);locm=find(eng_resx1>-1);eng_resx1(locm)=1;
%% display all the results


figure,plot(snrrange,enegy_propval,'r-.s')
hold on,plot(snrrange,match_filt_resx,'g-^')
hold on,plot(snrrange,cyclo_stat_resx,'k:<')
axis([min(snrrange) max(snrrange) -0.2 1.2]);
legend('Energy Detector','Matched Filter','Cyclostationary Feature Detector');
xlabel('SNR(dB)');
ylabel('Detection Probability of Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Absent');
pause(2);
%%

figure,plot(snrrange,enegy_propval,'r-.s')
hold on,plot(snrrange,match_filt_resx,'g-^')
hold on,plot(snrrange,cyclo_stat_resx,'k:<')
hold on,plot(snrrange,fuzzyres,'m-.o')
axis([min(snrrange) max(snrrange) -0.2 1.2]);
legend('Energy Detector','Matched Filter',...
    'Cyclostationary Feature Detector','Fuzzy Logic Based Detection');
xlabel('SNR(dB)');
ylabel('Detection Probability of Presence of Primary User');
title('Detection Using Minimizing Sensing Time Algorithm,Primary User Absent');

pause(2);


figure,plot(snrrange,eng_resx1,'r-.s')
hold on,plot(snrrange,match_filt_resx1,'g-^')
hold on,plot(snrrange,cyclo_stat_resxm,'k:<')
axis([min(snrrange) max(snrrange) -0.2 1.2]);
legend('Energy Detector','Matched Filter','Cyclostationary Feature Detector');
xlabel('SNR(dB)');
ylabel('Detection Probability of Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Present');
pause(2);

figure,plot(snrrange,eng_resx1,'r-.s')
hold on,plot(snrrange,match_filt_resx1,'g-^')
hold on,plot(snrrange,cyclo_stat_resxm,'k:<')
hold on,plot(snrrange,fuzzyresx,'m-.o')
axis([min(snrrange) max(snrrange) -0.2 1.2]);
legend('Energy Detector','Matched filter',...
    'Cyclostationary Feature detector','Fuzzy logic based detection');
xlabel('SNR(dB)');
ylabel('Detection Probability of Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Present');

pause(2);

figure,plot(snrrange,round(fuzzyres),'m-.o')
hold on,plot(snrrange,finalrespro3,'r-<');
axis([min(snrrange) max(snrrange) -1 1.5]);
legend('Fuzzy Logic Based detection','Minimizing Sensing Time Algorithm');
xlabel('SNR(dB)');
ylabel('Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Absent');

pause(2);
figure,plot(snrrange,finalrespor4f,'m-.o')
hold on,plot(snrrange,finalrespor4s,'r-<');
axis([min(snrrange) max(snrrange) -1 1.5]);
legend('Fuzzy Logic Based detection','Minimizing Sensing Time Algorithm');
xlabel('SNR(dB)');
ylabel('Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Present')
pause(2);
%%
figure,plot(snrrange,cs1,'r-.s')
hold on,plot(snrrange,cs2,'g-^')
hold on,plot(snrrange,cs3,'k:<')
hold on,plot(snrrange,cs4,'m-.o')
hold on,plot(snrrange,round(fuzzyres),'c-.o')
axis([min(snrrange) max(snrrange) -1 1.5]);
legend('Energy Detector','Matched Filter','Cyclostationary Feature Detector','Minimizing Sensing Time Algorithm',...
    'Location','SouthEast');
xlabel('SNR(dB)');
xlabel('SNR(dB)');
ylabel('Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Absent');
pause(2);
%%
figure,plot(snrrange,css1,'r-.s')
hold on,plot(snrrange,css2,'g-^')
hold on,plot(snrrange,css3,'k:<')
hold on,plot(snrrange,css4,'m-.o')
hold on,plot(snrrange,finalrespor4f,'c-.>');
axis([min(snrrange) max(snrrange) -1 1.5]);
legend('Energy Detector','Matched Filter','Cyclostationary Feature Detector','Minimizing Sensing Time Algorithm',...
    'Location','SouthEast');
xlabel('SNR(dB)');
xlabel('SNR(dB)');
ylabel('Presence of Primary User');
title('Detection using Minimizing Sensing Time Algorithm,Primary User Present');

%%