Я пытаюсь создать случайное поле Гаусса, создав сетку в пространстве Фурье, а затем обратный Фурье, переводя ее, чтобы получить случайное поле. Для этого обратное преобразованное Фурье изображение должно быть вещественным. Кажется, я получаю остатки в мнимой части сетки порядка 10^-18 - -22, поэтому я ожидал, что это будут числовые ошибки в БПФ. Реальная часть изображения отображает странный шаблон шахматной доски на пикселях, хотя пиксели переходят с положительного на отрицательный. Чтобы проверить, правильно ли функционирует FFT, я попытался преобразовать гауссовский язык, который должен вернуть другой гауссов, и снова изображение шахматной доски присутствует в изображении. Принимая абсолютное значение изображения, оно выглядит отлично, но мне также нужно, чтобы он учитывал отрицательные значения для моего гауссовского случайного поля.Рисунок шахматной доски после FFT
Для преобразования Фурье гауссовой я использую следующий код:
#! /usr/bin/env python
import numpy as n
import math as m
import pyfits
def fourierplane(a):
deltakx = 2*a.kxmax/a.dimkx #stepsize in k_x
deltaky = 2*a.kymax/a.dimky #stepsize in k_y
plane = n.zeros([a.dimkx,a.dimky]) #empty matrix to be filled in for the Fourier grid
for y in range(n.shape(plane)[0]):
for x in range(n.shape(plane)[1]):
#Defining coordinates centred at x = N/2, y = N/2
i1 = x - a.dimkx/2
j1 = y - a.dimky/2
#creating values to fill in in the grid:
kx = deltakx*i1 #determining value of k_x at gridpoint
ky = deltaky*j1 #determining value of k_y at gridpoint
k = m.sqrt(kx**2 + ky**2) #magnitude of k-vector
plane[y][x] = m.e**(-(k**2)/(2*a.sigma_k**2)) #gaussian
return plane
def substruct():
class fougrid:
pass
grid = fougrid()
grid.kxmax = 2.00 #maximum value k_x
grid.kymax = 2.00 #maximum value k_y
grid.sigma_k = (1./20.)*grid.kxmax #width of gaussian
grid.dimkx = 1024
grid.dimky= 1024
fplane = fourierplane(grid) #creating the Fourier grid
implane = n.fft.ifftshift(n.fft.ifft2(fplane)) #inverse Fourier transformation of the grid to get final image
##################################################################
#seperating real and imaginary part of the Fourier transformed grid
##################################################################
realimplane = implane.real
imagimplane = implane.imag
#taking the absolute value:
absimplane = n.zeros(n.shape(implane))
for a in range(n.shape(implane)[0]):
for b in range(n.shape(implane)[1]):
absimplane[a][b] = m.sqrt(implane[a][b].real**2 + implane[a][b].imag**2)
#saving images to files:
pyfits.writeto('randomfield.fits',realimplane) #real part of the image grid
pyfits.writeto('fplane.fits',fplane) #grid in fourier space
pyfits.writeto('imranfield.fits',imagimplane) #imaginary part of the image grid
pyfits.writeto('absranfield.fits',absimplane) #real part of the image grid
substruct() #running the script
Кто-нибудь есть идеи, как создается эта картина, и как решить эту проблему?
Решил. В коде, до того, как ifft2 выполняется на fplane, fplane должна быть смещена, так что: implane = n.fft.ifftshift (n.fft.ifft2 (n.fft.fftshift (fplane))) – Mizuti