【数字图像处理】Parametric Wiener Filter

实验要求：

(a) 编写一个给图像中添加高斯噪声的程序，程序的输入参数为噪声的均值与方差。
(b) 编写程序实现公式(5.6-11)所示的污损滤波；
(c) 如图 5.26(b)所示，对图像 5.26(a) 进行+45o 方向，T = 1 的污损滤波；
(d) 对污损后的图像加入均值为０，方差为 10 的高斯噪声；
(d) 编写程序使用公式(5.8-6)

技术论述：

实验代码：

% ----------------------Parametric Wiener Filter---------------------------%
function main_test
clc;clear all;close all;
 
img = imread('Fig5.26(a).jpg');                             % 原图像
img_fouling = fouling_filter(img,0.1,0.1,1);                % 污损滤波
img_fouling1 = fouling_filter(img,0.01,0.01,1);             % 污损滤波
img_gaussian_10 = gaussian_noise(img_fouling,0,10);         % 添加高斯噪声
img_gaussian_10_1 = gaussian_noise(img_fouling1,0,10);      % 添加高斯噪声
 
figure;
subplot(2,3,1);imshow(img);title('图1(a)：原图像');
subplot(2,3,2);imshow(img_fouling);title('图1(b)：污损滤波: +45°,T=1，a=b=0.1');
subplot(2,3,3);imshow(img_fouling1);title('图1(c)：污损滤波: +45°,T=1，a=b=0.01');
subplot(2,3,4);imshow(img_gaussian_10);title('图1(d)：对图(b)加入均值0,方差10的高斯噪声');
subplot(2,3,5);imshow(img_gaussian_10_1);title('图1(e)：对图(c)加入均值0,方差10的高斯噪声');
 
%% -------------------------维纳滤波函数-------------------------%
img_wiener = wiener_filter(img,img_gaussian_10,0.1);        % 维纳滤波
img_wiener1 = wiener_filter(img,img_gaussian_10_1,0.1);        % 维纳滤波
figure;
subplot(1,2,1);imshow(img_wiener);title('图2(a)：对图1(d)进行维纳滤波图像');
subplot(1,2,2);imshow(img_wiener1);title('图2(b)：对图1(e)进行维纳滤波图像');
 
%% ------------------------MATLAB自带函数-----------------------%
psf = fspecial('motion',25,45);
mf = imfilter(img,psf,'circular','conv');
noise = imnoise(img,'gaussian',0,10);
mfn = mf + noise;
nsr = sum(noise(:).^2)/sum(mf(:).^2);
img_out = deconvwnr(mfn,psf,nsr);
 
figure;
subplot(1,3,1);imshow(mf);title('图3(a)：Matlab自带污损滤波');
subplot(1,3,2);imshow(mfn);title('图3(b)：Matlab添加高斯噪声');
subplot(1,3,3);imshow(img_out);title('图3(c)：Matlab自带维纳滤波');
% ---------------------------------END---------------------------------%

(a) 编写一个给图像中添加高斯噪声的程序，程序的输入参数为噪声的均值与方差。 
% --------------------------给图像添加高斯噪声---------------------------%
function img_noise = gaussian_noise(img,mean,var) 
% 程序的输入为噪声的均值mean和方差var，输出为噪声图像
 
[M,N] = size(img);
add_noise = mean + randn(M,N)*sqrt(var);
img_noise = img + uint8(add_noise);
end

(b) 编写程序实现公式(5.6-11)所示的污损滤波； 
% ------------------------污损滤波函数程序------------------------%
function img_fouling = fouling_filter(img,a,b,T) 
% 输入参数：img为输入图像，a,b,T为别为技术论述中的公式(9)的三个参数
 
[M,N] = size(img);
F = fft2(img);
 
for u = 1:M
    for v = 1:N
        K(u,v) = pi * (u * a + v * b);
        H(u,v) = T * sin(K(u,v)) * exp(-1j * K(u,v))/K(u,v);
        G(u,v) = H(u,v) * F(u,v);
    end
end
 
img_fouling = ifft2(G);
img_fouling = uint8(abs(img_fouling));
end

(e) 编写程序使用公式(5.8-6)所示的参数维纳滤波对图像进行恢复。 
% ----------------------------维纳滤波程序-----------------------------%
function img_wiener = wiener_filter(img_src,img_degradation,K) 
 
[M,N] = size(img_src);
S = fft2(img_src);
G = fft2(img_degradation);
 
for u = 1:M
    for v = 1:N
        H(u,v) = G(u,v)/S(u,v);
        F(u,v) = 1/H(u,v)*(abs(H(u,v)))^2/((abs(H(u,v)))^2+K)*G(u,v);
    end
end
 
img_wiener = ifft2((F));
img_wiener = uint8(abs(img_wiener));
end