1三、OpenCV實現圖像的空間濾波——圖像平滑卷積在圖像處理中的應用（轉自https://medium.com/@irhumshafkat/intuitively-understanding-co

時間 2019-12-14

標籤 opencv 實現圖像空間濾波平滑圖像處理應用 https medium.com medium com irhumshafkat intuitively understanding 简体版

原文原文鏈接

一、空間濾波基礎概念

一、空間濾波基礎

空間濾波一詞中濾波取自數字信號處理，指接受或拒絕必定的頻率成分，可是空間濾波學習內容實際上和經過傅里葉變換實現的頻域的濾波是等效的，故而也稱爲濾波。空間濾波主要直接基於領域（空間域）對圖像中的像素執行計算，用濾波器（也成爲空間掩膜、核、模板和窗口）直接做用於圖像自己完成相似的平滑。html

二、空間濾波機理

對空間域中的每一點(x,y)，重複以下操做：ios

對預先定義的以（x,y)爲中心的領域內的像素進行預約義運算。
將（1）中運算的結果做爲(x,y)點新的響應。

上述過程就稱爲鄰域處理或空間域濾波。算法

三、空間濾波分類

根據預約義的操做，能夠將濾波器分爲：數組

線性濾波器
非線性濾波器

而根據濾波器最終對圖像形成的影響，能夠將濾波器分爲：dom

平滑濾波器，一般用於模糊圖像或者去除圖像中的噪聲
銳化濾波器，突出圖像中的邊緣細節部分

四、空間相關和卷積

　在執行線性空間vlbo時，有兩個相近的概念須要理解，相關和卷積。相關是濾波器模板移過圖像並計算每一個位置乘積之和的處理。卷積機理相似，可是濾波器首先須要旋轉180°。關於卷積的詳細說明，能夠參考文章函數

卷積在圖像處理中的應用（轉自https://medium.com/@irhumshafkat/intuitively-understanding-convolutions-for-deep-learning-1f6f42faee1）

一些算子及不一樣類型的卷積

這裏再也不作重複說明。post

二、爲圖片添加噪聲

使用攝像頭拍攝圖片，噪聲是難以免的，這裏重點提下兩種噪聲：學習

椒鹽噪聲：噪聲的幅值基本上相同，可是噪聲出現的位置是隨機的；（中值濾波效果好）
高斯噪聲：每一點都存在噪聲，但噪聲的幅值是隨機分佈的。

可是，咱們網上找到的圖片大可能是通過處理的，或者說噪聲不是明顯的，這裏爲了研究，爲圖片添加噪聲。優化

椒鹽噪聲表如今圖片上是圖片中離散分佈的黑點和白點ui

// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
    // 隨機數產生器
    std::random_device rd; //種子
    std::mt19937 gen(rd()); // 隨機數引擎

    auto cols = m.cols * m.channels();

    for (int i = 0; i < num; i++)
    {
        auto row = static_cast<int>(gen() % m.rows);
        auto col = static_cast<int>(gen() % cols);

        auto p = m.ptr<uchar>(row);
        p[col++] = 255;
        p[col++] = 255;
        p[col] = 255;
    }
}

高斯噪聲是一種加性噪聲，爲圖像添加高斯噪聲的代碼以下：

// 添加Gussia噪聲
// 使用指針訪問
void addGaussianNoise(Mat &m, int mu, int sigma)
{
    // 產生高斯分佈隨機數發生器
    std::random_device rd;
    std::mt19937 gen(rd());

    std::normal_distribution<> d(mu, sigma);

    auto rows = m.rows; // 行數
    auto cols = m.cols * m.channels(); // 列數

    for (int i = 0; i < rows; i++)
    {
        auto p = m.ptr<uchar>(i); // 取得行首指針
        for (int j = 0; j < cols; j++)
        {
            auto tmp = p[j] + d(gen);
            tmp = tmp > 255 ? 255 : tmp;
            tmp = tmp < 0 ? 0 : tmp;
            p[j] = tmp;
        }
    }
}

示例：

　　　爲圖像添加兩種噪聲的程序以下：

//實現直方圖的反向投影
#include "stdafx.h"
#include <math.h>
#include <random>
#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>
using namespace cv;
using namespace std;
// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
	// 隨機數產生器
	std::random_device rd; //種子
	std::mt19937 gen(rd()); // 隨機數引擎

	auto cols = m.cols * m.channels();

	for (int i = 0; i < num; i++)
	{
		auto row = static_cast<int>(gen() % m.rows);
		auto col = static_cast<int>(gen() % cols);

		auto p = m.ptr<uchar>(row);
		p[col++] = 255;
		p[col++] = 255;
		p[col] = 255;
	}
}
void addGaussianNoise(Mat &m, int mu, int sigma)
{
	// 產生高斯分佈隨機數發生器
	std::random_device rd;
	std::mt19937 gen(rd());

	std::normal_distribution<> d(mu, sigma);

	auto rows = m.rows; // 行數
	auto cols = m.cols * m.channels(); // 列數

	for (int i = 0; i < rows; i++)
	{
		auto p = m.ptr<uchar>(i); // 取得行首指針
		for (int j = 0; j < cols; j++)
		{
			auto tmp = p[j] + d(gen);
			tmp = tmp > 255 ? 255 : tmp;
			tmp = tmp < 0 ? 0 : tmp;
			p[j] = tmp;
		}
	}
}

int main()
{
	Mat srcImage = imread("111.jpg");
	Mat src_salt, src_gaussian;
	src_salt = imread("111.jpg");
	src_gaussian = imread("111.jpg");
	if (!srcImage.data)
	{
		cout << "圖像打開失敗！" << endl;
		return -1;
	}
	addSaltNoise(src_salt, 100);
	addGaussianNoise(src_gaussian,100,1);
	imshow("原圖", srcImage);
	imshow("添加椒鹽噪聲後", src_salt);
	imshow("添加高斯噪聲後", src_gaussian);
	waitKey();
	return 0;
}

添加完噪聲後的圖像以下：

三、線性平滑空間濾波

平滑濾波器主要用於模糊處理和下降噪音。

一、平滑線性濾波器(又稱均值濾波）

　　首先須要說明的是圖像平滑是一種減小和抑制圖像噪聲的實用數字圖像處理技術，通常來講，圖像具備局部連續性質，即相鄰像素的數值相近，而噪聲的存在使得在噪聲點處產生灰度跳躍，但通常能夠合理的假設偶爾出現的噪聲並無改變圖像局部連續的性質。其中，線性濾波是現實將空域模板對應鄰域內像素的灰度值加權之和做爲鄰域內中心像素的相應輸出。線性平滑模板的權係數全爲正值並且係數之和等於1，以下圖所示：

　　所以這種方法不會增長圖像中整體的灰度程度。也就是說在灰度一致的區域，線性平滑濾波的相應輸出不變。此外，能夠知道，線性平滑濾波等效於低通濾波。可是一樣的咱們能夠看出，線性平滑模板在模糊和降噪的同時，圖像中的邊緣和細節的銳度也都丟失了，也就是說，在平滑的過程當中，使得圖像的一部分細節信息丟失。

在OpenCV中，函數blur表示使用該模板的均值濾波器，其聲明以下：

void blur( InputArray src, OutputArray dst, Size ksize, Point anchor = Point(-1,-1), int borderType = BORDER_DEFAULT );

src是輸入圖像，dst爲輸出圖像；ksize是濾波器模板窗口的大小；後兩個參數分別表示，待處理像素在模板窗口的位置，默認值是窗口的中心位置，因此窗口的大小通常爲奇數，最後一個參數表示對編解類型的處理，使用默認值便可。其調用示例blur(src,dst,Size(5,5)，模板窗口的大小爲 $5 \times 5$

示例：

下面經過程序對兩種含噪聲的圖片進行線性平滑濾波處理

//實現直方圖的反向投影
#include "stdafx.h"
#include <math.h>
#include <random>
#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;

// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
	// 隨機數產生器
	std::random_device rd; //種子
	std::mt19937 gen(rd()); // 隨機數引擎

	auto cols = m.cols * m.channels();

	for (int i = 0; i < num; i++)
	{
		auto row = static_cast<int>(gen() % m.rows);
		auto col = static_cast<int>(gen() % cols);

		auto p = m.ptr<uchar>(row);
		p[col++] = 255;
		p[col++] = 255;
		p[col] = 255;
	}
}
void addGaussianNoise(Mat &m, int mu, int sigma)
{
	// 產生高斯分佈隨機數發生器
	std::random_device rd;
	std::mt19937 gen(rd());

	std::normal_distribution<> d(mu, sigma);

	auto rows = m.rows; // 行數
	auto cols = m.cols * m.channels(); // 列數

	for (int i = 0; i < rows; i++)
	{
		auto p = m.ptr<uchar>(i); // 取得行首指針
		for (int j = 0; j < cols; j++)
		{
			auto tmp = p[j] + d(gen);
			tmp = tmp > 255 ? 255 : tmp;
			tmp = tmp < 0 ? 0 : tmp;
			p[j] = tmp;
		}
	}
}

int main()
{
	Mat srcImage = imread("111.jpg");
	Mat src_salt, src_gaussian;
	src_salt = imread("111.jpg");
	src_gaussian = imread("111.jpg");
	if (!srcImage.data)
	{
		cout << "圖像打開失敗！" << endl;
		return -1;
	}
	addSaltNoise(src_salt, 100);
	addGaussianNoise(src_gaussian,100,1);

    //平滑濾波
	Mat aft_mean_filter_salt, aft_mean_filter_gaussian;
	blur(src_salt, aft_mean_filter_salt, Size(5, 5), Point(-1, -1));
	blur(src_gaussian, aft_mean_filter_gaussian, Size(5, 5), Point(-1, -1));


	//加強圖像對比度
	//Mat kernel = (Mat_<int>(3, 3) << 0, -1, 0, -1, 5, -1, 0, -1, 0);
	//filter2D(src_salt, src_salt, src_salt.depth(), kernel, Point(-1, -1));

	imshow("原圖", srcImage);
	imshow("添加椒鹽噪聲後", src_salt);
	imshow("添加高斯噪聲後", src_gaussian);
	imshow("通過平滑濾波後含椒鹽噪聲的的結果", aft_mean_filter_salt);
	imshow("通過平滑濾波後含高斯噪聲的的結果", aft_mean_filter_gaussian);
	waitKey();
	return 0;
}

二、高斯平滑

　　均值平滑對於領域內的像素一視同仁，爲了減小平滑處理中的模糊，獲得更天然的平滑效果，咱們能夠很天然的想到適當加大模板中心點的權重，隨着距離中心點的距離增長，權重迅速減少，從而能夠確保中心點看起來更接近於與它距離更近的點，基於這樣的考慮獲得的模板即爲高斯模板。經常使用的3X3高斯模板以下所示：

　　高斯模板名字的由來是高斯函數，即二維正態分佈密度函數，一個二維的高斯函數以下,

　　$$
　　h(x, y)=e^{-\frac{x^{2}+y^{2}}{2 \sigma^{2}}}
　　$$

公式中(x,y)是點座標 $(x, y)$

對於窗口模板的大小爲 (2K+1)X(2K+1)的矩陣，其(i，j）位置的元素值可以下肯定：

高斯模板的生成方法以下所示：

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;
            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;
        }
    }
}

　　須要一個二維數組，存放生成的係數（這裏假設模板的最大尺寸不會超過11）；第二個參數是模板的大小（不要超過11）；第三個參數就比較重要了，是高斯分佈的標準差。

生成的過程，首先根據模板的大小，找到模板的中心位置ksize/2。而後就是遍歷，根據高斯分佈的函數，計算模板中每一個係數的值。

須要注意的是，最後歸一化的過程，使用模板左上角的係數的倒數做爲歸一化的係數（左上角的係數值被歸一化爲1），模板中的每一個係數都乘以該值（左上角係數的倒數），而後將獲得的值取整，就獲得了整數型的高斯濾波器模板。

至於小數形式的生成也比較簡單，去掉歸一化的過程，而且在求解過程後，模板的每一個係數要除以全部係數的和。具體代碼以下：

void generateGaussianTemplate(double window[][11], int ksize, double sigma)
{
    static const double pi = 3.1415926;
    int center = ksize / 2; // 模板的中心位置，也就是座標的原點
    double x2, y2;
    double sum = 0;
    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;
            sum += g;
            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] /= sum;
        }
    }
}

$σ$

　　當標準差 $σ$

$σ$
$σ$

$σ$

void GaussianBlur( InputArray src, OutputArray dst, Size ksize, double sigmaX, double sigmaY=0, int borderType=BORDER_DEFAULT );

　這裏因爲高斯函數的可分離性，尺寸較大的高斯濾波器能夠分紅兩步進行：首先將圖像在水平（豎直）方向與一維高斯函數進行卷積；而後將卷積後的結果在豎直（水平）方向使用相同的一維高斯函數獲得的模板進行卷積運算。

你也能夠用以下源碼進行實現：

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;
}

示例：

//實現直方圖的反向投影
#include "stdafx.h"
#include <math.h>
#include <random>
#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;

// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
	// 隨機數產生器
	std::random_device rd; //種子
	std::mt19937 gen(rd()); // 隨機數引擎

	auto cols = m.cols * m.channels();

	for (int i = 0; i < num; i++)
	{
		auto row = static_cast<int>(gen() % m.rows);
		auto col = static_cast<int>(gen() % cols);

		auto p = m.ptr<uchar>(row);
		p[col++] = 255;
		p[col++] = 255;
		p[col] = 255;
	}
}
void addGaussianNoise(Mat &m, int mu, int sigma)
{
	// 產生高斯分佈隨機數發生器
	std::random_device rd;
	std::mt19937 gen(rd());

	std::normal_distribution<> d(mu, sigma);

	auto rows = m.rows; // 行數
	auto cols = m.cols * m.channels(); // 列數

	for (int i = 0; i < rows; i++)
	{
		auto p = m.ptr<uchar>(i); // 取得行首指針
		for (int j = 0; j < cols; j++)
		{
			auto tmp = p[j] + d(gen);
			tmp = tmp > 255 ? 255 : tmp;
			tmp = tmp < 0 ? 0 : tmp;
			p[j] = tmp;
		}
	}
}

int main()
{
	Mat srcImage = imread("111.jpg");
	Mat src_salt, src_gaussian;
	src_salt = imread("111.jpg");
	src_gaussian = imread("111.jpg");
	if (!srcImage.data)
	{
		cout << "圖像打開失敗！" << endl;
		return -1;
	}
	addSaltNoise(src_salt, 100);
	addGaussianNoise(src_gaussian,100,1);

    //平滑濾波
	Mat aft_mean_filter_salt, aft_mean_filter_gaussian;
	blur(src_salt, aft_mean_filter_salt, Size(5, 5), Point(-1, -1));
	blur(src_gaussian, aft_mean_filter_gaussian, Size(5, 5), Point(-1, -1));


	//加強圖像對比度
	//Mat kernel = (Mat_<int>(3, 3) << 0, -1, 0, -1, 5, -1, 0, -1, 0);
	//filter2D(src_salt, src_salt, src_salt.depth(), kernel, Point(-1, -1));

	imshow("原圖", srcImage);
	imshow("添加椒鹽噪聲後", src_salt);
	imshow("添加高斯噪聲後", src_gaussian);
	imshow("通過平滑濾波後含椒鹽噪聲的的結果", aft_mean_filter_salt);
	imshow("通過平滑濾波後含高斯噪聲的的結果", aft_mean_filter_gaussian);
	waitKey();
	return 0;
}

四、非線性平滑空間濾波

一、中值濾波

　　中值濾波本質上是一種統計排序濾波器，對於原圖像中某點(i,j)，中值濾波以該點爲中心的領域內的全部像素的統計排序中值做爲(i,j)的響應，即它將每一像素點的灰度值設置爲該點某鄰域窗口內的全部像素點灰度值的中值，也就是將中心像素的值用全部像素值的中間值（不是平均值）替換。中值濾波不一樣於均值，是指排序隊列中位於中間位置的元素的值。

　　中值濾波對於某些類型的隨機噪聲具備很是理想的降噪能力，對於線性平滑濾波而言，在處理的像素領域以內包含噪聲點時，噪聲的存在總會或多或少的影響該點的像素值的計算(對於高斯平滑，影響程度同噪聲點到中點的距離成正比），但在中值濾波中，噪聲點則經常是直接忽略掉的；並且在同線性平滑濾波器相比，中值濾波在降噪的同時引發的模糊效應較低。

　　OpenCV提供了中值濾波的API函數以下：

void medianBlur(InputArray src,OutputArray dst,int ksize)

　　參數解釋：原圖像，目標圖像，模板尺寸。

示例：

//實現直方圖的反向投影
#include "stdafx.h"
#include <math.h>
#include <random>
#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;

// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
	// 隨機數產生器
	std::random_device rd; //種子
	std::mt19937 gen(rd()); // 隨機數引擎

	auto cols = m.cols * m.channels();

	for (int i = 0; i < num; i++)
	{
		auto row = static_cast<int>(gen() % m.rows);
		auto col = static_cast<int>(gen() % cols);

		auto p = m.ptr<uchar>(row);
		p[col++] = 255;
		p[col++] = 255;
		p[col] = 255;
	}
}
void addGaussianNoise(Mat &m, int mu, int sigma)
{
	// 產生高斯分佈隨機數發生器
	std::random_device rd;
	std::mt19937 gen(rd());

	std::normal_distribution<> d(mu, sigma);

	auto rows = m.rows; // 行數
	auto cols = m.cols * m.channels(); // 列數

	for (int i = 0; i < rows; i++)
	{
		auto p = m.ptr<uchar>(i); // 取得行首指針
		for (int j = 0; j < cols; j++)
		{
			auto tmp = p[j] + d(gen);
			tmp = tmp > 255 ? 255 : tmp;
			tmp = tmp < 0 ? 0 : tmp;
			p[j] = tmp;
		}
	}
}

int main()
{
	Mat srcImage = imread("111.jpg");
	Mat src_salt, src_gaussian;
	src_salt = imread("111.jpg");
	src_gaussian = imread("111.jpg");
	if (!srcImage.data)
	{
		cout << "圖像打開失敗！" << endl;
		return -1;
	}
	addSaltNoise(src_salt, 100);
	addGaussianNoise(src_gaussian,100,1);

    //平滑濾波
	Mat aft_median_filter_salt, aft_median_filter_gaussian;
	medianBlur(src_salt, aft_median_filter_salt, 3);
	medianBlur(src_gaussian, aft_median_filter_gaussian, 3);
	//加強圖像對比度
	//Mat kernel = (Mat_<int>(3, 3) << 0, -1, 0, -1, 5, -1, 0, -1, 0);
	//filter2D(src_salt, src_salt, src_salt.depth(), kernel, Point(-1, -1));

	imshow("原圖", srcImage);
	imshow("添加椒鹽噪聲後", src_salt);
	imshow("添加高斯噪聲後", src_gaussian);
	imshow("通過平滑濾波後含椒鹽噪聲的的結果", aft_median_filter_salt);
	imshow("通過平滑濾波後含高斯噪聲的的結果", aft_median_filter_gaussian);
	waitKey();
	return 0;
}

二、自適應中值濾波

　　中值濾波的效果依賴於濾波窗口的大小，太大會使圖像邊緣模糊，過小則去噪效果很差。由於噪聲點和邊緣點一樣是灰度變化較爲劇烈的像素，普通中值濾波在改變噪聲點灰度值的時候，會必定程度上改變邊緣像素的灰度值。可是噪聲點幾乎就是鄰近像素的極值，而邊緣每每不是，咱們能夠經過這個特性來限制中值濾波。

　　具體的改進方法以下：追行掃描圖像，當處理每個像素時，判斷該像素是不是濾波窗口所覆蓋下領域像素的極大或者極小值，若是是，則採用正常的中值濾波處理該像素，若是不是，則不予處理。這種方法可以很是有效的去除突發噪聲點，尤爲是椒鹽噪聲，並且幾乎不影響邊緣。這種方法又稱爲自適應中值濾波。

　　算法實現以下所示：

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 adaptiveProcess(im, row, col, kernelSize, maxSize); // 增大窗口尺寸，繼續A過程。
		else
			return med;
	}
}

　示例：

//實現直方圖的反向投影
#include "stdafx.h"
#include <math.h>
#include <random>
#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;
uchar adaptiveProcess(const Mat &im, int row, int col, int kernelSize, int maxSize);
// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
	// 隨機數產生器
	std::random_device rd; //種子
	std::mt19937 gen(rd()); // 隨機數引擎

	auto cols = m.cols * m.channels();

	for (int i = 0; i < num; i++)
	{
		auto row = static_cast<int>(gen() % m.rows);
		auto col = static_cast<int>(gen() % cols);

		auto p = m.ptr<uchar>(row);
		p[col++] = 255;
		p[col++] = 255;
		p[col] = 255;
	}
}
void addGaussianNoise(Mat &m, int mu, int sigma)
{
	// 產生高斯分佈隨機數發生器
	std::random_device rd;
	std::mt19937 gen(rd());

	std::normal_distribution<> d(mu, sigma);

	auto rows = m.rows; // 行數
	auto cols = m.cols * m.channels(); // 列數

	for (int i = 0; i < rows; i++)
	{
		auto p = m.ptr<uchar>(i); // 取得行首指針
		for (int j = 0; j < cols; j++)
		{
			auto tmp = p[j] + d(gen);
			tmp = tmp > 255 ? 255 : tmp;
			tmp = tmp < 0 ? 0 : tmp;
			p[j] = tmp;
		}
	}
}

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 adaptiveProcess(im, row, col, kernelSize, maxSize); // 增大窗口尺寸，繼續A過程。
		else
			return med;
	}
}

int main()
{
	Mat srcImage = imread("111.jpg",0);
	Mat src_salt, src_gaussian;
	src_salt = imread("111.jpg",0);
	src_gaussian = imread("111.jpg",0);
	if (!srcImage.data)
	{
		cout << "圖像打開失敗！" << endl;
		return -1;
	}
	addSaltNoise(src_salt, 100);
	addGaussianNoise(src_gaussian,100,1);

	int minSize = 3; // 濾波器窗口的起始尺寸
	int maxSize = 7; // 濾波器窗口的最大尺寸
	//平滑濾波
	Mat aft_admedian_filter_salt, aft_admedian_filter_gaussian;
	// 擴展圖像的邊界
	copyMakeBorder(src_salt, aft_admedian_filter_salt, maxSize / 2, maxSize / 2, maxSize / 2, maxSize / 2, BorderTypes::BORDER_REFLECT);
	for (int j = maxSize / 2; j < aft_admedian_filter_salt.rows - maxSize / 2; j++)
	{
		for (int i = maxSize / 2; i < aft_admedian_filter_salt.cols * aft_admedian_filter_salt.channels() - maxSize / 2; i++)
		{
			aft_admedian_filter_salt.at<uchar>(j, i) = adaptiveProcess(aft_admedian_filter_salt, j, i, minSize, maxSize);
		}
	}

	copyMakeBorder(src_salt, aft_admedian_filter_gaussian, maxSize / 2, maxSize / 2, maxSize / 2, maxSize / 2, BorderTypes::BORDER_REFLECT);
	for (int j = maxSize / 2; j < aft_admedian_filter_gaussian.rows - maxSize / 2; j++)
	{
		for (int i = maxSize / 2; i < aft_admedian_filter_gaussian.cols * aft_admedian_filter_gaussian.channels() - maxSize / 2; i++)
		{
			aft_admedian_filter_gaussian.at<uchar>(j, i) = adaptiveProcess(aft_admedian_filter_gaussian, j, i, minSize, maxSize);
		}
	}

	medianBlur(src_gaussian, aft_admedian_filter_gaussian, 3);
	//加強圖像對比度
	//Mat kernel = (Mat_<int>(3, 3) << 0, -1, 0, -1, 5, -1, 0, -1, 0);
	//filter2D(src_salt, src_salt, src_salt.depth(), kernel, Point(-1, -1));

	imshow("原圖", srcImage);
	imshow("添加椒鹽噪聲後", src_salt);
	imshow("添加高斯噪聲後", src_gaussian);
	imshow("通過平滑濾波後含椒鹽噪聲的的結果", aft_admedian_filter_salt);
	imshow("通過平滑濾波後含高斯噪聲的的結果", aft_admedian_filter_gaussian);
	waitKey();
	return 0;
}

三、雙邊濾波

　　雙邊濾波（Bilateral filter）是一種非線性的濾波方法，是結合圖像的空間鄰近度和像素值類似度的一種折衷處理，同時考慮空域信息和灰度類似性，達到保邊去噪的目的。雙邊濾波可以提供一種不會將邊緣平滑掉的方法，但做爲代價，須要更多的處理時間。與高斯濾波相似，雙邊濾波會依據每一個像素及其領域構造一個加權平均值，加權計算包括兩個部分，其中第一部分加權方式與高斯平滑中相同，第二部分也屬於高斯加權，但不是基於中心像素點與其餘像素點的空間距離之上的加權，而是基於其餘像素與中心像素的亮度差值的加權。能夠將雙邊濾波視爲高斯平滑，對類似的像素賦予較高的權重，不類似的像素賦予較小的權重，也可用於圖像分割之中。

　　雙邊濾波器之因此可以作到在平滑去噪的同時還可以很好的保存邊緣（Edge Preserve），是因爲其濾波器的核由兩個函數生成：

一個函數由像素歐式距離決定濾波器模板的係數
另外一個函數由像素的灰度差值決定濾波器的係數

　　其綜合了高斯濾波器（Gaussian Filter）和 $α$ $α$ $α$ $α$

雙邊濾波器使用二維高斯函數生成距離模板，使用一維高斯函數生成值域模板。
距離模板係數的生成公式以下：

$$
d(i, j, k, l)=\exp \left(-\frac{(i-k)^{2}+(j-l)^{2}}{2 \sigma_{d}^{2}}\right)
$$

其中,(k,l)爲模板窗口的中心座標;(i,j)爲模板窗口的其餘係數的座標；σ_d爲高斯函數的標準差。使用該公式生成的濾波器模板和高斯濾波器使用的模板是沒有區別的。 $σ_{d}$

值域模板係數的生成公式以下：

$$
r(i, j, k, l)=\exp \left(-\frac{\|f(i, j)-f(k, l)\|^{2}}{2 \sigma_{r}^{2}}\right)
$$

其中，函數f(x,y)表示要處理的圖像，f(x,y)表示圖像在點(x,y)處的像素值；(k,l)爲模板窗口的中心座標；(i,j)爲模板窗口的其餘係數的座標；σ_r爲高斯函數的標準差。

將上述兩個模板相乘就獲得了雙邊濾波器的模板

$$
w(i, j, k, l)=d(i, j, k, l) * r(i, j, k, l)=\exp \left(-\frac{(i-k)^{2}+(j-l)^{2}}{2 \sigma_{d}^{2}}-\frac{\|f(i, j)-f(k, l)\|^{2}}{2 \sigma_{r}^{2}}\right)
$$

其中定義域濾波和值域濾波能夠用下圖形象表示：

　　代碼實現以下：

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提供了對應雙邊濾波的API接口

void bilateralFilter(InputArray src,OutputArray dst,int d,double sigmaColor,double sigmaSpace,int borderType=BORDER_DEFAULT)

第一個參數，InputArray類型的src，輸入圖像，即源圖像，須要爲8爲或者浮點型單通道、三通道的圖像；
第二個參數，OutputArray類型的dst，即目標圖像，須要和源圖像有同樣的尺寸和類型；
第三個參數，int類型的d，表示在過濾過程當中每一個像素鄰域的直徑。若是這個參數被設置爲負值，那麼OpenCV會從第五個參數sigmaSpace來計算出它；
第四個參數，double類型的sigmaColor，顏色空間濾波器的sigma值，這個參數的值越大，就代表該像素鄰域內有越寬廣的顏色會被混合到一塊兒，產生較大的半相等顏色區域；
第五個參數，double類型的 sigmaSpace，座標空間中的sigma值，座標空間的標註方差。它的數值越大，意味着越遠的像素會相互影響，從而使更大區域中足夠類似的顏色獲取相同的顏色，當d>0時，d制定了鄰域大小且與sigmaSpace無關不然，d正比於sigmaSpace；
第六個參數，int類型的borderType，用於推斷圖像外部像素的某種邊界模式，有默認值BORDER_DEFAULT。

示例：

//實現直方圖的反向投影
#include "stdafx.h"
#include <math.h>
#include <random>
#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;
// 添加椒鹽噪聲
void addSaltNoise(Mat &m, int num)
{
	// 隨機數產生器
	std::random_device rd; //種子
	std::mt19937 gen(rd()); // 隨機數引擎

	auto cols = m.cols * m.channels();

	for (int i = 0; i < num; i++)
	{
		auto row = static_cast<int>(gen() % m.rows);
		auto col = static_cast<int>(gen() % cols);

		auto p = m.ptr<uchar>(row);
		p[col++] = 255;
		p[col++] = 255;
		p[col] = 255;
	}
}
void addGaussianNoise(Mat &m, int mu, int sigma)
{
	// 產生高斯分佈隨機數發生器
	std::random_device rd;
	std::mt19937 gen(rd());

	std::normal_distribution<> d(mu, sigma);

	auto rows = m.rows; // 行數
	auto cols = m.cols * m.channels(); // 列數

	for (int i = 0; i < rows; i++)
	{
		auto p = m.ptr<uchar>(i); // 取得行首指針
		for (int j = 0; j < cols; j++)
		{
			auto tmp = p[j] + d(gen);
			tmp = tmp > 255 ? 255 : tmp;
			tmp = tmp < 0 ? 0 : tmp;
			p[j] = tmp;
		}
	}
}
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;
			}
		}
	}
}
int main()
{
	Mat srcImage = imread("111.jpg");
	Mat src_salt, src_gaussian;
	src_salt = imread("111.jpg");
	src_gaussian = imread("111.jpg");
	if (!srcImage.data)
	{
		cout << "圖像打開失敗！" << endl;
		return -1;
	}
	addSaltNoise(src_salt, 100);
	addGaussianNoise(src_gaussian,100,1);

	//平滑濾波
	Mat aft_bilateral_filter_salt, aft_bilateral_filter_gaussian, srcImage1;

	bilateralFilter(srcImage, srcImage1, 5, 100, 3);
	bilateralFilter(src_salt, aft_bilateral_filter_salt, 5, 100, 3);
	bilateralFilter(src_gaussian, aft_bilateral_filter_gaussian, 5, 100, 3);
	//加強圖像對比度
	Mat kernel = (Mat_<int>(3, 3) << 0, -1, 0, -1, 5, -1, 0, -1, 0);
	filter2D(srcImage1, srcImage1, srcImage1.depth(), kernel, Point(-1, -1));

	imshow("原圖", srcImage);
	imshow("雙邊濾波處理事後的原圖", srcImage1);
	imshow("添加椒鹽噪聲後", src_salt);
	imshow("添加高斯噪聲後", src_gaussian);
	imshow("通過平滑濾波後含椒鹽噪聲的的結果", aft_bilateral_filter_salt);
	imshow("通過平滑濾波後含高斯噪聲的的結果", aft_bilateral_filter_gaussian);
	waitKey();
	return 0;
}

　　程序運行結果以下，雙邊濾波很好的對原圖像進行了優化，可是對應椒鹽噪聲之類的高頻噪聲沒有濾波效果。