为什么 FFTW 输入非零，输出返回全零（C++）

Question

好的，所以我一直在研究以下问题一段时间，但似乎无法找到答案。简而言之，我正在为雷达信号实现断层成像，似乎无法让 FFTW 给我除零以外的输出。我正在通过 C++ 实现并使用 FFTW3 库。

现在的输入是场景原点的模拟点散射体，它具有所有真实的相位响应（参见rxphase 变量）。由于使用的带宽等原因，过滤后的接收信号的 IFFT 应该相当大（我还没有对其进行缩放），但每个脉冲我都得到了零。

摘自我的代码：

void ifftshift(std::vector<phdata> &argin, std::size_t count){
  int k = 0;
  int c = (int)floor((float)count/2);

  if (count % 2 == 0)
  {
    for (k=0; k<c;k++)
    {
      complex<double> tmp = argin[k];
      argin[k] = argin[k+c];
      argin[k+c] = tmp;
    }
  }
  else
  {
    complex<double> tmp = argin[count-1];
    for (k=c-1; k>=0; k--)
    {
      argin[c+k+1] = argin[k];
      argin[k] = argin[c+k];
    }
    argin[c] = tmp;
  }
};

void main(){

  std::vector<complex<double> > filteredData;
  // Number of pulses across the aperture
  std::int pulses = 11;
  // Define the frequency vector
  std::vector<double> freq;
  std::double freqmin = 9 * pow(10,9);
  std::double freqmax = 11 * pow(10,9);
  std::vector<complex<double> > outData;
  std::vector<complex<double> > rxphase;

  for  (int i = 0; i<64; i++)
  {
    freq.push_back(freqmin + (((freqmax-freqmin)/64)*i));
  }

  // Create a (samples X pulses) array of complex doubles
  rxphase.assign(64, std::vector<complex<double> > (11, {1,0}) );

  fftw_complex *in, *out;
  fftw_plan plan;

  in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * image.Nfft );
  out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * image.Nfft );
  plan = fftw_plan_dft_1d(image.Nfft, in, out, FFTW_BACKWARD, FFTW_MEASURE);

  // iterate through pulses to back project across the aperture
  for ( int i=0; i<=pulses; i++ )
  {
    // Reset vectors for each pulse
    filteredData.clear();
    outData.clear();

    for ( int ii=0; ii<freq.size(); ii++ )
    {
      // Apply simple ramp filter
      filteredData.push_back(freq[ii]*rxphase[ii][i]);
    }

    // Shift the data to put the center frequency at DC position (zero index)
    ifftshift(filteredData, filteredData.size());

    for (int ii=0; ii<filteredData.size(); ii++)
    {
      std::cout << filteredData[ii] ;
    }
    std::cout << std::endl;
    // filteredData is what I expect it to be.

    // Recast the vector to the type needed for FFTW and compute FFT
    in = reinterpret_cast<fftw_complex*>(&filteredData[0]);

    for (int ii=0; ii<filteredData.size(); ii++)
    {
      std::cout << "(" << (in[ii])[0] << "," << (in[ii])[1] << ");";
    }
    std::cout << std::endl;
    // values for in are the same as for filteredData, as expected.

    fftw_execute(plan);

    for (int ii=0; ii<filteredData.size(); ii++)
    {
      std::cout << "(" << (out[ii])[0] << "," << (out[ii])[1] << ");";
    }
    std::cout << std::endl;
    // The values for out are all (0,0)
  }

  fftw_destroy_plan(plan);
  //fftw_free(in); error here
  //fftw_free(out); error here
};

任何帮助将不胜感激。

编辑：代码应该更完整一些。

Answer 1

混乱可能源于这条线

in = reinterpret_cast<fftw_complex*>(&filteredData[0]);

虽然这会更新局部变量 in 以指向内存中存储数据的区域（filteredData 向量的存储缓冲区），但它不会更新 plan 内部保存的指针。因此，FFTW 知道并用作输入缓冲区的内存区域仍然是您最初分配的内存区域

in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * image.Nfft );

并指定为

的参数

fftw_plan_dft_1d(image.Nfft, in, out, FFTW_BACKWARD, FFTW_MEASURE);
//                           ^^

请注意，此内存区域在整个程序中保持未初始化（可能恰好为零）。

您应该改为直接写入in 输入缓冲区：

for ( int ii=0; ii<freq.size(); ii++ )
{
    // Apply simple ramp filter
    std::complex<double> value = freq[ii]*rxphase[ii][i];
    in[ii][0] = value.real();
    in[ii][1] = value.imag();
}

不用说，在这种情况下，in = reinterpret_cast<fftw_complex*>(&filteredData[0]) 行也应该被删除。

或者，如果您使用的是fftw_complex 和std::complex<double> 二进制兼容的编译器（请参阅FFTW documentation of Complex numbers type），您可能更喜欢：

std::complex<double>* filteredData = reinterpret_cast<std::complex<double>*>(in);
for ( int ii=0; ii<freq.size(); ii++ )
{
    // Apply simple ramp filter
    filteredData[ii] = freq[ii]*rxphase[ii][i];
}

Answer 2

为了快速而肮脏地思考这个问题，您使用fftw_plan_dft_1d创建fftw_plan首先，将有效输入数组的内容修改为fftw_plan_dft_1d，正如SleuthEye提到的那样，然后然后致电fftw_execute。