www.gusucode.com > 压缩感知重构算法 压缩感知源码程序 > shijun.m

    clc;
clear all;
%n=2^10;
%t = (1:n) ./n;
t=0:1:2047;
%pos = [ .1 .13 .15 .23 .25 .40 .44 .65  .76];
% pos=[50 150 220 280 350 550 620 720 800 850 980];
% hgt=[5.7 -3 -2.6 5.7 -5.8 3.8 -3.65 0.1  4.5  -4.6 0];
pos=[300 500 700 900 1100 1300 1500 1700 1900 2100];
hgt=[5.8 -5.8 3.6 -3.6 5.7 -5.7  5.8 -5.80  0 0];
sig = zeros(size(t));
for j=1:length(pos)
    sig = 0.01+sig + (1 + sign(t-pos(j))).*(hgt(j)/2) ;
end

sig0=sig';
sig0=sig0/max(sig0);
xfs=40;
points=2048;
f = xfs/2*(1:points)/points;
M = 2^11;
N = 750;

%sig0 = MakeBlocks(M)';
qmf = MakeONFilter('Haar',1);

SNR =15;
%randn('seed', 2.055615866000000e+09);
%[sig0 y]= NoiseMaker(sig0, SNR);
% y=awgn(sig0,SNR);
% y=y/max(y);
y=sig0;

figure;
y1=ifft(y,M);
plot([0:length(y1)-1],y1);
grid on;

f1 = xfs/2*(1:points)/(points);
figure;
y2=fft(y1,M);
plot(f1,y2);
axis([0 20 -0.1 1.1]);
grid on;


L = 3;
alpha0 = FWT_PO(y, L, qmf);

figure;
plot(f,alpha0);
grid on;

% Generate Random Dictionary
Phi = MatrixEnsemble(N,M);

% generate the vector S
S = Phi * alpha0;

% Solve the BP problem
alpha = SolveBP(Phi, S, M, 20);
sig = IWT_PO(alpha,L,qmf);

% Solve the BPDN problem
alpha2 = SolveOMP(Phi, S, M, 20);
sig2 = IWT_PO(alpha2,L,qmf);

figure;
subplot(2,1,1);
plot(f,y); 
axis([0 20 -0.1 1.1]);
% axis([0 M -1 7]); 
title('(a) Noisy Signal,SNR=15');
grid on;
%subplot(3,2,2); 
%PlotWaveCoeff(alpha0, L, 0);  
%title('(b) Noisy wavelet coefficients');

subplot(2,1,2);
plot(f,sig,'r'); 
axis([0 20 -0.1 1.1]);
% axis([0 M -1 7]); 
title('(b) CS Reconstruction by BP, n = 512');
grid on;
%subplot(3,2,4); 
%PlotWaveCoeff(alpha, L, 0);  
%title('(d) CS Reconstruction');

%subplot(3,1,3);
%plot(sig2); 
%axis([0 M -1 7]); 
%title('(e) CSDN Reconstruction, n = 512');
%grid on;
%subplot(3,2,6); 
%PlotWaveCoeff(alpha2, L, 0);  
%title('(f) CSDN Reconstruction');