bpfig24.m

％ Atomizer Main Directory, Version .802 里面信号含有分解去噪合成过程的代码 ％---------------------------------------

% bpfig24: BP Figure 2.4 -- Super-resolving two close frequencies
%		within the Rayleigh distance
%----------------------------------------------------------------------
%
%       BP is able to achieve super-resolution, while
%	MOF and MP sometimes cannot
%
% Signal:        TwinSine (two frequencies not in the dictionary)
%                cos(pi * w1 * t) + cos(pi * w2 * t)
% Signal Length: 256
% Dictionary:    factor 4-overcomplete discrete cosine dictionary
% Observations:
%    (a) Signal: TwinSine.
%    (b) MOF results in a Dirichlet-like oscillatory structure
%        consisting of not two but many frequencies.
%    (c) MP fails to superresolve. MP initially selects the frequency
%        between the two frequencies making up the signal. Because of
%        this mistake, MP is forced to make a series of alternating
%        corrections that suggest a highly complex and organized
%        structure.
%    (d) BP result resolves the doublet structure.
%
% Use:
%    bpfig24                   uses the current solver.
%    ATOMIZER_ENGINE = 1998;   selects original solver.
%    ATOMIZER_ENGINE = 2001;   selects later    solver.
%    Default is most recent solver.
%----------------------------------------------------------------------

%----------------------------------------------------------------------
%        1998: (S. Chen) Original script for BP paper.
% 09 Apr 2001: (M. Saunders) Choice of solvers implemented.
%----------------------------------------------------------------------

help bpfig24
bpengine;

n        = 256;
fineness =   4;
m        = n * fineness;
t        = ((1:n)' - .5) / n;
const    = sqrt(2/n);

x1 = const * cos(pi * (125.2  - 1) / fineness * t);
x2 = const * cos(pi * (127.2  - 1) / fineness * t);
x  = x1 + x2;

%------------------
% MOF
%------------------
time  = cputime;
cMOF  = MOF(x, 'DCT', fineness, 0, 0);
time1 = cputime - time;

%------------------
% MP
%------------------
natom = n;
frac  = 1e-2;
time  = cputime;
cMP   = MP(x, 'DCT', fineness, 0, 0, n, frac);
time3 = cputime - time;

%------------------
% BP
%------------------
time  = cputime;
switch ATOMIZER_ENGINE
  case 1998
    cBP = BP_Interior ( x, 'DCT', fineness, 0, 0, 1e-2, 1e-2, 1e-2 );
  case 2001
    cBP = BP_Interior2( x, 'DCT', fineness, 0, 0 );
end
time4 = cputime - time;

fprintf('\n')
fprintf('CPU Running Time of MOF = %8.4e\n', time1);
fprintf('CPU Running Time of  MP = %8.4e\n', time3);
fprintf('CPU Running Time of  BP = %8.4e\n', time4);

%------------------
%True frequencies
%------------------
P11 = [(125.2-1)/(fineness * n) (125.2-1)/(fineness * n)];
P12 = [-.5 2];
P21 = [(127.2-1)/(fineness * n) (127.2-1)/(fineness * n)];
P22 = P12;

%--------------------
% Plots
%--------------------
fprintf('\nConstructing figure(%1g) ...\n', FIGURE)
figure(FIGURE);   clf reset;

subplot(2,2,1);   plot(t, x);
                  title('(a) Signal: TwinSine')

subplot(2,2,2);   PlotCoeff(cMOF, n, 'DCT', fineness, 0, 0);
hold on;
axis([(101-1)/(fineness*n) (150-1)/(fineness*n) -.5 2])
plot(P11, P12, ':')
plot(P21, P22, ':')
hold off
xlabel('Frequency/Nyquist')
ylabel('Amplitude')
                  title('(b) MOF Coefs')

subplot(2,2,3);   PlotCoeff(cMP, n, 'DCT', fineness, 0, 0);
                  title('(c) MP Coefs');
hold on;
axis([(101-1)/(fineness*n) (150-1)/(fineness*n) -.5 2])
plot(P11, P12, ':')
plot(P21, P22, ':')
hold off
xlabel('Frequency/Nyquist')
ylabel('Amplitude')

subplot(2,2,4);   PlotCoeff(cBP, n, 'DCT', fineness, 0, 0);
                  title('(d) BP Coefs');
hold on;
axis([(101-1)/(fineness*n) (150-1)/(fineness*n) -.5 2])
plot(P11, P12, ':')
plot(P21, P22, ':')
hold off
xlabel('Frequency/Nyquist')
ylabel('Amplitude')