mpcosin.m

利用系数分解算法从噪声中检测出微弱信号

% MATLAB script for estimation the parameter of the cosinoidal signal by
% Matching Pursuit(MP)
%
%
% *******************************Preparation part*******************
%
% *******************************set parameter**********************
f=9.251;                       % freqency of the signal
fai=0.512*pi;                  % phase of the signal
snr=-20;                       % SNR in dB(logarithmic)
t=0:0.001:32;                  % time vector
s=cos(2*pi*f*t+fai);           % signal with the Amplify of 1
%
% *******************************nosie generation**********************
% 
tr=randn(size(s));
Amplifytr=max(abs(tr));
ampsig=Amplifytr*10^(snr/20);
s=ampsig*s;
%
% *******************************signal with nosie generation***********
%
s_o=s+tr;
%
% *******************************dictionary generation*****************
%
error=zeros(101,201);
for fatom=0:1000                          % freqency of the atom
    for faiatom1=0:1:200
        faiatom=faiatom1*0.01*pi;         % phase of the atom
        fatom1=fatom*0.01;
        s1=cos(2*pi*fatom1*t+faiatom);     % generation the atom
        s_s1=sqrt(sum(s1.*s1));
        s1=s1/s_s1;                       % normalizes to the norm of the atom
% ** calculate the decomposition cofficients of all the atoms and the
% signal**
        error(fatom+1,faiatom1+1)=sum(s_o.*s1);
    end
end



yabs=abs(error);            % save all the expansion coefficients
[m,n]=max(yabs);            % choose the best atom
[k,l]=max(m);
result_f=(n(l)-1)/100;      % save the freqency parameter of the best atom
result_fai=(l-1)*0.01;
result_a=error(n(l),1);


if result_a<0
    result_a=-result_a;
    result_fai=result_fai-1;
end
result_a=result_a/s_s1
result_fai;
%
% ******************************* renew the cosinoidal signal ************
%
sr=result_a*cos(2*pi*result_f*t+result_fai*pi);
%
% ******************************* the plot of the results ***************
%
tt=0:1:1000;
subplot(2,2,1);
plot(tt,s(1:1001));    % original signal 
subplot(2,2,2);
plot(tt,s_o(1:1001));    % signal with nosie
subplot(2,2,3);
plot(tt,sr(1:1001));      % renewed signal
subplot(2,2,4); 
plot(tt,yabs(1:1001));        % decomposed cofficients of freqency
% ******************************* end of file ****************************