• 自适应中值滤波
• 高斯滤波
• 双边滤波
• 导向滤波

### 自适应中值滤波

#### 自适应中值滤波器

• 滤除椒盐噪声
• 平滑其他非脉冲噪声
• 尽可能的保护图像中细节信息，避免图像边缘的细化或者粗化。

• 窗口中的最小灰度值
• 窗口中的最大灰度值
• 窗口中的灰度值的中值
• 表示坐标 处的灰度值
• 允许的最大窗口尺寸

A :

B:

#### 算法实现

``````uchar adaptiveProcess(const Mat &im, int row,int col,int kernelSize,int maxSize)
{
vector<uchar> pixels;
for (int a = -kernelSize / 2; a <= kernelSize / 2; a++)
for (int b = -kernelSize / 2; b <= kernelSize / 2; b++)
{
pixels.push_back(im.at<uchar>(row + a, col + b));
}
sort(pixels.begin(), pixels.end());
auto min = pixels[0];
auto max = pixels[kernelSize * kernelSize - 1];
auto med = pixels[kernelSize * kernelSize / 2];
auto zxy = im.at<uchar>(row, col);
if (med > min && med < max)
{
// to B
if (zxy > min && zxy < max)
return zxy;
else
return med;
}
else
{
kernelSize += 2;
if (kernelSize <= maxSize)
return adpativeProcess(im, row, col, kernelSize, maxSize); // 增大窗口尺寸，继续A过程。
else
return med;
}
}``````

https://github.com/zhangqizky/common-image-filteringgithub.com

#### 高斯滤波

``````//O(m * n * ksize^2)
void GaussianFilter(const Mat &src, Mat &dst, int ksize, double sigma)
{
CV_Assert(src.channels() || src.channels() == 3); //只处理3通道或单通道的图片
double **GaussianTemplate = new double *[ksize];
for(int i = 0; i < ksize; i++){
GaussianTemplate[i] = new double [ksize];
}
generateGaussianTemplate(GaussianTemplate, ksize, sigma);
int border = ksize / 2;
copyMakeBorder(src, dst, border, border, border, border, BORDER_CONSTANT);
int channels = dst.channels();
int rows = dst.rows - border;
int cols = dst.cols - border;
for(int i = border; i < rows; i++){
for(int j = border; j< cols; j++){
double sum[3] = {0};
for(int a = -border; a <= border; a++){
for(int b = -border; b <= border; b++){
if(channels == 1){
sum[0] += GaussianTemplate[border+a][border+b] * dst.at<uchar>(i+a, j+b);
}else if(channels == 3){
Vec3b rgb = dst.at<Vec3b>(i+a, j+b);
auto k = GaussianTemplate[border+a][border+b];
sum[0] += k * rgb[0];
sum[1] += k * rgb[1];
sum[2] += k * rgb[2];
}
}
}
for(int k = 0; k < channels; k++){
if(sum[k] < 0) sum[k] = 0;
else if(sum[k] > 255) sum[k] = 255;
}
if(channels == 1){
dst.at<uchar >(i, j) = static_cast<uchar >(sum[0]);
}else if(channels == 3){
Vec3b rgb = {static_cast<uchar>(sum[0]), static_cast<uchar>(sum[1]), static_cast<uchar>(sum[2])};
dst.at<Vec3b>(i, j) = rgb;
}
}
}
for(int i = 0; i < ksize; i++)
delete[] GaussianTemplate[i];
delete[] GaussianTemplate;
}``````

``````//分离实现高斯滤波
//O(m*n*k)
void separateGaussianFilter(const Mat &src, Mat &dst, int ksize, double sigma){
assert(src.channels()==1 || src.channels() == 3); //只处理单通道或者三通道图像
//生成一维的
double *matrix = new double[ksize];
double sum = 0;
int origin = ksize / 2;
for(int i = 0; i < ksize; i++){
double g = exp(-(i-origin) * (i-origin) / (2 * sigma * sigma));
sum += g;
matrix[i] = g;
}
for(int i = 0; i < ksize; i++) matrix[i] /= sum;
int border = ksize / 2;
copyMakeBorder(src, dst, border, border, border, border, BORDER_CONSTANT);
int channels = dst.channels();
int rows = dst.rows - border;
int cols = dst.cols - border;
//水平方向
for(int i = border; i < rows; i++){
for(int j = border; j < cols; j++){
double sum[3] = {0};
for(int k = -border; k<=border; k++){
if(channels == 1){
sum[0] += matrix[border + k] * dst.at<uchar>(i, j+k);
}else if(channels == 3){
Vec3b rgb = dst.at<Vec3b>(i, j+k);
sum[0] += matrix[border+k] * rgb[0];
sum[1] += matrix[border+k] * rgb[1];
sum[2] += matrix[border+k] * rgb[2];
}
}
for(int k = 0; k < channels; k++){
if(sum[k] < 0) sum[k] = 0;
else if(sum[k] > 255) sum[k] = 255;
}
if(channels == 1)
dst.at<Vec3b>(i, j) = static_cast<uchar>(sum[0]);
else if(channels == 3){
Vec3b rgb = {static_cast<uchar>(sum[0]), static_cast<uchar>(sum[1]), static_cast<uchar>(sum[2])};
dst.at<Vec3b>(i, j) = rgb;
}
}
}
//竖直方向
for(int i = border; i < rows; i++){
for(int j = border; j < cols; j++){
double sum[3] = {0};
for(int k = -border; k<=border; k++){
if(channels == 1){
sum[0] += matrix[border + k] * dst.at<uchar>(i+k, j);
}else if(channels == 3){
Vec3b rgb = dst.at<Vec3b>(i+k, j);
sum[0] += matrix[border+k] * rgb[0];
sum[1] += matrix[border+k] * rgb[1];
sum[2] += matrix[border+k] * rgb[2];
}
}
for(int k = 0; k < channels; k++){
if(sum[k] < 0) sum[k] = 0;
else if(sum[k] > 255) sum[k] = 255;
}
if(channels == 1)
dst.at<Vec3b>(i, j) = static_cast<uchar>(sum[0]);
else if(channels == 3){
Vec3b rgb = {static_cast<uchar>(sum[0]), static_cast<uchar>(sum[1]), static_cast<uchar>(sum[2])};
dst.at<Vec3b>(i, j) = rgb;
}
}
}
delete [] matrix;
}``````

https://github.com/zhangqizky/common-image-filtering/tree/maingithub.com

#### 双边滤波

opencv中提供了bilateralFilter()函数来实现双边滤波操作，其原型如下：

``````void cv::bilateralFilter(InputArray src,
OutputArray 	dst,
int 	d,
double 	sigmaColor,
double 	sigmaSpace,
int 	borderType = BORDER_DEFAULT
)``````
• InputArray src: 输入图像，可以是Mat类型，图像必须是8位整型或浮点型单通道、三通道的图像。
• OutputArray dst: 输出图像，和原图像有相同的尺寸和类型。
• int d: 表示在过滤过程中每个像素邻域的直径范围。如果这个值是非正数，则函数会从第五个参数sigmaSpace计算该值。
• double sigmaColor: 颜色空间过滤器的值，这个参数的值月大，表明该像素邻域内有越宽广的颜色会被混合到一起，产生较大的半相等颜色区域。（这个参数可以理解为值域核的 和 )
• double sigmaSpace: 坐标空间中滤波器的sigma值，如果该值较大，则意味着越远的像素将相互影响，从而使更大的区域中足够相似的颜色获取相同的颜色。当d>0时，d指定了邻域大小且与sigmaSpace无关，否则d正比于sigmaSpace. （这个参数可以理解为空间域核的 和 ）
• int borderType=BORDER_DEFAULT: 用于推断图像外部像素的某种边界模式，有默认值BORDER_DEFAULT.

``````#include <iostream>
#include <opencv2/opencv.hpp>

using namespace std;
using namespace cv;

//定义全局变量
const int g_ndMaxValue = 100;
const int g_nsigmaColorMaxValue = 200;
const int g_nsigmaSpaceMaxValue = 200;
int g_ndValue;
int g_nsigmaColorValue;
int g_nsigmaSpaceValue;

Mat g_srcImage;
Mat g_dstImage;

//定义回调函数
void on_bilateralFilterTrackbar(int, void*);

int main()
{

//判断图像是否加载成功
if(g_srcImage.empty())
{
cout << "图像加载失败!" << endl;
return -1;
}
else
cout << "图像加载成功!" << endl << endl;

namedWindow("原图像", WINDOW_AUTOSIZE);
imshow("原图像", g_srcImage);

//定义输出图像窗口属性和轨迹条属性
namedWindow("双边滤波图像", WINDOW_AUTOSIZE);
g_ndValue = 10;
g_nsigmaColorValue = 10;
g_nsigmaSpaceValue = 10;

char dName[20];
sprintf(dName, "邻域直径 %d", g_ndMaxValue);

char sigmaColorName[20];
sprintf(sigmaColorName, "sigmaColor %d", g_nsigmaColorMaxValue);

char sigmaSpaceName[20];
sprintf(sigmaSpaceName, "sigmaSpace %d", g_nsigmaSpaceMaxValue);

//创建轨迹条
createTrackbar(dName, "双边滤波图像", &g_ndValue, g_ndMaxValue, on_bilateralFilterTrackbar);
on_bilateralFilterTrackbar(g_ndValue, 0);

createTrackbar(sigmaColorName, "双边滤波图像", &g_nsigmaColorValue,
g_nsigmaColorMaxValue, on_bilateralFilterTrackbar);
on_bilateralFilterTrackbar(g_nsigmaColorValue, 0);

createTrackbar(sigmaSpaceName, "双边滤波图像", &g_nsigmaSpaceValue,
g_nsigmaSpaceMaxValue, on_bilateralFilterTrackbar);
on_bilateralFilterTrackbar(g_nsigmaSpaceValue, 0);

waitKey(0);

return 0;
}

void on_bilateralFilterTrackbar(int, void*)
{
bilateralFilter(g_srcImage, g_dstImage, g_ndValue, g_nsigmaColorValue, g_nsigmaSpaceValue);
imshow("双边滤波图像", g_dstImage);
}``````

### 导向滤波(Guide Filter)

``````void cv::ximgproc::guidedFilter	(	InputArray 	guide,
InputArray 	src,
OutputArray 	dst,
double 	eps,
int 	dDepth = -1
)``````
 guide 引导图像，3通道，如果大于3通道则只有前三个会被用到 src 待滤波图像 dst 输出图像 radius 导向滤波的窗口 eps 正则化参数 dDepth 可选，图像的深度参数

``````import numpy as np
import scipy as sp
import scipy.ndimage

def box(img, r):
""" O(1) box filter
img - >= 2d image
r   - radius of box filter
"""
(rows, cols) = img.shape[:2]
imDst = np.zeros_like(img)

tile = [1] * img.ndim
tile[0] = r
imCum = np.cumsum(img, 0)
imDst[0:r+1, :, ...] = imCum[r:2*r+1, :, ...]
imDst[r+1:rows-r, :, ...] = imCum[2*r+1:rows, :, ...] - imCum[0:rows-2*r-1, :, ...]
imDst[rows-r:rows, :, ...] = np.tile(imCum[rows-1:rows, :, ...], tile) - imCum[rows-2*r-1:rows-r-1, :, ...]

tile = [1] * img.ndim
tile[1] = r
imCum = np.cumsum(imDst, 1)
imDst[:, 0:r+1, ...] = imCum[:, r:2*r+1, ...]
imDst[:, r+1:cols-r, ...] = imCum[:, 2*r+1 : cols, ...] - imCum[:, 0 : cols-2*r-1, ...]
imDst[:, cols-r: cols, ...] = np.tile(imCum[:, cols-1:cols, ...], tile) - imCum[:, cols-2*r-1 : cols-r-1, ...]

return imDst

def _gf_color(I, p, r, eps, s=None):
""" Color guided filter
I - guide image (rgb)
p - filtering input (single channel)
eps - regularization (roughly, variance of non-edge noise)
s - subsampling factor for fast guided filter
"""
fullI = I
fullP = p
if s is not None:
I = sp.ndimage.zoom(fullI, [1/s, 1/s, 1], order=1)
p = sp.ndimage.zoom(fullP, [1/s, 1/s], order=1)
r = round(r / s)

h, w = p.shape[:2]
N = box(np.ones((h, w)), r)

mI_r = box(I[:,:,0], r) / N
mI_g = box(I[:,:,1], r) / N
mI_b = box(I[:,:,2], r) / N

mP = box(p, r) / N

# mean of I * p
mIp_r = box(I[:,:,0]*p, r) / N
mIp_g = box(I[:,:,1]*p, r) / N
mIp_b = box(I[:,:,2]*p, r) / N

# per-patch covariance of (I, p)
covIp_r = mIp_r - mI_r * mP
covIp_g = mIp_g - mI_g * mP
covIp_b = mIp_b - mI_b * mP

# symmetric covariance matrix of I in each patch:
#       rr rg rb
#       rg gg gb
#       rb gb bb
var_I_rr = box(I[:,:,0] * I[:,:,0], r) / N - mI_r * mI_r;
var_I_rg = box(I[:,:,0] * I[:,:,1], r) / N - mI_r * mI_g;
var_I_rb = box(I[:,:,0] * I[:,:,2], r) / N - mI_r * mI_b;

var_I_gg = box(I[:,:,1] * I[:,:,1], r) / N - mI_g * mI_g;
var_I_gb = box(I[:,:,1] * I[:,:,2], r) / N - mI_g * mI_b;

var_I_bb = box(I[:,:,2] * I[:,:,2], r) / N - mI_b * mI_b;

a = np.zeros((h, w, 3))
for i in range(h):
for j in range(w):
sig = np.array([
[var_I_rr[i,j], var_I_rg[i,j], var_I_rb[i,j]],
[var_I_rg[i,j], var_I_gg[i,j], var_I_gb[i,j]],
[var_I_rb[i,j], var_I_gb[i,j], var_I_bb[i,j]]
])
covIp = np.array([covIp_r[i,j], covIp_g[i,j], covIp_b[i,j]])
a[i,j,:] = np.linalg.solve(sig + eps * np.eye(3), covIp)

b = mP - a[:,:,0] * mI_r - a[:,:,1] * mI_g - a[:,:,2] * mI_b

meanA = box(a, r) / N[...,np.newaxis]
meanB = box(b, r) / N

if s is not None:
meanA = sp.ndimage.zoom(meanA, [s, s, 1], order=1)
meanB = sp.ndimage.zoom(meanB, [s, s], order=1)

q = np.sum(meanA * fullI, axis=2) + meanB

return q

def _gf_gray(I, p, r, eps, s=None):
""" grayscale (fast) guided filter
I - guide image (1 channel)
p - filter input (1 channel)
r - window raidus
eps - regularization (roughly, allowable variance of non-edge noise)
s - subsampling factor for fast guided filter
"""
if s is not None:
Isub = sp.ndimage.zoom(I, 1/s, order=1)
Psub = sp.ndimage.zoom(p, 1/s, order=1)
r = round(r / s)
else:
Isub = I
Psub = p

(rows, cols) = Isub.shape

N = box(np.ones([rows, cols]), r)

meanI = box(Isub, r) / N
meanP = box(Psub, r) / N
corrI = box(Isub * Isub, r) / N
corrIp = box(Isub * Psub, r) / N
varI = corrI - meanI * meanI
covIp = corrIp - meanI * meanP

a = covIp / (varI + eps)
b = meanP - a * meanI

meanA = box(a, r) / N
meanB = box(b, r) / N

if s is not None:
meanA = sp.ndimage.zoom(meanA, s, order=1)
meanB = sp.ndimage.zoom(meanB, s, order=1)

q = meanA * I + meanB
return q

def _gf_colorgray(I, p, r, eps, s=None):
""" automatically choose color or gray guided filter based on I's shape """
if I.ndim == 2 or I.shape[2] == 1:
return _gf_gray(I, p, r, eps, s)
elif I.ndim == 3 and I.shape[2] == 3:
return _gf_color(I, p, r, eps, s)
else:
print("Invalid guide dimensions:", I.shape)

def guided_filter(I, p, r, eps, s=None):
""" run a guided filter per-channel on filtering input p
I - guide image (1 or 3 channel)
p - filter input (n channel)
r - window raidus
eps - regularization (roughly, allowable variance of non-edge noise)
s - subsampling factor for fast guided filter
"""
if p.ndim == 2:
p3 = p[:,:,np.newaxis]

out = np.zeros_like(p3)
for ch in range(p3.shape[2]):
out[:,:,ch] = _gf_colorgray(I, p3[:,:,ch], r, eps, s)
return np.squeeze(out) if p.ndim == 2 else out

def test_gf():
import imageio

r = 8
eps = 0.05

cat_smoothed = guided_filter(cat, cat, r, eps)
cat_detail = cat / cat_smoothed
print(cat_detail.shape)
cat_smoothed_s4 = guided_filter(cat, cat, r, eps, s=4)

imageio.imwrite('cat_smoothed.png', cat_smoothed)
imageio.imwrite('cat_smoothed_s4.png', cat_smoothed_s4)
imageio.imwrite('cat_smoothed_detailed.png',cat_detail)

tulips_smoothed4s = np.zeros_like(tulips)
tulips_detailed = np.zeros_like(tulips)
for i in range(3):
tulips_smoothed4s[:,:,i] = guided_filter(tulips, tulips[:,:,i], r, eps, s=4)

tulips_detailed = tulips / tulips_smoothed4s
imageio.imwrite('tulips_detailed.png',tulips_detailed)
imageio.imwrite('tulips_smoothed4s.png', tulips_smoothed4s)

tulips_smoothed = np.zeros_like(tulips)
for i in range(3):
tulips_smoothed[:,:,i] = guided_filter(tulips, tulips[:,:,i], r, eps)
imageio.imwrite('tulips_smoothed.png', tulips_smoothed)

if __name__ == '__main__':
test_gf()``````

https://github.com/atilimcetin/guided-filter