使用python进行二维FFT会使频率稍微偏移 [英] Two dimensional FFT using python results in slightly shifted frequency
问题描述
我知道在python中使用快速傅立叶变换(FFT)方法存在一些问题,但不幸的是,这些问题都无法帮助我解决问题:
I know there have been several questions about using the Fast Fourier Transform (FFT) method in python, but unfortunately none of them could help me with my problem:
我想使用python计算给定二维信号f(即f(x,y))的快速傅立叶变换. Python文档提供了很多帮助,解决了FFT带来的一些问题,但是与我希望它显示的频率相比,我最终仍然会稍微改变一下频率.这是我的python代码:
I want to use python to calculate the Fast Fourier Transform of a given two dimensional signal f, i.e. f(x,y). Pythons documentation helps a lot, solving a few issues, which the FFT brings with it, but i still end up with a slightly shifted frequency compared to the frequency i expect it to show. Here is my python code:
from scipy.fftpack import fft, fftfreq, fftshift
import matplotlib.pyplot as plt
import numpy as np
import math
fq = 3.0 # frequency of signal to be sampled
N = 100.0 # Number of sample points within interval, on which signal is considered
x = np.linspace(0, 2.0 * np.pi, N) # creating equally spaced vector from 0 to 2pi, with spacing 2pi/N
y = x
xx, yy = np.meshgrid(x, y) # create 2D meshgrid
fnc = np.sin(2 * np.pi * fq * xx) # create a signal, which is simply a sine function with frequency fq = 3.0, modulating the x(!) direction
ft = np.fft.fft2(fnc) # calculating the fft coefficients
dx = x[1] - x[0] # spacing in x (and also y) direction (real space)
sampleFrequency = 2.0 * np.pi / dx
nyquisitFrequency = sampleFrequency / 2.0
freq_x = np.fft.fftfreq(ft.shape[0], d = dx) # return the DFT sample frequencies
freq_y = np.fft.fftfreq(ft.shape[1], d = dx)
freq_x = np.fft.fftshift(freq_x) # order sample frequencies, such that 0-th frequency is at center of spectrum
freq_y = np.fft.fftshift(freq_y)
half = len(ft) / 2 + 1 # calculate half of spectrum length, in order to only show positive frequencies
plt.imshow(
2 * abs(ft[:half,:half]) / half,
aspect = 'auto',
extent = (0, freq_x.max(), 0, freq_y.max()),
origin = 'lower',
interpolation = 'nearest',
)
plt.grid()
plt.colorbar()
plt.show()
运行它时我从中得到的是:
And what i get out of this when running it, is:
现在您看到x方向的频率不完全在fq = 3
处,而是稍微向左移动.为什么是这样?
我认为这与事实有关,FFT是一种使用对称参数和
Now you see that the frequency in x direction is not exactly at fq = 3
, but slightly shifted to the left. Why is this?
I would assume that is has to do with the fact, that FFT is an algorithm using symmetry arguments and
half = len(ft) / 2 + 1
用于显示正确位置的频率.但是我不太清楚确切的问题是什么以及如何解决.
is used to show the frequencies at the proper place. But I don't quite understand what the exact problem is and how to fix it.
我也尝试过使用更高的采样频率(N = 10000.0),这不能解决问题,而是将频率稍微向右移了一点.因此,我很确定问题不在于采样频率.
I have also tried using a higher sampling frequency (N = 10000.0), which did not solve the issue, but instead shifted the frequency slightly too far to the right. So i am pretty sure that the problem is not the sampling frequency.
注意:我知道一个事实,泄漏效应会导致此处出现非自然的振幅,但是在这篇文章中,我主要对正确的频率感兴趣.
Note: I'm aware of the fact, that the leakage effect leads to unphysical amplitudes here, but in this post I am primarily interested in the correct frequencies.
推荐答案
我发现了很多问题
您使用2 * np.pi
两次,如果您想要一个整数倍的周期,则应该选择linspace或arg中的一个作为弧度正弦值.
you use 2 * np.pi
twice, you should choose one of either linspace or the arg to sine as radians if you want a nice integer number of cycles
另外np.linspace
默认为endpoint=True
,为您提供101而不是100的加分
additionally np.linspace
defaults to endpoint=True
, giving you an extra point for 101 instead of 100
fq = 3.0 # frequency of signal to be sampled
N = 100 # Number of sample points within interval, on which signal is considered
x = np.linspace(0, 1, N, endpoint=False) # creating equally spaced vector from 0 to 2pi, with spacing 2pi/N
y = x
xx, yy = np.meshgrid(x, y) # create 2D meshgrid
fnc = np.sin(2 * np.pi * fq * xx) # create a signal, which is simply a sine function with frequency fq = 3.0, modulating the x(!) direction
您可以检查以下问题:
len(x)
Out[228]: 100
plt.plot(fnc[0])
现在修复linspace端点意味着您有偶数个fft箱,因此您将+ 1
放在half
计算中
fixing the linspace endpoint now means you have an even number of fft bins so you drop the + 1
in the half
calc
matshow()
似乎具有更好的默认值,您在imshow
中的extent = (0, freq_x.max(), 0, freq_y.max()),
似乎无法使用fft bin编号
matshow()
appears to have better defaults, your extent = (0, freq_x.max(), 0, freq_y.max()),
in imshow
appears to fubar the fft bin numbering
from scipy.fftpack import fft, fftfreq, fftshift
import matplotlib.pyplot as plt
import numpy as np
import math
fq = 3.0 # frequency of signal to be sampled
N = 100 # Number of sample points within interval, on which signal is considered
x = np.linspace(0, 1, N, endpoint=False) # creating equally spaced vector from 0 to 2pi, with spacing 2pi/N
y = x
xx, yy = np.meshgrid(x, y) # create 2D meshgrid
fnc = np.sin(2 * np.pi * fq * xx) # create a signal, which is simply a sine function with frequency fq = 3.0, modulating the x(!) direction
plt.plot(fnc[0])
ft = np.fft.fft2(fnc) # calculating the fft coefficients
#dx = x[1] - x[0] # spacing in x (and also y) direction (real space)
#sampleFrequency = 2.0 * np.pi / dx
#nyquisitFrequency = sampleFrequency / 2.0
#
#freq_x = np.fft.fftfreq(ft.shape[0], d=dx) # return the DFT sample frequencies
#freq_y = np.fft.fftfreq(ft.shape[1], d=dx)
#
#freq_x = np.fft.fftshift(freq_x) # order sample frequencies, such that 0-th frequency is at center of spectrum
#freq_y = np.fft.fftshift(freq_y)
half = len(ft) // 2 # calculate half of spectrum length, in order to only show positive frequencies
plt.matshow(
2 * abs(ft[:half, :half]) / half,
aspect='auto',
origin='lower'
)
plt.grid()
plt.colorbar()
plt.show()
放大了剧情:
这篇关于使用python进行二维FFT会使频率稍微偏移的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!