Le filtrage des signaux numériques¶

Nous allons voir quelques moyens de filtrer les signaux numériques :

Lissage d'un signal par une moyenne mobile exponentielle
Filtre du premier ordre (passe-bas, passe-haut)
Filtre numérique (passe-bas, passe-bande, coupe-bande, passe-haut)

Info

Ceci n'est pas un cours sur le filtrage, juste quelques méthodes pour filtrer un signal après acquisition.

Génération du signal pour les tests¶

import numpy as np
from numpy.random import randn
import matplotlib.pyplot as plt

sig = np.cumsum(randn(800))  # Brownian noise
plt.plot(sig, color='silver', label='Signal')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Signal d'origine")
plt.show()

png

Lissage d'un signal par une moyenne mobile exponentielle¶

La moyenne mobile exponentielle permet de lisser une courbe pour mettre en évidence son allure en supprimer les écarts ponctuels et le bruit.

L'équation permettant de calculer la moyenne mobile $s_n$ quelque soit $n$ est :

$s_n=\alpha.e_{n} +(1-\alpha)s_{n-1}$

Avec :

$s_0$ une moyenne des premières valeurs ;
$\alpha = \frac{2}{N+1}$ la constante de lissage où $N$ est un entier à choisir en fonction de la précision voulue (tester plusieurs valeurs pour obtenir un résultat cohérent).

import numpy as np
import matplotlib.pyplot as plt

# Constante du système
N = 10
alpha = 2/(N+1)

# Préparation de la liste de sortie
s_m = []
s_m.append(0)

# Application du filtre
for e in sig:
    s_m.append(alpha*e+(1-alpha)*s_m[-1])

# Affichage du signal filtré
plt.plot(sig, color='silver', label='Signal')
plt.plot(s_m, color='#cc0000', label='Signal filtré')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre par moyenne glissante exponetielle")
plt.show()

png

Filtre passe-bas¶

Filtre passe-bas du 1er ordre¶

La fonction de transfert d'un filtre passe-bas du premier ordre est : $H(p)=\frac{1}{1+\tau.p}\quad$
Avec $\frac{1}{\tau}$ la pulsation de coupure du filtre.

L'équation différentielle associée à ce filtre est : $\tau\frac{ds(t)}{dt}+s(t)=e(t)$

L'équation discrétisée devient alors : $\tau\frac{s_{n+1}-s_n}{T_e}+s_n=e_n$

L'équation permettant de calculer $s_n$ quelque soit $n$ devient :

$s_n=s_{n-1}+\frac{T_e}{\tau}(e_{n-1}-s_{n-1}) \quad \text{Avec}\, s_0=0$

import numpy as np
import matplotlib.pyplot as plt

# Fréquence de coupure
fc = 0.1  # Hz
tau = 1/(2*np.pi*fc)

# Période d'échantillonnage
Te = 1  # s

# Préparation de la liste de sortie
s_pb = []
s_pb.append(0)

# Application du filtre
for i in range(1, len(sig)):
    s_pb.append(s_pb[i-1]+Te/tau*(sig[i-1]-s_pb[i-1]))

# Affichage du signal filtré
plt.plot(sig, color='silver', label='Signal')
plt.plot(s_pb, color='#cc0000', label='Signal filtré')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre passe-bas du 1er ordre")
plt.show()

png

Passe bas avec un filtre de Butterworh¶

On va ici utiliser une méthode numérique un filtre de Butterworth.

Le filtre numérique utilise la fréquence de Nyquist. Cette fréquence correspond à la fréquence d'échantillonnage divisée par deux.

Attention

Ce filtre ne peut fonctionner que si l'on a un jeu de données complet, car il a un fonctionnement acausal. Il ne peut pas être utilisé pour du temps réel.
L'avantage c'est qu'il ne créait pas de déphasage, si on utilise la fonction filfilt pour appliquer le filtre.
Pour appliquer le filtre en mode causal, vous pouvez vous orienter vers la fonction lfilter.

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

# Fréquence d'échantillonnage
fe = 1  # Hz

# Fréquence de nyquist
f_nyq = fe / 2.  # Hz

# Fréquence de coupure
fc = 0.1  # Hz

# Préparation du filtre de Butterworth en passe-bas
b, a = signal.butter(4, fc/f_nyq, 'low', analog=False)

# Application du filtre
s_but = signal.filtfilt(b, a, sig)


# Affichage du signal filtré
plt.plot(sig, color='silver', label='Signal')
plt.plot(s_but, color='#cc0000', label='Signal filtré')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre passe-bas numérique")
plt.show()

png

Comparaison filtre du 1er ordre et filtre de Butterworth¶

import matplotlib.pyplot as plt

# Affichage du signal filtré
plt.plot(sig, color='silver', label='Signal')

plt.plot(s_pb, label='1er ordre')
plt.plot(s_but, label='Butterworth')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.xlim(200, 240)
plt.ylim((-24, -19))
plt.title("Comparaison")
plt.show()

png

On peut observer au moins deux choses :

Le filtre de Butterworth filtre beaucoup plus le signal que le 1er ordre ;
Le filtre basé sur un premier ordre présente un retard, contrairement au filtre de Butterworth.

Filtre passe-haut¶

On construit un signal qui est la somme de deux signaux sinusoïdaux de fréquence différente.

L'objectif est de retrouver le signal de haute fréquence.

import numpy as np
import matplotlib.pyplot as plt

# Génération d'un signal à filtrer
x = np.linspace(0, 4*np.pi, 1000)
e = np.sin(x)+0.2*np.sin(10*x)

plt.plot(x, e, color='silver', label='Signal')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre passe-bas du 1er ordre")
plt.show()

png

Filtre passe-haut du 1er ordre¶

La fonction de transfert d'un filtre passe-haut du premier ordre est : $H(p)=\frac{\tau.p}{1+\tau.p}\quad$
Avec $\frac{1}{\tau}$ la pulsation de coupure du filtre.

L'équation différentielle associée à ce filtre est : $\tau\frac{ds(t)}{dt}+s(t)=\tau\frac{de(t)}{dt}$

L'équation discrétisée devient alors : $\tau\frac{s_{n+1}-s_n}{T_e}+s_n=\tau\frac{e_{n+1}-e_n}{T_e}$

L'équation permettant de calculer $s_n$ quelque soit $n$ devient :

$s_n=s_{n-1}\left(1-\frac{T_e}{\tau}\right)+e_{n}-e_{n-1} \quad \text{Avec}\, s_0=0$

import numpy as np
import matplotlib.pyplot as plt

# Fréquence de coupure
fc = 1  # Hz
tau = 1/(2*np.pi*fc)

# Période d'échantillonnage
Te = 4*np.pi/1000  # s

# Préparation de la liste de sortie
s_ph = []
s_ph.append(0)

# Application du filtre
for i in range(1, len(e)):
    s_ph.append(s_ph[i-1]*(1-Te/tau)+e[i]-e[i-1])

# Affichage du signal filtré
plt.plot(x, e, color='silver', label='Signal')
plt.plot(x, s_ph, color='#cc0000', label='Signal filtré')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre passe-haut du 1er ordre")
plt.show()

png

Passe haut avec un filtre de Butterworth¶

On applique la même méthode que pour le passe bas avec le filtre de Butterworth.

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

# Fréquence d'échantillonnage
fe = 1/Te  # Hz

# Fréquence de nyquist
f_nyq = fe / 2.  # Hz

# Fréquence de coupure
fc = 1  # Hz

# Préparation du filtre de Butterworth en passe-haut
b, a = signal.butter(4, fc/f_nyq, 'high', analog=False)

# Application du filtre
s_but2 = signal.filtfilt(b, a, e)


# Affichage du signal filtré
plt.plot(x, e, color='silver', label='Signal')
plt.plot(x, s_but2, color='#cc0000', label='Signal filtré')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre passe-haut numérique")
plt.show()

png

Comparaison des deux méthodes¶

import matplotlib.pyplot as plt
import numpy as np

# Signal à retrouver
s_ori = 0.2*np.sin(10*x)
plt.plot(x, s_ori, label='Signal à retrouver')
# Affichage du signal filtré
plt.plot(x, s_ph, label='1er ordre')
plt.plot(x, s_but2,  label='Filtre de Butterworth', color='red')

# Affichage d'une barre verticale pour le sommet
s_max = np.where(e == max(e))
plt.axvline(x=x[s_max], color='red')

plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre passe-haut numérique")
plt.xlim(1, 2)
plt.show()

png

Le constat est facile, le filtre de Butterworth est le plus efficace. Le signal filtré se superpose quasiment sur le signal d'origine.

Configuration du filtre de Butterworth pour d'autres type de filtre¶

la fonction signal.butter(N, Wn, btype='low', analog=False, fs=None) peut prendre les paramètres suivants :

Paramètre	Détails
`N`	Nombre entier donnant l'ordre du filtre.
`Wn`	La ou les fréquences de coupures du filtre. Pour un filtre numérique la valeur est comprise entre 0 et 1 avec 1 la fréquence de Nyquist. Pour un filtre analogique, `Wn` correspond à la pulsation de coupure en rad/s.
`btype`	Type de filtre : ‘lowpass’, ‘highpass’, ‘bandpass’, ‘bandstop’
`analog`	Modèle analogique si vrai, sinon numérique.
`fs`	Fréquence d'échantillonnage pour un système numérique.

Passe bande¶

import numpy as np
from scipy import signal

# Période d'échantillonnage
Te = 4*np.pi/1000  # s

# Génération du signal à filtrer
x = np.linspace(0, 4*np.pi, 1000)
e = np.sin(x)+0.4*np.sin(10*x)+0.1*np.sin(100*x)

# Fréquence d'échantillonnage
fe = 1/Te  # Hz

# Fréquence de nyquist
f_nyq = fe / 2.  # Hz

# Fréquence de coupure
fc1 = 1  # Hz
fc2 = 10  # Hz

# Préparation du filtre de Butterworth en passe bande
b, a = signal.butter(4, (fc1/f_nyq, fc2/f_nyq), 'bandpass', analog=False)
s_but_bp = signal.filtfilt(b, a, e)


# Préparation du filtre de Butterworth en passe-bas
b, a = signal.butter(4, fc1/f_nyq, 'low', analog=False)
s_but_pb = signal.filtfilt(b, a, e)

# Préparation du filtre de Butterworth en passe-haut
b, a = signal.butter(4, fc2/f_nyq, 'high', analog=False)
s_but_ph = signal.filtfilt(b, a, e)

import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1.inset_locator import zoomed_inset_axes, mark_inset

fig, ax = plt.subplots()

# Affichage du signal filtré
plt.plot(x, e, color='silver', label='Signal')
plt.plot(x, s_but_bp, label='Passe-bandes')
plt.plot(x, s_but_pb, label='Passe bas')
plt.plot(x, s_but_ph, label='Passe haut')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre numérique")

# Zoom sur le passe-haut
# zoom-factor: 2.5, location: upper-left
axins = zoomed_inset_axes(ax, 4, loc=8)
axins.plot(x, s_but_ph, label='Passe haut', color='green')
axins.set_xlim(2.5, 4)  # apply the x-limits
axins.set_ylim(0.1, -0.1)  # apply the y-limits
plt.yticks(visible=False)
plt.xticks(visible=False)
mark_inset(ax, axins, loc1=2, loc2=1, fc="none", ec="0.5")

plt.show()

png

Coupe-bandes¶

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

# Période d'échantillonnage
Te = 4*np.pi/1000  # s

# Génération du signal à filtrer
x = np.linspace(0, 4*np.pi, 1000)
e = np.sin(x)+0.4*np.sin(10*x)+0.1*np.sin(100*x)

# Fréquence d'échantillonnage
fe = 1/Te  # Hz

# Fréquence de nyquist
f_nyq = fe / 2.  # Hz

# Fréquence de coupure
fc1 = 1  # Hz
fc2 = 10  # Hz

# Préparation du filtre de Butterworth en coupe bande
b, a = signal.butter(4, (fc1/f_nyq, fc2/f_nyq), 'bandstop', analog=False)
s_but_bs = signal.filtfilt(b, a, e)

# Affichage du signal filtré
plt.plot(x, e, color='silver', label='Signal')
plt.plot(x, s_but_bs, label='Coupe bande')
plt.grid(True, which='both')
plt.legend(loc="best")
plt.title("Filtre numérique")
plt.show()

png