clc;
clear;
tic;
im=phantom（‘Modified?Shepp-Logan‘256）;
G=abs（fftshift（fft2（im）））;
figure;
imshow（G[]）;
colorbar
title（‘二維Fourier變換頻譜圖‘）;
im=im2uint8（mat2gray（im））;
figure;
imshow（im[]）;
title（‘256*256原圖像‘）;

[m?n]=size（im）;
im2=zeros（mn）;
bm=12;
bn=12;	%分塊大小
p=0.7;	%采樣率
d=bm*bn;	%信號長度
N=round（d*p）;	%測量值個數
bx=ceil（m/bm）;
by=ceil（n/bn）;	%分塊數
imc=zeros（bx*bmby*bn）;
for?i=1:bx*bm	%分塊沒有覆蓋到的地方補0
	for?j=1:by*bn
		if?i<=m&&j<=n
			imc（ij）=im（ij）;
		else?imc（ij）=0;
		end
	end
end
for?x=1:bx	%對每個小塊處理
	for?y=1:by
		b=imc（（1+（x-1）*bm）:（x*bm）（1+（y-1）*bn）:（y*bn））;
		xin=reshape（bd1）;
		xin=double（xin）;	%輸入信號
		A=dctmtx（d）;		%DCT矩陣
		Phi=randn（Nd）;	%高斯隨機矩陣
		T=Phi*inv（A）;		%傳感矩陣
		s=Phi*xin;			%測量值
		hat_y=zeros（1d）;	%OMP算法
		Aug_t=[];
		aug_y=[];
		pos_array=[];
		r_n=s;
		times=1;
		while?norm（r_n）>0.01
			for?col=1:d;
				product（col）=abs（T（:col）‘*r_n）;
			end
			[valpos]=max（product）;
			Aug_t=[Aug_tT（:pos）];
			T（:pos）=zeros（N1）;
			aug_y=pinv（Aug_t‘*Aug_t）*Aug_t‘*s;
			r_n=s-Aug_t*aug_y;
			pos_array（times）=pos;
			times=times+1;
		end
		hat_y（pos_array）=aug_y;
		hat_x=hat_y;
		hat_x2=inv（A）*hat_x‘;	%復原信號
		b2=reshape（hat_x2bmbn）;
		%將計算好的分塊組合
		imc2（（1+（x-1）*bm）:（x*bm）（1+（y-1）*bn）:（y*bn））=b2;
	end
end
imc2=uint8（imc2）;
im2=imc2（1:m1:n）;
mse=sum（sum（abs（im-im2）.^2））/（m*n）;	%計算均方誤差
toc;		%計算時間
figure;	%顯示圖像
imshow（im2[0?255]）;
title（strcat（‘OMP重構結果:‘‘采樣率=‘num2str（p）‘分塊=‘...
num2str（bm）‘בnum2str（bn）‘計算時間=‘...
num2str（round（toc））‘s?MSE=‘?...
num2str（round（mse*10）/10）））;
colormap?gray;

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----

?????文件???????1730??2014-12-11?19:12??BCS.m

-----------?---------??----------?-----??----

?????????????????1730????????????????????1