求高手指教！！！压缩感知重构算法GPSR，分块实现，无法重建。

光脚丫丫 2013-05-15 04:17:59

程序如下：
clc
clear
totalframes = seq_frames('foreman_cif.yuv',[352 288],'YUV420_8');

frame1stimage = save_yuvframe('foreman_cif.yuv',[352 288],1,'foreman_1stframe.jpg'); %提取第一帧,保存图片的格式可以为.jpg .bmp .png
%saveY(;,;,frame) = save_eframes( 'foreman_cif.yuv',[352 288],totalframes)

% 读文件
t0=cputime;
%im=uint8(im);
%seg = frame_block(im,8,1); %可改变，块大小8*8
X=imread('foreman_1stframe.jpg');
if ndims(X)>2
X=rgb2gray(X); %转化成灰度图像
end
X=imresize(X,[256 256]); %大小调整为256*256
X=double(X);
[a,b]=size(X);
%size_block=8*8;
n=8*8;
X2=zeros(n);
X3=zeros(a,b); % 恢复矩阵

% 小波变换矩阵生成
ww=DWT(n);
% 随机矩阵生成
M=48;
%M=4;
R=randn(M,n);

R = orth(R')';

hR = @(x1) (R)*x1; %定义转置操作
hRt = @(x1) (R)'*x1;

tic
for i_x=1:ceil(a/n)
for i_y=1:ceil(b/n)
XX=X((i_x-1)*n+1:i_x*n,(i_y-1)*n+1:i_y*n);
% 小波变换让图像稀疏化（注意该步骤会耗费时间，但是会增大稀疏度）
X1=ww*sparse(XX)*ww';
X1=full(X1);
% 测量
%Y=R*X1;
%x1=X1(:,i_y);
y = hR(x1);

% GPSR 算法
%for i=1:n % 列循环

% regularization parameter
tau =0.01*max(abs(R'*y)); %tau[0.01-0.02] 与X同维数大小
debias = 0;
stopCri = 3;
tolA = 0.01; %默认参数值的设置，调用函数式用
[x_BB_notmono,x_debias_BB_notmono,obj_BB_notmono,...
times_BB_notmono,debias_start_BB_notmono,mse]= ...
GPSR_BB(y,hR,tau,...
'Debias',debias,...
'AT',hRt,...
'Monotone',0,...
'Initialization',0,...
'StopCriterion',stopCri,...
'ToleranceA',tolA,...
'ToleranceD',0.0001);
t_BB_notmono = times_BB_notmono(end);
X2(:,i_y)=x_BB_notmono;
% end

X3((i_x-1)*n+1:i_x*n,(i_y-1)*n+1:i_y*n)=ww'*sparse(X2)*ww; % 小波反变换
end

end
X3=full(X3);
%use_time=toc
%errorx=sum(sum(abs(X3-X).^2)); % MSE误差
%psnr(I)=10*log10(255*255/(errorx/a/b)) % PSNR
%save GPSR_50_0.5_PSNR.mat psnr time;
%psnr=s+psnr(I);
t1=cputime-t0
errorx=sum(sum(abs(X3-X).^2)) % MSE误差
psnr=10*log10(255*255/(errorx/a/b)) % PSNR

%load GPSR_50_0.5_PSNR.mat

%figure(1)
%plot(time,psnr)
%title('抽样比例为0.5的psnr')
%grid on
%xlabel('times')
%ylabel('psnr')

% 原始图像
figure(1);
imshow(uint8(X));
title('原始图像');
% 压缩传感恢复的图像
figure(2);
imshow(uint8(X3));
title('恢复的图像');

能跑通但是出来的图形，不是重建图形。
如下：
原图