Преобразование Фурье в Python 2D - PullRequest
0 голосов
/ 13 мая 2018

Я хочу выполнить численно Fourier transform из Gaussian function, используя fft2.При этом преобразовании функция сохраняется до постоянной.

Я создаю 2 сетки: одна для real space, вторая для frequency (импульс, k и т. Д.).(Частоты смещены в ноль).Я оцениваю функции и в конечном итоге отображаю результаты.

Вот мой код

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

Итак, графики для gaussian, fourier(gaussian), inverse_fourier(fourier(gaussian)) следующие: Initial , Фурье , Обратный Фурье

Используя plt.imshow(), я дополнительно нанесу фурье на гауссиану:

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

Результат выглядит следующим образом: imshow

Это не имеет смысла.Я ожидаю увидеть то же самое gaussian function, что и начальное, вплоть до некоторого постоянного порядка единства.

Я был бы очень рад, если бы кто-то мог прояснить это для меня.

1 Ответ

0 голосов
/ 13 мая 2018

Я думаю, вы немного озадачены формой вашего вывода F.В частности, вы можете спросить, почему вы видите такой резкий пик, а не широко распространенный гауссиан.

Я немного изменил ваш код:

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

Обратите внимание, что я ввел sigma параметр для контроля ширины гауссиана.Теперь я приглашаю вас поиграть со следующими параметрами: N_x и N_y, d_x и d_y и sigma.

Затем вы должны увидеть обратное поведение гауссиана в реальном пространстве и в пространстве Фурье: чем больше гауссиан в реальном пространстве, тем уже в пространстве Фурье и наоборот.

Таким образом, с текущими параметрами, установленными в моем коде, вы получаете следующие графики:

Реальное пространство: enter image description here

Пространство Фурье: enter image description here

...