使用 GNU Octave FFT 函数

发布于 2024-08-31 14:23:46 字数 528 浏览 6 评论 0原文

我正在使用八度的 fft 函数，但我无法真正弄清楚如何缩放它们的输出：我使用以下（非常短的）代码来近似函数：

function y = f(x)
    y = x .^ 2;
endfunction;

X=[-4096:4095]/64;
Y = f(X);
# plot(X, Y);

F = fft(Y);
S = [0:2047]/2048;

function points = approximate(input, count)
    size    = size(input)(2);
    fourier = [fft(input)(1:count) zeros(1, size-count)];
    points  = ifft(fourier);
endfunction;

Y = f(X); plot(X, Y, X, approximate(Y, 10));

基本上，它的作用是采用函数，计算图像一个音程，对其进行 fft，然后保留一些谐波，并对结果进行 ifft。然而我得到了一个垂直压缩的图（输出的垂直比例是错误的）。有什么想法吗？

原文

I'm playing with octave's fft functions, and I can't really figure out how to scale their output: I use the following (very short) code to approximate a function:

function y = f(x)
    y = x .^ 2;
endfunction;

X=[-4096:4095]/64;
Y = f(X);
# plot(X, Y);

F = fft(Y);
S = [0:2047]/2048;

function points = approximate(input, count)
    size    = size(input)(2);
    fourier = [fft(input)(1:count) zeros(1, size-count)];
    points  = ifft(fourier);
endfunction;

Y = f(X); plot(X, Y, X, approximate(Y, 10));

Basically, what it does is take a function, compute the image of an interval, fft-it, then keep a few harmonics, and ifft the result. Yet I get a plot that is vertically compressed (the vertical scale of the output is wrong). Any ideas?

分享到QQ

分享到微博

如果你对这篇内容有疑问，欢迎到本站社区发帖提问参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

发布评论

需要登录才能够评论，你可以免费注册一个本站的账号。

反目相谮 2024-09-07 14:23:46

您正在放弃转换的后半部分。对于实值输入，变换是埃尔米特对称的，您必须保留这些线。试试这个：

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;

由于数值计算误差，逆变换总是会产生一些微小的虚部，因此是实数提取。

注意，input和size是Octave中的关键字；用你自己的变量破坏它们是一个很好的方法，可以在路上得到非常奇怪的错误！

You are throwing out the second half of the transform. The transform is Hermitian symmetric for real-valued inputs and you have to keep those lines. Try this:

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;

The inverse transform will invariably have some tiny imaginary part due to numerical computation error, hence the real extraction.

Note that input and size are keywords in Octave; clobbering them with your own variables is a good way to get really weird bugs down the road!

回复收藏 0 原文

两相知 2024-09-07 14:23:46

你可能做错了。您删除代码中的所有“负”频率。您应该保留正低频和负低频。这是 python 代码和结果。该情节具有正确的规模。

替代文本 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))

You are probably doing it wrong. You remove all the "negative" frequencies in your code. You should keep both positive and negative low frequencies. Here is a code in python and the result. The plot has the right scale.

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

The code:

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))

回复收藏 0 原文

~没有更多了~