Python-FFT und PSD Berechnung

Question

Ich habe eine Reihe von Daten, die aus 36002 Elemente besteht und ich möchte FFT und PSD davon zu wissen, welche Frequenz es enthält und entsprechende Leistungsdichte der Frequenz.

Mein Code ist:

from __future__ import division 
import numpy 
import matplotlib.pyplot as plt 

#read in the pressure p_dot and time t, they are [36002,] vector 
nSteps=36002 
p_dot=numpy.genfromtxt((r'E:\p_dot.dat'), delimiter=' ')[:,2] 
t=numpy.genfromtxt((r'E:\t.dat'), delimiter=' ')[:,0] 

T=(t[-1]-t[0])/nSteps # the interval between two data points 
N=len(p_dot)//2+1 # FFT is symmetrical, so plot one half 
Y=numpy.fft.fft(p_dot) # to compare with Yhann 
hann=numpy.hanning(len(p_dot)) 
Yhann=numpy.fft.fft(hann*p_dot) 
fa=1.0/T # scan frequency 
X=numpy.linspace(0, fa/2, N, endpoint=True) # Nyquest frequency=fa/2 
plt.close() 
plt.subplot(2,1,1) 
plt.plot(x,y) 
plt.subplot(2,1,2) 
plt.plot(X, 2.0*abs(Yhann[:N])/N) 
plt.tight_layout() 
plt.show() 

Aber die Ergebnisse erwies sich als falsch, zumindest glaube ich so. Anscheinend sind meine Daten periodisch, aber die FFT hat nur eine Spitze bei 0 Hz. Sehen Sie das Ergebnis der ersten 1000 Artikel.

result 1000 Und das Ergebnis von 36002 result 36002 Also, was ist das Problem? Danke vielmals.

Btw, kann jemand erklären, wie man Überlappungsfenster in Python verwendet? Und wann sollten wir es benutzen? Seit wann ich die Lösung meines Problems suche, sehe ich immer diese Methoden und keine Ahnung, wie man sie benutzt. Vielen Dank!

Answer 1

Der Grund für die Spitze bei 0Hz ist, dass Ihr Eingangssignal einen von Null verschiedenen Mittelwert hat. Dies ist die Gleichstromkomponente des Signals. Es gibt einen anderen, kleineren Peak bei etwa 1Hz oder so, der in etwa die Hauptfrequenz Ihres Signals ist. Die Fensterfunktionen werden normalerweise angewendet, wenn man die Leistung in kleinen Segmenten der Daten berechnet, die dann kombiniert oder gemittelt werden, um das Rauschen bei der Schätzung des Spektrums zu reduzieren. Es ist in Ordnung, es auf das volle Signal anzuwenden, aber es wird weniger häufig auf diese Weise verwendet. Sie können die gesamte Fensterung und Mittelwertbildung manuell durchführen, aber ich würde wirklich empfehlen, die scipy.signal.welch-Funktion (oder eine ähnliche Funktion) zu verwenden, um die PSD zu schätzen. Es kümmert sich um alle Buchhaltung, Fensterung, Mittelwertbildung usw. für Sie.