在OpenCV環境下寫的灰度影象二維傅立葉換,幅值計算,頻譜平移和將數值歸一化到0到255區間的四個函式
阿新 • • 發佈:2019-01-10
影象處理開發資料、影象處理開發需求、影象處理接私活掙零花錢,可以搜尋公眾號"qxsf321",並關注!
影象處理開發資料、影象處理開發需求、影象處理接私活掙零花錢,可以搜尋公眾號"qxsf321",並關注!
影象處理開發資料、影象處理開發需求、影象處理接私活掙零花錢,可以搜尋公眾號"qxsf321",並關注!
灰度影象的二維傅立葉變換(cv_gray_fft2[我用C語言寫成]和DFT[用C++寫成]),二維傅立葉變換結果的幅值計算(cv_abs),頻譜平移(cv_gray_fft2shift)【頻譜平移的作用詳見我的博文http://blog.csdn.net/wenhao_ir/article/details/51689960
傅立葉變換的意義這裡就不多說了,大家在高等數學、訊號與系統中應該都學過,不清楚的可以去查閱相關資料!這裡我簡單說兩點!①影象的頻率是表徵影象中灰度變化劇烈程式的指標;②影象低頻分量表示背景和緩慢變化區域,高頻分量表示影象的邊緣、細節及相關噪聲。
原始碼如下:
#include <opencv2/opencv.hpp> #include <opencv2/legacy/compat.hpp> #include <fstream> using namespace std; #pragma comment(linker, "/subsystem:\"windows\" /entry:\"mainCRTStartup\"") void cv_gray_fft2(IplImage *src, IplImage *dst) //注意:此函式僅適用灰度影象,src因為是灰度影象,所以通道數要求為1,dst的屬性為:IPL_DEPTH_64F,2(兩通道) { IplImage *image_Re = 0, *image_Im = 0, *Fourier = 0; image_Re = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1); //實部 image_Im = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1); //虛部 Fourier = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 2);//2 channels,分別儲存image_Re、image_Im // Real part conversion from u8 to 64f (double) cvConvertScale(src, image_Re); // Imaginary part (zeros) cvZero(image_Im); // Join real and imaginary parts and stock them in Fourier image cvMerge(image_Re, image_Im, 0, 0, Fourier); // Application of the forward Fourier transform cvDFT(Fourier, dst, CV_DXT_FORWARD); cvReleaseImage(&image_Re); cvReleaseImage(&image_Im); cvReleaseImage(&Fourier); } void cv_abs(IplImage *src,IplImage *fft2_out, IplImage *dst)//src就是灰度影象原始資料,fft2_out要求是cv_gray_fft2的dst資料型別,即"64F, 2";dst為"64F, 1" { IplImage *image_Re = 0, *image_Im = 0; image_Re = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1); //實部 image_Im = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1); //虛部 cvSplit(fft2_out,image_Re,image_Im,0,0); cvPow(image_Re,image_Re,2); cvPow(image_Im,image_Im,2); cvAdd(image_Re,image_Im,image_Re,NULL); cvPow(image_Re,dst,0.5); cvReleaseImage(&image_Re); cvReleaseImage(&image_Im); } void cv_gray_fft2shift(IplImage *src)//引數是"64F, 2" { int nRow, nCol, i, j, cy, cx; double tmp13, tmp24; //Rearrange the quadrants of Fourier image so that the origin is at the image center nRow = src->height; nCol = src->width; cy = nRow/2; // image center cx = nCol/2; //CV_IMAGE_ELEM為OpenCV定義的巨集,用來讀取影象的畫素值,這一部分就是進行中心變換 for( j = 0; j < cy; j++ ){ for( i = 0; i < cx; i++ ){ //中心化,將整體份成四塊進行對角交換 tmp13 = CV_IMAGE_ELEM( src, double, j, i); CV_IMAGE_ELEM( src, double, j, i) = CV_IMAGE_ELEM( src, double, j+cy, i+cx); CV_IMAGE_ELEM( src, double, j+cy, i+cx) = tmp13; tmp24 = CV_IMAGE_ELEM( src, double, j, i+cx); CV_IMAGE_ELEM( src, double, j, i+cx) = CV_IMAGE_ELEM( src, double, j+cy, i); CV_IMAGE_ELEM( src, double, j+cy, i) = tmp24; } } } void cv_range_0to255(IplImage *src, IplImage *dst)//兩個引數都是"64F, 2" { double m,M; double scale; double shift; cvMinMaxLoc(src,&m,&M,NULL,NULL); scale = 255/(M - m); shift = -m * scale; cvConvertScale(src,dst,scale,shift); } int main() { int i=0;//迴圈變數 //從檔案中載入原圖 IplImage *pSrcImage = cvLoadImage("coins.png", CV_LOAD_IMAGE_UNCHANGED); //建立輸出的影象 IplImage *pOutImage = cvCreateImage(cvGetSize(pSrcImage), IPL_DEPTH_64F,1); //建立中間結果影象 IplImage *cv_gray_fft2_out = cvCreateImage(cvGetSize(pSrcImage), IPL_DEPTH_64F,2);//注意這裡是二通道! IplImage *cv_abs_out = cvCreateImage(cvGetSize(pSrcImage), IPL_DEPTH_64F, 1); cv_gray_fft2(pSrcImage,cv_gray_fft2_out); double watch_cv_gray_fft2_out_Re[100];//利用除錯時的區域性變數視窗觀察cvDFTOut實部的前100個值,看與MATLAB執行的結果否相同 double watch_cv_gray_fft2_out_Im[100];//利用除錯時的區域性變數視窗觀察cvDFTOut虛部的前100個值,看與MATLAB執行的結果否相同 for(i=0;i<10;i++)//一般情況下觀察10個值就夠了 { watch_cv_gray_fft2_out_Re[i]=cvGet2D(cv_gray_fft2_out,0,i).val[0]; } for(i=0;i<10;i++)//一般情況下觀察10個值就夠了 { watch_cv_gray_fft2_out_Im[i]=cvGet2D(cv_gray_fft2_out,0,i).val[1]; } cv_abs(pSrcImage,cv_gray_fft2_out,cv_abs_out); double watch_cv_abs_out[100];//利用除錯時的區域性變數視窗觀察cv_abs_out的前100個值,看與MATLAB執行的結果否相同 for(i=0;i<10;i++)//一般情況下觀察10個值就夠了 { watch_cv_abs_out[i]=cvGet2D(cv_abs_out,0,i).val[0]; } cv_gray_fft2shift(cv_abs_out); double watch_cv_gray_fft2shift_out[100];//利用除錯時的區域性變數視窗觀察cv_gray_fft2shift_out的前100個值,看與MATLAB執行的結果否相同 for(i=0;i<10;i++)//一般情況下觀察10個值就夠了 { watch_cv_gray_fft2shift_out[i]=cvGet2D(cv_abs_out,0,i).val[0]; } cv_range_0to255(cv_abs_out,pOutImage); double watch_cv_range_0to255_out[100];//利用除錯時的區域性變數視窗觀察cv_range_0to255_out的前100個值,看與MATLAB執行的結果否相同 for(i=0;i<10;i++)//一般情況下觀察10個值就夠了 { watch_cv_range_0to255_out[i]=cvGet2D(pOutImage,0,i).val[0]; } return 0; }
上面四個函式與MATLAB執行的結果一致,MATLAB程式如下:
clear all; clc; I=imread('coins.png'); F=fft2(I); MatlabFFt2Out=F; MatlabAbsOut=abs(F); MatlabFftshiftOut=fftshift(MatlabAbsOut); T=MatlabFftshiftOut; max_T=max(max(T)); min_T=min(min(T)); shift_T=-min_T*255/(max_T-min_T); MatlabRange0to255=T.*255/(max_T-min_T)+shift_T;
執行結果截圖如下:
cv_gray_fft2
cv_abs
cv_gray_fft2shift
cv_range_0to255
下面再給出OpenCV下的C++程式碼實現,程式碼中的DFT函式已經包含了上述函式的功能:
//OpenCV版本2.4.9
//交流QQ2487872782
# include <opencv2/core/core.hpp>
# include <opencv2/highgui/highgui.hpp>
# include <opencv2/imgproc/imgproc.hpp>
#include <iostream>
using namespace cv;
using namespace std;
cv::Mat DFT(cv::Mat srcImage)
{
cv::Mat srcGray;
cvtColor(srcImage,srcGray,CV_RGB2GRAY);
// 將輸入影象延擴到最佳的尺寸(快速傅立葉變換是基於影象尺寸是2、3或5的倍數完成的,因此對於輸入源影象,首先應將其轉換成DFTsize,OpenCV中提供了函式getOptimalDFTSize來實現尺寸轉換)
int nRows = getOptimalDFTSize(srcGray.rows);
int nCols = getOptimalDFTSize(srcGray.cols);
cv::Mat resultImage;
// 把灰度影象放在左上角,在右邊和下邊擴充套件影象,
// 新增的畫素初始化為0
copyMakeBorder(srcGray, resultImage, 0,
nRows - srcGray.rows,
0, nCols - srcGray.cols,
BORDER_CONSTANT, Scalar::all(0));
// 為傅立葉變換的結果(實部和虛部)分配儲存空間
cv::Mat planes[] = { cv::Mat_<float>(resultImage),
cv::Mat::zeros(resultImage.size(), CV_32F)};//這裡實際上是建立MAT陣列,陣列有兩個成員:
//第一個就是sizeConvMat這個物件(只是資料型別轉換成了float型別)
//第二個是全0的型別為32F的物件
Mat completeI;
// 為延擴後的影象增添一個初始化為0的通道
merge(planes,2,completeI); //把groupMats的第0和第1個物件合併到mergeMat,通過這個操作mergeMat是雙通道的資料陣列了
// 進行離散傅立葉變換
dft(completeI,completeI);
// 將複數轉換為幅度
split(completeI,planes);
magnitude(planes[0],planes[1],planes[0]);//0中存的是實部值,1中存的是虛部值
cv::Mat dftResultImage = planes[0];
// 對數尺度(logarithmic scale)縮放
dftResultImage += 1;//陣列加1作對數變換,以擴大頻域動態顯示範圍
log(dftResultImage,dftResultImage);//作對數變換
// 剪下和重分佈幅度圖象限
dftResultImage= dftResultImage(Rect(0,
0,srcGray.cols,srcGray.rows));
// 歸一化影象
normalize(dftResultImage,dftResultImage,
0,1,CV_MINMAX);
int cx = dftResultImage.cols/2;
int cy = dftResultImage.rows/2;
Mat tmp;
// Top-Left - 為每一個象限建立ROI
Mat q0(dftResultImage,Rect(0,0,cx,cy));
// Top-Right
Mat q1(dftResultImage,Rect(cx,0,cx,cy));
// Bottom-Left
Mat q2(dftResultImage,Rect(0,cy,cx,cy));
// Bottom-Right
Mat q3(dftResultImage,Rect(cx,cy,cx,cy));
// 交換象限 (Top-Left with Bottom-Right)
q0.copyTo(tmp);
q3.copyTo(q0);
tmp.copyTo(q3);
// 交換象限 (Top-Right with Bottom-Left)
q1.copyTo(tmp);
q2.copyTo(q1);
tmp.copyTo(q2);
return dftResultImage;
}
int main()
{
cv::Mat srcImage = imread("coins.png");
if(srcImage.empty())
return-1;
imshow("srcImage", srcImage);
cv::Mat resultImage = DFT(srcImage);
imshow("resultImage", resultImage);
cv::waitKey(0);
return 0;
}
執行結果如下圖所示: