opencv-9-影像雜訊以及評估指標 PSNR 與SSIM

2020 年 5 月 1 日
筆記

開始之前

我們在將 opencv 的影像顯示在了 qt 的label 上, 我們能夠將圖顯示在label 上, 用於顯示我們的演算法,
我們在 opencv 上一篇文章中介紹了 opencv 的核操作, 我們這裡就要進入一個很重要的章節了,影像濾波操作, 也是影像核操作應用的一個很重要的章節,

那我們就從降噪的角度完整的講一下, 並通過 opencv 核的方式進行影像演算法操作, 【技術綜述】一文道盡傳統影像降噪方法這篇文章寫的還算比較完整, 也是傳統的演算法的一個綜述過程,

前言

數字成像過程中由於電雜訊以及其他因素, 導致我們獲取到的影像存在雜訊,雜訊出現在輸入部分, 在後續的每個步驟都會受到影響, 所以在數字影像處理的前面必須要進行的一個步驟就是影像降噪

每個做訊號處理的都會接觸到一類問題 , 訊號降噪, 讓人最頭疼的一門課,真是感謝老師給過, 但是後面自己用到的時候反而感覺真的好用, 原來是這樣, 然後就慢慢學會了怎麼使用吧..(感覺還是弱雞)
知乎可復現的影像降噪演算法總結這篇文章列出了一個能夠復現的影像降噪操作演算法列表, 近年來實現了的演算法可以見reproducible-image-denoising-state-of-the-art, 之後使用相應的文章進行演算法實現吧.( 又立了一個 flag )

影像降噪理論基礎

影像降噪主要的目的是在進行去除影像雜訊的同時保留儘可能多的主要特徵, 對於人眼來說, 區分雜訊還算比較容易, 但是對於電腦來說,輸入的都是數據, 我怎麼區分哪個是雜訊, 哪個不是雜訊呢, 這裡就要引入雜訊的理論基礎了

影像雜訊的產生

我們在之前的章節介紹了影像的程式系統, 實際上在成像過程中可能由於點雜訊, 量化過程等造成雜訊,
實際上的雜訊主要分為三種:

加性雜訊: 與輸入無關, \(f(x,y) = g(x,y) +n(x,y)\)
乘性雜訊: 與輸入訊號有關, \(f(x,y) = g(x,y) + n(x,y) \cdot g(x,y)\)
量化雜訊: 與輸入無關, 在影像量化過程雜訊的量化誤差導致的雜訊,

實際上後兩種很難解決, 目前處理的都是以加性雜訊為主, 屬於隨機的雜訊訊號, 根據統計學的觀點, 雜訊在無限長時間窗內的雜訊和為0, 在第一類中的 \(n(x,y)\) 隨著時間存在正負訊號的不確定變化.

上圖所示虛線代表真實訊號，紅藍線表示的就是隨機雜訊訊號，所有的隨機雜訊訊號求和後結果為0。

這裡關於雜訊的說明可以參考影像去噪演算法簡介

雜訊在理論上可以定義為「不可預測，只能用概率統計方法來認識的隨機誤差」。因此將影像雜訊看成是多維隨機過程是合適的，因而描述雜訊的方法完全可以借用隨機過程的描述，即用其概率分布函數和概率密度分布函數。但在很多情況下，這樣的描述方法是很複雜的，甚至是不可能的。而實際應用往往也不必要。通常是用其數字特徵，即均值方差，相關函數等。因為這些數字特徵都可以從某些方面反映出雜訊的特徵。

我認為影像雜訊的成因分類與常見影像去噪演算法簡介這篇文章關於雜訊的分類部分講的還比較細, 可以參考

影像雜訊的模型

由於我們認為雜訊在時間尺度的隨機性, 但是我們可以使用雜訊的概率分布與概率密度函數進行描述, 那麼我們就能將雜訊根據其分布特點進行分類,
我們稍微介紹一下常見的雜訊模型吧

雜訊模型主要可以分為:

高斯雜訊，高斯雜訊模型經常被用於實踐中。
脈衝雜訊（椒鹽雜訊），影像上一個個點，也可稱為散粒和尖峰雜訊。
伽馬雜訊
瑞利雜訊
指數分布雜訊
均勻分布雜訊

這裡能查到的資料很多, 可以看我們的參考部分, 內容都一樣, 再寫只是浪費時間和精力, 有興趣的可以自己翻閱

影像降噪操作

其實吧, 我就不應該講那麼多, 直接開始影像處理部分就行了, 為了開始進行影像處理, 我們要先進行一點小工作, 我們要按造以下步驟進行降噪演算法的比較,

選擇標準影像— lena.png
添加雜訊
量化雜訊
降噪操作
量化結果值
比較結果

在我們進行演算法比對之前, 我們選擇的是 lena 的影像, 加入隨機雜訊, 然後計算出來一個雜訊的比例, 進行降噪操作, 再次計算以下雜訊參數, 看下效果值.

如果是進行演算法比較的時候, 最好選擇現有的降噪的數據集進行比較, 比如, Kodak , BSD

雜訊添加

我們認為雜訊是隨機的, 我們生成隨機數加在原始影像上便能夠得到雜訊影像, opencv 沒有提供相應的實現, 但是知道原理了, 寫起來都比較簡單, 我比較喜歡
影像處理基礎(1)：雜訊的添加和過濾使用的方法, 他使用的是梅森旋轉演算法來實現的偽隨機演算法,

其實吧這裡我也不懂, 但是隨機數能用就行了, 我又不是數學家, 然後看到了談談梅森旋轉：演算法及其爆破

這裡就不重複造輪子了, 直接複製他給出的程式碼就好,

// 添加椒鹽雜訊 // 生成 隨機 num 個 白點
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;
		}
	}
}

雜訊量化方法

這裡其實涉及到影像品質評估的領域,可以參考影像品質評價概述（評估指標、傳統檢測方法）介紹的方法, 存在太多的計算方式,

我們必須選擇一個量化雜訊的方式進行影像品質的評估, 一般進行雜訊評估手段就是雜訊比(Signal to Noise Ratio,SNR),峰值信噪比(Peak Signal to Noise Ratio, PSNR) , 均方差值(Mean Square Error, MSE), 結構相似性(Structural SIMilarity, SSIM),

我們一個一個來看, 均方差值是用於比較兩幅影像 K, I 的均方差值

\[MSE=\frac{1}{mn}\sum_{i=0}^{n-1}\sum_{j=0}^{m-1}\parallel K(i,j)-I(i,j)\parallel^{2}
\]

峰值信噪比PSNR衡量影像失真或是雜訊水平的客觀標準。2個影像之間PSNR值越大，則越相似。普遍基準為30dB，30dB以下的影像劣化較為明顯。定義為，

\[PSNR=10log_{10}(\frac{MAX^{2}}{MSE})
\]

其中\(MAX^2\) 為圖片可能的最大像素值。如果每個像素都由 8 位二進位來表示，那麼就為 255。

SNR用於描述訊號與雜訊的比值

\[SNR (dB)=10 log_{10}
\left[
\frac{\sum_{x=0}^{m-1} \sum_{y=0}^{n-1}(f(x, y))^{2}}{\sum_{x=0}^{m-1} \sum_{y=0}^{n-1}(f(x, y)-\hat{f}(x, y))^{2}}\right]
\]

SSIM 描述兩個影像的相似性, 通過三個進行比較, 亮度,對比度和結構, 參考影像品質評價指標之 PSNR 和 SSIM

\[l(x, y)=\frac{2 \mu_{x} \mu_{y}+c_{1}}{\mu_{x}^{2}+\mu_{y}^{2}+c_{1}} c(x, y)=\frac{2 \sigma_{x} \sigma_{y}+c_{2}}{\sigma_{x}^{2}+\sigma_{y}^{2}+c_{2}} s(x, y)=\frac{\sigma_{x y}+c_{3}}{\sigma_{x} \sigma_{y}+c_{3}}
\]

\[\operatorname{SSIM}(x, y)=\frac{\left(2 \mu_{x} \mu_{y}+c_{1}\right)\left(2 \sigma_{x y}+c_{2}\right)}{\left(\mu_{x}^{2}+\mu_{y}^{2}+c_{1}\right)\left(\sigma_{x}^{2}+\sigma_{y}^{2}+c_{2}\right)}
\]

一般取\(c_3 = \frac{c_2}{2}\)。
\(u_x\) 為 \(x\) 的均值
\(u_y\) 為 \(y\) 的均值
\(\sigma_x^2\) 為\(x\) 的方差
\(\sigma_y^2\) 為\(y\) 的方差
\(\sigma_{xy}\) 為\(x\) 和\(y\) 的協方差
\(c_1 = (k_1 L)^2, c_2=(k_2 L)^2\) 為兩個常數，避免除零
\(L\) 為像素值的範圍，\((0,255)\)
\(k_1 = 0.01, k_2 = 0.03\) 為默認值
默認參數\(\alpha = 1, \beta = 1, \gamma = 1\)

opencv 計算 PSNR 和 SSIM

本來不想寫這麼多的, 但是 opencv 給出了一個常式Similarity check (PNSR and SSIM) on the GPU, 提供了計算的方法, 自己不用去寫了, 豈不是很爽, 所以上面就詳細介紹了各個方法的使用.
官方給出了普通版本以及 GPU 加速的版本, 我們暫時只使用基礎的版本就好,
PSNR返回一個浮點數，如果兩個輸入在30到50之間相似（越高越好）。
SSIM返回影像的MSSIM。這也是一個介於零和一之間的浮點數（越高越好），但是每個通道都有一個浮點數。因此，我們返回一個Scalar OpenCV數據結構：

double getPSNR(const Mat& I1, const Mat& I2)
{
    Mat s1;
    absdiff(I1, I2, s1);       // |I1 - I2|
    s1.convertTo(s1, CV_32F);  // cannot make a square on 8 bits
    s1 = s1.mul(s1);           // |I1 - I2|^2
    Scalar s = sum(s1);         // sum elements per channel
    double sse = s.val[0] + s.val[1] + s.val[2]; // sum channels
    if( sse <= 1e-10) // for small values return zero
        return 0;
    else
    {
        double  mse =sse /(double)(I1.channels() * I1.total());
        double psnr = 10.0*log10((255*255)/mse);
        return psnr;
    }
}
Scalar getMSSIM( const Mat& i1, const Mat& i2)
{
    const double C1 = 6.5025, C2 = 58.5225;
    /***************************** INITS **********************************/
    int d     = CV_32F;
    Mat I1, I2;
    i1.convertTo(I1, d);           // cannot calculate on one byte large values
    i2.convertTo(I2, d);
    Mat I2_2   = I2.mul(I2);        // I2^2
    Mat I1_2   = I1.mul(I1);        // I1^2
    Mat I1_I2  = I1.mul(I2);        // I1 * I2
    /*************************** END INITS **********************************/
    Mat mu1, mu2;   // PRELIMINARY COMPUTING
    GaussianBlur(I1, mu1, Size(11, 11), 1.5);
    GaussianBlur(I2, mu2, Size(11, 11), 1.5);
    Mat mu1_2   =   mu1.mul(mu1);
    Mat mu2_2   =   mu2.mul(mu2);
    Mat mu1_mu2 =   mu1.mul(mu2);
    Mat sigma1_2, sigma2_2, sigma12;
    GaussianBlur(I1_2, sigma1_2, Size(11, 11), 1.5);
    sigma1_2 -= mu1_2;
    GaussianBlur(I2_2, sigma2_2, Size(11, 11), 1.5);
    sigma2_2 -= mu2_2;
    GaussianBlur(I1_I2, sigma12, Size(11, 11), 1.5);
    sigma12 -= mu1_mu2;
    Mat t1, t2, t3;
    t1 = 2 * mu1_mu2 + C1;
    t2 = 2 * sigma12 + C2;
    t3 = t1.mul(t2);              // t3 = ((2*mu1_mu2 + C1).*(2*sigma12 + C2))
    t1 = mu1_2 + mu2_2 + C1;
    t2 = sigma1_2 + sigma2_2 + C2;
    t1 = t1.mul(t2);               // t1 =((mu1_2 + mu2_2 + C1).*(sigma1_2 + sigma2_2 + C2))
    Mat ssim_map;
    divide(t3, t1, ssim_map);      // ssim_map =  t3./t1;
    Scalar mssim = mean( ssim_map ); // mssim = average of ssim map
    return mssim;
}

演算法雜訊數據

我們完成了雜訊添加以及雜訊的量化, 我們來試一下, 給影像隨機添加一定的雜訊, 然後看下相應的參數變化情況對比來看就好

椒鹽雜訊測試

我們先來測試椒鹽雜訊分別計算沒有雜訊的圖, 以及添加了 1000個和10000個雜訊的數據結果, 並將後面兩個顯示出來


void MainWindow::testFunc1(void)
{
    // 添加椒鹽雜訊 並計算 PSNR和 SSIM
    cv::Mat salt_img;
    double psnr = 0;
    cv::Scalar mssim;

    QString res_temp = "Salt-%1 : psnr:%2, mssim: B:%3 G:%4 R:%5 ";
    QString res_str;

    // 計算三組影像的參數 0, 1000, 10000

    // 複製原始影像, 添加雜訊, 計算 psnr和ssim  顯示在 ui上
    salt_img = gSrcImg.clone();
    addSaltNoise(salt_img,0);

    psnr = getPSNR(gSrcImg, salt_img);
    mssim = getMSSIM(gSrcImg,salt_img);
    res_str = res_temp.arg(0)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);

    salt_img = gSrcImg.clone();
    addSaltNoise(salt_img,1000);

    psnr = getPSNR(gSrcImg, salt_img);
    mssim = getMSSIM(gSrcImg,salt_img);
    res_str = res_temp.arg(1000)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);

    // 左側顯示 1000 雜訊 右側顯示 10000 雜訊
    ShowMatOnQtLabel(salt_img,ui->lb_src);

    salt_img = gSrcImg.clone();
    addSaltNoise(salt_img,10000);

    psnr = getPSNR(gSrcImg, salt_img);
    mssim = getMSSIM(gSrcImg,salt_img);
    res_str = res_temp.arg(10000)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);

    ShowMatOnQtLabel(salt_img,ui->lb_dst);
}

我們可以直接計算得到椒鹽雜訊 psnr 和 ssim 都是越大越好的, 可以明顯的看到影像品質退化

Salt-0 : psnr:0, mssim: B:1 G:1 R:1 
Salt-1000 : psnr:27.7528, mssim: B:0.865341 G:0.870555 R:0.914122 
Salt-10000 : psnr:17.8062, mssim: B:0.311999 G:0.327485 R:0.493874

高斯雜訊測試

高斯雜訊我們測試了四組分別使用參數(0,1) (0,10)(10,1)(10,10) 作為高斯參數, 最終得到後面的圖, 然後計算得到的結果, 我們做的結果比較簡單, 可以參考數字影像處理——添加高斯雜訊&椒鹽雜訊, 給出了很多的圖, 可以參考學


void MainWindow::testFunc2(void)
{
    // 添加高斯雜訊 並計算 PSNR和 SSIM
    cv::Mat guass_img;
    double psnr = 0;
    cv::Scalar mssim;

    QString res_temp = "gauss-%1- %2 : psnr:%3, mssim: B:%4 G:%5 R:%6 ";
    QString res_str;

    // 計算三組影像的參數 (0,1) (0,10), (10,1), (10,10)

    // 複製原始影像, 添加雜訊, 計算 psnr和ssim  顯示在 ui上
    guass_img = gSrcImg.clone();
    addGaussianNoise(guass_img,0,1);

    psnr = getPSNR(gSrcImg, guass_img);
    mssim = getMSSIM(gSrcImg,guass_img);
    res_str = res_temp.arg(0)
                        .arg(1)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);
    guass_img = gSrcImg.clone();
    addGaussianNoise(guass_img,0,10);

    psnr = getPSNR(gSrcImg, guass_img);
    mssim = getMSSIM(gSrcImg,guass_img);
    res_str = res_temp.arg(0)
                        .arg(10)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);

    guass_img = gSrcImg.clone();
    addGaussianNoise(guass_img,10,1);

    psnr = getPSNR(gSrcImg, guass_img);
    mssim = getMSSIM(gSrcImg,guass_img);
    res_str = res_temp.arg(10)
                        .arg(1)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);

    guass_img = gSrcImg.clone();
    addGaussianNoise(guass_img,10,10);

    psnr = getPSNR(gSrcImg, guass_img);
    mssim = getMSSIM(gSrcImg,guass_img);
    res_str = res_temp.arg(10)
                        .arg(10)
                        .arg(psnr)
                        .arg(mssim.val[0])
                        .arg(mssim.val[1])
                        .arg(mssim.val[2]);
    ui->pt_log->appendPlainText(res_str);
}

gauss-0- 1 : psnr:46.8791, mssim: B:0.991811 G:0.991622 R:0.992751 
gauss-0- 10 : psnr:28.1229, mssim: B:0.614219 G:0.608773 R:0.648285 
gauss-10- 1 : psnr:28.5293, mssim: B:0.978448 G:0.980308 R:0.987926 
gauss-10- 10 : psnr:25.3511, mssim: B:0.605665 G:0.600491 R:0.646768

總結

原本想把濾波一起做了的, 但是越寫越, 就不做太多的處理了, 我們算是介紹了雜訊的來源, 雜訊的模型, 以及個雜訊的量化方式,
然後介紹了影像添加雜訊的方法我們分別給影像添加椒鹽雜訊與高斯雜訊, 然後分別量化了雜訊的結果值, 進行對比展示,

示例的圖不是很多, 程式是在程式碼庫裡面的, 可以直接去自己實現, 然後進行進行更多圖的展示

參考

《高斯雜訊_百度百科》. 見於 2020年4月30日. //baike.baidu.com/item/高斯雜訊.
知乎專欄. 《【技術綜述】一文道盡傳統影像降噪方法》. 見於 2020年4月29日. //zhuanlan.zhihu.com/p/51403693.
知乎專欄. 《可復現的影像降噪演算法總結》. 見於 2020年4月29日. //zhuanlan.zhihu.com/p/32502816.
《梅森旋轉演算法》. 收入維基百科，自由的百科全書, 2019年11月4日. //zh.wikipedia.org/w/index.php?title=梅森旋轉演算法&oldid=56745942.
《實現灰度影像峰值信噪比計算_人工智慧_松子茶的專欄-CSDN部落格》. 見於 2020年4月30日. //blog.csdn.net/songzitea/article/details/17529445.
《數字影像處理-雜訊 – Mohanson》. 見於 2020年4月30日. //accu.cc/content/pil/noise/.
《影像處理基礎(1)：雜訊的添加和過濾 – Brook_icv – 部落格園》. 見於 2020年4月30日. //www.cnblogs.com/wangguchangqing/p/6372025.html.
《影像處理PSNR及其計算（OpenCV和matlab實現）_人工智慧_無機器不學習-加大碼的分享-CSDN部落格》. 見於 2020年4月30日. //blog.csdn.net/laoxuan2011/article/details/51519062.
《影像的 SNR 和 PSNR 的計算 – rldts – 部落格園》. 見於 2020年4月30日. //www.cnblogs.com/qrlozte/p/5340216.html.
《影像去噪演算法簡介 – InfantSorrow – 部落格園》. 見於 2020年4月29日. //www.cnblogs.com/CCBB/archive/2011/01/06/1929033.html.
《影像雜訊的成因分類與常見影像去噪演算法簡介_Java_qq_27606639的部落格-CSDN部落格》. 見於 2020年4月30日. //blog.csdn.net/qq_27606639/article/details/80912071.
《影像品質評估指標 SSIM / PSNR / MSE_人工智慧_兔角與禪-CSDN部落格》. 見於 2020年4月30日. //blog.csdn.net/edogawachia/article/details/78756680.
《影像品質評價概述（評估指標、傳統檢測方法）_人工智慧_qq_23304241的部落格-CSDN部落格》. 見於 2020年4月30日. //blog.csdn.net/qq_23304241/article/details/80953613.
《影像降噪》. 收入維基百科，自由的百科全書, 2018年9月20日. //zh.wikipedia.org/w/index.php?title=影像降噪&oldid=51354600.