分塊壓縮感知程序matlab

1、%本程序?qū)崿F(xiàn)圖像LENA的壓縮傳感%程序作者：沙威，香港大學(xué)電氣電子工程學(xué)系，wshaeee.hku.hk%算法采用正交匹配法，參考文獻(xiàn)Joel A. Tropp and Anna C. Gilbert% Signal Recovery From Random Measurements Via Orthogonal Matching% Pursuit, IEEE TRANSACTIONS ON INFORMATION THEORY , VOL. 53, NO. 12,% DECEMBER 2007.%該程序沒有經(jīng)過任何優(yōu)化%function Wavelet_OMPclcclear%讀文件X=i

2、mread(lena256.bmp);X=double(X);a,b=size(X);size_kuai=16*4;X2=zeros(size_kuai); %恢復(fù)矩陣X3=zeros(a,b); %恢復(fù)矩陣%小波變換矩陣生成ww=DWT(size_kuai);%隨機(jī)矩陣生成M=12*4;R=randn(M,size_kuai);ticfor i_x=1:ceil(a/size_kuai)for i_y=1:ceil(b/size_kuai)XX=X(i_x-1)*size_kuai+1:i_x*size_kuai,(i_y-1)*size_kuai+1:i_y*size_kuai);%小波變

3、換讓圖像稀疏化(注意該步驟會耗費(fèi)時間，但是會增大稀疏度)X1=ww*sparse(XX)*ww;X1=full(X1);%測量Y=R*X1;% OMP算法for i=1:size_kuai %列循環(huán)rec=omp_fenkuai(Y(:,i),R,size_kuai);X2(:,i)=rec;endX3(i_x-1)*size_kuai+1:i_x*size_kuai,(i_y-1)*size_kuai+1:i_y*size_kuai)=ww*sparse(X2)*ww;%小波反變換end endX3=full(X3);use_time=toc%原始圖像figure(1);imshow(uin

4、t8(X);title(原始圖像);%壓縮傳感恢復(fù)的圖像figure(2);imshow(uint8(X3);title(分塊恢復(fù)的圖像);%誤差(PSNR) errorx=sum(sum(abs(X3-X).A2); psnr=10*log10(255*255/(errorx/a/b)endhat_y(pos_array)=aug_y;% MSE誤差% PSNR% OMP的函數(shù)% s-測量；T-觀測矩陣；N-向量大小function hat_y=omp_fenkuai(s,T,N)Size=size(T);M=Size(1);hat_y=zeros(1,N);Aug_t=;r_n=s;for

5、 times=1:M/4for col=1:Nproduct(col)=abs(T(:,col)*r_n);數(shù)(內(nèi)積值)endval,pos=max(product);Aug_t=Aug_t,T(:,pos);T(:,pos)=zeros(M,1);為了簡單我把它置零)aug_y=(Aug_t*Aug_t)A(-1)*Aug_t*s;r_n=s-Aug_t*aug_y;pos_array(times)=pos;if (norm(r_n)40)break;end%觀測矩陣大小%測量% 待重構(gòu)的譜域（變換域）向量%增量矩陣（初始值為空矩陣）%殘差值%迭代次數(shù)（稀疏度是測量的1/4）%恢復(fù)矩陣的所有

6、列向量%恢復(fù)矩陣的列向量和殘差的投影系%最大投影系數(shù)對應(yīng)的位置%矩陣擴(kuò)充%選中的列置零（實(shí)質(zhì)上應(yīng)該去掉,%最小二乘，使殘差最小%殘差%紀(jì)錄最大投影系數(shù)的位置%殘差足夠小%重構(gòu)的向量%程序作者：沙威，香港大學(xué)電氣電子工程學(xué)系，wshaeee.hku.hk%參考文獻(xiàn)：小波分析理論與MA TLAB R2007實(shí)現(xiàn)，葛哲學(xué)，沙威，第20章 小波變換在矩陣方程求解中的應(yīng)用（沙威、陳明生編寫）%構(gòu)造正交小波變換矩陣，圖像大小N*N , N=2AP , P是整數(shù)。function ww=DWT(N)%矩陣構(gòu)造for jj=rank_min:rank_maxnn=2Ajj；%構(gòu)造向量p1_0=sparse(h

7、,zeros(1,nn-L);p2_0=sparse(g,zeros(1,nn-L);%向量圓周移位for ii=1:nn/2p1(ii,:)=circshift(p1_0,2*(ii-1);p2(ii,:)=circshift(p2_0,2*(ii-1);end%構(gòu)造正交矩陣w1=p1;p2;mm=2Arank_max-length(w1);w=sparse(w1,zeros(length(w1),mm);zeros(mm,length(w1),eye(mm,mm);ww=ww*w;clear p1;clear p2;endh,g= wfilters(sym8,d);%分解低通和高通濾波器%