【问题标题】:fftw3 inverse transform not workfftw3逆变换不起作用
【发布时间】:2015-05-22 17:06:40
【问题描述】:

我正在使用 fftw3 库在 c 中为复杂到复杂的 DFT 编写一个简单的代码。 我已经编写了一个带有输入数组双精度数据的文件,因此我可以与 matlab fft 函数进行比较。 我尝试从变换数组执行反向变换,但结果和第一个输入数组不同。这是我的结果:

FFTW3 TRANSFORM WELCOME <<<<<

enter the number (integer) N of samples (Bit: 64) (preferably power of 2):8

SAMPLE INPUT
in[0][0] = -216448918.015237        in[0][1] = 0.000000 
in[1][0] = 948904790.062151         in[1][1] = 0.000000
in[2][0] = 826811206.185300         in[2][1] = 0.000000
in[3][0] = 1868763250.342451        in[3][1] = 0.000000
in[4][0] = 703135606.077152         in[4][1] = 0.000000
in[5][0] = -1989016445.622210       in[5][1] = 0.000000
in[6][0] = 1912963650.704585        in[6][1] = 0.000000
in[7][0] = 811527262.805480         in[7][1] = 0.000000

 Hit enter to continue ... 


 FORWARD TRANSFORM COEFFICIENTS

out[0][0] = 4866640402.539672       out[0][1] = 0.000000
out[1][0] = 410260768.150135        out[1][1] = -1738850319.926936
out[2][0] = -2253088168.827970      out[2][1] = 3720402168.707990
out[3][0] = -2249429816.334913      out[3][1] = -3911155208.965507
out[4][0] = 1586282687.363928       out[4][1] = 0.000000
out[5][0] = -2249429816.334913      out[5][1] = 3911155208.965507
out[6][0] = -2253088168.827970      out[6][1] = -3720402168.707990
out[7][0] = 410260768.150135        out[7][1] = 1738850319.926936
do you want to calculate the inverse-transform? (y/n) 
y


INVERSE TRANSFORM COEFFICIENTS
rev[0][0] = -1731591344.121896      rev[0][1] = 0.000000
rev[1][0] = 7591238320.497208       rev[1][1] = 0.000000
rev[2][0] = 6614489649.482399       rev[2][1] = 0.000000
rev[3][0] = 14950106002.739609      rev[3][1] = 0.000000
rev[4][0] = 5625084848.617215       rev[4][1] = 0.000000
rev[5][0] = -15912131564.977680        rev[5][1] = 0.000000
rev[6][0] = 15303709205.636681      rev[6][1] = 0.000000
rev[7][0] = 6492218102.443840       rev[7][1] = 0.000000

正如您所见,'in' 和 'rev' 数组是不同的,但直接转换是正确的。我已经将它与matlab进行了比较,结果是相同的。 当我用 matlab 执行逆变换时,我得到了输入数组。 我能做什么?

这是我的 c 代码:

#include <fftw3.h>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#define PI 3.141592653589

int main()
{

  fftw_complex *in, *out, *rev;
  int i,f0,A,N;
  char no;  
  FILE *fp;
  fftw_plan p;

  printf("\n\n>>>>> FFTW3 TRANSFORM WELCOME <<<<<");
  printf("\n\n enter the number (integer) N of samples (bit: %ld) (preferably power of 2):",(sizeof(fftw_complex)/2)*8);
  scanf("%d",&N);

  //f0 = 50;
  //A = 1;


  //allocating memory for input & output arrays
  in = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*N);
  out = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*N);
  rev = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*N);



  //Opening the data file
  if((fp=fopen("lista_numeri_double.dat","rb"))==NULL)
  {
    printf("\nError reading file\n");
    exit(1);
  }


  printf("\nSAMPLE INPUT");
  //assigning samples read from the file
  for(i=0;i<N;i++)
  {
    //in[i][0] = A * cos(2*PI*f0*i/N);
    fread(&in[i][0],sizeof(double),1,fp);   
    in[i][1]=0;

    printf("\nin[%d][0] = %f \t\tin[%d][1] = %f",i,in[i][0],i,in[i][1]);
  }


  //plan and execute transform
  p = fftw_plan_dft_1d(N,in,out,FFTW_FORWARD,FFTW_ESTIMATE);
  fftw_execute(p);

  printf("\n\n Hit enter to continue ... \n");
  scanf("%c",&no);


  //print output values 
  printf("\n\nFORWARD TRANSFORM COEFFICIENTS\n");
  for(i=0;i<N;i++)
  {
    printf("\nout[%d][0] = %f \t\tout[%d][1] = %f",i,out[i][0],i,out[i][1]);
  }

  fftw_destroy_plan(p); 

  printf("\n do you want to calculate the inverse-transform? (y/n)  \n");
  scanf ("%c",&no);

  if(no=='y')
  {

    //plan and execute inverse transform
    p = fftw_plan_dft_1d(N,out,rev,FFTW_BACKWARD,FFTW_ESTIMATE);
    fftw_execute(p);

    printf("\n\nINVERSE TRANSFORM COEFFICIENTS\n");
    for(i=0;i<N;i++)
    {
      printf("rev[%d][0] = %f \t\trev[%d][1] = %f\n",i,rev[i][0],i,rev[i][1]);
    }

    fftw_destroy_plan(p);           
  }

  return 0;
}

【问题讨论】:

    标签: c matlab transform fft fftw


    【解决方案1】:

    Matlab 和 FFTW 的区别在于应用于变换的缩放因子。

    虽然Matlab's FFT 已标准化,但in FFTW's documentation 所描述的FFTW 使用的算法未标准化。换句话说,使用 FFTW 的全圆变换(向前后向后)将结果缩放了一个因子N

    相应地,比较inrev 数组表明rev 缩放了一致的8 倍(您的示例中使用的变换的大小N)。

    【讨论】:

    • 我认为这些天 matlab 窃取了 fftw 并将其合并为自己的。 [3] FFTW (fftw.org) [4] Frigo, M. 和 S. G. Johnson,“FFTW:FFT 的自适应软件架构”,国际声学、语音和信号处理会议论文集,卷。 3,1998 年,第 1381-1384 页。
    • @GAlexander afaik matlab 正在使用 fftw 的修改版本,为了向后兼容方便地对其进行了缩放。
    • 在 pld 的日子里,fftw 会粉碎 matlab fft,以至于我知道没有人使用 matlabs 实现。
    猜你喜欢
    • 1970-01-01
    • 2021-01-01
    • 1970-01-01
    • 2020-11-28
    • 1970-01-01
    • 2015-08-07
    • 2019-01-20
    • 2020-12-07
    • 1970-01-01
    相关资源
    最近更新 更多