使用 DFT 进行卷积答案

【问题标题】：convolution using DFT使用 DFT 进行卷积
【发布时间】：2015-03-09 20:51:07
【问题描述】：

我使用以下代码计算图像与指定内核的卷积（在我的情况下为高斯）。每次我得到不同的结果，结果图像甚至不接近我在空间域中通过卷积获得的图像。首先，我认为问题出在图像的数据类型上。我将它们更改为 32 和 64，但结果仍然相同。谁能告诉我可能出了什么问题？

http://opencv.willowgarage.com/documentation/cpp/core_operations_on_arrays.html#dft 上面的这个函数给了我一个黑色的图像。我有灰度输入。

void convol_fft(const Mat& A,const vector<vector<float>>& kernel2d,Mat& result)
{

    Mat B = Mat(3,3,CV_64F);
    for (int row = 0; row < kernel2d.size(); row++)
        for (int col = 0; col < kernel2d[row].size(); col++){
            B.at<uchar>(row,col) = (uchar)kernel2d[row][col];
        }

    int dft_M = getOptimalDFTSize( A.rows+B.rows-1 );
    int dft_N = getOptimalDFTSize( A.cols+B.cols-1 );
    Mat dft_A = Mat::zeros(dft_M, dft_N, CV_64F);
    Mat dft_B = Mat::zeros(dft_M, dft_N, CV_64F);

    Mat dft_A_part = dft_A(Rect(0, 0, A.cols,A.rows));
    A.convertTo(dft_A_part, dft_A_part.type(), 1, -mean(A)[0]);
    Mat dft_B_part = dft_B(Rect(0, 0, B.cols,B.rows));
    B.convertTo(dft_B_part, dft_B_part.type(), 1, -mean(B)[0]);

    dft(dft_A, dft_A, 0, A.rows);
    dft(dft_B, dft_B, 0, B.rows);

    // set the last parameter to false to compute convolution instead of correlation
    mulSpectrums( dft_A, dft_B, dft_A, 0, false );
    idft(dft_A, dft_A, DFT_SCALE, A.rows + B.rows - 1 );

    result = dft_A(Rect(0, 0, A.cols + B.cols - 1, A.rows + B.rows - 1));
    normalize(result, result, 0, 1, NORM_MINMAX, result.type());
    pow(result, 3., result);

  //  B ^= Scalar::all(255);

}

【问题讨论】：

如果你总是得到不同的结果，那么你可能会错过一些初始化，因此使用随机值......
我没有使用任何随机值来初始化函数
我的意思是 kernel2d（乍一看是唯一可能的罪魁祸首）
是的，我正在使用函数制作 kernel2d。但即使是来自 OpenCV 的代码也给了我一个全黑的图像。我不得不改变这一行 Mat tempA(dftSize, CV_64FC1, Scalar::all(0));并在此之后使该代码工作，因为它不是多频谱或 dft 只接受 4 种输入数据类型。最终的转换不会对结果产生任何影响。

标签： c++ opencv fft convolution dft

【解决方案1】：

以下基于 openCV 的 phaseCorrelateRes() 的代码将在二维中进行关联。

static void fftShift(InputOutputArray _out)
{
    Mat out = _out.getMat();

    if(out.rows == 1 && out.cols == 1)
    {
        // trivially shifted.
        return;
    }

    vector<Mat> planes;
    split(out, planes);

    int xMid = out.cols >> 1;
    int yMid = out.rows >> 1;

    bool is_1d = xMid == 0 || yMid == 0;

    if(is_1d)
    {
        xMid = xMid + yMid;

        for(size_t i = 0; i < planes.size(); i++)
        {
            Mat tmp;
            Mat half0(planes[i], Rect(0, 0, xMid, 1));
            Mat half1(planes[i], Rect(xMid, 0, xMid, 1));

            half0.copyTo(tmp);
            half1.copyTo(half0);
            tmp.copyTo(half1);
        }
    }
    else
    {
        for(size_t i = 0; i < planes.size(); i++)
        {
            // perform quadrant swaps...
            Mat tmp;
            Mat q0(planes[i], Rect(0,    0,    xMid, yMid));
            Mat q1(planes[i], Rect(xMid, 0,    xMid, yMid));
            Mat q2(planes[i], Rect(0,    yMid, xMid, yMid));
            Mat q3(planes[i], Rect(xMid, yMid, xMid, yMid));

            q0.copyTo(tmp);
            q3.copyTo(q0);
            tmp.copyTo(q3);

            q1.copyTo(tmp);
            q2.copyTo(q1);
            tmp.copyTo(q2);
        }
    }

    merge(planes, out);
}

void Correlate2d(
    const cv::Mat& src1, 
    const cv::Mat& src2, 
    cv::Mat& dst,
    double* response)
{

    CV_Assert( src1.type() == src2.type());
    CV_Assert( src1.type() == CV_32FC1 || src1.type() == CV_64FC1 );
    CV_Assert( src1.size == src2.size);

    int M = getOptimalDFTSize(src1.rows);
    int N = getOptimalDFTSize(src1.cols);

    Mat padded1, padded2, paddedWin;

    if(M != src1.rows || N != src1.cols)
    {
        copyMakeBorder(src1, padded1, 0, M - src1.rows, 0, N - src1.cols, BORDER_CONSTANT, Scalar::all(0));
        copyMakeBorder(src2, padded2, 0, M - src2.rows, 0, N - src2.cols, BORDER_CONSTANT, Scalar::all(0));
    }
    else
    {
        padded1 = src1;
        padded2 = src2;
    }

    Mat FFT1, FFT2, P, Pm, C;

    // correlation equation
    // Reference: http://en.wikipedia.org/wiki/Phase_correlation
    dft(padded1, FFT1, DFT_REAL_OUTPUT);
    dft(padded2, FFT2, DFT_REAL_OUTPUT);

    mulSpectrums(FFT1, FFT2, dst, 0, true);
    idft(dst, dst, DFT_SCALE); // gives us the correlation result...
    fftShift(dst); // shift the energy to the center of the frame.

    // locate the highest peak
    Point peakLoc;
    minMaxLoc(dst, NULL, NULL, NULL, &peakLoc);

    // max response is scaled
    if( response )
        *response = dst.at<float>(peakLoc);
}

你可以在\opencv\sources\modules\imgproc\src\phasecorr.cpp中找到代码

为了将代码更改为卷积，只需更改此行：

mulSpectrums(FFT1, FFT2, dst, 0, true);

到

mulSpectrums(FFT1, FFT2, dst, 0, false);

这相当于在matlab中做：

dst = fftshift(ifft2(fft2(src1).*conj(fft2(src2))))

【讨论】：

【解决方案2】：

我不确定 OpenCV...但这看起来很可疑。

for (int row = 0; row < kernel2d.size(); row++)
    for (int col = 0; col < kernel2d[row].size(); col++){
        B.at<uchar>(row,col) = (uchar)kernel2d[row][col];
 }

如果您正在填充 B 内核，那么该行应该是 kernel2d[col].size()。看起来您正在超出 B 内核。 kernel2d.size() 的值是多少？

为什么不直接加载值？保存所有函数调用。

对于高斯内核，它应该类似于 {1,2,1,2,3,2,1,2,1}。

【讨论】：