Tiefpassfilter und FFT für Anfänger mit Python

23

Ich bin neu in der Signalverarbeitung und besonders in der FFT. Daher bin ich mir nicht sicher, ob ich hier das Richtige tue und bin ein bisschen verwirrt mit dem Ergebnis.

Ich habe eine diskrete reelle Funktion (Messdaten) und möchte darauf einen Tiefpassfilter aufbauen. Das Werkzeug der Wahl ist Python mit dem Numpy-Paket. Ich folge dieser Prozedur:

  • berechne die fft meiner Funktion
  • Hochfrequenzen abschneiden
  • Führe das inverse fft aus

Hier ist der Code, den ich verwende:

import numpy as np
sampling_length = 15.0*60.0 # measured every 15 minutes
Fs = 1.0/sampling_length
ls = range(len(data)) # data contains the function
freq = np.fft.fftfreq(len(data), d = sampling_length)
fft = np.fft.fft(data)
x = freq[:len(data)/2] 
for i in range(len(x)):
if x[i] > 0.005: # cut off all frequencies higher than 0.005
    fft[i] = 0.0
    fft[len(data)/2 + i] = 0.0
inverse = np.fft.ifft(fft)

Ist das die richtige Vorgehensweise? Das Ergebnis inverseenthält komplexe Werte, die mich verwirren.

Bis B
quelle
1
Als ich FFT lernte, fand ich diesen Blogeintrag sehr hilfreich. glowingpython.blogspot.com/2011/08/…
David Poole

Antworten:

23

Dass das Ergebnis komplex ist, ist zu erwarten. Ich möchte auf ein paar Dinge hinweisen:

Sie wenden einen Brick-Wall-Frequenzbereich-Filter auf die Daten an, versuchen, alle FFT-Ausgänge auf Null zu setzen, die einer Frequenz von mehr als 0,005 Hz entsprechen, und transformieren sie dann rückwärts, um wieder ein Zeitbereichssignal zu erhalten. Damit das Ergebnis real ist, muss die Eingabe in die inverse FFT konjugiert symmetrisch sein . Dies bedeutet, dass für eine Länge FFTN

X[k]=X[N-k],k=1,2,,N2-1(Neven)

X[k]=X[N-k],k=1,2,,N2(NOdd)
  • Beachten Sie, dass für gerade X [ 0 ] und X [ NNX[0]sind im Allgemeinen nicht gleich, aber beide sind real. Für ungeradesNmussX[0]reell sein.X[N2]NX[0]

Ich sehe, dass Sie in Ihrem obigen Code versucht haben, so etwas zu tun, aber es ist nicht ganz richtig. Wenn Sie die obige Bedingung für das Signal erzwingen, das Sie an die inverse FFT übergeben, sollten Sie ein echtes Signal ausgeben.

sichnc(x)sichnc

sichnc

Handlung der sinc-Funktion

sichncsichnc

Es gibt praktischere Möglichkeiten, Tiefpassfilter sowohl im Zeit- als auch im Frequenzbereich anzuwenden. Finite-Impulsantwort- und Infinite-Impulsantwort- Filter können direkt unter Verwendung ihrer Differenzgleichungsdarstellung angewendet werden . Oder, wenn Sie einen Filter , eine ausreichend lange Impulsantwort hat, können Sie oft Performance - Vorteile mit erhalten schnelle Faltungstechniken auf der Basis der FFT (Anwendung den Filters durch im Frequenzbereich multipliziert wird anstelle der Faltung in der Zeitdomäne), wie die überlagernde Methoden speichern und überlappen .

Jason R
quelle
Die sinc-Funktion ist jedoch die ideale Filterung, nicht wahr? Das ist es, was alle anderen Filter anstreben, aber nicht erreichen. Es ist schlecht für die Bildverarbeitung, weil Bilder nicht zuerst antialiasing sind und daher ein schreckliches Klingeln erzeugen. Aber für Audio- oder andere Signale, die vor dem Sampling mit Antialias gefiltert wurden, ist es nicht der beste Filter, den Sie bekommen können?
Endolith
1
Ja, mein Ergebnis war nicht konjugiert symmetrisch. Ich habe den Code korrigiert, jetzt funktioniert alles gut. Vielen Dank!
Bis
3
@endolith - ein Sinc ist ein idealer Interpolator für bestimmte Arten von Interpolation, kann jedoch bei weitem nicht ideal als Filter für die meisten gängigen Filteranforderungen sein, z. B. für die Flachheit der
Passbandantwort, die Sperrbandunterdrückung
+1 für die nette Erklärung zu "Warum die Leute den Filter nicht implementieren, wie es die PO tut"
Sibbs Gambling
Sie müssen eine Sinc mit Fenster verwenden. Wenn Sie nicht an die Zeit gebunden sind, ist dies der optimale Filter, weitaus besser als Chebichev.