基於opencv的理想低通濾波器和巴特沃斯低通濾波器
阿新 • • 發佈:2019-02-17
首先看個圖瞭解下什麼是理想低通濾波器公式和圖是轉自Rolin的專欄
低通濾波器
1.理想的低通濾波器
其中,D0表示通帶的半徑。D(u,v)的計算方式也就是兩點間的距離,很簡單就能得到。
使用低通濾波器所得到的結果如下所示。低通濾波器濾除了高頻成分,所以使得影象模糊。由於理想低通濾波器的過度特性過於急峻,所以會產生了振鈴現象。
2.巴特沃斯低通濾波器
同樣的,D0表示通帶的半徑,n表示的是巴特沃斯濾波器的次數。隨著次數的增加,振鈴現象會越來越明顯。
void ideal_Low_Pass_Filter(Mat src){ Mat img; cvtColor(src, img, CV_BGR2GRAY); imshow("img",img); //調整影象加速傅立葉變換 int M = getOptimalDFTSize(img.rows); int N = getOptimalDFTSize(img.cols); Mat padded; copyMakeBorder(img, padded, 0, M - img.rows, 0, N - img.cols, BORDER_CONSTANT, Scalar::all(0)); //記錄傅立葉變換的實部和虛部 Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) }; Mat complexImg; merge(planes, 2, complexImg); //進行傅立葉變換 dft(complexImg, complexImg); //獲取影象 Mat mag = complexImg; mag = mag(Rect(0, 0, mag.cols & -2, mag.rows & -2));//這裡為什麼&上-2具體檢視opencv文件 //其實是為了把行和列變成偶數 -2的二進位制是11111111.......10 最後一位是0 //獲取中心點座標 int cx = mag.cols / 2; int cy = mag.rows / 2; //調整頻域 Mat tmp; Mat q0(mag, Rect(0, 0, cx, cy)); Mat q1(mag, Rect(cx, 0, cx, cy)); Mat q2(mag, Rect(0, cy, cx, cy)); Mat q3(mag, Rect(cx, cy, cx, cy)); q0.copyTo(tmp); q3.copyTo(q0); tmp.copyTo(q3); q1.copyTo(tmp); q2.copyTo(q1); tmp.copyTo(q2); //Do為自己設定的閥值具體看公式 double D0 = 60; //處理按公式保留中心部分 for (int y = 0; y < mag.rows; y++){ double* data = mag.ptr<double>(y); for (int x = 0; x < mag.cols; x++){ double d = sqrt(pow((y - cy),2) + pow((x - cx),2)); if (d <= D0){ } else{ data[x] = 0; } } } //再調整頻域 q0.copyTo(tmp); q3.copyTo(q0); tmp.copyTo(q3); q1.copyTo(tmp); q2.copyTo(q1); tmp.copyTo(q2); //逆變換 Mat invDFT, invDFTcvt; idft(mag, invDFT, DFT_SCALE | DFT_REAL_OUTPUT); // Applying IDFT invDFT.convertTo(invDFTcvt, CV_8U); imshow("理想低通濾波器", invDFTcvt); } void Butterworth_Low_Paass_Filter(Mat src){ int n = 1;//表示巴特沃斯濾波器的次數 //H = 1 / (1+(D/D0)^2n) Mat img; cvtColor(src, img, CV_BGR2GRAY); imshow("img", img); //調整影象加速傅立葉變換 int M = getOptimalDFTSize(img.rows); int N = getOptimalDFTSize(img.cols); Mat padded; copyMakeBorder(img, padded, 0, M - img.rows, 0, N - img.cols, BORDER_CONSTANT, Scalar::all(0)); Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) }; Mat complexImg; merge(planes, 2, complexImg); dft(complexImg, complexImg); Mat mag = complexImg; mag = mag(Rect(0, 0, mag.cols & -2, mag.rows & -2)); int cx = mag.cols / 2; int cy = mag.rows / 2; Mat tmp; Mat q0(mag, Rect(0, 0, cx, cy)); Mat q1(mag, Rect(cx, 0, cx, cy)); Mat q2(mag, Rect(0, cy, cx, cy)); Mat q3(mag, Rect(cx, cy, cx, cy)); q0.copyTo(tmp); q3.copyTo(q0); tmp.copyTo(q3); q1.copyTo(tmp); q2.copyTo(q1); tmp.copyTo(q2); double D0 = 100; for (int y = 0; y < mag.rows; y++){ double* data = mag.ptr<double>(y); for (int x = 0; x < mag.cols; x++){ //cout << data[x] << endl; double d = sqrt(pow((y - cy), 2) + pow((x - cx), 2)); //cout << d << endl; double h = 1.0 / (1 + pow(d / D0, 2 * n)); if (h <= 0.5){ data[x] = 0; } else{ //data[x] = data[x]*0.5; //cout << h << endl; } //cout << data[x] << endl; } } q0.copyTo(tmp); q3.copyTo(q0); tmp.copyTo(q3); q1.copyTo(tmp); q2.copyTo(q1); tmp.copyTo(q2); //逆變換 Mat invDFT, invDFTcvt; idft(complexImg, invDFT, DFT_SCALE | DFT_REAL_OUTPUT); // Applying IDFT invDFT.convertTo(invDFTcvt, CV_8U); imshow("巴特沃斯低通濾波器", invDFTcvt); }