Softwareentwicklung
›
Thread
Einfache Lock-In Verstärkung in Python - wo liegt der Fehler?
T
Tom Gast
09.05.2022 20:42 #7059838
Hallo zusammen
Ich möchte mittels Python eine einfache Lock-In Verstärkung simulieren.
Dazu generiere ich ein 50Hz Signal welches mit der Zeit abnimmt und ein
100Hz Signal welches mit der Zeit zunimmt, siehe angehängtes Bild der
STFT.
Im 3ten Plot stelle ich den Ausgang der Lock-In Verstärkung dar, die
allerdings keinen Sinn ergibt. Der Ausgang der Lock-In Verstärkung
scheint nicht von der Frequenz "flockin" abzuhängen, wenn ich die
Variable "flockin" (Zeile 28) ändere kommt immer das selbe raus.
Sieht vielleicht jemand den Fehler?
1 from scipy.signal import stft
2 import numpy as np
3 from numpy import sin,pi,exp,sqrt
4 import matplotlib.pyplot as plt
5 from scipy import signal
6
7 dt=1e-3
8 T=10
9 ts=dt*np.arange(0,T/dt)
10 fs=1/dt
11
12 f0=50
13 f1=100
14 tau=4
15 sigs=3*exp(-ts/tau)*sin(2*pi*f0*ts)+0.3*exp(ts/tau)*sin(2*pi*f1*ts)
16
17 n_stft=256
18 fs_,ts_,Cs_=stft(sigs,1/dt,nperseg=n_stft)
19
20 #perform lockin amplification
21 flockin=50
22 s_lockin_d=np.sin(2*pi*ts*flockin)
23 s_lockin_q=np.cos(2*pi*ts*flockin)
24 #compute the measured signal with the d and q component
25 s_prod=np.sqrt((s_lockin_d*sigs)**2+(s_lockin_q*sigs)**2)
26 #filter the signal using a low pass filter
27 flowpass=10
28 flowpass_norm=flowpass/(fs/2)
29 b,a=signal.butter(3,flowpass_norm,'low')
30 lockin_output=signal.filtfilt(b,a,s_prod)
31
32 plt.figure(1)
33 plt.pcolormesh(ts_,fs_,np.abs(Cs_))
34 plt.title("$N_{window}=%d$"%(n_stft))
35 plt.xlabel("Time t [s]")
36 plt.ylabel("DTFT[f](t)")
37
38 plt.figure(2)
39 plt.plot(ts,sigs)
40 plt.xlabel("Time t [s]")
41 plt.ylabel("Signal s(t)")
42
43 plt.figure(3)
44 plt.plot(lockin_output)
T
Tom Gast
09.05.2022 21:26 #7059878
Hmm, das Problem scheint zu sein, dass ich zuerst tiefpassfiltern und
erst naher die Norm der d/q Komponenten berechnen muss, der folgende
Code funktioniert richtig:
1 # -*- coding: utf-8 -*-
2
3 from scipy.signal import stft
4 import numpy as np
5 from numpy import sin,pi,exp,sqrt,square
6 import matplotlib.pyplot as plt
7 from scipy import signal
8
9 dt=1e-3
10 T=10
11 ts=dt*np.arange(0,T/dt)
12 fs=1/dt
13
14 f0=50
15 f1=100
16 tau=4
17 sigs=3*exp(-ts/tau)*sin(2*pi*f0*ts)+0.3*exp(ts/tau)*sin(2*pi*f1*ts)
18
19 n_stft=256
20 fs_,ts_,Cs_=stft(sigs,1/dt,nperseg=n_stft)
21
22 #perform lockin amplification
23 flockin=float(100)
24 s_lockin_d=np.sin(2*pi*ts*flockin)
25 s_lockin_q=np.cos(2*pi*ts*flockin)
26 #compute the measured signal with the d and q component
27 s_prod_d=s_lockin_d*sigs
28 s_prod_q=s_lockin_q*sigs
29 #filter the signal using a low pass filter
30 flowpass=10
31 flowpass_norm=flowpass/(fs/2)
32 b,a=signal.butter(3,flowpass_norm,'low')
33 lockin_output_d=signal.filtfilt(b,a,s_prod_d)
34 lockin_output_q=signal.filtfilt(b,a,s_prod_q)
35
36 plt.figure(1)
37 plt.pcolormesh(ts_,fs_,np.abs(Cs_))
38 plt.title("$N_{window}=%d$"%(n_stft))
39 plt.xlabel("Time t [s]")
40 plt.ylabel("DTFT[f](t)")
41
42 plt.figure(2)
43 plt.plot(ts,sigs)
44 plt.xlabel("Time t [s]")
45 plt.ylabel("Signal s(t)")
46
47 plt.figure(3)
48 plt.plot(ts,np.sqrt(np.square(lockin_output_d))+np.square(lockin_output_q))
49 plt.ylabel("Lockin Output")
50 plt.xlabel("Time t [s]")
↑ Nach oben← Zurück zu Softwareentwicklung
Antwort schreiben
Eine eigene Frage?
Dieser Thread ist über 4 Jahre alt.
Stell sie in einem eigenen Thread – dort lesen mehr Leute mit, und sie wird öfter und schneller beantwortet.
Neuen Thread starten
Eine Lösung, Ergänzung oder Rückmeldung zu diesem Thread? Dann antworte hier.
stft_plot.png
(12 KB)
