cs2_2sig.m

压缩传感理论仿真的一个实际例子

tic

clear
format short
c=3*10.^8;
L=1;
N=256;
M=80;

fm1=0.35;%信号1起始频率
Km1=0.05;%信号2调频斜率
fm2=0.25;%信号2起始频率
Km2=0.12;%信号2调频斜率
f0=fm1+Km1/2;
lamta=c/f0;
len=lamta/2;
doature=10;
P1=doature*pi/180;


t=1:N;

%%%%%%%%%%%%%%%%%%
%CS的方法

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%构造原子库
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
atom=zeros(51*201,N);
for l=1:51
    fa=(l-1)*0.01;%0~0.5
    for k=1:201
        Ka=k*0.001;    %-0.05~0.05
         atom((l-1)*201+k,:)=1.*exp(j.*2*pi*(Ka/2*t.^2/N+fa*t));
    end
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%设置信噪比个数及每个信噪比试验次数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

time=1;
A=(1/sqrt(M))*randn(M,N);       %构造观测矩阵

sig1=exp(j*2*pi*(Km1/2*(t.^2)/N+fm1*t));
sig2=exp(j*2*pi*(Km2/2*(t.^2)/N+fm2*t)); 

for ii=1:6;
  snr=(ii-1)*4-10;

 Amp=sqrt(2*10^(snr/10));


%生成信号

  

%控制循环次数 
 for times=1:time
     
 noise=randn(1,N)+j*randn(1,N);
z=Amp*sig1+Amp*sig2+noise; 

%求出观测值y
y1=A*z';


%解优化问题

cvx_begin
 variable x(51*201)
 minimize(sum_square(real(y1)-A*real(atom')*x)+0.01*norm(x,1))
cvx_end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%求解起始频率和调频斜率
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[m,n]=max(abs(x));
e_fm1=floor(n/201)/100;
e_Km1=(n-floor(n/201)*201)/1000;
x(n)=0;
[m,n]=max(abs(x));
e_fm2=floor(n/201)/100;
e_Km2=(n-floor(n/201)*201)/1000;
if e_fm1>e_fm2
ef1(times)=e_fm1;
eK1(times)=e_Km1;
ef2(times)=e_fm2;
eK2(times)=e_Km2;
else
ef1(times)=e_fm2;
eK1(times)=e_Km2;
ef2(times)=e_fm1;
eK2(times)=e_Km1;
end
end
  
f1_mean(ii)=mean(ef1);
K1_mean(ii)=mean(eK1);
f2_mean(ii)=mean(ef2);
K2_mean(ii)=mean(eK2);

f1_resm(ii)=rmse(ef1,fm1);
K1_resm(ii)=rmse(eK1,Km1);
f2_resm(ii)=rmse(ef2,fm2);
K2_resm(ii)=rmse(eK2,Km2);

end







%%%%%%%%%%%%%%%%%显示结果
 i=1:6;
 snr=(i-1)*4-10;

figure(1)
plot(snr,f1_resm,'r-o',snr,f2_resm,'b-*');
xlabel('SNR/dB','fontsize',12,'fontweight','b');
ylabel('RMSE','fontsize',12,'fontweight','b');
title('起始频率');
legend('信号1-起始频率','信号2-起始频率')

figure(2)
plot(snr,K1_resm,'r-o',snr,K2_resm,'b-*');
xlabel('SNR/dB','fontsize',12,'fontweight','b');
ylabel('RMSE','fontsize',12,'fontweight','b');
title('调频斜率');

legend('信号1-调频斜率','信号2-调频斜率')

toc