盒子(方框)滤波(BoxFilter)原理及C++及Matlab实现

本文主要是介绍盒子(方框)滤波(BoxFilter)原理及C++及Matlab实现，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

写在前面

盒子滤波是一种非常有用的线性滤波，也叫方框滤波，最简单的均值滤波就是盒子滤波归一化的情况。

应用：可以说，一切需要求某个邻域内像素之和的场合，都有盒子滤波的用武之地，比如：均值滤波、引导滤波、计算Haar特征等等。

优势：就一个字：快！它可以使复杂度为O(MN)的求和，求方差等运算降低到O(1)或近似于O(1)的复杂度，也就是说与邻域尺寸无关了，有点类似积分图吧，但是貌似比积分图更快（与它的实现方式有关）。

opencv函数：

void boxFilter( InputArray src, OutputArray dst, int ddepth,Size ksize, Point anchor = Point(-1,-1),bool normalize = true,int borderType = BORDER_DEFAULT );

原理

在原理上，和均值滤波一样，用一个内核和图像进行卷积：

$\texttt{K} = \alpha \begin{bmatrix} 1 & 1 & 1 & \cdots & 1 & 1 \\ 1 & 1 & 1 & \cdots & 1 & 1 \\ \hdotsfor{6} \\ 1 & 1 & 1 & \cdots & 1 & 1 \end{bmatrix}$

其中：

$\alpha = \fork{\frac{1}{\texttt{ksize.width*ksize.height}}}{when \texttt{normalize=true}}{1}{otherwise}$

可见，归一化了就是均值滤波；不归一化则可以计算每个像素邻域上的各种积分特性，方差、协方差，平方和等等。

实现

c++实现

Note:

1、我这里用的积分图思想实现的，虽然效果一样，但速度慢一些，所以算不上真正意义上的盒子滤波实现形式，若要看真正的实现方式，可以参考：https://www.cnblogs.com/lwl2015/p/4460711.html。

2、这个c++程序只是实验，仅仅为了学习盒子滤波的原理。若真正的去应用，例如用到引导滤波中，这个程序还不够稳健，或许会出问题，因为没有考虑多个通道以及多种数据类型的情况。建议可以进一步看看OpenCV关于boxfitler的源码。

#include <iostream>
#include <opencv2/core.hpp>
#include <opencv2/highgui.hpp>
#include <opencv2/imgproc.hpp>/
//求积分图-优化方法
//由上方negral(i-1,j)加上当前行的和即可
//对于W*H图像：2*(W-1)*(H-1)次加减法
//比常规方法快1.5倍左右
/
void Fast_integral(cv::Mat& src, cv::Mat& dst){int nr = src.rows;int nc = src.cols;int sum_r = 0;dst = cv::Mat::zeros(nr + 1, nc + 1, CV_64F);for (int i = 1; i < dst.rows; ++i){for (int j = 1, sum_r = 0; j < dst.cols; ++j){//行累加，因为积分图相当于在原图上方加一行，左边加一列，所以积分图的(1,1)对应原图(0,0),(i,j)对应(i-1,j-1)sum_r = src.at<uchar>(i - 1, j - 1) + sum_r; //行累加dst.at<double>(i, j) = dst.at<double>(i - 1, j) + sum_r;}}
}//
//盒子滤波-均值滤波是其特殊情况
/
void BoxFilter(cv::Mat& src, cv::Mat& dst, cv::Size wsize, bool normalize){//图像边界扩充if (wsize.height % 2 == 0 || wsize.width % 2 == 0){fprintf(stderr, "Please enter odd size!");exit(-1);}int hh = (wsize.height - 1) / 2;int hw = (wsize.width - 1) / 2;cv::Mat Newsrc;cv::copyMakeBorder(src, Newsrc, hh, hh, hw, hw, cv::BORDER_REFLECT);//以边缘为轴，对称src.copyTo(dst);//计算积分图cv::Mat inte;Fast_integral(Newsrc, inte);//BoxFilterdouble mean = 0;for (int i = hh + 1; i < src.rows + hh + 1; ++i){  //积分图图像比原图（边界扩充后的）多一行和一列 for (int j = hw + 1; j < src.cols + hw + 1; ++j){double top_left = inte.at<double>(i - hh - 1, j - hw - 1);double top_right = inte.at<double>(i - hh - 1, j + hw);double buttom_left = inte.at<double>(i + hh, j - hw - 1);double buttom_right = inte.at<double>(i + hh, j + hw);if (normalize == true)mean = (buttom_right - top_right - buttom_left + top_left) / wsize.area();elsemean = buttom_right - top_right - buttom_left + top_left;//一定要进行判断和数据类型转换if (mean < 0)mean = 0;else if (mean>255)mean = 255;dst.at<uchar>(i - hh - 1, j - hw - 1) = static_cast<uchar>(mean);}}
}int main(){cv::Mat src = cv::imread("I:\\Learning-and-Practice\\2019Change\\Image process algorithm\\Img\\woman2.jpeg");if (src.empty()){return -1;}//自编BoxFilter测试cv::Mat dst1;double t2 = (double)cv::getTickCount(); //测时间if (src.channels() > 1){std::vector<cv::Mat> channel;cv::split(src, channel);BoxFilter(channel[0], channel[0], cv::Size(7, 7), true);//盒子滤波BoxFilter(channel[1], channel[1], cv::Size(7, 7), true);//盒子滤波BoxFilter(channel[2], channel[2], cv::Size(7, 7), true);//盒子滤波cv::merge(channel,dst1);}elseBoxFilter(src, dst1, cv::Size(7, 7), true);//盒子滤波t2 = (double)cv::getTickCount() - t2;double time2 = (t2 *1000.) / ((double)cv::getTickFrequency());std::cout << "FASTmy_process=" << time2 << " ms. " << std::endl << std::endl;//opencv自带BoxFilter测试cv::Mat dst2;double t1 = (double)cv::getTickCount(); //测时间cv::boxFilter(src, dst2, -1, cv::Size(7, 7), cv::Point(-1, -1), true, cv::BORDER_CONSTANT);//盒子滤波t1 = (double)cv::getTickCount() - t1;double time1 = (t1 *1000.) / ((double)cv::getTickFrequency());std::cout << "Opencvbox_process=" << time1 << " ms. " << std::endl << std::endl;cv::namedWindow("src");cv::imshow("src", src);cv::namedWindow("ourdst",CV_WINDOW_NORMAL);cv::imshow("ourdst", dst1);cv::namedWindow("opencvdst", CV_WINDOW_NORMAL);cv::imshow("opencvdst", dst2);cv::waitKey(0);}

Matlab实现

Note: 来自何恺明大神主页引导滤波代码 http://kaiminghe.com/

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 axis
imCum = cumsum(imSrc, 1);
%difference over Y axis
imDst(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 axis
imCum = cumsum(imDst, 2);
%difference over Y axis
imDst(:, 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