图像去模糊（维纳滤波）

在数学应用上，对于运动引起的图像模糊，最简单的方法是直接做逆滤波，但是逆滤波对加性噪声特别敏感，使得回复的图像几乎不可用。最小均方差（维纳）滤波用来去除含有噪声的模糊图像，其目标是找到未污染图像的一个估计，使它们之间的均方差最小，可以去除噪声，同时清晰化模糊图像。

定义

给定一个系统

$$y(t)=h(t)*x(t)+n(t)$$
这里，$*$是卷积符号

• $x(t)$是在时间$t$刻输入的信号（未知）
• $h(t)$是一个线性时间不变系统的脉冲响应（已知）
• $n(t)$是加性噪声，与$x(t)$不相关（未知）
• $y(t)$是我们观察到的信号
  我们的目标是找出这样的卷积函数$g(t)$，这样我们可以如下得到估计的$x(t)$：
  $$\hat x(t) = g(t)*y(t)$$
  这里$\hat x(t)$是$x(t)$的最小均方差估计。
  基于这种误差度量， 滤波器可以在频率域如下描述
  $$\begin{split} G(f) &= {{H^*(f)S(f)} \over {|H(f){|^2}S(f) + N(f)}}\\ &= {{H^*(f)} \over {|H(f){|^2} + N(f)/S(f)}} \end{split}$$
  这里：
• $G(f)$和$H(f)$是$g$和$h$在频率域$f$的傅里叶变换。
• $S(f)$是输入信号$x(t)$的功率谱。
• $N(f)$是噪声的$n(t)$的功率谱。
• 上标$*$代表复数共轭。
  滤波过程可以在频率域完成：
  $$\hat X(f) = G(f)*Y(f)$$
  这里 $\hat X(f)$是 $\hat x(t)$的傅里叶变换，通过逆傅里叶变化可以得到去卷积后的结果$\hat x(t)$。

解释

上面的式子可以改写成更为清晰的形式

$$\begin{split} G(f) &= {1 \over {H(f)}}\left[ {{{|H(f){|^2}} \over {|H(f){|^2} + {{N(f)} \over {S(f)}}}}} \right] \\ &={1 \over {H(f)}}\left[ {{{|H(f){|^2}} \over {|H(f){|^2} + {1 \over {SNR(f)}}}}} \right] \end{split}$$
这里$H(f)$是$h$在频率域$f$的傅里叶变换。$SNR(f)=S(f)/N(f)$是信号噪声比。当噪声为零时（即信噪比趋近于无穷），方括号内各项也就等于1，意味着此时刻维纳滤波也就简化成逆滤波过程。但是当噪声增加时，信噪比降低，方括号里面值也跟着降低。这说明，维纳滤波的带通频率依赖于信噪比。

推导

上面直接给出了维纳滤波的表达式，接下来介绍推导过程。
上面提到，维纳滤波是建立在最小均方差，可以如下表示：

$$e(f) = E|X(f) - \hat X(f){|^2}$$
这里$E$是期望
如果我们替换表达式中的$\hat X(f)$，上面可以重新组合成
$$\begin{split} e(f) &= E|X(f) - G(f)Y(f){|^2} \ &= E|X(f) - G(f)[H(f)X(f) + V(f)]{|^2}\ &= E|[1 - G(f)H(f)]X(f) - G(f)V(f){|^2} \end{split}$$
展开二次方，得到下式：
$$\eqalign{ e(f)&= [1 - G(f)H(f)]{[1 - G(f)H(f)]^*}E{\rm{|X(f)}}{{\rm{|}}^2} \cr &- [1 - G(f)H(f)]{G^*}(f)E\{ X(f){V^*}(f)\} \cr &- G(f){[1 - G(f)H(f)]^*}E\{ V(f){X^*}(f)\} \cr &+ G(f){G^{\rm{*}}}(f)E|V(f){|^2} \cr}$$
然而，我们假设噪声与信号独立无关，这样有
$$E\{ X(f){V^*}(f)\} =E\{ V(f){X^*}(f)\}=0$$
并且我们如下定义功率谱
$$S(f) =E|X(f)|^2\N(f) =E|V(f)|^2$$
这样我们有
$$\eqalign{ & e(f) = [1 - G(f)H(f)]{[1 - G(f)H(f)]^*}S(f) + G(f){G^*}(f)N(f)}$$
为了得到最小值，我们对$G(f)$求导，令方程等于零。
$${{{\rm d}(f)} \over {{\rm d}G(f)}} = {G^*}(f)N(f) - H(f){[1 - G(f)H(f)]^*}S(f) = 0$$
由此最终推出维纳滤波器。

测试

Matlab自带了示例程序，如下

%Read image
I = im2double(imread(‘cameraman.tif‘));
figure,subplot(2,3,1),imshow(I);
title(‘Original Image (courtesy of MIT)‘);

%Simulate a motion blur
LEN = 21;
THETA = 11;
PSF = fspecial(‘motion‘, LEN, THETA);
blurred = imfilter(I, PSF, ‘conv‘, ‘circular‘);
subplot(2,3,2),imshow(blurred);
title(‘Blurred Image‘);

%Restore the blurred image
wnr1 = deconvwnr(blurred, PSF, 0);
subplot(2,3,3),imshow(wnr1);
title(‘Restored Image‘);

%Simulate blur and noise
noise_mean = 0;
noise_var = 0.0001;
blurred_noisy = imnoise(blurred, ‘gaussian‘, ...
                        noise_mean, noise_var);
subplot(2,3,4),imshow(blurred_noisy)
title(‘Simulate Blur and Noise‘)

%Restore the blurred and noisy image:First attempt
wnr2 = deconvwnr(blurred_noisy, PSF, 0);
subplot(2,3,5);imshow(wnr2);title(‘Restoration of Blurred, Noisy Image Using NSR = 0‘)

%Restore the Blurred and Noisy Image: Second Attempt
signal_var = var(I(:));
wnr3 = deconvwnr(blurred_noisy, PSF, noise_var / signal_var);
subplot(2,3,6),imshow(wnr3)
title(‘Restoration of Blurred, Noisy Image Using Estimated NSR‘);

维纳滤波需要估计图像的信噪比(SNR)或者噪信比(NSR)，信号的功率谱使用图像的方差近似估计，噪声分布是已知的。从第一排中可以看出，若无噪声，此时维纳滤波相当于逆滤波，恢复运动模糊效果是极好的。从第二排可以看出噪信比估计的准确性对图像影响比较大的，二排中间效果几乎不可用。

参考阅读

http://en.wikipedia.org/wiki/Wiener_deconvolution 英文维基百科
http://www.owlnet.rice.edu/~elec539/Projects99/BACH/proj2/wiener.html 莱斯大学的项目资料

转载保留声明

作者	日期	联系方式
风吹夏天	2015年5月29日	wincoder@qq.com