2011-07-31 72 views
2

我正在尝试创建一个程序,该程序将绘制给定图像的二维灰度光谱。我正在使用OpenCV和FFTW库。通过使用来自互联网的提示和代码并修改它们,我设法加载一张图片,计算这张图片的fft并从fft重新创建图片(这是相同的)。我无法做的是画出傅立叶频谱本身。你可以帮我吗? 下面的代码(除去不太重要的行):在C++中绘制图像的频谱(fftw,OpenCV)

/* Copy input image */ 

/* Create output image */ 

/* Allocate input data for FFTW */ 
in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N); 
dft = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N); 

/* Create plans */ 
plan_f = fftw_plan_dft_2d(w, h, in, dft, FFTW_FORWARD, FFTW_ESTIMATE); 

/* Populate input data in row-major order */ 
for (i = 0, k = 0; i < h; i++) 
{ 
    for (j = 0; j < w; j++, k++) 
    { 
     in[k][0] = ((uchar*)(img1->imageData + i * img1->widthStep))[j]; 
     in[k][1] = 0.; 
    } 
} 

/* forward DFT */ 
fftw_execute(plan_f); 

/* spectrum */ 
for (i = 0, k = 0; i < h; i++) 
{ 
    for (j = 0; j < w; j++, k++) 
     ((uchar*)(img2->imageData + i * img2->widthStep))[j] = sqrt(pow(dft[k][0],2) + pow(dft[k][1],2)); 
}  

cvShowImage("iplimage_dft(): original", img1); 
cvShowImage("iplimage_dft(): result", img2); 
cvWaitKey(0); 

/* Free memory */ 

}

的问题是 “频谱” 部分。我得到一些噪音,而不是谱。我究竟做错了什么?我会很感激你的帮助。

+0

听起来像一个缩放问题 - 你应该检查你的FFT输出幅度值的范围。 –

+0

你能建议我应该怎么做?从我在不同论坛上阅读的内容来看,如果是关于幅度的话,我会得到一个黑色的图像(中心有一个白点)。谢谢您的回答。 – Narren

+1

既然你没有对幅度进行任何范围检查,那么如果它很大,当你将它分配给一个8位像素(它看起来像噪声)时,它会将模2包裹。这就是为什么我说你应该检查范围 - 例如添加另一个循环来查找最大和最小幅度,然后相应地缩放您的值,以便您知道它们将适合8位范围。 –

回答

0

你可以尝试做IFFT步骤,看看你是否恢复原始图像?那么,你可以逐步检查你的问题在哪里。找到问题的另一个解决方案是用你预先定义好的一个小矩阵来完成这个过程,然后在MATLAB中计算出FFT,然后逐步检查,它对我有用!

1

您需要绘制频谱的大小。这里是代码。

void ForwardFFT(Mat &Src, Mat *FImg) 
{ 
    int M = getOptimalDFTSize(Src.rows); 
    int N = getOptimalDFTSize(Src.cols); 
    Mat padded;  
    copyMakeBorder(Src, padded, 0, M - Src.rows, 0, N - Src.cols, BORDER_CONSTANT, Scalar::all(0)); 
    // Создаем комплексное представление изображения 
    // planes[0] содержит само изображение, planes[1] его мнимую часть (заполнено нулями) 
    Mat planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)}; 
    Mat complexImg; 
    merge(planes, 2, complexImg); 
    dft(complexImg, complexImg);  
    // После преобразования результат так-же состоит из действительной и мнимой части 
    split(complexImg, planes); 

    // обрежем спектр, если у него нечетное количество строк или столбцов 
    planes[0] = planes[0](Rect(0, 0, planes[0].cols & -2, planes[0].rows & -2)); 
    planes[1] = planes[1](Rect(0, 0, planes[1].cols & -2, planes[1].rows & -2)); 

    Recomb(planes[0],planes[0]); 
    Recomb(planes[1],planes[1]); 
    // Нормализуем спектр 
    planes[0]/=float(M*N); 
    planes[1]/=float(M*N); 
    FImg[0]=planes[0].clone(); 
    FImg[1]=planes[1].clone(); 
} 
void ForwardFFT_Mag_Phase(Mat &src, Mat &Mag,Mat &Phase) 
{ 
    Mat planes[2]; 
    ForwardFFT(src,planes); 
    Mag.zeros(planes[0].rows,planes[0].cols,CV_32F); 
    Phase.zeros(planes[0].rows,planes[0].cols,CV_32F); 
    cv::cartToPolar(planes[0],planes[1],Mag,Phase); 
} 
Mat LogMag; 
    LogMag.zeros(Mag.rows,Mag.cols,CV_32F); 
    LogMag=(Mag+1); 
    cv::log(LogMag,LogMag); 
    //--------------------------------------------------- 
    imshow("Логарифм амплитуды", LogMag); 
    imshow("Фаза", Phase); 
    imshow("Результат фильтрации", img);