Filterung eines gestörten Signals

[Toben A. (toben_a)]

Hallo,

ich habe vorliegendes (habe ich angehängt) Signal gemessen. Man sieht 
einen leicht nach oben wandernden Sinus mit ziemlich vielen 
unerwünschten Peaks nach oben und unten. Der wandernde Sinus ist mein 
Nutzsignal und hat eine ziemlich konstante Amplitude und Frequenz bzw. 
wie man sieht ändern sich beide nur sehr langsam. Das Störsignal jedoch 
ist wie ersichtlich meistens sehr hochfrequent mit großer amplitude.

Kann mir jemand hinweise geben mit was für einem Filter ich die 
unerwünschten Störsignale herausfiltern kann? Klar, man könnte einen 
Tiefpassfilter nehmen, ich finde bei diesem Fall nur interessant, dass 
Amplitude -und- Frequenz des Nutzsignales annähernd konstant sind. 
Machen sich dies komplexere Filter zu nutze? hier bräuchte ich Hinweise 
und Schlagworte, wonach ich suchen muss.

In Schritt 1 verwende ich Ipython bzw. Matlab um einen Filter zu 
implementieren und die Daten zu analysieren, hier kann ich also leicht 
komplexere filter realisieren. Für Schritt 2 verwende ich TwinCat, hier 
bin ich nach meinem Verständnis auf simple Mittelwertfilter o.ä. 
beschränkt. Hinweise hierzu sind ebenfalls willkommen

Ich danke vielmals

Toben

[A. F. (chefdesigner)]

Ich könnte Dir eine SW anbieten, die das löst. Wäre in C und prozessiert 
Amplitude und Phase des Sinus mit einer Latenz von 2 Perioden. Ist nach 
etwa 4-6 Perioden stabil.

Permanent appliziert (also auf immer eine ausreichende Anzahl von 
Samples, damit auch wirklich 2 volle Perioden der Frequenz drin sind) 
ist das gegen jede Form von Störung iummun. Ist ein diskretes 
Filtermodell.

Den echten Sinus kann man daraus auch regenerieren, wenn man möchte. Er 
passt sich den Ergebnissen an und geht dynamisch mit Frequenzänderungen 
mit, wenn sie nicht zu sprunghaft sind.

So ca 5Hz/s sind kein Problem.

[Joe F. (easylife)]

Zunächst würde ich es mit einem Medianfilter (Rangordnungsfilter) 
versuchen.
Tiefpass höherer Ordnung könnte auch sehr vielversprechend sein.

Ich hätte schon längst ein paar Filter ausprobiert und mir das Resultat 
angeguckt. Ist ja in der Regel eine Sache von ein paar wenigen Minuten.

Eine weitere Frage wäre natürlich, warum sieht das Signal überhaupt so 
aus?
Evtl. könnte man da bereits auf Hardwareebene die Störungen 
beseitigen...

[Toben A. (toben_a)]

erstmal danke für euren Input.

@An Fi: ja das klingt doch interessant. An soetwas hatte ich auch schon 
gedacht: einen Sinus zB per Least-Squares über das echte signal legen 
und alles was zuweit abweicht rauslöschen bzw. der langsamen 
Signaländerung folgen. Könntest du mir dienen Code zukommenlassen? PN ? 
gibt es einen gesonderten Begriff für diesen Filter, mit dem ich 
weiterrecherchieren kann?

@Joe: wir tun schon unser möglichstes, die Störungen gar nicht erst auf 
das Signal kommen zu lassen, ganz unterdrücken lassen wird es sich aber 
nicht. Ich habe bereits die von dir angesprochenen Filter verwendet und 
es funktioniert soweit ok. einen Medianfilter könnte man auch mit 
TwinCat gut implementieren.

Wie gesagt interessiert mich vor allem, ob sich aus den meiner Meinung 
nach besonderen Eigenschaften des Nutzsignales noch Vorteile ziehen 
lassen bzw. ob es Ansätze gibt, diese zu nutzen.

[Joe F. (easylife)]

Toben A. schrieb:
> Wie gesagt interessiert mich vor allem, ob sich aus den meiner Meinung
> nach besonderen Eigenschaften des Nutzsignales noch Vorteile ziehen
> lassen

Das kommt halt darauf an, wie diese Eigenschaften genau sind, und welche 
Daten/Eigenschaften du aus dem Nutzsignal herausziehen möchtest.

Wenn du sagst, Amplitude und Freqzenz ist immer annähernd gleich, worin 
besteht dann überhaupt die Information, die dieses Signal in sich trägt?

[Toben A. (toben_a)]

man sieht das der Mittelwert des "Sinus" sich langsam nach oben bewegt. 
diese Änderung des Mittelwertes sowie die Amplitude mit der der Sinus 
schwingt sind meine Informationen. Wie man im Plot sieht ändert sich 
beides, langsam.

[Joe F. (easylife)]

Ich würde es mal klassisch aus einer Kombination aus Medianfilter und 
steilem Tiefpass/Bandpass versuchen.

Der Aussetzer, bzw. zeitlich längerer Offset des Messwertes am Ende 
deiner Messreihe könnte allerdings zum Problem werden (siehe roter 
Pfeil).
Lade doch mal die Rohdaten einer Messung hoch, dann könnte man sich das 
im Detail genauer angucken.

[Toben A. (toben_a)]

voila.

ich bekomme nun für folgendes sehr schöne Ergebnisse:

1)Signal fenstern und Mittelwert bilden
1.1)Wenn das delta zwischen Laser(i) und Laser(i+1) größer als 0.009 
ist, überschreibe Laser(i+1) mit diesem Mittelwert

Hier nutze ich ganz gezielt, dass Amplitude und Frequenz des 
Nutzsignales bekannt sind. Diese Lösung passt mir sehr gut, weil ich sie 
auch einfach mit TwinCat umsetzen kann.

Weitere (komplexere :) ) Vorschläge sind gern willkommen.

[Joe F. (easylife)]

Au weia sieht das Signal schlimm aus.
Die Amplitude der Störungen übersteigt die Amplitude des Nutzsignals ja 
mindestens um Faktor 10, und die Periodendauer ist recht nahe an der des 
Nutzsignals (ca. Faktor 2).
Bin ehrlich gesagt etwas ratlos, wie man das anständig filtern könnte.
Bist du sicher, dass alle hardwareseitigen Möglichkeiten (Schirmung, 
symmetrisches Messignal, Pegelerhöhung (aktive Verstärkung am Sensor)) 
ausgeschöpft wurden?

[Raymund H. (raymund_h753)]

Es sind also der DC und die Amplitude eines Sinus mit konstanter 
Frequenz gefragt.

Der DC kann wohl am besten über Median/Tiefpass extrahiert werden.

Dann könnte man eine PLL mit sehr kleiner Schleifenbandbreite auf die 
konstante Signalfrequenz locken und damit das Signal synchron IQ 
(Quadratur) demodulieren.

Auf DC des IQ liegt dann die Amplitude des Signals.

Tiefpass + Median Filterung, Absolutwertbildung des IQ gibt die 
Amplitude.

[Raymund H. (raymund_h753)]

Ich wollte etwas mit Jupyter Notebook und Scipy/Numpy arbeiten und habe 
mich an der Datenauswertung versucht:

An den Daten ist auffällig das sie auf 1/3 dezimiert werden können ohne 
Informationsverlust.

Die Frequenz zu Synchrondemodulation habe ich von Hand ermittelt, wenn 
diese leicht abweicht/driftet ist das bei hier ein vernachlässigbarer 
Amplitudenfehler.

%matplotlib nbagg
import numpy
import scipy.signal
import matplotlib.pyplot as plt
my_data = numpy.genfromtxt('/cygdrive/c/Users/Visitor/Downloads/out.csv', delimiter=',')
data=(my_data[:,3])[1:]
decimate=3
freqest=797/len(data)*decimate

ilog=[]
qlog=[]

# decimate data
data=data[numpy.arange(0,len(data)-1,decimate, dtype='int32')]
# differentiate and clip
data=numpy.clip(0.01,-0.01,numpy.diff(data))

#synchronous IQ demodulation
vcophase=0
for sample in data:
    ilog.append(scipy.cos(vcophase*2*scipy.pi)*sample)
    qlog.append(scipy.sin(vcophase*2*scipy.pi)*sample)
    vcophase+=freqest

# gaussian low pass
ilog=scipy.signal.convolve(ilog,scipy.signal.gaussian(51,51/scipy.pi/2))
qlog=scipy.signal.convolve(qlog,scipy.signal.gaussian(51,51/scipy.pi/2))
# median filter
ilog=scipy.signal.medfilt(ilog,333)
qlog=scipy.signal.medfilt(qlog,333)
#plt.plot (scipy.signal.convolve(ilog,scipy.ones(1)))
#plt.plot (scipy.signal.convolve(qlog,scipy.ones(1)))
#plot demodulated amplitude
plt.plot (scipy.sqrt(ilog**2+qlog**2))
plt.grid()
plt.show()

[Toben A. (toben_a)]

Hi,

sorry, ich komme da nicht hinterher. Könntest du etwas klarer formuliert 
erläutern was dein Vorschlag ist?

Vielen Dank!

Beste Grüße und guten Rutsch :)

[Raymund H. (raymund_h753)]

Hier noch mal ein schönerer Python Script, mit Kommentaren.

Ich lasse den in Ipython Jupyter Notebook laufen:

https://ipython.org/notebook.html

Scheint brauchbar zu sein, bin von GNU Octave abgekommen.

%matplotlib nbagg
import numpy
import math
import scipy.signal as signal
import matplotlib.pyplot as plt

# read data
my_data = numpy.genfromtxt('/cygdrive/c/Users/Visitor/Downloads/out.csv', delimiter=',')

# extract data column
data=(my_data[:,3])[1:]

# decimate data, looks like 3x sample replication
decimate=3
data=data[numpy.arange(0,len(data)-1,decimate, dtype='int32')]

# differentiate and clip, plot data
data=numpy.clip(0.012,-0.012,numpy.diff(data))
plt.plot(data)
plt.grid()
plt.savefig("data.svg")
plt.show()

# plot fft of windowed data
plt.figure()
plt.grid()
plt.plot(20*numpy.log10(numpy.absolute(scipy.fft(data*signal.blackman(len(data)))))[0:int(len(data)/2)])
plt.savefig("freq.svg")
plt.show()

#synchronous IQ demodulation
periods=797 # frequency in periods in the samples looked up from fft plot
lastphase=periods*(len(data)-1)/len(data)*2*math.pi

ilog=scipy.cos(numpy.linspace(0,lastphase, len(data)))*data
qlog=scipy.sin(numpy.linspace(0,lastphase, len(data)))*data

# gaussian low pass of IQ
ilog=scipy.signal.convolve(ilog,signal.gaussian(51,51/scipy.pi/2))
qlog=scipy.signal.convolve(qlog,signal.gaussian(51,51/scipy.pi/2))

# median filter of IQ
ilog=scipy.signal.medfilt(ilog,333)
qlog=scipy.signal.medfilt(qlog,333)

#plot demodulated amplitude
plt.figure()
plt.plot (scipy.sqrt(ilog**2+qlog**2))

#plot DC
plt.plot(100*scipy.signal.medfilt(scipy.convolve(data,scipy.signal.gaussian(101,101/scipy.pi/2)),2001))

plt.grid()
plt.savefig("rslt.svg")
plt.show()

[Kritiker (Gast)]

yeah, Python rulez!