使用GNU八度FFT函数

你可能做错了。删除代码中的所有“负”频率。你应该保持正负低频。下面是一个用python编写的代码和结果。这个图的比例尺是正确的。

alt text http://files.droplr.com/files/35740123/XUl90.fft.png

代码:

from __future__ import division

from scipy.signal import fft, ifft
import numpy as np

def approximate(signal, cutoff):
    fourier = fft(signal)
    size = len(signal)
    # remove all frequencies except ground + offset positive, and offset negative:
    fourier[1+cutoff:-cutoff] = 0
    return ifft(fourier)

def quad(x):
    return x**2

from pylab import plot

X = np.arange(-4096,4096)/64
Y = quad(X)

plot(X,Y)
plot(X,approximate(Y,3))

你正在扔掉转换的后半部分。对于实值输入,转换是Hermitian对称的,必须保留这些线。试试这个:

function points = approximate(inp, count)
    fourier = fft(inp);
    fourier((count+1):(length(fourier)-count+1)) = 0;
    points  = real(ifft(fourier)); %# max(imag(ifft(fourier))) should be around eps(real(...))
endfunction;

由于数值计算的误差,逆变换总是会有一些微小的虚部,因此 real 提取。

注意 input 和 size 是八度音阶的关键词;用你自己的变数来打击它们是一个很好的方法,让真正奇怪的虫子上路!