Description of bpdq_demo

0001 % bpdq_demo_1d : Demo for BPDQ reconstruction in 1-d
0002 %
0003 % This demo first creates a random, sparse 1-d signal, creates a gaussian
0004 % random sensing matrix, and computes the quantized measurements.
0005 %
0006 % Then the BPDQ algorithm is run to reconstruct the signal from the quantized
0007 % measurements.
0008 %
0009 % For Comparison, BPDN (corresponding to p=2) is also run, then the results
0010 % from both reconstruction programs are compared.
0011 %
0012 % Plots of relative errors ( ||x_n - x_(n-1)||/||x_n|| )
0013 % as a function of iteration number are given to demonstrate empirical
0014 % convergence of the douglas-rachford algorithm used to compute the
0015 % BPDQ program.
0016 %
0017 % This demonstrates the method described in
0018 % "Dequantizing Compressed Sensing : When Oversampling and Non-Gaussian
0019 % Constraints Combine"
0020 % by Laurent Jacques, David Hammond, M. Jalal Fadili
0021 % ( submitted to IEEE Transactions on Signal Processing)
0022 %
0023 % This file is part of BPDQ Toolbox (Basis Pursuit DeQuantizer)
0024 % Copyright (C) 2009, the BPDQ Team (see the file AUTHORS distributed with
0025 % this library) (See the notice at the end of the file.)
0026 
0027 
0028 
0029 % With parameters as set in this demo, execution time was 18 seconds
0030 % on 2.4 GHz Intel Core 2 Duo processor.
0031 
0032 
0033 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0034 % Demo Parameters :
0035 % N - signal dimension
0036 % K - signal sparsity
0037 % M - # of measurements
0038 % p - BPDQ moment (p=2 corresponds to standard BPDN)
0039 % nbins - desired # of quantization bins
0040 % (Quantization threshold alpha is computed depending on nbins and on
0041 %  the signal instance)
0042 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0043 N=1024;
0044 K=16;
0045 M=K*20;
0046 p=10;
0047 nbins=10;
0048 
0049 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0050 % Parameters for controlling bpdq algorithm.
0051 % The Douglas-Rachford algorithm parameters
0052 % dr_lambda must be between 0 and 2, dr_gamma must be positive.
0053 % dr_gamma and dr_lambda are described in Jacques et. al. section
0054 % IV-A; dr_lambda corresponds to alpha_t in equation 10
0055 % Note : a priori these parameters could vary at each iteration, in this code they are
0056 % kept constant.
0057 %
0058 % dr_verbosity and fca_verbosity control output of information about
0059 % speed of convergence at each iteration, set to 0 to supress output, to 1
0060 % to show output.
0061 % See bpdq_d_r_iter.m and bpdq_prox_f_comp_A.m for details.
0062 %
0063 % Other options
0064 % dr_maxiter, dr_relerr_thresh, fca_maxiter, fca_relerr_thresh
0065 % are defined in in control_params in bpdq_1d.m, see bpdq_1d.m for details.
0066 
0067 bpdq_opts={'dr_gamma',.05,'dr_lambda',1.5,'dr_verbosity',0,'fca_verbosity',0};
0068 
0069 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0070 % Create signal and sensing matrix
0071 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0072 fprintf('Generating signal of length N=%g and sparsity K=%g\n',N,K);
0073 x=bpdq_generate_1d_signal(N,K);
0074 fprintf('Generating sensing matrix for M=%g random measurements\n',M);
0075 A=bpdq_generate_sensing_matrix(M,N);
0076 
0077 % compute measurements and quantize
0078 fprintf('Computing quantized measurements with %g quantization levels\n',nbins);
0079 y=A*x;
0080 alpha = 2.01*norm(y, Inf)/nbins;
0081 yq = bpdq_quantize(y, alpha);
0082 
0083 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0084 % Call BPDQ decoder
0085 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0086 % Compute apropriate epsilon (ref Jacques et. al section III-C)
0087 
0088 epsilon=bpdq_err_p(p,alpha,M);
0089 fprintf('Running BPDQ reconstruction with p=%g\n',p);
0090 [x_bpdq,D_bpdq] = bpdq_1d(yq,A,epsilon,p,bpdq_opts{:});
0091 
0092 % compare to BPDN
0093 epsilon2=bpdq_err_p(2,alpha,M);
0094 fprintf('Running BPDN (p=2) reconstruction\n');
0095 [x_bpdn,D_bpdn]  = bpdq_1d(yq,A,epsilon2,2,bpdq_opts{:});
0096 
0097 % examine results
0098 SNR_bpdq = bpdq_compute_snr(x,x_bpdq);
0099 SNR_bpdn = bpdq_compute_snr(x,x_bpdn);
0100 
0101 fprintf('Displaying results and convergence curves\n');
0102 % plot results and convergence curves
0103 
0104 dmax=1.05*max(abs([x(:);x_bpdq(:);x_bpdn(:)]));
0105 fig1ylim=[-dmax,dmax];
0106 figure(1)
0107 subplot(3,1,1);
0108 tstr=sprintf('original signal, N=%g, K=%g',N,K);
0109 plot(x);
0110 ylim(fig1ylim); 
0111 title(tstr);
0112 
0113 subplot(3,1,2);
0114 tstr=sprintf('bpdq (p=%g) reconstruction, SNR %g dB',p,SNR_bpdq);
0115 plot(x_bpdq);
0116 ylim(fig1ylim); 
0117 title(tstr);
0118 
0119 subplot(3,1,3);
0120 tstr=sprintf('bpdn (p=2) reconstruction, SNR %g dB',SNR_bpdn);
0121 plot(x_bpdn); 
0122 ylim(fig1ylim); 
0123 title(tstr);
0124 
0125 figure(2)
0126 subplot(2,1,1)
0127 plot(log10(D_bpdq.dr_relerr_save));
0128 title('DR convergence (for bpdq)');
0129 xlabel('iteration');
0130 ylabel('log relative error');
0131 
0132 subplot(2,1,2)
0133 plot(log10(D_bpdn.dr_relerr_save));
0134 title('DR convergence (for bpdn)');
0135 xlabel('iteration')
0136 ylabel('log relative error');
0137 
0138 figure(3)
0139 diff_bpdq=x-x_bpdq;
0140 diff_bpdn=x-x_bpdn;
0141 dmax=1.05*max([abs(diff_bpdq);abs(diff_bpdn)]);
0142 fig3ylim=[-dmax,dmax];
0143 
0144 subplot(2,1,1);
0145 tstr=sprintf('Error of bpdq (p=%g) reconstruction, SNR %g',p,SNR_bpdq);
0146 plot(diff_bpdq);
0147 ylim(fig3ylim);
0148 title(tstr);
0149 
0150 subplot(2,1,2);
0151 tstr=sprintf('Error of bpdq (p=2) reconstruction, SNR %g',SNR_bpdn);
0152 plot(diff_bpdn);
0153 ylim(fig3ylim);
0154 title(tstr);
0155 
0156 % The BPDQ Toolbox is free software: you can redistribute it and/or modify
0157 % it under the terms of the GNU General Public License as published by
0158 % the Free Software Foundation, either version 3 of the License, or
0159 % (at your option) any later version.
0160 %
0161 % The BPDQ Toolbox is distributed in the hope that it will be useful,
0162 % but WITHOUT ANY WARRANTY; without even the implied warranty of
0163 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0164 % GNU General Public License for more details.
0165 %
0166 % You should have received a copy of the GNU General Public License
0167 % along with The BPDQ Toolbox.  If not, see <http://www.gnu.org/licenses/>.
bpdq_demo_1d

PURPOSE

SYNOPSIS

DESCRIPTION

CROSS-REFERENCE INFORMATION

SOURCE CODE
