自适应滤波-最小均方误差滤波

1.最小均方误差滤波原理

        低通滤波不能像中值滤波那样很好的滤除冲激噪声。因为低通滤波的最终结果混合了图像信号无关的噪声和信号本身。相反，中值滤波能够在保护图像边缘不受损失的情况下，滤除与图像信号无关的噪声。但是当噪声不完全和图像信号无关，比如被混合了图像本身信号和一定噪声的加性噪声或乘性噪声污染的图像，我们该如何消除此种噪声呢？或许需要构造一种更加智能的滤波。

       自适应滤波的提出就是为了解决上述问题。自适应滤波利用图像局部特性和结构自适应选择合适的方法滤除噪声。根据图像局部邻域统计信息自适应滤波的时候，如果我们知道被处理图像的先验知识将更高效的滤除噪声。比如图像的噪声统计指标-噪声方差。

    自适应滤波被广泛应用于受斑点噪声污染的图像，比如超声波检查中所拍摄的图像和合成孔径雷达所获取的图像。斑点噪声的模型一般用如下公式表示：

    为受噪声污染的图像， 为给定标准差，均值为0的高斯噪声。由公式可知，斑点噪声取决于输入信号的强度。斑点噪声由于影响了图像的对比图和空间分辨率，阻碍了对图像结果的解析。如下图分别是中值滤波，MATLAB工具箱中维纳滤波以及自适应滤波的实际效果图：

        绝大部分自适应滤波利用图像邻域局部均值、局部方差等图像局部统计信息自适应调整滤波强度。局部均值，即邻域内所有像素求其平均值；局部方差分两步计算可得到：第一步计算邻域内像素平方和均值；第二步，从第一步中的结果减去邻域像素平均值的平方。经过这两步就可以得到邻域内像素的方差。具体推导步骤可参看之前的博客《图像比较之模板匹配》。

其中，为待处理图像，为邻域大小。为邻域均值，为邻域方差。

最小均方误差滤波一般用来滤除加性白噪声和斑点噪声。我们假设待处理图像为，邻域大小为，图像中噪声的方差为，局部的均值为，局部方差为。则最小均方误差滤波可以表示为：

 

 

其中噪声的方差我们可以对邻域内所有像素所求的方差的平均值进行估算。这种方法只适合针对加性噪声，并不适合乘性噪。

我们可以从上述最小均方误差数学表达式可知：

1. 当邻域局部方差远大于噪声的方差，或者说噪声方差比较小或者接近于0，则经过滤波之后输出像素为图像本身对应的像素值。当邻域局部方差远大于噪声方差，在原始图像有可能表示的是边缘区域，我们最好保留这种边缘区域信息；
2. 当噪声方差数值超过局部邻域方差，则返回为图像局部均值；
3. 当噪声方差和局部邻域方差相差不大的时候，返回原始像素本身和局部均值加权后的结果；

 

小均方误差滤波的保边效果，主要有由性质1决定。性质1决定了最小均方误差滤波能够保留图像绝大部分分边缘信息，尽管边缘部分的噪声不能消除。该滤波的伪代码可表示为：

 

       上图是受斑点噪声污染图像，经过3*3邻域中值滤波和最小均方误差滤波的结果对比。可以看到，最小均方误差滤波比中值滤波稍好，滤除噪声的同时保留了图像清晰度。从图像峰值信噪比可知，最小均方误差滤波比中值滤波提高了1/3db左右。

     最小均方误差滤波适合处理斑点噪声污染的图像，对其他噪声污染的图像并不合适。下图是盐椒噪声污染的图像经过中值滤波和最小均方误差滤波的结果对比。此种情况下，中值滤波效果要好于最小均方误差滤波效果。

2.最小均方误差滤波实现

 

这是一个定点优化的版本。其中的定点优化方面有几处可以学习。

1. 在上一节最小均方误差滤波的伪代码中，对局部均值和方差求解过程中，针对邻域内每一次局部均值求解需要9次加法操作，每一次局部方差计算时，分别需要9次加法和9次乘法操作。由于相邻两个像素的邻域有重叠的部分，其计算局部均值和方差的时候我们可以保留其重叠的部分像素，然后随着像素从左到右移动，计算相邻像素邻域均值和方差的时候，可以减去邻域最左部分像素的灰度值和像素灰度平方和，然后再加上最右边像素的灰度值和像素灰度平方和。这样，只需要计算最左边和最右边灰度值和像素灰度平方和，即求均值的时候只需要计算6次加法，计算方差时只需要6次加法和6次乘法运算。其计算数量比原来减少了1/3.

其具体方法类似于下图所示操作：

1. 为了避免比较耗时的除法操作，当局部邻域像素和求取之后，乘以了用Q12格式表示小数1/9的定点化所对应的数值455，然后再进行右移12计算其局部均值结果；
2. 同样的，在最小均方误差滤波最后结果需要对局部均值和原始输入值，根据噪声方差和局部方差比值进行线性插值。而噪声方差和局部方差的比值需要涉及到除法运算。在此，运用了定点格式化的牛顿-拉夫逊方法对除法进行估算计算。针对a/b的除法操作，改写为a*(1/b)，针对Q15格式所表示的小数1/x定点表示为,即32768/x。
3. 牛顿-拉夫逊方法中对其进行倒数估算的时候，使用了查找表的方式指定其估算的初始值。牛顿-拉夫逊方法对1/x求其倒数的时候，其迭代关系式为：

   

   即Q15格式浮点数定点化之后的牛顿-拉夫逊方法迭代关系式用C语言表示为：

   

   当时，或跳出循环时认为找到x的倒数。对于上述最小均方误差滤波定点格式化代码实现中，只需要经过六次就可以找到其对应的倒数。这个应该和构建牛顿-拉夫逊方法初始值的查找表有关。为什么如此构建初始值的查找表不太明白。

参考资料：Embedded Image Processing on the TMS320C6000 dsp