Hallo allerseits,
ich hatte folgenden Text bei Stackoverflow versucht zu platzieren (deshalb auf Englisch), weil ich dieses Forum noch nicht kannte. Antworten gerne auf Deutsch. Auf Wunsch kann ich den Text auch noch übersetzen.
We need to analyse voltages (UL1L2, UL2L3, UL3L1) and currents (IL1L2, IL2L3, IL3L1) which will be importet frome a static file. We want to know things like frequencies of the fundamental, amplitudes, phase differences between voltages an currents, active power, reactive power, complex powers, frequency analysis (Fourier spectrum) and all the typcal things you learn at university.
Of course, you can look up the formulas in a textbook and programm all the functions you need, but I thought this should really be a stand problem and hence there might be ready to use functions in some Python library.
However, using Google didn't result in the desired result. For example I ended at PyPSA, but this seems more like something to analyze electrical networks or something like that and not a simple tool to evalute Measured data.
So, does anyone know a library as described above?
Gruß
Florian
Python Bibliothek um gemessene Spannungen und Ströme von 3 Phasigen Systemen zu analysieren
Ja, habe gesehen, das scipy einiges mitbringt wie beispielsweise curve fitting, aber bis man dann bei dem ist was ich suche muss man doch noch Arbeit reinstecken und da es wirklich ein Standardproblem ist, dachte ich das gibts vielleicht schon.
Also nach längerem rumsuchen, habe ich diverse Ansätze gefunden, als trivial würde ich das jedoch nicht bezeichnen.
Einfaches Beispiel:
Dies ergibt genau das Frequenzspektrum mit Peaks bei 50 Hz, 100Hz und 150 Hz.
Aber wie erhalte ich nun die Phase, sprich wie müssten die Codezeilen aussehen, die mir den Wert 0.15 zumindest näherungsweise liefern?
Einfaches Beispiel:
Code: Alles auswählen
from scipy.fft import fft, fftfreq
import matplotlib.pyplot
import numpy
# Messwerte simulieren
N=10001 # Anzahl Messpunkte
tMax=5 # Maximale Zeit in s
f=50 # Netzfrequenz in Hz
omega=2*numpy.pi*f # Umrechnung in Omega
phase = 0.15 # Phasenverschiebung
t = numpy.linspace(0, tMax, N) # Array mit "Zeitstützstellen" der Messung
U0L1 = 325 * numpy.sin(omega*t+phase) + 5 * numpy.sin(2*omega * t) + 3 * numpy.sin(3*omega * t) #Spannungswerte mit ein paar oberwellen erzeugen
# Fouriertransformation berechnen
U0L1f = fft(U0L1)
tf = fftfreq(N, tMax/N)[:N//2]
Aber wie erhalte ich nun die Phase, sprich wie müssten die Codezeilen aussehen, die mir den Wert 0.15 zumindest näherungsweise liefern?
Die Phaseninformation steckt in den komplexen Werten der FFT. Mit fftfreq holst du dir da ja nur die Frequenzen raus, nicht die Phasen von denen. Hier zB ein Beispiel, welches anschaulich ist, aber leider in matlab: https://www.gaussianwaves.com/2015/11/i ... formation/ - doch das sollte auch ergoogelbar sein mit numpy/matplotlib.
So, weil ich meine eigene Intuition um Signalverarbeitung verbessern will, habe ich nochmal selbst damit gespielt. Und ich muss mich da korrigieren, soooo simpel ist es tatsaechlich nicht.
Hier erstmal ein Stueck Code, dass das ganze illusrtiert:
Ist nicht vorbildlich, aber erstmal ein bisschen sauberer. Im Kommentar ist der Link, der erklaert, was hier eigentlich abgeht: die Phaseninformation wird negativ beeinflusst durch die "offenen" Enden des Signals. Die exakte Phasenlage bekommt man in meinem Script, wenn man mit USE_PERIODIC = True ein Signal erzeugt, welches periodengenau ist. Mit USE_WINDOW kann man noch ein Hanning-Window dazu aktivieren.
Wenn ich das laufen lasse, bekomme ich 1a-Werte fuer den "synthetischen" Fall, aber dein Signal liegt um 5 Grad daneben. Ob das fuer den gewuenschten Anwendungsfall ausreichend ist, ober nochmal genauer muss, musst du entscheiden. Ich werde auch noch mal versuchen, das Ergebnis zu verbessern.
Hier erstmal ein Stueck Code, dass das ganze illusrtiert:
Code: Alles auswählen
# https://stackoverflow.com/questions/54454723/scipy-fft-how-to-get-phase-angle
from scipy.fft import fft, fftfreq
from scipy.signal.windows import hann
import matplotlib.pyplot as plt
import numpy as np
USE_PERIODIC = False
USE_WINDOW = True
def plot_signal(signal, tf, fft_plot, phase_plot, c="blue"):
if USE_WINDOW:
signal = signal * hann(len(signal))
signal_fft = fft(signal)
highest_frequency_index = np.argmax(np.abs(signal_fft))
highest_frequency = tf[highest_frequency_index]
threshold = max(np.abs(signal_fft)) / 10000
clean_fft = np.array(signal_fft)
clean_fft[abs(signal_fft) < threshold] = 0
phase = np.arctan2(np.imag(clean_fft), np.real(clean_fft)) * 180 / np.pi
highest_phase = phase[highest_frequency_index]
print(highest_frequency, highest_phase)
fft_plot.plot(tf, 2.0/len(signal_fft) * np.abs(signal_fft), c=c)
phase_plot.set(ylim=(-180, 180))
phase_plot.scatter(tf, phase, s=1, c=c)
def generate_signals(use_periodic):
phase = np.pi / 4
if use_periodic:
n = 200
# Exakt periodische Signale
t = np.linspace(0, 10, num=n, endpoint=False)
signal = np.cos(2 * np.pi * t)
shifted_signal = np.cos(2 * np.pi * t + phase)
tf = fftfreq(n)
else:
# Deine simulierten Messwerte, allerdings
# ohne die falschen (weil keine Phase) Oberwellen
# und mit cos, weil sonst der phase-angle schon - 90 ist
n = 10001 # Anzahl Messpunkte
tMax = 5 # Maximale Zeit in s
f = 50 # Netzfrequenz in Hz
omega = 2*np.pi*f # Umrechnung in Omega
t = np.linspace(0, tMax, n) # Array mit "Zeitstützstellen" der Messung
shifted_signal = 325 * np.cos(omega * t + phase)
signal = 325 * np.cos(omega * t)
tf = fftfreq(n, tMax/n)
return signal, shifted_signal, tf
def main():
signal, shifted_signal, tf = generate_signals(USE_PERIODIC)
fig, ax = plt.subplots(nrows=2, ncols=2)
plot_signal(signal, tf, ax[0, 0], ax[0, 1], c="blue")
plot_signal(shifted_signal, tf, ax[1, 0], ax[1, 1], c="red")
plt.show()
if __name__ == '__main__':
main()
Wenn ich das laufen lasse, bekomme ich 1a-Werte fuer den "synthetischen" Fall, aber dein Signal liegt um 5 Grad daneben. Ob das fuer den gewuenschten Anwendungsfall ausreichend ist, ober nochmal genauer muss, musst du entscheiden. Ich werde auch noch mal versuchen, das Ergebnis zu verbessern.