Matlab图像退化与复原

退化：

% checkerboard产生测试板图像，第一个参数是每个正方形一边的像素数，第二个参数行数，第三为列数（缺省则等于行数）
f = checkerboard(8);    % 产生一个一面为8个正方形的测试板
PSF = fspecial('motion',7,45);      % 运动模糊，PSF刚好为空间滤波器
gb = imfilter(f,PSF,'circular');   % 减少边界效应
noise = imnoise(zeros(size(f)),'gaussian',0,0.001);     % 高斯噪声
g = gb + noise;     % 添加高斯噪声构造退化的图像模型
subplot(2,2,1);imshow(pixeldup(f,8),[ ]);title('原图像'); % 大图像运算过慢，故选用小图像来节省时间，   
subplot(2,2,2);imshow(gb);title('运动模糊图像');    % 以显示为目的，可通过像素赋值来放大图像。      
subplot(2,2,3);imshow(noise,[ ]);title('高斯噪声图像');
subplot(2,2,4);imshow(g);title('运动模糊+高斯噪声');

复原：

f = checkerboard(8);                                     % 产生一个一面为8个正方形的测试板
PSF = fspecial('motion',7,45);            % 运动模糊，PSF刚好为空间滤波器
gb = imfilter(f,PSF,'circular');     % 减少边界效应
noise = imnoise(zeros(size(f)),'gaussian',0,0.001);      % 高斯噪声
g = gb + noise;                    % 添加高斯噪声构造退化的图像模型
fr1 = deconvwnr(g,PSF);          % 直接逆滤波
Sn = abs(fft(noise)).^2;               % 噪声功率谱
nA = sum(Sn(:))/prod(size(noise));       % 平均噪声功率，prod计算数组元素的连乘积。
Sf = abs(fft2(f)).^2;                    % 图像功率谱
fA = sum(Sf(:))/prod(size(f));       % 平均图像功率
R = nA/fA;                     % 噪信功率比
fr2 = deconvwnr(g,PSF,R);       % 参数维纳滤波器
NCORR = fftshift(real(ifft2(Sn)));      % 自相关函数
ICORR = fftshift(real(ifft2(Sf)));  
fr3 = deconvwnr(g,PSF,NCORR,ICORR);
subplot(2,2,1);imshow(g,[]);title('噪声图像'); 
subplot(2,2,2);imshow(fr1);title('逆滤波结果');          
subplot(2,2,3);imshow(fr2);title('使用常数比率的维纳滤波的结果');
subplot(2,2,4);imshow(fr3);title('使用自相关函数的维纳滤波的结果');

 

正则滤波：

f = checkerboard(8);     % 产生一个一面为8个正方形的测试板
PSF = fspecial('motion',7,45);       % 运动模糊，PSF刚好为空间滤波器
gb = imfilter(f,PSF,'circular');                         % 减少边界效应
noise = imnoise(zeros(size(f)),'gaussian',0,0.001);      % 高斯噪声
g = gb + noise;               % 添加高斯噪声构造退化的图像模型
fr1 = deconvreg(g,PSF,4);            % 正则滤波器
fr2 = deconvreg(g,PSF,0.4,[1e-7 1e7]);
subplot(1,3,1);imshow(g);title('运动模糊+高斯噪声(g)');
subplot(1,3,2);imshow(fr1);title('正则滤波器');
subplot(1,3,3);imshow(fr2);title('正则滤波器(缩小范围)');