1.算法描述
图像去噪是用于解决图像由于噪声干扰而导致其质量下降的问题，通过去噪技术可以有效地提高图像质量，增大信噪比，更好的体现原来图像所携带的信息。在我们的图像中常见的噪声主要有以下4种：加性噪声、乘性噪声、量化噪声、椒盐噪声。根据不同的噪声特点，我们可以采用不同的去噪算法，按照数学运算主要分为两大类，一类是通过滤波（相当于积分的过程），又可以在空域（和频率域（傅立叶变换和小波变换）中分别采用此操作，比如说空域中值滤波对于椒盐噪声有很好的处理效果，另一类是通过偏微分方程，具有各向异性的特点，具有平滑图像和将边缘尖锐化的能力，在低噪声密度的图像处理中取得了较好的效果，但是在处理高噪声密度图像时去噪效果不好。
普通的时空域的低通滤波器，在像素空间完成滤波以后，导致图像的边缘部分也变得不那么明显，整张图像都变得同样的模糊，图像边缘细节丢失。双边滤波器（ABilateral Filter）可以很好的保留边缘的同时消除噪声。双边滤波器能做到这些原因在于它不像普通的高斯/卷积低通滤波，只考虑了位置对中心像素的影响，它还考虑了卷积核中像素与中心像素之间相似程度的影响，根据位置影响与像素值之间的相似程度生成两个不同的权重表(WeightTable)，在计算中心像素的时候加以考虑这两个权重，从而实现双边低通滤波。据说AdobePhotoshop的高斯磨皮功能就是应用了双边低通滤波算法实现。
双边滤波（Bilateral filter）是一种非线性的滤波方法，是结合图像的空间邻近度和像素值相似度的一种折中处理，同时考虑空域信息和灰度相似性，达到保边去噪的目的。具有简单、非迭代、局部的特点。双边滤波器的好处是可以做边缘保存（edge preserving），一般过去用的维纳滤波或者高斯滤波去降噪，都会较明显地模糊边缘，对于高频细节的保护效果并不明显。
双边滤波（Bilateral filter）是一种非线性的滤波方法，是结合图像的空间邻近度和像素值相似度的一种折衷处理，同时考虑空域信息和灰度相似性，达到保边去噪的目的。具有简单、非迭代、局部的特点。双边滤波器的好处是可以做边缘保存（edge preserving），一般过去用的维纳滤波或者高斯滤波去降噪，都会较明显地模糊边缘，对于高频细节的保护效果并不明显。双边滤波器顾名思义比高斯滤波多了一个高斯方差sigma－d，它是基于空间分布的高斯滤波函数，所以在边缘附近，离的较远的像素不会太多影响到边缘上的像素值，这样就保证了边缘附近像素值的保存。但是由于保存了过多的高频信息，对于彩色图像里的高频噪声，双边滤波器不能够干净的滤掉，只能够对于低频信息进行较好的滤波。
双边滤波之所以能够既作平滑处理又保留边界，是因为它综合了高斯滤波器和α-截尾均值滤波器的特点，同时考虑了空间域与值域，即其核是由空间域核和值域核相乘得到。滤波算法中，目标点上的像素值通常是由其所在位置上的周围的一个小局部邻居像素的值所决定。在2D高斯滤波中的具体实现就是对周围的一定范围内的像素值分别赋以不同的高斯权重值，并在加权平均后得到当前点的最终结果。而这里的高斯权重因子是利用两个像素之间的空间距离（在图像中为2D）关系来生成。通过高斯分布的曲线可以发现，离目标像素越近的点对最终结果的贡献越大，反之则越小。其公式化的描述一般如下所述：
高斯滤波在低通滤波算法中有不错的表现，但是其却有另外一个问题，那就是只考虑了像素间的空间位置上的关系，因此滤波的结果会丢失边缘的信息。这里的边缘主要是指图像中主要的不同颜色区域（比如蓝色的天空，黑色的头发等），而Bilateral就是在Gaussian blur中加入了另外的一个权重分部来解决这一问题。Bilateral滤波中对于边缘的保持通过下述表达式来实现：
双边滤波器中，输出像素的值依赖于邻域像素的值的加权组合，
权重系数w(i,j,k,l)取决于定义域核
和值域核
的乘积同时考虑了空间域与值域的差别。一般过去用的维纳滤波或者高斯滤波去降噪，只考虑了空间域差别，都会较明显地模糊边缘，对于高频细节的保护效果并不明显；α-截尾均值滤波器，去掉百分率为α的最小值和最大之后剩下像素的均值作为滤波器，只考虑了值域差别。
2.仿真效果预览
matlab2022a仿真结果如下：
3.MATLAB核心程序
`Img = double( imread('./images/ckb.jpg') );
[m, n] = size(Img);
% create noisy image (additive Gaussian noise)
sigma = 20;
inImg = Img + sigma * randn(m, n);
% filter parameters
sigma1 = 4;
sigma2 = 30;
tol = 0.01;
% Set window for spatial Gaussian
w = 6*sigma1;
if (mod(w,2) == 0)
w = w+1;
end
% call bilateral filter
tic;
[outImg, param] = shiftableBF(inImg, sigma1, sigma2, w, tol);
toc;
% plot results
T = param.T;
N = param.N;
M = param.M;
gamma = 1 / (sqrt(N) * sigma2);
twoN = 2^N;
warning('off'); %#ok
s = linspace(-T, T, 200);
g = exp( -s.^2 / (2 * sigma2 *sigma2) );
gApprox = cos(gamma * s).^N;
if M == 0
gTrunc = gApprox;
else
gTrunc = zeros( 1, length(s) );
for k = M : N - M
gTrunc = gTrunc + (nchoosek(N, k) / twoN) * ...
cos( (2*k - N) * gamma * s );
end
end
figure('Units','normalized','Position',[0 0.5 1 0.5]);
plot(s, g, 'b');
hold on,
plot(s, gApprox, 'm'),
hold on,
plot(s, gTrunc, 'r');
axis('tight'), grid('on'),
legend('Gassian','Raised cosine','Truncated raised cosine','FontSize', 10);
title('Comparison of the range kernels', 'FontSize', 10),
peak = 255;
PSNR0 = 10 * log10(m * n * peak^2 / sum(sum( (inImg - Img).^2)) );
PSNR1 = 10 * log10(m * n * peak^2 / sum(sum((outImg - Img).^2)) );`
