执行FFT计划后,FFTW输出矢量大小错误

FFTW output vector size is wrong after executing the FFT plan

本文关键字:错误 输出 FFTW FFT 计划 执行      更新时间:2023-10-16

我在使用fftw库执行FFT计划时遇到问题。我声明大小为 2048 的fftInfftOut向量,并用它们定义一个fftwf_plan_dft_1d

但是当我执行计划时,突然间fftOut向量都错了。在我的代码中,我可以在调试时看到它的大小从 2048 变为 0。当我在下面编写小型可执行文件示例时,我只是得到了诸如17293858104588369909之类的错误大小。

当然,当我尝试访问向量的第一项时,会发生SIGSEGV

我的代码:

#include <complex>
#include <iostream>
#include <fftw3.h>
#include <vector>
using namespace std;
typedef std::vector<std::complex<float>> fft_vector;
int main() {
const unsigned int VECTOR_SIZE = 2048;
fft_vector* fftIn = new fft_vector(VECTOR_SIZE);
fft_vector* fftOut = new fft_vector(VECTOR_SIZE);
fftwf_plan fftPlan = fftwf_plan_dft_1d(VECTOR_SIZE, reinterpret_cast<fftwf_complex*>(fftIn), reinterpret_cast<fftwf_complex*>(fftOut), FFTW_FORWARD, FFTW_ESTIMATE);
fftwf_execute(fftPlan);
std::cout << fftOut->size() << std::endl;
std::cout << fftOut->at(0).real() << fftOut->at(0).imag() << std::endl;
return 0;
}

当然,我知道在这个例子中fftIn向量是空的,但是当它不是时,输出就会被破坏。在这种情况下,SIGSEGV 发生在第二个 cout 中,如前所述。

我的完整代码有线程(但 FFT 都发生在同一个线程中,所以竞争条件不应该适用),这是试图在这个小示例中隔离代码的原因之一,以防万一,但它似乎有问题无论如何。

知道吗?

主要问题是你正在向它传递一个向量,这是行不通的;你需要传递向量内容:&((*fftIn)[0])和&((*fftOut)[0])或类似的东西。事实上,你告诉fftw踩踏矢量对象的元数据(包括长度,这解释了为什么它有时是0,有时是胡言乱语)。它实际上是写入向量结构的开头,因为这就是指针指向的内容。此外,它还使用 fftIn 的元数据作为 fft 输入的一部分,这也不是您想要的。

您可以考虑改用fftw_complex和fftw_malloc,这将确保您的数据存储符合 fftw 需要的方式:http://www.fftw.org/doc/SIMD-alignment-and-fftw_005fmalloc.html 如果你真的需要,你可以在之后把fftOut放在一个向量中。这将确保您获得 SIMD 的优势(如果可用),并避免任何特定于编译器或平台的复杂类型和内存分配行为。使用 fftw 的类型和分配器意味着您的代码将始终以最佳方式工作并按预期工作,无论您在哪个平台上运行它以及使用哪个编译器。

以下是从 fftw 的文档 (http://www.fftw.org/doc/Complex-One_002dDimensional-DFTs.html) 转换的示例:

#include <fftw3.h>
...
{
fftw_complex *in, *out;
fftw_plan p;
...
in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
p = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
...
fftw_execute(p); /* repeat as needed */
...
fftw_destroy_plan(p);
fftw_free(in); fftw_free(out);
}

尝试匹配此示例,您的代码应该会执行您的期望。

编辑: 如果您必须使用C++复数和向量,这应该有效:

const unsigned int VECTOR_SIZE = 2048;
fft_vector* fftIn = new fft_vector(VECTOR_SIZE);
fft_vector* fftOut = new fft_vector(VECTOR_SIZE);
fftwf_plan fftPlan = fftwf_plan_dft_1d(VECTOR_SIZE, 
reinterpret_cast<fftwf_complex*>(&(*fftIn)[0]),
reinterpret_cast<fftwf_complex*>(&(*fftOut)[0]),
FFTW_FORWARD, FFTW_ESTIMATE);
fftwf_execute(fftPlan);
std::cout << fftOut->size() << std::endl;
std::cout << fftOut->at(0).real() << fftOut->at(0).imag() << std::endl;

请注意,您必须取消引用矢量指针才能使 [] 运算符工作;获取索引 0 会为您提供实际矢量数据的第一个元素,因此该元素的地址是您需要提供给fftw_plan的地址(指针)。

你不能reinterpret_cast向量到复杂元素上的指针:你访问的是向量内部变量(如大小),不一定是向量原始数据,这可能解释了为什么对象在输出中被丢弃,以及接下来发生的所有坏事。

另一方面,转换矢量数据的地址是有效的。

只需做(也不需要new向量,使用标准声明):

fft_vector fftIn(VECTOR_SIZE);
fft_vector fftOut(VECTOR_SIZE);
// load your data in fftIn
fftwf_plan_dft_1d(VECTOR_SIZE, static_cast<const fftwf_complex*>(&fftIn[0]), static_cast<fftwf_complex*>(&fftOut[0]), FFTW_FORWARD, FFTW_ESTIMATE);