Teorema de Parseval en Python

Estoy tratando de controlar la funcionalidad fft de Python, y una de las cosas extrañas con las que me he topado es que el teorema de Parseval no parece aplicarse, ya que da una diferencia de aproximadamente 50 ahora, mientras que debería serlo. 0.

import numpy as np import matplotlib.pyplot as plt import scipy.fftpack as fftpack pi = np.pi tdata = np.arange(5999.)/300 dt = tdata[1]-tdata[0] datay = np.sin(pi*tdata)+2*np.sin(pi*2*tdata) N = len(datay) fouriery = abs(fftpack.rfft(datay))/N freqs = fftpack.rfftfreq(len(datay), d=(tdata[1]-tdata[0])) df = freqs[1] - freqs[0] parceval = sum(datay**2)*dt - sum(fouriery**2)*df print parceval plt.plot(freqs, fouriery, 'b-') plt.xlim(0,3) plt.show() 

Estoy bastante seguro de que es un factor de normalización, pero no parece poder encontrarlo, ya que toda la información que puedo encontrar sobre esta función es la documentación de scipy.fftpack.rfft .

Su factor de normalización proviene de intentar aplicar el teorema de Parseval para la transformada de Fourier de una señal continua a una secuencia discreta. En el panel lateral del artículo de wikipedia sobre la transformada discreta de Fourier hay una discusión sobre la relación de la transformada de Fourier, la serie de Fourier, la transformada discreta de Fourier y el muestreo con peines de Dirac .

Para resumir , el teorema de Parseval, cuando se aplica a DFT , no requiere integración, sino sum: un 2*pi que está creando multiplicando por dt y df sus sums.

Tenga en cuenta también que, debido a que está utilizando scipy.fftpack.rfft , lo que está obteniendo no es correctamente la DFT de sus datos, sino solo la mitad positiva de los mismos, ya que la negativa sería simétrica. Entonces, como solo estás agregando la mitad de los datos, más el 0 en el término de DC, faltan los 2 para llegar al 4*pi que encontró @unutbu.

En cualquier caso, si datay mantiene su secuencia, puede verificar el teorema de Parseval de la siguiente manera:

 fouriery = fftpack.rfft(datay) N = len(datay) parseval_1 = np.sum(datay**2) parseval_2 = (fouriery[0]**2 + 2 * np.sum(fouriery[1:]**2)) / N print parseval_1 - parseval_2 

Usando scipy.fftpack.fft o numpy.fft.fft la segunda sum no necesita tomar una forma tan extraña:

 fouriery_1 = fftpack.fft(datay) fouriery_2 = np.fft.fft(datay) N = len(datay) parseval_1 = np.sum(datay**2) parseval_2_1 = np.sum(np.abs(fouriery_1)**2) / N parseval_2_2 = np.sum(np.abs(fouriery_2)**2) / N print parseval_1 - parseval_2_1 print parseval_1 - parseval_2_2