选择维度以在 C 中使用 FFTW 执行 FFT答案

【问题标题】：Select dimension to perform FFT with FFTW in C选择维度以在 C 中使用 FFTW 执行 FFT
【发布时间】：2018-07-05 13:25:47
【问题描述】：

我有一些存储为D[x][y][z][t] 的数据，我想使用fftw 库执行fft，但仅在t 维度上，这意味着我想获得该矩阵的频率响应。当然D[x][y][z][t] 只是一个向量化数组的包装器（因为这种数组可以用作fftw 中的参数）。在 MATLAB 中可以直接计算它，但在 c 中我根本找不到任何类型的参考，即使在 fftw 参考中也是如此。有什么想法吗？

数组存储为D[l + k*Nt + j*Nt*Nz + i*Nt*Nz*Ny]，其中i,j,k,l 是x,y,z,t 的索引，对应的数组大小为Nx, Ny, Nz, Nt。因此，我想要的是在 t 维度上计算 FFT，并带有 Nt 点。

相应的计划应该类似于myplan = fftw_plan_dft_r2c_1d(Nt, &D[0], &out_array, FFTW_ESTIMATE);，但我不确定这是否可行。我希望输出数组是一个类似大小的数组，包含 FFT 转换 x,y,z，其方式与原始 D 数组相同。请注意，output_array 是一个复杂数组，其分配方式与双精度数组 D 相同。

【问题讨论】：

请发帖Minimal, Complete, and Verifiable example，这样人们至少可以看到你在这方面的尝试
我已经发布了我想做的事情。在这一点上一个例子是不可行的，因为我不知道从哪里开始。我熟悉一维傅里叶变换，但我想知道如何仅在特定维度上执行一维

标签： c fftw

【解决方案1】：

您不能使用 FFTW 对上述数据的第 n 维进行 fft。您可以做的是沿相关维度复制“列”并执行许多 ffts 并将数据复制回来。换句话说，您会将 FFT 计划为一个或多个正确长度的一维 fft。将数据“column”逐个“column”复制到您为操作准备的内存中，沿着该数据复制回它所属的位置。

另一个更昂贵的操作是置换数据，以便第一个维度是您想要 fft 并遍历所有列的维度。

在攻读博士学位期间，我编写了这个用于医学图像重建的库 (https://github.com/kvahed/codeare)。 DFT 实现允许沿任何维度进行 fft，但它是按照我描述的方式完成的： https://github.com/kvahed/codeare/blob/master/src/matrix/ft/DFT.hpp。这是 C++，但你会明白的。寻找功能

template<class T> inline static Matrix<T> fft (
  const Matrix<T>& in, size_t dim, bool shift, bool fwd) noexcept;

具体的精度代码在这里： https://github.com/kvahed/codeare/blob/master/src/matrix/ft/FFTWTraits.hpp

【讨论】：

【解决方案2】：

这对于您的具体情况实际上很简单，并且可以推广到任意维度。您的目标是在“t”维度上进行一系列 FFT，因此您只需告诉 FFTW 数据的内存布局。如果您只想在 t 维度中进行 FFT，则“异地”实现将类似于以下内容

int rank=1;
const int * N={Nx,Ny,Nz,Nt};
int howmany=Nx*Ny*Nz;
double *in = D;
const int *inembed=NULL;
int istride=1; 
int idist=Nt;
fftw_complex *out = Df;
const int *onembed = NULL;
int ostride=1; 
int odist = (Nt+2)/2;
unsigned flags = FFTW_ESTIMATE;

fftw_plan plan = fftw_plan_many_dft_r2c(rank, &N[3], howmany, 
                                        D, inembed,istride, idist, 
                                        DF, onembed, ostride, odist, flags);

有关详细信息，请参阅 http://www.fftw.org/fftw3_doc/Advanced-Real_002ddata-DFTs.html，您需要遵守输出数组中 r2c 转换的特定内存要求。

【讨论】：