Introfiltres PDF
Introfiltres PDF
Introfiltres PDF
1. Introduction
Le ltrage du signal numérique intervient dans de nombreux domaines : transmission
du signal, traitement du son, des images, traitement des données scientiques, etc. Ce
document explique le principe de fonctionnement d'un ltre numérique et donne quelques
exemples de ltres.
On se limite au cas des signaux unidimensionnels. On considère principalement des
ltres à réponse impulsionnelle nie (ltres RIF). Une étude plus détaillée des ltres RIF
est présentée dans Filtres à réponse impulsionnelle nie. Voir aussi Exemples de ltres
RIF.
On verra aussi un exemple de ltre à réponse impulsionnelle innie, le ltre intégra-
teur.
2. Signal échantillonné
Lorsqu'un signal analogique x(t) est échantillonné avec une période d'échantillonnage
Te , on obtient la suite suivante :
xn = x(nTe ) (1)
où n est un entier. La suite xn constitue un signal numérique.
Lorsque le spectre de x(t) a une fréquence maximale fmax , la condition de Nyquist-
Shannon permet de s'assurer que toute l'information du signal analogique est présente
dans le signal échantillonné : la fréquence d'échantillonnage fe = 1/Te doit être supérieure
à 2fmax .
La valeur de la sortie à l'instant n est donc une combinaison linéaire des valeurs de
l'entrée à l'instant n et aux N − 1 instants antérieurs. Si les valeurs de l'entrée ne sont
disponibles qu'à partir de l'indice 0, la première valeur de la sortie est yN −1 , puisqu'il
faut au moins N éléments en entrée pour faire le calcul. Pour voir la signication des
coecients hk , considérons en entrée une impulsion, dénie par :
x0 = 1 (3)
xn = 0 (n 6= 0) (4)
Frédéric Legrand Licence Creative Commons 2
yn = 0 (n < 0) (5)
y 0 = b0 x 0 = b0 (6)
y 1 = b1 x 0 = b1 (7)
··· (8)
yN −1 = bN −1 x0 = bN −1 (9)
yn = 0 (n ≥ N ) (10)
Les coecients b0 , b1 ...bN −1 constituent donc la réponse impulsionnelle du ltre.
L'opération dénie par la relation (2) est un produit de convolution. La sortie est
donc le produit de convolution de l'entrée par la réponse impulsionnelle. Ce type de
ltrage est aussi appelé ltrage par convolution.
Le produit de convolution peut être calculé avec la fonction scipy.signal.convolve,
en suivant la syntaxe suivante :
y = scipy.signal.convolve(x,b,mode='valid')
N −1
(12)
X
yn = bk xn−k
k=0
N −1
(13)
X
= exp(i2π n f Te ) bk exp(−i2π f k Te )
k=0
=xn H(Z) (14)
avec :
N −1
(15)
X
H(Z) = bk Z −k
k=0
Z = exp(i2π f Te ) (16)
Frédéric Legrand Licence Creative Commons 3
On voit donc que la sortie est dans ce cas proportionnelle à l'entrée. H(z) est la fonction
de transfert en Z , analogue à la fonction de transfert des signaux continus. La réponse
fréquentielle est nalement :
N −1
(17)
X
Hf (f ) = bk exp(−i2π k f Te )
k=0
On remarque que cette réponse fréquentielle est en fait une fonction de f Te , qui est
le rapport de la fréquence du signal sinusoïdal sur la fréquence d'échantillonnage fe =
1/Te . D'après le théorème de Shannon, ce rapport doit être inférieur à 1/2. La réponse
fréquentielle est de période fe mais en pratique on la trace sur l'intervalle [0, fe /2].
La fonction suivante renvoit le module et l'argument de la réponse fréquentielle. La
réponse impulsionnelle et le nombre de points de la courbe sont fournis en argument.
import math
import cmath
import numpy
def reponseFreq(b,nf):
f = numpy.arange(0.0,0.5,0.5/nf)
g = numpy.zeros(f.size)
phi = numpy.zeros(f.size)
for m in range(f.size):
H = 0.0
for k in range(b.size):
H += b[k]*cmath.exp(-1j*2*math.pi*k*f[m])
g[m] = abs(H)
phi[m] = cmath.phase(H)
return (f,g,numpy.unwrap(phi))
4. Exemples
4.a. Synthèse d'un signal bruité
Pour tester les diérents ltres, nous allons utiliser un signal périodique auquel on
ajoute un bruit gaussien :
On commence par dénir le signal par ses harmoniques :
def signal(t):
return 1.0*math.cos(2*math.pi*t) \
+0.5*math.cos(3*2*math.pi*t+math.pi/3) \
+0.2*math.cos(5*2*math.pi*t+math.pi/5) \
+0.2*math.cos(10*2*math.pi*t)
Frédéric Legrand Licence Creative Commons 4
import numpy
import random
from matplotlib.pyplot import *
fe=500
T=2.0
N=int(fe*T)
te = T/N
x = numpy.zeros(N)
t = numpy.zeros(N)
sigma = 0.05
for k in range(N):
t[k] = k*te
x[k] = signal(t[k])+random.gauss(0,sigma)
figure(figsize=(10,4))
plot(t,x)
On repère sur ce spectre les raies correspondant aux harmoniques dénies dans le signal.
Le bruit se manifeste à toute fréquence (c'est un bruit blanc). La fréquence d'échan-
tillonnage a été choisie très largement supérieure aux fréquences utiles du signal. On
parle dans ce cas de sur-échantillonnage. Cette condition est nécessaire pour eectuer un
ltrage ecace du bruit.
b = numpy.array([0.5,0.5])
(f,g,phi)=reponseFreq(b,100)
figure(figsize=(6,4))
plot(f,g)
xlabel('f/fe')
ylabel('|H|')
axis([0,0.5,0,1])
grid()
Frédéric Legrand Licence Creative Commons 6
figure(figsize=(6,4))
plot(f,phi)
xlabel('f/fe')
ylabel('Phase')
axis([0,0.5,-3,3])
grid()
Le ltre moyenneur est un ltre passe-bas (peu sélectif). Le déphasage varie linéairement
avec la fréquence.
Frédéric Legrand Licence Creative Commons 7
Pour obtenir le signal discret ltré, il sut d'eectuer la convolution avec la réponse
impulsionnelle. Le premier point calculé correspond à l'instant Te , c'est pourquoi l'échelle
de temps doit être modiée.
import scipy.signal
b = np.array([0.5,0.5])
y = scipy.signal.convolve(x,b,mode='valid')
ny = y.size
ty = te+np.arange(ny)
figure(figsize=(10,4))
plot(ty,y)
xlabel('t')
ylabel('y')
grid()
On voit que le bruit est légèrement réduit par ce ltre passe-bas. Pour augmenter l'ef-
cacité du ltrage, il faut augmenter l'ordre N du ltre (la longueur de la réponse
impulsionnelle).
Une première idée consiste à calculer une moyenne sur un plus grand nombre de
points. Voici la réponse fréquentielle pour une moyenne sur 10 points :
b = numpy.ones(10)/10.0
(f,g,phi)=reponseFreq(b,100)
figure(figsize=(6,4))
plot(f,g)
xlabel('f/fe')
ylabel('|H|')
axis([0,0.5,0,1])
grid()
Frédéric Legrand Licence Creative Commons 8
Il y a des rebonds dans la bande atténuée qui en font un très mauvais ltre. Nous allons
voir comment modier les coecients de la réponse impulsionnelle pour obtenir une
meilleure réponse.
bk = f (k − P ) (20)
Si P est donné, on choisit σ pour avoir une valeur très faible sur les bords :
P
σ=√ (21)
−2 ln
La réponse impulsionnelle est normalisée pour que la somme de ses coecients soit égale
à 1. Voici un exemple :
P_gauss=10
b_gauss = numpy.zeros(2*P_gauss+1)
epsilon=0.01
sigma=P_gauss/math.sqrt(-2.0*math.log(epsilon))
som = 0.0
for k in range(2*P_gauss+1):
Frédéric Legrand Licence Creative Commons 9
b_gauss[k] = math.exp(-(k-P_gauss)**2/(2*sigma**2))
som += b_gauss[k]
b_gauss = b_gauss/som
(f,g,phi)=reponseFreq(b_gauss,500)
figure(figsize=(6,4))
plot(f,g)
xlabel('f/fe')
ylabel('|H|')
axis([0,0.5,0,1])
grid()
La réponse fréquentielle d'un ltre gaussien est aussi une gausienne. Voyons la phase :
figure(figsize=(6,4))
plot(f,phi)
xlabel('f/fe')
ylabel('Phase')
grid()
Frédéric Legrand Licence Creative Commons 10
Dans la bande passante (et même au delà), la phase varie linéairement avec la fréquence.
Il y a donc un décalage temporel constant entre l'entrée et la sortie, égal à P Te .
Voyons l'eet sur le signal bruité déni plus haut. L'échelle de temps est calculée de
manière à annuler le décalage entre l'entrée et la sortie (cela n'est pas faisable dans un
ltre fonctionnant en temps réel).
y = scipy.signal.convolve(x,b_gauss,mode='valid')
ny = y.size
ty = np.arange(ny)+P_gauss*te
figure(figsize=(10,4))
plot(t,x,'b')
plot(ty,y,'r')
xlabel('t')
ylabel('y')
grid()
Frédéric Legrand Licence Creative Commons 11
Le bruit est éliminé sans que la forme du signal ne soit modiée. On remarque que P
points au début et à la n du signal ne sont pas traités. Dans un ltre temps-réel, cela
ne pose pas de problème car les 2P points non ltrés sont au début du signal.
P=10
b = numpy.zeros(2*P+1)
def sinc(u):
if u==0:
return 1.0
else:
return math.sin(u)/u
a=0.15
for k in range(2*P+1):
b[k] = 2*a*sinc(2*math.pi*(k-P)*a)
(f,g,phi)=reponseFreq(b,500)
figure(figsize=(6,4))
plot(f,g)
xlabel('f/fe')
ylabel('|H|')
axis([0,0.5,0,1.2])
grid()
Frédéric Legrand Licence Creative Commons 12
Ce ltre est beaucoup plus sélectif que le ltre gaussien, mais a l'inconvénient de com-
porter des ondulations dans la bande passante et dans la bande atténuée.
figure(figsize=(6,4))
plot(f,phi)
xlabel('f/fe')
ylabel('Phase')
grid()
Frédéric Legrand Licence Creative Commons 13
La phase est bien linéaire dans la bande passante. Voyons l'eet sur le signal :
y = scipy.signal.convolve(x,b,mode='valid')
ny = y.size
ty = np.arange(ny)+P*te
figure(figsize=(10,4))
plot(t,x,'b')
plot(ty,y,'r')
xlabel('t')
ylabel('y')
grid()
Ce ltre est moins ecace que le ltre gaussien pour réduire le bruit, en raison de la
présence de rebonds dans la bande atténuée.
b = numpy.array([1.0,-1.0])
(f,g,phi)=reponseFreq(b,500)
figure(figsize=(6,4))
plot(f,g)
xlabel('f/fe')
ylabel('|H|')
grid()
Frédéric Legrand Licence Creative Commons 14
figure(figsize=(6,4))
plot(f,phi)
xlabel('f/fe')
ylabel('Phase')
grid()
Le ltre dérivateur est passe-haut. Voyons son eet sur le signal bruité :
y_deriv = scipy.signal.convolve(x,b,mode='valid')
Frédéric Legrand Licence Creative Commons 15
ny = y_deriv.size
ty = te+np.arange(ny)
figure(figsize=(10,4))
plot(ty,y_deriv,'r')
xlabel('t')
ylabel('y')
grid()
Le bruit est tellement amplié par le ltre dérivateur qu'il couvre complètement le signal
utile. Cela est dû au caractère passe-haut du dérivateur. Pour dériver un signal bruité, il
faut impérativement eectuer un ltrage passe-bas au préalable, par exemple un ltrage
gaussien :
y = scipy.signal.convolve(x,b_gauss,mode='valid')
y_deriv = scipy.signal.convolve(y,b,mode='valid')
ny = y_deriv.size
ty = np.arange(ny)+(P_gauss+1)*te
figure(figsize=(10,4))
plot(ty,y_deriv,'r')
xlabel('t')
ylabel('y')
grid()
Frédéric Legrand Licence Creative Commons 16
5. Filtre intégrateur
Un ltre intégrateur peut être obtenu par la relation suivante :
xn + xn−1
yn = yn−1 + Te (25)
2
Voir Filtres intégrateur et dérivateur pour la justication de cette relation.
Cette relation est diérente de la relation (2) d'un ltre à réponse impulsionnelle
nie : la sortie à l'instant n dépend de l'état de la sortie à l'instant n − 1. Ce type de
ltre est appelé ltre récursif. Pour l'appliquer à un signal, il faut xer une condition
initiale, ici la valeur y0 .
Appliquons ce ltre à l'impulsion unité, avec la condition initiale y0 = 0. On obtient :
x0 1
y1 = y0 + = (26)
Te Te
1
y2 = y1 = (27)
Te
1
y3 = y2 = (28)
Te
··· (29)
La réponse impulsionnelle est un échelon. C'est un exemple de réponse impulsionnelle
innie. Un ltre est stable si sa réponse impulsionnelle tend vers zéro. Le ltre intégrateur
est donc instable, mais sa réponse impulsionnelle est tout de même bornée (stabilité
marginale).
La réponse fréquentielle de ce ltre peut être obtenue en suivant la même méthode
qu'en 3.b. On obtient ainsi la fonction de transfert en Z :
Te 1 + Z −1
H(Z) = (30)
2 1 − Z −1
et la réponse fréquentielle se déduit en posant :
Z = exp(i2π f Te ) (31)
D'une manière générale, la fonction de transfert en Z d'un ltre linéaire est une frac-
tion rationnelle en Z −1 . La fonction scipy.signal.freqz permet d'obtenir la réponse
fréquentielle d'un ltre déni par sa fonction de transfert en Z :
b=[0.5,0.5]
a=[1,-1]
[w,h] = scipy.signal.freqz(b,a)
figure(figsize=(6,4))
plot(w/(2*math.pi),numpy.abs(h))
xlabel('f/fe')
ylabel('|H|')
grid()
Frédéric Legrand Licence Creative Commons 17
y = scipy.signal.lfilter(b,a,x)
figure(figsize=(10,4))
plot(t,y)
xlabel('t')
ylabel('y')
grid()
Références
[1] M. Bellanger, Traitement numérique du signal, (Dunod, 1998)