Преобразование Фурье гауссиана не является гауссовским, но это неправильно! - питон

Question

Я пытаюсь использовать функцию fft от Numpy, однако, когда я даю функции простую гауссовскую функцию, fft этой гауссовой функции не является гауссовой, она близка, но делится пополам, так что каждая половина находится на любом конце оси x .

Гауссовская функция, которую я вычисляю, y = exp (-x ^ 2)

Вот мой код:

from cmath import *
from numpy import multiply
from numpy.fft import fft
from pylab import plot, show

""" Basically the standard range() function but with float support """
def frange (min_value, max_value, step):
    value = float(min_value)
    array = []
    while value < float(max_value):
        array.append(value)
        value += float(step)
    return array


N = 256.0 # number of steps
y = []
x = frange(-5, 5, 10/N)

# fill array y with values of the Gaussian function   
cache = -multiply(x, x)
for i in cache: y.append(exp(i))

Y = fft(y)

# plot the fft of the gausian function
plot(x, abs(Y))
show()

Результат не совсем верный, потому что БПФ гауссовой функции должно быть самой гауссовой функцией ...

Answer 1

np.fft.fft возвращает результат в так называемом «стандартном порядке»: ( из документов )

Если A = fft(a, n), то A[0] содержит термин с нулевой частотой ( среднее значение сигнала), который всегда чисто реальный для реальных входов. затем A[1:n/2] содержит положительно-частотные термины и A[n/2+1:] содержит термины с отрицательной частотой, в порядке убывающая отрицательная частота.

Функция np.fft.fftshift упорядочивает результат в порядке, ожидаемом большинством людей (и который хорош для построения графиков):

Рутина np.fft.fftshift(A) сдвиги преобразований и их частоты поставить нулевую частоту компоненты в середине ...

Итак, используя np.fft.fftshift:

import matplotlib.pyplot as plt
import numpy as np

N = 128
x = np.arange(-5, 5, 10./(2 * N))
y = np.exp(-x * x)
y_fft = np.fft.fftshift(np.abs(np.fft.fft(y))) / np.sqrt(len(y))
plt.plot(x,y)
plt.plot(x,y_fft)
plt.show()

Answer 2

Ваш результат даже не близок к гауссовскому, даже одно деление на две половины.

Чтобы получить ожидаемый результат, вам нужно расположить свой собственный гауссиан с центром в индексе 0, и результат также будет позиционироваться таким образом. Попробуйте следующий код:

from pylab import *
N = 128
x = r_[arange(0, 5, 5./N), arange(-5, 0, 5./N)]
y = exp(-x*x)
y_fft = fft(y) / sqrt(2 * N)
plot(r_[y[N:], y[:N]])
plot(r_[y_fft[N:], y_fft[:N]])
show()

Команды построения разбивают массивы на две половины и меняют их местами, чтобы получить более красивую картинку.

Answer 3

Отображается с центром (то есть средним значением) с нулевым индексом коэффициента.Вот почему кажется, что правая половина слева, и наоборот.

РЕДАКТИРОВАТЬ: Исследуйте следующий код:

import scipy
import scipy.signal as sig
import pylab
x = sig.gaussian(2048, 10)
X = scipy.absolute(scipy.fft(x))
pylab.plot(x)
pylab.plot(X)
pylab.plot(X[range(1024, 2048)+range(0, 1024)])

Последняя строка будет отображаться X, начиная сцентр вектора, затем переходите к началу.

Answer 4

Исходя из ответа Свена Марнаха, более простая версия будет такой:

from pylab import *
N = 128

x = ifftshift(arange(-5,5,5./N))

y = exp(-x*x)
y_fft = fft(y) / sqrt(2 * N)

plot(fftshift(y))
plot(fftshift(y_fft))

show()

Это дает график, идентичный приведенному выше.

Ключ (и мне это кажется странным) заключается в том, что предполагаемый порядок данных NumPy - во обеих частотных и временных областях - должен иметь "ноль" значение в первую очередь. Это не то, что я ожидаю от других реализаций FFT, таких как библиотеки FFTW3 в C.

Это было слегка обманчиво в ответах от unutbu и Стива Тджоа выше, потому что они берут абсолютное значение БПФ перед его построением, таким образом стирая фазовые проблемы, возникающие из-за неиспользования «стандартного порядка» во времени.

Answer 5

Преобразование Фурье неявно повторяется бесконечно, так как это преобразование сигнала, которое неявно повторяется бесконечно. Обратите внимание, что когда вы передаете y для преобразования, значения x не передаются, поэтому на самом деле преобразованный гауссиан равен центру от медианного значения от 0 до 256, поэтому 128.

Помните также, что перевод f (x) - это изменение фазы F (x).