Python 2D中的傅立叶变换 [英] Fourier Transform in Python 2D

问题描述

我想使用 fft2 对 Gaussian函数的傅里叶变换进行数字化处理.在这种转换下，该函数将保留到一个常数.

我创建2个网格:一个用于真实空间，第二个用于 frequency (动量，k等).(频率变为零).我评估函数并最终绘制结果.

这是我的代码

 将numpy导入为np导入matplotlib.pyplot作为plt从scipy.fftpack导入fft2，ifft2从mpl_toolkits.mplot3d导入Axes3D"创建真实和动量空间网格""N_x，N_y = 2 ** 11，2 ** 11range_x，range_y = np.arange(N_x)，np.arange(N_y)dx，dy = 0.005，0.005#真实空间网格向量xv，yv = dx *(range_x-0.5 * N_x)，dy *(range_y-0.5 * N_y)dk_x，dk_y = np.pi/np.max(xv)，np.pi/np.max(yv)#动量空间网格矢量，移到零频率的中心k_xv，k_yv = dk_x * np.append(range_x [:N_x//2]，-range_x [N_x//2:0:-1])，\dk_y * np.append(range_y [:N_y//2]，-range_y [N_y//2:0:-1])#创建真实和动量空间网格x，y = np.meshgrid(xv，yv，sparse = False，indexing ='ij')kx，ky = np.meshgrid(k_xv，k_yv，sparse = False，indexing ='ij')功能"""f = np.exp(-0.5 *(x ** 2 + y ** 2))F = fft2(f)f2 = ifft2(F)"绘图""无花果= plt.figure()ax = Axes3D(图)surf = ax.plot_surface(x，y，np.abs(f)，cmap ='viridis')#我更改为的其他地块#surf = ax.plot_surface(kx，ky，np.abs(F)，cmap ='viridis')#surf = ax.plot_surface(x，y，np.abs(f2)，cmap ='viridis')plt.show() 

因此， gaussian，fourier(gaussian)，inverse_fourier(fourier(gaussian))的图如下:

傅立叶空间:

I want to perform numerically Fourier transform of Gaussian function using fft2. Under this transformation the function is preserved up to a constant.

I create 2 grids: one for real space, the second for frequency (momentum, k, etc.). (Frequencies are shifted to zero). I evaluate functions and eventually plot the results.

Here is my code

import numpy as np
import matplotlib.pyplot as plt
from scipy.fftpack import fft2, ifft2
from mpl_toolkits.mplot3d import Axes3D

"""CREATING REAL AND MOMENTUM SPACES GRIDS"""
N_x, N_y = 2 ** 11, 2 ** 11
range_x, range_y = np.arange(N_x), np.arange(N_y)
dx, dy = 0.005, 0.005
# real space grid vectors
xv, yv = dx * (range_x - 0.5 * N_x), dy * (range_y - 0.5 * N_y)
dk_x, dk_y = np.pi / np.max(xv), np.pi / np.max(yv)
# momentum space grid vectors, shifted to center for zero frequency
k_xv, k_yv = dk_x * np.append(range_x[:N_x//2], -range_x[N_x//2:0:-1]), \
            dk_y * np.append(range_y[:N_y//2], -range_y[N_y//2:0:-1])

# create real and momentum spaces grids
x, y = np.meshgrid(xv, yv, sparse=False, indexing='ij')
kx, ky = np.meshgrid(k_xv, k_yv, sparse=False, indexing='ij')

"""FUNCTION"""
f = np.exp(-0.5 * (x ** 2 + y ** 2))
F = fft2(f)
f2 = ifft2(F)
"""PLOTTING"""
fig = plt.figure()
ax = Axes3D(fig)
surf = ax.plot_surface(x, y, np.abs(f), cmap='viridis')
# for other plots I changed to
# surf = ax.plot_surface(kx, ky, np.abs(F), cmap='viridis')
# surf = ax.plot_surface(x, y, np.abs(f2), cmap='viridis')
plt.show()

So, the plots for gaussian, fourier(gaussian), inverse_fourier(fourier(gaussian)) are the following:Initial, Fourier, Inverse Fourier

Using plt.imshow(), I additionally plot fourier of gaussian:

   plt.imshow(F)
   plt.colorbar()
   plt.show()

The result is as follows: imshow

That doesn't make sense. I expect see the same gaussian function as the initial up to some constant order of unity.

I would be very glad if someone could clarify this for me.

解决方案

I think you are a bit puzzled by the shape of your output F. Especially, you might wonder why you see such a sharp peak and not a wide-spread gaussian.

I changed your code a little bit:

 import numpy as np
 import matplotlib.pyplot as plt
 from scipy.fftpack import fft2, ifft2
 from mpl_toolkits.mplot3d import Axes3D

 """CREATING REAL AND MOMENTUM SPACES GRIDS"""
 N_x, N_y = 2 ** 10, 2 ** 10
 range_x, range_y = np.arange(N_x), np.arange(N_y)
 dx, dy = 0.005, 0.005
 # real space grid vectors
 xv, yv = dx * (range_x - 0.5 * N_x), dy * (range_y - 0.5 * N_y)
 dk_x, dk_y = np.pi / np.max(xv), np.pi / np.max(yv)
 # momentum space grid vectors, shifted to center for zero frequency
 k_xv, k_yv = dk_x * np.append(range_x[:N_x//2], -range_x[N_x//2:0:-1]), \
             dk_y * np.append(range_y[:N_y//2], -range_y[N_y//2:0:-1])

 # create real and momentum spaces grids
 x, y = np.meshgrid(xv, yv, sparse=False, indexing='ij')
 kx, ky = np.meshgrid(k_xv, k_yv, sparse=False, indexing='ij')

 """FUNCTION"""
 sigma=0.05
 f = 1/(2*np.pi*sigma**2) * np.exp(-0.5 * (x ** 2 + y ** 2)/sigma**2)
 F = fft2(f)
 """PLOTTING"""
 fig = plt.figure()
 ax = Axes3D(fig)
 surf = ax.plot_surface(x, y, np.abs(f), cmap='viridis')
 # for other plots I changed to
 fig2 = plt.figure()
 ax2 =Axes3D(fig2)
 surf = ax2.plot_surface(kx, ky, np.abs(F)*dx*dy, cmap='viridis')
 plt.show()

Notice that I introduced a sigma parameter to control the width of the gaussian. I now invite you to play with the following parameters: N_x and N_y, d_x and d_y and sigma.

You should then see the inverse behaviour of gaussian in real-space and in fourier space: The larger the gaussian in real-space, the narrower in fourier-space and vice-versa.

So with the currently set parameters in my code, you get the following plots:

Real space:

Fourier Space:

这篇关于Python 2D中的傅立叶变换的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！