BoxFilter包滤波器的Matlab代码实现分析（基础）

最近在使用引导滤波器，速度相对于双边滤波器快了很多，主要是算法中涉及到大量的窗口内的求和、求均值、求方差等，利用了BoxFilter滤波器，使得上述过程得到了极大的加速，所以效率很高，这里我查了相关的资料，结合经典的暗通道去雾算法中给出的引导滤波的Matlab代码，分析其在Matlab中实现的过程，比较基础。

BoxFilter包滤波器的原理，在这篇文章中有介绍，http://www.cnblogs.com/easymind223/archive/2012/11/13/2768680.html

总的来说，BoxFilter包的功能可以大致概括如下：
输入： 有 2个 对象：
第一、待处理的矩阵（input）（在此讨论的是二维矩阵，即针对灰度图像的处理中的应用）；
第二、窗口半径（r）；
输出：out 和 input 维数大小相等，每个元素out(i,j)表示以该input(i,j)为中心的半径为r的窗口内的所有元素的和，当处于边缘时，窗口自动根据边框缩小。

首先给出BoxFilter包滤波器的Matlab代码如下：

function imDst = boxfilter(imSrc, r)%   BOXFILTER   O(1) time box filtering using cumulative sum%%   - Definition imDst(x, y)=sum(sum(imSrc(x-r:x+r,y-r:y+r)));%   - Running time independent of r; %   - Equivalent to the function: colfilt(imSrc, [2*r+1, 2*r+1], 'sliding', @sum);%   - But much faster.[hei, wid] = size(imSrc);imDst = zeros(size(imSrc));%cumulative sum over Y axisimCum = cumsum(imSrc, 1);%列向量主次向上求和，输出同大小的矩阵%difference over Y axisimDst(1:r+1, :) = imCum(1+r:2*r+1, :);imDst(r+2:hei-r, :) = imCum(2*r+2:hei, :) - imCum(1:hei-2*r-1, :);imDst(hei-r+1:hei, :) = repmat(imCum(hei, :), [r, 1]) - imCum(hei-2*r:hei-r-1, :);%cumulative sum over X axisimCum = cumsum(imDst, 2);%difference over Y axisimDst(:, 1:r+1) = imCum(:, 1+r:2*r+1);imDst(:, r+2:wid-r) = imCum(:, 2*r+2:wid) - imCum(:, 1:wid-2*r-1);imDst(:, wid-r+1:wid) = repmat(imCum(:, wid), [1, r]) - imCum(:, wid-2*r:wid-r-1);end

在理解了上面所给出的那篇文章对于BoxFilter包滤波器原理的介绍之后，通过一个输入例子，给出Matlab的实现过程。
设定输入矩阵，A（10，10）=imput；如下图：

整个过程分两个部分，分别对行和列对矩阵操作，其原理和方法相同，在此只介绍处理矩阵列的方法。步骤如下：
1）、利用Matlab自身工具箱函数 comsum ，得到输入矩阵的每个元素对应的列向上求和矩阵。
comsum 函数：逐项求和。调用格式详见matlab帮助

imCum = cumsum(imSrc, 1);%列向量主次向上求和，输出同大小的矩阵

此处，imSrc 为输入的缓存，imSrc=imput；参数‘1’表示沿元素所在列对 向上所有行求和，当参数为‘2’时刚好相反。
得到的输出矩阵 imCum 的每一个元素 imCum(i,j) 对应于输入的 imput(i,j) 所在 j 列，i 行以前 的imput(1 : i ,j)的和。如图所示输出imCum。

2）、对矩阵顶部边缘，求每个元素中心半径为 r 的不完整窗口内的和。

```imDst(1:r+1, :) = imCum(1+r:2*r+1, :);```此处，从第1行 到 第 r+1 行，该元素中心半径为 r 的窗口，包含矩阵中的行元素的和 即是 该元素 沿列向上求和，，也就是imCum对应的值。 结果如图![顶部边缘求和](http://img.blog.csdn.net/20151215103248245)

3）、对矩阵中间部分，每个元素的半径为 r 的窗口都能完整的包含在矩阵中，则窗口内的和为 2*r+1 个元素的和。

```imDst(r+2:hei-r, :) = imCum(2*r+2:hei, :) - imCum(1:hei-2*r-1, :);```此处，利用imCum的特性，相隔 n 行的 imCum(i+n,:) 与 imCum(i,:) 的差即是中间 n 个行元素的和。也即是BoxFilter的巧妙的利用这种方法，降低了循环累加求和的时间复杂度。在上述代码中 n=2*r+1 ，即对于任意窗口大小 r ，计算元素半径r窗口大小的元素和只需一次减法就能完成。![对应相减](http://img.blog.csdn.net/20151215105755962)![中间部分求和](http://img.blog.csdn.net/20151215110812356)当imCum(2*r+2,;) - imCum(1,:)时，得到中间2*r+1个元素的和，这些元素对应的中心元素 为imput(r+2,:),同理直元素imput(hei-r,:),这些元素半径为r的窗口包含元素的和都是按上述方法计算。

4）、对矩阵底部边缘部分，仿照对顶端不完整窗口求和方法，稍加改动即可。

```imDst(hei-r+1:hei, :) = repmat(imCum(hei, :), [r, 1]) - imCum(hei-2*r:hei-r-1, :);```此处，原理同样是imCum后面的元素减去前面的元素，得到相隔的元素的和。repmat 函数，是matlab自身工具箱函数，对矩阵整理，具体功能见帮助，不做赘述。![用最后一行元素做差](http://img.blog.csdn.net/20151215111907642)![底部求和](http://img.blog.csdn.net/20151215143648643)

通过以上步骤，就可以计算出输入imput的每个元素以r为半径的中心窗口涵盖的行元素的和，同理依照此方法，对输出imDst，求元素中心窗口涵盖的列元素的和，即得到了整个窗口元素的和。代码如下：

```imCum = cumsum(imDst, 2);%difference over Y axisimDst(:, 1:r+1) = imCum(:, 1+r:2*r+1);imDst(:, r+2:wid-r) = imCum(:, 2*r+2:wid) - imCum(:, 1:wid-2*r-1);imDst(:, wid-r+1:wid) = repmat(imCum(:, wid), [1, r]) - imCum(:, wid-2*r:wid-r-1);```

最终，得到输入input的每个元素的以r为中心的窗口所涵盖的所有元素的和imDst。如图

以上是建立在对之前提的那篇关于BoxFilter原理的文章的理解的基础上，对Matlab平台的实现过程的分析。关于在其他平台上的实现还没有研究，有问题和疑问，欢迎大家一起交流。