IIR-Filter

Bei der Verarbeitung von Sensorwerten mit einem Mikrocontroller treten häufig Probleme auf, da die Sensorwerte verrauscht oder mit Störimpulsen überlagert sind. In der Analogtechnik würde man ein analoges Filter, beispielsweise ein RC-Glied als Tiefpass, einsetzen. In der Digitaltechnik gibt es entsprechend digitale Filter. Ein RC-Glied in der Analogtechnik entspricht einem digitalen IIR-Filter (Infinite Impulse Response). Ein solches Filter entsteht, indem ein neuer Sensorwert FneuF_{neu} mit einem aus den vergangenen Sensorwerten gespeicherten Wert FaltF_{alt} verrechnet wird und der berechnete Wert für den nächsten Sensorwert gespeichert wird

Theorie

Die Formel für einen IIR-Tiefpass 1. Ordnung lautet:
Falt=(Fneu⋅k2)+(Falt⋅k1)F_{alt}=\left(F_{neu}\cdot k_2\right)+\left(F_{alt}\cdot k_1\right)
Da im Durchlassbereich keine Veränderung der Amplitude erfolgen soll, muss k1+k2=1k_1+k_2=1 sein. Damit gilt dann
Falt=Falt+(Fneu−Falt)⋅k2F_{alt}=F_{alt}+\left(F_{neu}-F_{alt}\right)\cdot k_2
Bei einem 8-Bit-Microcontroller wie dem Arduino UNO ist eine Berechnung mit Gleitkommazahlen sehr zeitaufwendig. Auch eine Multiplikation oder Division mit ganzzahligen Werten ist „teuer”. Beides sollte vermieden werden. Um auf ganzzahlige Werte zu kommen, wird die obige Gleichung mit einem Faktor n multipliziert, sodass k2⋅n=1k_2\cdot n=1 wird:
Falt⋅n=Falt⋅n+(Fneu−Falt)⋅k2⋅n=Falt⋅n+(Fneu−Falt)F_{alt}\cdot n=F_{alt}\cdot n+\left(F_{neu}-F_{alt}\right)\cdot k_2\cdot n=F_{alt}\cdot n+\left(F_{neu}-F_{alt}\right)
oder
Falt=(Falt⋅n+(Fneu−Falt))/nF_{alt}=\left(F_{alt}\cdot n+\left(F_{neu}-F_{alt}\right)\right)/n
Mit etwas geringerer Genauigkeit können wir auch mit der folgenden Gleichung rechnen:
Falt=Falt+(Fneu−Falt)/nF_{alt}=F_{alt}+\left(F_{neu}-F_{alt}\right)/n
Wenn der Faktor n nun eine 2er-Potenzen ist (n=2sn=2^{s}), dann kann die Multiplikation oder Division auch vermieden und durch Shift-Operationen um s Bits ersetzt werden. Dies schränkt allerdings die Wahl der Grenzfrequenz des Filters sehr stark ein. Diese Grenzfrequenz FgrenzF_{grenz} des Filters berechnet sich in Anhängigkeit des Sampling Frequenz fsamplf_{sampl} nach der Formel:
k2=1−e−2π⋅(fgrenz/fsampl)k_2=1-e^{-2\pi\cdot(f_{grenz}/f_{sampl})} Hier einige praxisnahe Werte:

k2k_2 n s fgrenz≈f_{grenz}\approx
0,5 2 1 fsampl/10f_{sampl}/10
0,25 4 2 fsampl/20f_{sampl}/20
0,125 8 3 fsampl/50f_{sampl}/50

Implementierung

Für die Implantierung betrachten wir die folgende Formel, bei der Shift-Operatoren anstelle von Multiplikation und Division verwendet werden:
(Falt⋅≪s)=(Falt⋅≪s)+(Fneu−Falt)(F_{alt}\cdot \ll s)=(F_{alt}\cdot \ll s)+\left(F_{neu}-F_{alt}\right)
Wie man sieht, berechnen wir den Wert Falt≪sF_{alt}\ll s, verwenden aber nachher in der Subtraktion den Wert FaltF_{alt}. Dadurch leidet die Genauigkeit. Es empfiehlt sich daher, die Subtraktion mit höherer Genauigkeit zu berechnen. Wir kommen dann zu folgender Formel:
(Falt≪s)=(Falt≪s)+(((Fneu≪s)−(Falt≪s))≫s)\left(F_{alt}\ll s\right)=\left(F_{alt}\ll s\right)+\left(\left(\left(F_{neu}\ll s\right)-\left(F_{alt}\ll s\right)\right)\gg s\right)
Es liegt jetzt nahe, statt dem Wert FaltF_{alt} für die nächste Berechnung den Wert Falt≪sF_{alt}\ll s zu speichern.

Überläufe bei der Berechnung müssen unbedingt vermieden werden. Man sieht sofort, dass die Breite der Eingabewerte (z.B. Sensorwerte) plus dem Shiftfaktor s die Wortbreite der Berechnung (beispielswese 16-Bit Arithmetik) nicht überschreiten darf. Kritisch ist aber wieder die Subtraktion. Hier kann es bei einem Sprung vom höchsten auf den niedrigsten Wert zu einem Überlauf kommen. Wir benötigen also noch ein Bit in Reserve. Bei einer angenommenen 16-Bit-Arithmetik können aber mit 12-Bit-Sensorwerten aber immer noch die Shiftfaktoren 1 bis 3 verwendet werden. Entsprechende Überlegungen können auch für andere Wertebereiche und Shiftfaktoren angestellt werden.

Damit kommen wir jetzt zu dieser sehr einfachen und schnellen Implementation eines Tiefpass Filters mit 6dB Flankensteilheit:

int16_t IIRfilter_6dB(int16_t input, int8_t s) {
    static int16_t alt = 0; // wir speichern hier schon den geschifteten Wert alt << s
    alt = alt + (((input << s) - alt) >> s);
    return alt >> s; // den Shitfaktor im gespeicherten Wert müssen wir bei der Rückgabe wieder abziehen
}

Ein Filter höherer Ordnung erhält man durch Hintereinanderschaltung zweier Filter:

int16_t IIRfilter_12dB(int16_t input, int8_t s) {
    static int16_t alt1 = 0; 
    static int16_t alt2 = 0; 
    alt1 = alt1 + (((input << s) - alt1) >> s);
    alt2 = alt2 + ((alt1 - alt2) >> s);
    return alt2 >> s;
}

Um in anderen Fällen (beispielsweise 16Bit ADC) Überläufe zu vermeiden, muss man mit geringerer Genauigkeit rechnen,.

int16_t IIRfilter(int16_t input, int8_t s) {
    static int16_t alt = 0; // hier wird ein 15-Bit Wert gespeichert
    alt = alt + (((input >> 1) - alt) >> s);
    return alt << 1; // intern wird ja mit 15-Bit gerechnet
}

Die Alternative wäre, in der Routine mit höherer Genauigkeit zu rechnen:

int16_t IIRfilter(int16_t input, int8_t s) {
    static int24_t alt = 0; // wir speichern hier schon den geschifteten Wert alt << s
    alt = alt + ((((int24_t)input << s) - alt) >> s);
    return int16_t(alt >> s);
}

Hier der Frequenzgang zweier IIR-Filter mit 6db und 12db Steilheit und jeweils den Grenzfrequenzen von ca 1/10 und 1/20 der Samplingfrequenz.


Hier der Vergleich zweier IIR-Filter mit 6dB und Grenzfrequenz ca. 1/10 der Sampling Frequenz und eine 12dB Filter mit der halben Grenzfrequenz.

Zu dem letzten Fall hier noch die Sprungantwort der beiden Filter

Das es sich bei der Simulation tatsächlich um ein digitale Filter handelt sieht man, wenn man den Frequenzgang etwas weiter ausdehnt. Die maximale Filterdämpfung ist bei der halben Samplingfrequenz, bei Vielfachen dieser Frequenz ist keine Filterwirkung zu erkennen.

Fazit

Diese Art der Filterung ist nur für bestimmte Anwendungsfälle geeignet.

  • Der genaue Wert der Grenzfrequenz ist nicht entscheidend.
  • Die Berechnungsgeschwindigkeit ist entscheidend.

Dies ist beispielsweise bei 8-Bit-Prozessoren wie dem Arduino Uno und der Verarbeitung von Sensorwerten der Fall. Wenn man dort mit einer Sampling-Frequenz von 10 kHz oder bis zu 50 kHz arbeitet, bleiben nur wenige µs zur Verarbeitung der Werte übrig.

Wenn genauere Werte für die Grenzfrequenz eingehalten werden müssen, kommt man um eine Rechnung mit höherer Wortbreite nicht herum. Wir müssen dann zur ursprünglichen Formel für das Filter zurückkehren und unseren Shiftfaktor so wählen, dass k2≪s=Kk_2\ll s=K ganzzahlig wird. K liegt dann im Wertebereich 0<K<2s0<K<2^s.

(Falt≪s)=(Falt≪s)+(((Fneu≪s)−(Falt≪s))⋅K≫s)\left(F_{alt}\ll s\right)=\left(F_{alt}\ll s\right)+\left(\left(\left(F_{neu}\ll s\right)-\left(F_{alt}\ll s\right)\right)\cdot K\gg s\right)
Diese Formel liegt dem Algorithmus in dem im Link genannten Artikel zugrunde. Wie viele Bit benötigen wir jetzt für die Berechnung? Bisher waren es die Breite der Eingabewerte plus der Shiftfaktor s plus 1, jetzt kommen durch die Multiplikation mit K noch einmal s Bits dazu. Dann sind wir schnell bei einer Wortbreite von 32 Bit für die Berechnung. Für einen 32-Bit-Prozessor ist das alles kein Problem.

  • Carsten Meyer: Digital filtern und oszillieren, Make-Magazin 6/2022: Heise Verlag Hannover 2022

#Arduino, #Analogeingang, #Filter