FFTs in Python mit UliEngineering einfach berechnen & visualisieren
UliEngineering ist eine gemischte Data-Analytics-Bibliothek in Python – eines der Werkzeuge, die sie bereitstellt, ist ein einfach zu verwendendes Paket zur Berechnung von FFTs. Im Gegensatz zu anderen Paketen ist diese Bibliothek auf praktische Anwendungsfälle ausgerichtet und ermöglicht dir, die FFT in nur einer Codezeile durchzuführen! Keine Mathematikkenntnisse erforderlich. Installiere zunächst UliEngineering.
Erste Schritte
Zuerst generieren wir einige Testdaten. Siehe diesen vorherigen Beitrag für weitere Details zur Generierung von Sinus-Testdaten:
from UliEngineering.SignalProcessing.Simulation import *
# Testdaten generieren: 100 Hz + 400 Hz Ton
data = sine_wave(frequency=100.0, samplerate=1000, amplitude=1.) \
+ sine_wave(frequency=400.0, samplerate=1000, amplitude=0.5)Die Testdaten bestehen aus einer 100 Hz-Sinusschwingung plus einer 400 Hz-Sinusschwingung (mit halber Amplitude). Dieses Signal wird mit einer Abtastrate von 1000 Hz abgetastet.
Nun können wir die FFT mit matplotlib berechnen & visualisieren:
# FFT berechnen. Gleiche Abtastrate verwenden
# HINWEIS: Fenster ist standardmäßig "blackman"!
from UliEngineering.SignalProcessing.FFT import compute_fft
fft = compute_fft(data, samplerate=1e3)
# Plot
from matplotlib import pyplot as plt
plt.style.use("ggplot")
plt.gcf().set_size_inches(10, 5) # (20, 10) für einen größeren Plot verwenden
plt.plot(fft.frequencies, fft.amplitudes)
plt.xlabel("Frequenz")
plt.ylabel("Amplitude")compute_fft(data, samplerate=1e3) gibt ein FFT-Objekt zurück, das Felder wie frequencies, amplitude und phase enthält. Es führt eine FFT der Größe des Inputs durch (d.h. da data ein Array der Länge 1000 ist, ist die FFT der Größe 1000).
fft.frequencies ist ein Array von Frequenzen (in Hz), entsprechend den Werten in fft.amplitudes. Du kannst auch fft.angles verwenden, um relative Winkel in Grad zu erhalten, aber dies wird in diesem Blogbeitrag nicht behandelt.
Wie im oben gezeigten Plot zu erkennen ist, ist die maximale Frequenz, die erkannt werden kann, immer die Hälfte der Abtastrate, d.h. für unsere Abtastrate von $f_s = 1000,\text{Hz}$ ist es $500,\text{Hz}$. Siehe diese mathematikorientiertere FFT-Erklärung wenn du weitere Details wissen möchtest.
Intern führt compute_fft() diese Berechnung durch:
- $2$ ist ein Korrekturfaktor, der berücksichtigt, dass wir die hintere Hälfte des rohen FFT-Ergebnisses verwerfen (da wir eine reelle FFT durchführen)
- $\frac{1}{\text{len(data)}}$ normalisiert die FFT-Ergebnisse, sodass sie unabhängig von der Datenlänge sind (d.h. wenn du ein längeres Sample derselben Sinusschwingung übergibst, erhältst du weiterhin dasselbe Ergebnis)
- $\text{abs}\left(\cdots\right)$ konvertiert das komplexe phasenbehaftete Ergebnis der FFT in ein Spektrum, das leichter zu lesen & visualisieren ist.
- Fenster (standardmäßig
blackman) ist das Fenster, das auf die Daten angewendet wird, um einige mathematische Effekte am Anfang und am Ende des Datensatzes zu mildern. Siehe Wikipedia zu Fensterfunktionen für weitere Details. UliEngineering bietet derzeit diese Liste von Fensterfunktionen:blackmanbartletthamminghanningkaiser(Parameter is fixed to2.0)none
Frequenzbereiche auswählen
Mit der UliEngineering-API ist die Auswahl eines Frequenzbereichs der FFT trivial einfach: Verwende einfach fft[lowfreq:highfreq]. Du kannst fft[lowfreq:] verwenden, um alles ab lowfreq auszuwählen, oder fft[:highfreq] verwenden, um alles bis highfreq auszuwählen.
from UliEngineering.SignalProcessing.FFT import compute_fft
fft = compute_fft(data, samplerate=1e3)
# Frequenzbereich auswählen: 50 bis 200 Hz
fft = fft[50.0:200.0]
# Plot
from matplotlib import pyplot as plt
plt.style.use("ggplot")
plt.gcf().set_size_inches(10, 5) # (20, 10) für einen größeren Plot verwenden
plt.plot(fft.frequencies, fft.amplitudes)
plt.xlabel("Frequenz")
plt.ylabel("Amplitude")
plt.savefig("/ram/fft-frequency-range.svg")
Amplitude & Winkel bei einer bestimmten Frequenz extrahieren
Durch die Verwendung von [frequency], d.h. dem Getitem-Operator mit einem einzelnen Wert, erhältst du ein FFTPoint()-Objekt, das die Frequenz, Amplitude und den relativen Winkel für eine gegebene Frequenz enthält. Die Bibliothek wählt automatisch den nächsten FFT-Bucket aus, sodass du auch dann sinnvolle Ergebnisse erhältst, wenn deine FFT keinen Bucket für diese spezifische Frequenz hat.
from UliEngineering.SignalProcessing.FFT import compute_fft
fft = compute_fft(data, samplerate=1e3)
# Wert bei einer bestimmten Frequenz anzeigen
print(fft[30]) # FFTPoint(frequency=30.0, value=2.91e-08, angle=0.0)
print(fft[100]) # FFTPoint(frequency=100.0, value=0.419, angle=0.0)Kurze FFTs für lange Daten
Bei Verwendung von compute_fft() bedeutet ein extrem langes Daten-Array, dass wir eine extrem lange FFT berechnen. In vielen Fällen ist dies nicht erwünscht und du möchtest eine FFT fester Größe berechnen (meistens eine FFT mit Zweierpotenz, z.B. 1024, 2048, 4096 etc).
Lass uns einige lange Testdaten generieren und annehmen, dass wir eine FFT der Größe 1024 darauf berechnen möchten
from UliEngineering.SignalProcessing.FFT import *
fft = simple_serial_fft_reduce(data, samplerate=1e3, fftsize=1024)
# fft.frequencies, fft.amplitudes etcDas Plotten der Daten wie oben ergibt
was fast genau so aussieht wie unser compute_fft()-Plot zuvor – genau wie erwartet.
simple_serial_fft_reduce() übernimmt für uns die gesamte Magie und Normalisierung, einschließlich der Partitionierung der Daten in überlappende Chunks, der Addition der FFTs und der korrekten Normalisierung des Ergebnisses.
Die Namenskonvention ist hier von Bedeutung:
simple_...._reducebedeutet, dass dies eine Variante mit sinnvollen Standardwerten (simple) für die Verwendung einerreduction-Funktion (Standard:sum) auf mehreren FFTs ist.serialbedeutet, dass die einzelnen FFTs
FFTs parallelisieren
Wenn du einen riesigen Datensatz hast, kannst du simple_parallel_fft_reduce() identisch zu simple_serial_fft_reduce() verwenden:
from UliEngineering.SignalProcessing.FFT import *
fft = simple_parallel_fft_reduce(data, samplerate=1e3, fftsize=1024)
# fft.frequencies, fft.amplitudes etc verwendenIn den meisten Fällen möchtest du den Executor jedoch manuell initialisieren, damit du ihn später wiederverwenden kannst:
from UliEngineering.SignalProcessing.FFT import *
from concurrent.futures import ThreadPoolExecutor
executor = ThreadPoolExecutor() # Kein Argument => verwende num_cpus Threads
fft = simple_parallel_fft_reduce(data, samplerate=1e3, fftsize=1024, executor=executor)Wir können einen ThreadPoolExecutor() verwenden, da scipy.fftpack (das UliEngineering für die schwere Mathematik verwendet) das Python-GIL freigibt.
Beachte, dass simple_parallel_fft_reduce() aufgrund der vielen Verwaltungsaufgaben viel langsamer ist als simple_serial_fft_reduce(), wenn dein Datensatz so klein ist, dass Parallelisierung nicht effektiv ist. Meine anfängliche Empfehlung ist, die parallele Variante in Betracht zu ziehen, wenn die gesamte Ausführungszeit der seriellen Variante größer als $0,5s$ ist.