GSL 快速傅立叶变换 - 变换高斯的非零虚数？

Question

作为我询问的 this question 的扩展。实高斯的傅里叶变换是实高斯。现在当然，一组仅类似于高斯的点的 DFT 并不总是完美的高斯，但它肯定应该接近。在下面的代码中，我使用 GSL 进行 [离散] 傅里叶变换。除了 returned/transformed 实分量的问题（在链接问题中概述）之外，我得到了虚分量（应该完全为零）的奇怪结果。当然，它的幅度很小，但它仍然很奇怪。 这种不对称和古怪的输出的原因是什么？

#include <gsl/gsl_fft_complex.h>
#include <gsl/gsl_errno.h>
#include <fstream>
#include <iostream>
#include <iomanip> 

#define REAL(z,i) ((z)[2*(i)]) //complex arrays stored as    [Re(z0),Im(z0),Re(z1),Im(z1),...]
#define IMAG(z,i) ((z)[2*(i)+1])
#define MODU(z,i) ((z)[2*(i)])*((z)[2*(i)])+((z)[2*(i)+1])*((z)[2*(i)+1])
#define PI 3.14159265359

using namespace std;

int main(){

    int n = pow(2,9);
    double data[2*n];
    double N = (double) n;

    ofstream file_out("out.txt");

    double xmin=-10.;
    double xmax=10.;
    double dx=(xmax-xmin)/N;
    double x=xmin;

    for (int i=0; i<n; ++i){
        REAL(data,i)=exp(-100.*x*x);
        IMAG(data,i)=0.;
        x+=dx;
    }

    gsl_fft_complex_radix2_forward(data, 1, n); 

    for (int i=0; i<n; ++i){
        file_out<<(i-n/2)<<"    "<<IMAG(data,((i+n/2)%n))<<'\n';
    }

    file_out.close();
}

Answer 1

你的虚部结果是正确的，符合预期。

与零 (10^-15) 的差值小于你给 pi 的精度（12 位数字，pi 在 FFT 中使用，但我不知道你是否覆盖了里面的 pi例程）。

实函数的 FFT 通常不是实函数。当您以分析方式进行数学计算时，您将对以下表达式进行积分：

f(t) e^{i w t} = f(t) cos wt  + i f(t) sin wt,

因此只有当函数 f(t) 为 real and even 时，虚部（否则为奇数）才会在积分过程中消失。不过这没什么意义，因为实部和虚部只有在特殊情况下才具有物理意义。

直接物理意义是在abs值（magnitude spectrum），abs。值的平方 (intensity spectrum) 和相位或角度 (phase spectrum)。

如果虚数部分不以时间向量的中心为中心，则会发生更显着的零偏移。尝试将 x 向量移动 dx 的一部分。

请参阅下文 dx/2（右列）的输入偏移如何影响虚部，但不影响幅度（用 Python、Numpy 编写的示例）。

from __future__ import division
import numpy as np
import matplotlib.pyplot as p
%matplotlib inline

n=512    # number of samples 2**9
x0,x1=-10,10
dx=(x1-x0)/n

x= np.arange(-10,10,dx)  # even number, asymmetric range [-10, 10-dx]   

#make signal
s1= np.exp(-100*x**2) 
s2= np.exp(-100*(x+dx/2 )**2)

#make ffts
f1=np.fft.fftshift(np.fft.fft(s1))
f2=np.fft.fftshift(np.fft.fft(s2))

#plots
p.figure(figsize=(16,12))
p.subplot(421)
p.title('gaussian (just ctr shown)')
p.plot(s1[250:262])
p.subplot(422)
p.title('same, shifted by dx/2')
p.plot(s2[250:262])

p.subplot(423)
p.plot(np.imag(f1))
p.title('imaginary part of FFT')
p.subplot(424)
p.plot(np.imag(f2))

p.subplot(425)
p.plot(np.real(f1))
p.title('real part of FFT')
p.subplot(426)
p.plot(np.real(f2))

p.subplot(427)
p.plot(np.abs(f1))
p.title('abs. value of FFT')
p.subplot(428)
p.plot(np.abs(f2))

GSL 快速傅立叶变换 - 变换高斯的非零虚数？

GSL Fast-Fourier Transform - Non-zero Imaginary for Transformed Gaussian?

fft

gsl

dft