fortran - Fortran 结果中的 FFTW 仅包含零

Question

我一直在尝试编写一个简单的程序来使用 fftw3 对一维输入数组执行 fft。这里我使用地震图作为输入。然而，输出数组仅包含零。

我知道输入是正确的，因为我也尝试在 MATLAB 中对同一输入文件进行 fft，这给出了正确的结果。没有编译错误。我使用 f95 来编译它，但是 gfortran 也给出了几乎相同的结果。这是我编写的代码:-

program fft

    use functions
    implicit none
    include 'fftw3.f90'
    integer nl,row,col
    double precision, allocatable :: data(:,:),time(:),amplitude(:)
    double complex, allocatable :: out(:)
    integer*8 plan


    open(1,file='test-seismogram.xy')
    nl=nlines(1,'test-seismogram.xy')
    allocate(data(nl,2))
    allocate(time(nl))
    allocate(amplitude(nl))
    allocate(out(nl/2+1))
    do row = 1,nl
        read(1,*,end=101) data(row,1),data(row,2)
        amplitude(row)=data(row,2)
    end do
    101 close(1)


    call dfftw_plan_dft_r2c_1d(plan,nl,amplitude,out,FFTW_R2HC,FFTW_PATIENT)
    call dfftw_execute_dft_r2c(plan, amplitude, out)
    call dfftw_destroy_plan(plan)


    do row=1,(nl/2+1)
        print *,out(row)
    end do


    deallocate(data)
    deallocate(amplitude)
    deallocate(time)
    deallocate(out)
end program fft

nlines() 函数是一个用于计算文件中行数的函数，并且它工作正常。它在名为“functions”的模块中定义。

该程序几乎尝试遵循 http://www.fftw.org/fftw3_doc/Fortran-Examples.html 中的示例

我可能只是犯了一个非常简单的逻辑错误，但我真的无法弄清楚这里出了什么问题。任何指示都会非常有帮助。

这几乎就是整个输出的样子:-

           .
           .
           .
           (0.0000000000000000,0.0000000000000000)
           (0.0000000000000000,0.0000000000000000)
           (0.0000000000000000,0.0000000000000000)
           (0.0000000000000000,0.0000000000000000)
           (0.0000000000000000,0.0000000000000000)
           .
           .
           .

我的疑问直接与 fftw 有关，因为 SO 上有 fftw 的标签，所以我希望这个问题不是偏离主题

Answer 1

正如 @roygvib 和 @Ross 首先在评论中所解释的，计划子例程会覆盖输入数组，因为它们使用不同的参数多次尝试转换。我会添加一些实际使用注意事项。

您声称您确实关心性能。那么有两种可能:

1. 正如代码中所示，仅进行一次转换。那么就没有必要使用 FFTW_MEASURE。计划子例程比实际计划执行子例程慢很多倍。使用 FFTW_ESTIMATE 会快得多。

FFTW_MEASURE tells FFTW to find an optimized plan by actually computing several FFTs and measuring their execution time. Depending on your machine, this can take some time (often a few seconds). FFTW_MEASURE is the default planning option.

FFTW_ESTIMATE specifies that, instead of actual measurements of different algorithms, a simple heuristic is used to pick a (probably sub-optimal) plan quickly. With this flag, the input/output arrays are not overwritten during planning.

http://www.fftw.org/fftw3_doc/Planner-Flags.html

您对不同的数据多次进行相同的转换。然后，您必须在第一次转换之前仅进行一次规划，然后重新使用该计划。只需先制定计划，然后再用第一个输入数据填充数组。在每次传输之前制定计划会使程序变得非常慢。