開始之前 這幾天由於自己的原因沒有寫, 一個是因為自己懶了, 一個是感覺這裡遇到點問題不想往下寫了, 我們先努力結束這個章節吧, 之前介紹了比較常用而且比較好理解的均值和中值濾波, 但是呢,在常式 "Smoothing Images" , 還有給出的其他的濾波方式, 主要是高斯濾波和雙邊濾波, 我們 ...
開始之前
這幾天由於自己的原因沒有寫, 一個是因為自己懶了, 一個是感覺這裡遇到點問題不想往下寫了, 我們先努力結束這個章節吧, 之前介紹了比較常用而且比較好理解的均值和中值濾波, 但是呢,在常式Smoothing Images, 還有給出的其他的濾波方式, 主要是高斯濾波和雙邊濾波,
我們這一次完結掉濾波與平滑的這個部分, 寫的有點多了,反而不想再寫了, 加油
目錄
本文目標
本文主要是介紹
- 高斯濾波
- 雙邊濾波
和之前介紹的一樣, 我們仍然還是 介紹一下原理, 給出一下具體的形式, 然後使用 opencv 進行一下實現的過程, 最後使用我們之前的圖像進行測試 進行演算法的分析與總結.
正文
高斯濾波(Gaussian Filter)
我們在之前介紹了中值濾波是統計排序的結果, 屬於非線性的結果, 均值濾波是使用模板核進行的操作, 我們在的文章中也提到了均值濾波在計算的過程中必須要考慮權重的問題, 進而提出了加權的均值濾波的操作, 比如最常見的加權均值濾波的操作核.
\[M = \frac{1}{16} \left [ \begin{array}{c} 1 & 2 & 1 \\ 2& 4 & 2 \\ 1 & 2 & 1 \end{array} \right ] \]
其實呢,這個核也就是高斯濾波器在 3*3
視窗的離散取整的結果值, 最明顯的特點就是模板的繫數隨著距離模板中心的距離而變換, 能夠有效的抑制雜訊,平滑圖像, 相比均值濾波能夠更好的平滑圖像, 保留圖像邊緣.
高斯濾波原理
由於我們的圖像是二維的, 但是高斯分佈是一維的, 那我們先考慮一維的高斯分佈, 就是我們常用的正太分佈曲線,
\[G(x) = \frac{1}{\sqrt{2\pi \sigma}} e^{-\frac{x^2}{2\sigma^2}} \]
對於二維的高斯分佈其實可以考慮成兩個方向的運算相疊加的得到的結果
\[G(x,y) = \frac{1}{2\pi \sigma^2} e^{-\frac{x^2+y^2}{2\sigma^2}} = G(x)*G(y) \]
考慮到圖像的計算實際上是離散的座標, 對於視窗大小為 \((2k + 1) \times (2k + 1)\) 模板, 我們可以表示成
\[G{i,j} = \frac{1}{2\pi \sigma ^ 2}e ^{-\frac{(i - k - 1)^2 + (j - k - 1)^2}{2 \sigma ^ 2}} \]
可以參考圖像處理基礎(4):高斯濾波器詳解
裡面給出的方法, 使用
void generateGaussianTemplate(double window[][11], int ksize, double sigma)
{
static const double pi = 3.1415926;
int center = ksize / 2; // 模板的中心位置,也就是坐標的原點
double x2, y2;
for (int i = 0; i < ksize; i++)
{
x2 = pow(i - center, 2);
for (int j = 0; j < ksize; j++)
{
y2 = pow(j - center, 2);
double g = exp(-(x2 + y2) / (2 * sigma * sigma));
g /= 2 * pi * sigma*sigma; //
window[i][j] = g;
}
}
double k = 1 / window[0][0]; // 將左上角的繫數歸一化為1
for (int i = 0; i < ksize; i++)
{
for (int j = 0; j < ksize; j++)
{
window[i][j] *= k;
}
}
}
生成了\(3 \times 3, \sigma = 0.8\) 的高斯模板, 對應的將其取整就得到了
\[M = \frac{1}{16} \left [ \begin{array}{c} 1 & 2 & 1 \\ 2& 4 & 2 \\ 1 & 2 & 1 \end{array} \right ] \]
上面給出的文章同樣的詳細介紹了 \(\sigma\) 在統計學中的意義, 可以去參考學習
不過根據高中的知識, 我們可以看到 正態分佈的曲線
C++ 實現
在我們之前提到的圖像處理基礎(4):高斯濾波器詳解 這裡給出了基於 opencv 的代碼實現, 這裡是\(O(m*n*k^2)\) 的演算法實現
// 來源鏈接: https://www.cnblogs.com/wangguchangqing/p/6407717.html
void GaussianFilter(const Mat &src, Mat &dst, int ksize, double sigma)
{
CV_Assert(src.channels() || src.channels() == 3); // 只處理單通道或者三通道圖像
const static double pi = 3.1415926;
// 根據視窗大小和sigma生成高斯濾波器模板
// 申請一個二維數組,存放生成的高斯模板矩陣
double **templateMatrix = new double*[ksize];
for (int i = 0; i < ksize; i++)
templateMatrix[i] = new double[ksize];
int origin = ksize / 2; // 以模板的中心為原點
double x2, y2;
double sum = 0;
for (int i = 0; i < ksize; i++)
{
x2 = pow(i - origin, 2);
for (int j = 0; j < ksize; j++)
{
y2 = pow(j - origin, 2);
// 高斯函數前的常數可以不用計算,會在歸一化的過程中給消去
double g = exp(-(x2 + y2) / (2 * sigma * sigma));
sum += g;
templateMatrix[i][j] = g;
}
}
for (int i = 0; i < ksize; i++)
{
for (int j = 0; j < ksize; j++)
{
templateMatrix[i][j] /= sum;
cout << templateMatrix[i][j] << " ";
}
cout << endl;
}
// 將模板應用到圖像中
int border = ksize / 2;
copyMakeBorder(src, dst, border, border, border, border, BorderTypes::BORDER_REFLECT);
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] += templateMatrix[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 = templateMatrix[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[] templateMatrix[i];
delete[] templateMatrix;
}
然後同樣的給出了分離的實現, 將圖像進行水平運算之後再進行豎直運算, 計算的時間上會有一定的速度提升
// 來源鏈接: https://www.cnblogs.com/wangguchangqing/p/6407717.html
// 分離的計算
void separateGaussianFilter(const Mat &src, Mat &dst, int ksize, double sigma)
{
CV_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, BorderTypes::BORDER_REFLECT);
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<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 = 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<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;
}
}
}
delete[] matrix;
}
這裡的演算法都是 上面提到的https://www.cnblogs.com/wangguchangqing/p/6407717.html 這篇文章, 具體可以去看內容
opencv 高斯濾波
其實這篇文章圖像處理--高斯濾波寫的很好
其實主要的結構也就是他給出的過程
其實整個高斯濾波的過程就是創建高斯核, 然後使用 filter2D 的方法進行的濾波操作, 具體要深入的話可以看函數的調用圖, 實現起來也是一樣的思路, 很簡單的操作, 我們之後測試一下效果..
// /modules\imgproc\src\smooth.dispatch.cpp:600
void GaussianBlur(InputArray _src, OutputArray _dst, Size ksize,
double sigma1, double sigma2,
int borderType)
- src 輸入圖像
- dst 輸出圖像
- ksize 核的尺寸 奇數
- sigmaX x 方向 的 sigma 值
- sigmaY y 方向 的 sigma 值
- borderType 邊界處理的方式
高斯濾波效果對比
我們還是使用之前的高椒鹽雜訊圖像, 然後直接進行演算法濾波, 計算結果就好, 跟之前的測試圖像很相似, 這裡
這裡的四張圖分別對應 高雜訊圖像, 直接高斯濾波的結果, 分離xy方向進行濾波結果,以及opencv 自帶的高斯濾波效果圖, 這裡是預覽圖像, 實際的檢測結果就是上面給出的參數值, 實際上效果只能說一般, 我們之後再進行演算法層面的對比.
image-noise: psnr:19.4727, mssim: B:0.353134 G:0.383638 R:0.629353
image-noise: psnr:26.3155, mssim: B:0.584585 G:0.617172 R:0.812303
image-noise: psnr:26.1721, mssim: B:0.574719 G:0.607494 R:0.809844
image-noise: psnr:26.4206, mssim: B:0.598176 G:0.630657 R:0.819658
雙邊濾波(Bilateral Filter)
雙邊濾波原理
我們在上面提出了高斯濾波的原理是對於距離模板中心 距離不同給予不同的權重, 而雙邊濾波則不僅考慮圖像的空間距離, 還要考慮其灰度距離, 對於越接近中間灰度值的點權重越高, 灰度值相差大的則權重更小.
雙邊濾波的原理可以參考雙邊濾波(Bilateral Filter)詳解,
可以參考Bilateral Filtering for Gray and Color Images
在文章圖像處理基礎(5):雙邊濾波器詳細介紹了雙邊濾波
其實跟上面給出的濾波演示一致, 都是在保證圖像邊緣信息的情況下進行雜訊的濾波..
可以參考bilateral filter雙邊濾波器的通俗理解 給出的雙邊濾波的數學表達
\[g(x,y) = \frac{\sum_{kl}f(k,l)w(i,j,k,l)}{\sum_{kl}w(i,j,k,l)} \]
對於不同的模板繫數又有兩個部分, 主要是 空間域模板權值 \(w_d\) 和 灰度域 模板權值 \(w_r\),
\[\begin{array}{rl} w_d(i,j,k,l) &= e^{-\frac{(i-k)^2 +(j-l)^2}{2\sigma_d^2}} \\ w_r(i,j,k,l) &= e^{-\frac{\left \| f(i,j) - f(k,l) \right \|} {2\sigma_r^2}} \\ w &= w_d * w_r \end{array} \]
其中,\(q(i,j)\) 為模板視窗的其他繫數的坐標,\(f(i,j)\) 表示圖像在點\(q(i,j)\) 處的像素值;\(p(k,l)\) 為模板視窗的中心坐標點,對應的像素值為\(f(k,l)\) ;\(\sigma_r\) 為高斯函數的標準差。
C++ 實現 雙邊濾波
感覺這裡寫的挺好的 圖像處理基礎(5):雙邊濾波器, 手動實現了雙邊濾波, 我們可以詳細的參考, 這裡
// 參考來源: https://www.cnblogs.com/wangguchangqing/p/6416401.html
void myBilateralFilter(const Mat &src, Mat &dst, int ksize, double space_sigma, double color_sigma)
{
int channels = src.channels();
CV_Assert(channels == 1 || channels == 3);
double space_coeff = -0.5 / (space_sigma * space_sigma);
double color_coeff = -0.5 / (color_sigma * color_sigma);
int radius = ksize / 2;
Mat temp;
copyMakeBorder(src, temp, radius, radius, radius, radius, BorderTypes::BORDER_REFLECT);
vector<double> _color_weight(channels * 256); // 存放差值的平方
vector<double> _space_weight(ksize * ksize); // 空間模板繫數
vector<int> _space_ofs(ksize * ksize); // 模板視窗的坐標
double *color_weight = &_color_weight[0];
double *space_weight = &_space_weight[0];
int *space_ofs = &_space_ofs[0];
for (int i = 0; i < channels * 256; i++)
color_weight[i] = exp(i * i * color_coeff);
// 生成空間模板
int maxk = 0;
for (int i = -radius; i <= radius; i++)
{
for (int j = -radius; j <= radius; j++)
{
double r = sqrt(i*i + j * j);
if (r > radius)
continue;
space_weight[maxk] = exp(r * r * space_coeff); // 存放模板繫數
space_ofs[maxk++] = i * temp.step + j * channels; // 存放模板的位置,和模板繫數相對應
}
}
// 濾波過程
for (int i = 0; i < src.rows; i++)
{
const uchar *sptr = temp.data + (i + radius) * temp.step + radius * channels;
uchar *dptr = dst.data + i * dst.step;
if (channels == 1)
{
for (int j = 0; j < src.cols; j++)
{
double sum = 0, wsum = 0;
int val0 = sptr[j]; // 模板中心位置的像素
for (int k = 0; k < maxk; k++)
{
int val = sptr[j + space_ofs[k]];
double w = space_weight[k] * color_weight[abs(val - val0)]; // 模板繫數 = 空間繫數 * 灰度值繫數
sum += val * w;
wsum += w;
}
dptr[j] = (uchar)cvRound(sum / wsum);
}
}
else if (channels == 3)
{
for (int j = 0; j < src.cols * 3; j+=3)
{
double sum_b = 0, sum_g = 0, sum_r = 0, wsum = 0;
int b0 = sptr[j];
int g0 = sptr[j + 1];
int r0 = sptr[j + 2];
for (int k = 0; k < maxk; k++)
{
const uchar *sptr_k = sptr + j + space_ofs[k];
int b = sptr_k[0];
int g = sptr_k[1];
int r = sptr_k[2];
double w = space_weight[k] * color_weight[abs(b - b0) + abs(g - g0) + abs(r - r0)];
sum_b += b * w;
sum_g += g * w;
sum_r += r * w;
wsum += w;
}
wsum = 1.0f / wsum;
b0 = cvRound(sum_b * wsum);
g0 = cvRound(sum_g * wsum);
r0 = cvRound(sum_r * wsum);
dptr[j] = (uchar)b0;
dptr[j + 1] = (uchar)g0;
dptr[j + 2] = (uchar)r0;
}
}
}
}
opencv 實現 雙邊濾波
void bilateralFilter( InputArray _src, OutputArray _dst, int d,
double sigmaColor, double sigmaSpace,
int borderType )
- InputArray src: 輸入圖像,可以是Mat類型,圖像必須是8位或浮點型單通道、三通道的圖像。
- OutputArray dst: 輸出圖像,和原圖像有相同的尺寸和類型。
- int d: 表示在過濾過程中每個像素鄰域的直徑範圍。如果這個值是非正數,則函數會從第五個參數sigmaSpace計算該值。
- double sigmaColor: 顏色空間過濾器的sigma值,這個參數的值月大,表明該像素鄰域內有越寬廣的顏色會被混合到一起,產生較大的半相等顏色區域。 (這個參數可以理解為值域核的)
- double sigmaSpace: 坐標空間中濾波器的sigma值,如果該值較大,則意味著越遠的像素將相互影響,從而使更大的區域中足夠相似的顏色獲取相同的顏色。當d>0時,d指定了鄰域大小且與sigmaSpace無關,否則d正比於sigmaSpace. (這個參數可以理解為空間域核的)
- int borderType=BORDER_DEFAULT: 用於推斷圖像外部像素的某種邊界模式,有預設值BORDER_DEFAULT.
雙邊濾波演算法對比
一開始的時候看雙邊濾波真的搞不懂, 也不知道這麼做有什麼目的, 最終的結果又代表什麼, 我們按照之前的方法去測試我們的圖像, 結果真的是幾種演算法中最差的了, 但是這隻是說不適用於我們的圖像結果, 在實際使用過程中還是要進行測試之後才能得出結論
測試結果如下: 對應原始圖和 手動實現的結果以及 opencv 的結果 都使用的 是3 的視窗, sigma 的值 為 255
, 這篇文章https://blog.csdn.net/Jfuck/article/details/8932978 講的很好, 介紹了參數對濾波的影響, 可以學習一下..
image-noise: psnr:19.4727, mssim: B:0.353134 G:0.383638 R:0.629353
image-noise: psnr:24.4502, mssim: B:0.538774 G:0.570666 R:0.776195
image-noise: psnr:24.4691, mssim: B:0.539177 G:0.571087 R:0.776461
總結
其實個人使用雙邊濾波真的不算很多, 在之前研究導向濾波的時候才瞭解過很多, 這裡寫的比較差吧, 只能說勉強能看, 強烈推薦 https://www.cnblogs.com/wangguchangqing/category/740760.html 這個系列, 將的很詳細, 很多都是博文裡面的內容, 可以參考學習, 高斯濾波就比較簡單了, 其實複雜的濾波過程主要是理解演算法, 然後根據演算法的思路進行代碼的實現過程, 最後做一定的程式上的優化就好, 理解第一, 實現其次.. 希望帶給讀者一點點啟發..
我這裡一開始不准備寫這麼多的, 結果越寫越多, 導致自己收不住了, 很多自己說不上很瞭解的地方, 這一次也是深入的瞭解了一下, 但是還是很僵硬, 只能說能用而已, 這裡還是推薦看我給出的鏈接或者自己去查閱相關的內容, 我這裡只是給出一個大略的介紹, 如果有錯誤還請指名, 十分感謝
參考鏈接
- 《快速高斯濾波、高斯模糊、高斯平滑(二維捲積分步為一維捲積)_人工智慧_青城山小和尚-CSDN博客》. 見於 2020年5月10日. https://blog.csdn.net/qq_36359022/article/details/80188873.
- 《雙邊濾波 - 旗亭涉的博客 | Qitingshe Blog》. 見於 2020年5月10日. https://qitingshe.github.io/2018/06/14/雙邊濾波/.
- 《雙邊濾波(Bilateral Filter)詳解_人工智慧_Jfuck的專欄-CSDN博客》. 見於 2020年5月10日. https://blog.csdn.net/Jfuck/article/details/8932978.
- 《雙邊濾波器》. 收入 維基百科,自由的百科全書, 2019年11月16日. https://zh.wikipedia.org/w/index.php?title=雙邊濾波器&oldid=56898678.
- 《圖像處理--高斯濾波_網路_L-inYi的專欄-CSDN博客》. 見於 2020年5月10日. https://blog.csdn.net/L_inYi/article/details/8915116.
- 《圖像處理基礎(4):高斯濾波器詳解 - Brook_icv - 博客園》. 見於 2020年5月10日. https://www.cnblogs.com/wangguchangqing/p/6407717.html.
- 《圖像處理基礎(5):雙邊濾波器 - Brook_icv - 博客園》. 見於 2020年5月10日. https://www.cnblogs.com/wangguchangqing/p/6416401.html.
- 《圖像處理-線性濾波-3 高斯濾波器 - Tony Ma - 博客園》. 見於 2020年5月10日. https://www.cnblogs.com/pegasus/archive/2011/05/20/2052031.html.
- 《【轉】高斯圖像濾波原理及其編程離散化實現方法_Smile_Gogo_新浪博客》. 見於 2020年5月10日. http://blog.sina.com.cn/s/blog_640577ed0100yz8v.html.
- 《bilateral filter雙邊濾波器的通俗理解_網路_pan_jinquan的博客-CSDN博客》. 見於 2020年5月10日. https://blog.csdn.net/guyuealian/article/details/82660826.
- 《Bilateral Filtering》. 見於 2020年5月10日. http://homepages.inf.ed.ac.uk/rbf/CVonline/LOCAL_COPIES/MANDUCHI1/Bilateral_Filtering.html.
- 《Cv圖像處理 - OpenCV China :圖像處理,電腦視覺庫,Image Processing, Computer Vision》. 見於 2020年5月10日. http://wiki.opencv.org.cn/index.php/Cv圖像處理.
- 《o(1)複雜度之雙邊濾波演算法的原理、流程、實現及效果。 - 雲+社區 - 騰訊雲》. 見於 2020年5月10日. https://cloud.tencent.com/developer/article/1011738.
本文由博客一文多發平臺 OpenWrite 發佈!