维纳滤波器(Wiener Filter)在图像处理中的应用（一）

　　维纳滤波器是一种自适应的滤波器，在数字信号处理中有着广泛的应用。ispforfun会在从今天开始定期给大家带来维纳滤波器在图像处理中应用。本节讲诉维纳滤波器在图像去噪中的简单应用。

       让我们从Matlab中的函数wiener2开始。Matlab的help中对wiener2的说明如下：

       wiener2估计图像中每个像素的局部均值和局部方差，

      其中，是图像A中每个像素的NXM邻域。维纳滤波器的去噪原理如下：

      其中是噪声的方差，如果没有给出来，那么会根据图像的局部方差进行估计。

     

     ispforfun根据这段原理，编写了如下的matlab代码，共享给大家。

function J = wiener_image_denoising_filter(I, Win, noise_var)
% J = wiener_image_denoising_filter(I, w)
% -----------------------------------------
% by: ispforfun
% I is the input image. Single channel image.
% Win is the window size for filtering.

[img_height, img_width] = size(I);
half_win_size = floor(Win/2.0);

[X, Y] = meshgrid(1:img_width, 1:img_height);
[dX, dY] = meshgrid(-half_win_size(2):half_win_size(2), -half_win_size(1):half_win_size(1));

dX = reshape(dX, [1 1 Win(1) Win(2)]);
dY = reshape(dY, [1 1 Win(1) Win(2)]);
X = repmat(X, [1 1 Win(1) Win(2)]) + repmat(dX, [img_height img_width 1 1]);
Y = repmat(Y, [1 1 Win(1) Win(2)]) + repmat(dY, [img_height img_width 1 1]);

X(X<1) = 2 - X(X<1);
X(X>img_width) = 2*img_width - X(X>img_width);
Y(Y<1) = 2 - Y(Y<1);
Y(Y>img_height) = 2*img_height - Y(Y>img_height);

patch = @(f)f(Y + (X-1)*img_height);

P = patch(I);

img_mean = local_image_mean(P);
img_var = local_image_var(P, img_mean, Win);

if nargin < 3
    estimated_noise_var = sum(sum(img_var))/(img_height * img_width);
else
    estimated_noise_var = noise_var;
end

gain = (estimated_noise_var == 0.0).*zeros(img_height, img_width) + ...
    (estimated_noise_var ~= 0.0).*((img_var - estimated_noise_var)./img_var);

J = img_mean + gain.*(double(I) - img_mean);

J = min(max(J, 0), 255);

J = uint8(J);

%gain = (img_var - estimated_noise_var)./img_var;



function img_mean = local_image_mean(P)
% img_mean = local_image_mean(P)

local_img_sum = sum(sum(double(P), 4), 3);
sum_size = size(P,3)*size(P,4);

img_mean = local_img_sum/sum_size;


function img_var = local_image_var(P, img_mean, Win)
% img_var = local_img_var(P, img_mean)

diff_img_patch = double(P) - repmat(img_mean, [1 1 Win(1) Win(2)]);
diff_img_patch_pow2 = diff_img_patch.^2;
local_diff_img_sum = sum(sum(double(diff_img_patch_pow2), 4), 3);
img_var = local_diff_img_sum/(Win(1)*Win(2) - 1);

实验如下，富含噪声的图像如下：

去噪后的图像如下：