2

由于 OpenCV 的许可更宽松，我正在尝试使用 OpenCV 而不是 fftw。在 FFTW fftwf_plan_r2r_2d() 和 FFTW_REDFT01 等价物的帮助下

我知道它目前是错误的，这只是我最初的想法。我需要一些指导才能通过 OpenCV 获得 fftw 的确切功能，我正在为可以帮助我做到这一点的人提供赏金。

我只需要那些特定的功能，我现在对转换其他 fftw 功能不感兴趣。

谢谢！

任何人？

H：

#pragma once

#define DO_NOT_USE_FFTW

#ifdef DO_NOT_USE_FFTW
enum fftwf_r2r_kind
{
    FFTW_R2HC,    //computes a real-input DFT with output in “halfcomplex” format, i.e. real and imaginary parts for a transform of size n stored as: r0, r1, r2, ..., rn/2, i(n+1)/2-1, ..., i2, i1 (Logical N=n, inverse is FFTW_HC2R.)
    FFTW_HC2R,    //computes the reverse of FFTW_R2HC, above. (Logical N=n, inverse is FFTW_R2HC.)
    FFTW_DHT,     //computes a discrete Hartley transform. (Logical N=n, inverse is FFTW_DHT.)
    FFTW_REDFT00, //computes an REDFT00 transform, i.e. a DCT-I. (Logical N=2*(n-1), inverse is FFTW_REDFT00.)
    FFTW_REDFT10, //computes an REDFT10 transform, i.e. a DCT-II (sometimes called “the” DCT). (Logical N=2*n, inverse is FFTW_REDFT01.)
    FFTW_REDFT01, //computes an REDFT01 transform, i.e. a DCT-III (sometimes called “the” IDCT, being the inverse of DCT-II). (Logical N=2*n, inverse is FFTW_REDFT=10.)
    FFTW_REDFT11, //computes an REDFT11 transform, i.e. a DCT-IV. (Logical N=2*n, inverse is FFTW_REDFT11.)
    FFTW_RODFT00, //computes an RODFT00 transform, i.e. a DST-I. (Logical N=2*(n+1), inverse is FFTW_RODFT00.)
    FFTW_RODFT10, //computes an RODFT10 transform, i.e. a DST-II. (Logical N=2*n, inverse is FFTW_RODFT01.)
    FFTW_RODFT01, //computes an RODFT01 transform, i.e. a DST-III. (Logical N=2*n, inverse is FFTW_RODFT=10.)
    FFTW_RODFT11, //computes an RODFT11 transform, i.e. a DST-IV. (Logical N=2*n, inverse is FFTW_RODFT11.)
};

struct fftwf_plan_imp
{
     int nx;
     int ny;
     float * in;
     float * out; 
     fftwf_r2r_kind kindx;
     fftwf_r2r_kind kindy;
     unsigned flags;
};

#define fftwf_plan         fftwf_plan_imp *
#define fftwf_malloc(x)    malloc(x)
#define fftwf_free(x)      free(x)
#define FFTW_DESTROY_INPUT 1
#define FFTW_ESTIMATE      1

fftwf_plan fftwf_plan_r2r_2d(int nx, int ny, float * in, float * out, fftwf_r2r_kind kindx, fftwf_r2r_kind kindy, unsigned flags);
void fftwf_execute(fftwf_plan x);
void fftwf_destroy_plan(fftwf_plan x);
void fftwf_cleanup();
#else
#include <fftw3.h>
#endif

源代码：

#include "KissInterface.h"
#include "kiss_fftr.h"
#include "kiss_fftndr.h"
#include "opencv.hpp"

//currently cv::dct supports even-size arrays (2, 4, 6 ...). For data analysis and approximation you can pad the array when necessary.

#ifdef DO_NOT_USE_FFTW
fftwf_plan fftwf_plan_r2r_2d(int nx, int ny, float * in, float * out, fftwf_r2r_kind kindx, fftwf_r2r_kind kindy, unsigned flags)
{
    fftwf_plan plan = new fftwf_plan_imp;

    plan->nx    = nx;
    plan->ny    = ny;
    plan->in    = in;
    plan->out   = out;
    plan->kindx = kindx;
    plan->kindy = kindy;
    plan->flags = flags;

    return plan;
}

void fftwf_execute(fftwf_plan x)
{
    int newx = x->nx;
    int newy = x->ny;
    if(newx % 2 != 0)++newx;
    if(newy % 2 != 0)++newy;

    cv::Mat Input (cv::Size(newx, newy), CV_32FC1, x->in); //AUTO_STEP or not?
    cv::Mat Output(cv::Size(newx, newy), CV_32FC1, x->out); //AUTO_STEP or not?

    if(x->kindx == FFTW_REDFT00 && x->kindy == FFTW_REDFT00)
    {
        cv::dct(Input, Output);
    }
    else if(x->kindx == FFTW_REDFT10 && x->kindy == FFTW_REDFT10)
    {
        cv::dct(Input, Output);

        Output *= (4 * sqrt(float(newx / 2)) * sqrt(float(newy / 2)));
        Output.row(0) *= 1.41421356237;
        Output.col(0) *= 1.41421356237;
    }
    else if(x->kindx == FFTW_REDFT01 && x->kindy == FFTW_REDFT01)
    {
        // First re-scale the data for idct():
        Input /= (4 * sqrt(float(newx / 2)) * sqrt(float(newy / 2)));
        Input.row(0) /= 1.41421356237;
        Input.col(0) /= 1.41421356237;

        cv::idct(Input, Output); // this will return the input exactly

        // However, the transforms computed by FFTW are unnormalized, exactly like the corresponding, 
        // so computing a transform followed by its inverse yields the original array scaled by N, where N is the logical DFT size. 
        // The logical DFT size: Logical N=2*n for each axis, this is th implicit symmetrization
        // of the image: reflect right and then reflect both halves down.
        int logicalSizeN = (2 * newx) * (2 * newy);
        Output *= logicalSizeN; // scale to be like FFTW result
    }

    memcpy(x->out, Output.ptr<float>(0), x->nx * x->ny);
}

void fftwf_destroy_plan(fftwf_plan x)
{
    delete x;
    x = NULL;
}

void fftwf_cleanup()
{

}
#endif
4

1 回答 1

0

OpenCV 的 FFT/DCT 的工作方式类似于 Matlab 的实现。了解这一点可以帮助您找到更好的参考资料和答案。
在我的问答中，您有几个链接可以帮助您入门。

最近我不得不回到这段代码并重新编写它，结果发现因为我想使用转换然后再转换回来，我可以避免所有额外的缩放，之前和之后，它就可以工作了。
当然，这意味着中间值（变换之间）与 fftw 不同，但逆变换后的结果却是。

于 2013-12-01T08:05:45.990 回答