Cours1 PDF

Télécharger au format pdf ou txt
Télécharger au format pdf ou txt
Vous êtes sur la page 1sur 79

ISSN xxxx-yyyy

Bases de Traitement du Signal


Notes de cours

Nicolas Moreau

Nicolas Moreau 1er octobre 2007

Licence de droits d’usage


en fin de document
Chapitre 1

Introduction

1.1 Aspect “Signal”


Signal = représentation physique d’une information envoyée d’une source vers un destinataire.
En pratique résultat d’une mesure par un capteur.
Exemple : signaux de parole (figures 1.1) et de musique (figures 1.2), signaux audiofréquences,
variation d’une pression en fonction du temps.

1.1.1 Signal de parole

0.8

0.6

0.4

0.2

−0.2

−0.4

−0.6

−0.8

−1
500 1000 1500 2000 2500 3000 3500 4000
60

0.5

0.4 50

0.3

0.2 40
Puissances [dB]

0.1

0 30

−0.1
20
−0.2

−0.3
10
−0.4

−0.5
0
0 0.5 1 1.5 2 2.5 3 3.5 4
1300 1350 1400 1450 1500 1550 1600 1650 1700 1750 1800 Frequences [kHz]

Fig. 1.1 – Exemple d’un signal de parole : début d’une phrase (“Des gens ...”) prononcée par un
locuteur féminin. Effet de zoom et représentation fréquentielle du signal correspondant.

– Trois types de sons : voisés (signaux périodiques sur une durée inférieure à 100 ms, existence
d’une fréquence fondamentale, d’harmoniques), non-voisés (variations plus rapides, signal

Nicolas Moreau 1er octobre 2007


page 1/78
Licence de droits d’usage
CHAPITRE 1. INTRODUCTION

plus “hautes fréquences”), plosives (support temporel très étroit).


En 64 ms, environ 12 “périodes” ⇒ 5 ms par période soit un premier pic à approximative-
ment 200 Hz ⇒ voix féminine.
On privilégie dans ce cours l’étude de la représentation fréquentielle des signaux. L’outil de
base est la transformée de Fourier. Le module au carré donne la répartition de la puissance
en fonction de la fréquence (spectre).
– Signal non-stationnaire : les propriétés statistiques évoluent au cours du temps. Signal
localement stationnaire sur des durées inférieures à une centaine de ms. Nécessité d’utiliser
des fenêtres de pondération.
La transformée de Fourier n’est pas utilisable directement. Nécessité de définir de nouveaux
outils.
– Mémorisation puis tracé par une machine d’un signal “à temps continu” (x(t) avec t ∈ R)
⇒ signal “à temps discret” (x(nT ) avec n ∈ Z).
On privilégie très fortement dans ce cours l’étude des signaux à temps discret. Existence
d’une transformée de Fourier à temps discret.
– Choix de la fréquence d’échantillonnage. Recherche de conditions de telle sorte qu’il n’y ait
pas perte d’information par échantillonnage (que l’on puisse revenir au signal de départ).
En pratique c’est le résultat d’un compromis. Le signal de parole est échantillonné soit à
8 kHz (réseau téléphonique commuté, bande téléphonique) soit à 16 kHz (téléconférences,
communications de groupe, bande élargie).
Dans ce cours, importance des ordres de grandeur. Ce n’est pas un cours abstrait.

1.1.2 Signal de musique

0.8

0.6

0.4

0.2

−0.2

−0.4

−0.6

−0.8

−1
2000 4000 6000 8000 10000 12000 14000 16000
100

0.5
90
0.4
80
0.3
70
0.2
Puissances [dB]

60
0.1

0 50

−0.1 40

−0.2 30

−0.3
20
−0.4
10
−0.5
0
0 2 4 6 8 10 12 14 16
7000 7200 7400 7600 7800 8000 8200 8400 8600 8800 9000 Frequences [kHz]

Fig. 1.2 – Exemple d’un signal de musique : violon. Effet de zoom et représentation fréquentielle
du signal correspondant.

Nicolas Moreau 1er octobre 2007


page 2/78
Licence de droits d’usage
1.2. ASPECT “SYSTÈME”

– Bande Hi-Fi : l’oreille est sensible à des fréquences comprises entre 20 Hz et 20 kHz ⇒ fe
= 32 kHz (bande FM, un peu juste), 44.1 kHz (CD), 48 kHz (studio production).
– Timbre d’un instrument caractérisé par les pics dans le spectre.
– Importance de la représentation fréquentielle. Interprétation généralement plus facile dans
le domaine fréquentiel que dans le domaine temporel. Exemple d’un signal dégradé par une
opération de compression.

1.1.3 Autres signaux


– Images. Luminance : fonction de deux variables spatiales et du temps l(x, y, t). Signal
vidéo : définition d’un balayage de façon à se ramener à une fonction d’une seule variable,
le temps.
Ordre de grandeur de la “largeur de bande”. 625 lignes, 25 fois par seconde ⇒ 64 µs
par ligne. Écran très bonne définition : 1024 points (pixels) par ligne ⇒ fréquence max :
512×(blanc, noir) ⇒ période max = 64/512 = 1/8 µs ⇒ 8 MHz. Largeur de bande 1000
fois plus importante que la parole. Problème du traitement temps réel dans une machine. A
complexité de traitement équivalente, il faudra un processeur 1000 fois plus puissant pour
traiter des images que pour traiter du signal de parole. Dans la pratique, généralement
traitements plus élaborés en parole qu’en image ...
– Signaux biomédicaux : EEG, ECG ...
– Radar, sonar ...
– Signaux sismiques (recherche pétrolière) ...

1.2 Aspect “système”


– Système = organe physique qui transforme un signal d’entrée en un signal de sortie. Re-
présentation par une boîte noire. Exemple : système phonatoire (entrée = vibrations des
cordes vocales, système = conduit vocal, sortie = pression acoustique).
– Systèmes à temps continu et à temps discret.
– Traitement de base : opération de filtrage.
– Pour caractériser un système, autres outils que la transformée de Fourier : transformée de
Laplace pour les systèmes à temps continu et la transformée en z pour les systèmes à temps
discret.

1.3 Exemples de traitement


– Parole : comprimer (transmission : téléphone mobile, GSM), synthétiser (serveurs vocaux),
reconnaître (commande vocale), améliorer la prise de son (filtrage d’antenne) ...
– Musique : comprimer (diffusion, stockage, MPEG-Audio), modéliser, simuler synthétiser
(instruments, salle), corriger des enregistrements dégradés ...
– Images : comprimer (archiver, MPEG), reconstruire (tomodensitométrie, scanner), recon-
naître (analyse de scènes, robotique), corriger ...
– Radar, sonar : détecter (applications militaires).
– Signaux sismiques : prospecter.
– Signaux biomédicaux : analyser, archiver.
– Autres applications : correction du canal de transmission (trajets multiples, interférences),
surveillance (non-destructive) de machines ...

Nicolas Moreau 1er octobre 2007


page 3/78
Licence de droits d’usage
CHAPITRE 1. INTRODUCTION

1.4 Plan du cours, idées directrices


– Etude des signaux et des systèmes monodimensionnels.
– Etude très superficielle des signaux et des systèmes à temps continu. Cette étude réclame
un développement important pour être traitée correctement (distribution). Plus de pro-
blèmes mathématiques délicats si on se restreint à l’étude des signaux et des systèmes à
temps discret. C’est un cours de base en traitement du signal et non en théorie du signal.
Impasse presque totale sur les filtres à temps continu (analogiques) : vu dans d’autres
cours. Uniquement quelques qualificatifs et notion de réponse impulsionnelle et de réponse
en fréquence. Impasse totale sur l’étude de la transformée de Laplace.
Etude des signaux et des systèmes à temps continu faite uniquement pour rappeler et
interpréter physiquement quelques propriétés de la transformée de Fourier (signaux de
module et carré sommables) et du développement en série de Fourier (signaux périodiques).
– Etude des signaux et des systèmes à temps discret en privilégiant l’interprétation fréquen-
tielle (ce n’est pas un cours de théorie des systèmes) :
– Comment passe-t-on de la transformée de Fourier à temps continu (TFTC) à la trans-
formée de Fourier à temps discret (TFTD), à la transformée de Fourier discrète (TFD,
FFT) et à la transformée de Fourier à court terme (introduction) ? Cf table 1.1.
Importance pratique de ces transformées.
– Comment filtre-t-on des signaux à temps discret ?
– Introduction aux processus aléatoires.
Densité spectrale de puissance. Filtrage.
Modélisation paramétrique.

TFTC d’une fonction x(t) ∈ L1 ∩ L2

Z +∞ Z +∞
x(t) = X(f )ej2πf t
df avec X(f ) = x(t)e−j2πf t dt
−∞ −∞

DSF d’une fonction x(t) périodique (T0 )

+∞ Z +T0 /2
X 1
x(t) = X(k)e j2πkt/T0
avec X(k) = x(t)e−j2πkt/T0 dt
T0 −T0 /2
k=−∞

TFTD d’une suite x(n) (série convergente)

Z +1/2 +∞
X
x(n) = X(f )ej2πf n df avec X(f ) = x(n)e−j2πf n
−1/2 n=−∞

TFD d’une suite x(n) périodique (N0 )

N −1 0 −1
NX
1 X
x(n) = X(k)ej2πkn/N0 avec X(k) = x(n)e−j2πkn/N0
N0
k=0 n=0

Tab. 1.1 – Résumé de quelques résultats (futurs) du cours

Nicolas Moreau 1er octobre 2007


page 4/78
Licence de droits d’usage
Chapitre 2

Signaux à temps continu : Analyse et


filtrage

2.1 Introduction
– Signaux déterministes à temps continu : modélisables par une fonction x(t) avec t ∈ R à
valeurs réelles ou complexes. Ce sont soit des signaux “test” (fonction rectangle, sinusoïde),
soit des “réponses impulsionnelles” de filtre.
– Classification de type continuité, dérivabilité peu utile. Distinction entre signaux d’énergie
finie et de puissance finie.
– Définitions :
Puissance instantanée
R +∞ = |x(t)|2
Energie = −∞ |x(t)|2 dt
R +T /2
Puissance moyenne = limT →∞ T1 −T /2 |x(t)|2 dt

2.2 Signaux d’énergie finie


2.2.1 Transformée de Fourier à temps continu (TFTC)
Existence de la transformée de Fourier. Notations
Z +∞ Z +∞
x(t) = X(f )ej2πf t
df
X(f ) = x(t)e−j2πf t dt
−∞ −∞

Interprétation : X(f ) est la projection du signal sur l’exponentielle complexe à la fréquence


f . Le signal est décomposé sur une base continue d’exponentielles complexes. Le paramètre f
s’interprète comme une “fréquence” au sens habituel (l’inverse d’une période) à un détail près :
existence de fréquences négatives. C’est dû à la projection sur des exponentielles complexes plutôt
que sur des sinusoïdes.
Problèmes mathématiques délicats : on projette sur une fonction qui n’est pas elle même
d’énergie finie.
Exemple : fenêtre rectangulaire x(t) = rectT (t).

sin(πf T )
X(f ) = = T sinc(f T ).
πf

Observation : x(t) de “durée” finie (support de x(t)), X(f ) de “bande” infinie (support de X(f )).
Existence de relations d’incertitude.

Nicolas Moreau 1er octobre 2007


page 5/78
Licence de droits d’usage
CHAPITRE 2. SIGNAUX À TEMPS CONTINU : ANALYSE ET FILTRAGE

2.2.2 Principales propriétés et commentaires

x(t)
X(f )
Linéarité
Homothétie temporelle : x(at)
X(f /a)/a
Exemple : magnétophone ralenti ⇒ perception plus graves des fréquences.
Translation temporelle : x(t − t0 )
X(f ) e−j2πf t0
Amplitude conservée, phase modifiée. La phase est directement reliée à l’origine des temps.
Modulation : x(t) ej2πf0 t
X(f − f0 )
Opération fondamentale en transmission. Création d’un “multiplex fréquentiel” :
K
X K
X
y(t) = xk (t) ej2πfk t
Y (f ) = Xk (f − fk )
k=1 k=1

Transmission de plusieurs signaux distincts sur un même support (fil de cuivre, liaison
hertzienne) à la condition que ces signaux soient à “bande limitée”.
R +∞
Pondération par une fenêtre : x(t) × y(t)
−∞ X(λ)Y (f − λ) dλ
Résultat très utile pour formaliser l’extraction d’une portion de signal. Cf exemple d’un
signal de parole en introduction. Utilisation de la fonction rectangle.

sin(πf T )
y(t) = rectT (t)
Y (f ) = T = T sinc(f T )
πf T
R +∞
Interprétation graphique de −∞ X(λ)Y (f − λ) dλ pas tout à fait évidente. Etude détaillée
plus tard.
R +∞
Convolution : x(t) ∗ y(t) = −∞ x(τ )y(t − τ ) dτ
X(f ) × Y (f )
Opération fondamentale dans l’étude des systèmes : opération de filtrage.
Exemple : Y (f ) = rectB (f ) ⇒ filtrage passe-bas “idéal”.
R +∞ R +∞
Parseval : −∞ |x(t)|2 dt = −∞ |X(f )|2 df
Conservation de l’énergie, isométrie.
Signal réel : X(−f ) = X ∗ (f ) ⇒ |X(f )|2 = X(f )X(−f )
Symétrie hermitienne.

2.2.3 Autres définitions


R +∞
– Energie = −∞ |X(f )|2 df
SX (f ) = |X(f )|2 s’appelle la “densité spectrale énergétique” ou “spectre”. C’est une fonc-
tion à valeur réelle, positive, paire. Interprétation physique fondamentale. On a perdu
toute information relative à la phase. Généralement la phase joue un rôle beaucoup moins
important surtout lorsque le récepteur est l’oreille (peu sensible à la phase).
– Le support de x(t) s’appelle la “durée” du signal.
– Le support de X(f ) s’appelle la “largeur de bande” (ou bande) du signal.
– Un signal vérifiant X(f ) = 0 pour f 6∈ [−B1 , B2 ] est “à bande limitée”. Si x(t) ∈ R alors
B1 = B2 .
– Signal “bande étroite” : X(f ) 6= 0 pour f0 − B < |f | < f0 + B avec B << f0 . Importance
des ordres de grandeur, exemple : pour une chaîne de radio FM, f0 ≈ 100 MHz et B ≈ 20
kHz.

Nicolas Moreau 1er octobre 2007


page 6/78
Licence de droits d’usage
2.3. SIGNAUX PÉRIODIQUES

2.2.4 Théorème de Bernstein


Existence de “relations d’incertitude” entre les supports temporels et fréquentiels. Exemple
du signal rectangulaire, figure ??. A un support fini dans le domaine temporel correspond un
support infini dans le domaine fréquentiel. Conclusion : un signal à bande limitée n’existe pas
(en théorie). En pratique, importance des ordres de grandeur.
Le théorème de Bernstein donne des relations entre la largeur de bande du signal et les
dérivées n-èmes du signal (plus la bande du signal est importante, plus il est susceptible de varier
rapidement).
Autre interprétation : plus un signal est impulsif (concentré dans le temps) plus son spectre
est large (et réciproquement).

2.3 Signaux périodiques


Nécessité de définir une “fonction” qui chiffre la distribution de la puissance en fonction de la
fréquence.

2.3.1 Développement en série de Fourier (DSF)


Energie infinie ⇒ pas de transformée de Fourier au sens des fonctions. Existence d’un déve-
loppement en série de Fourier.
Si x(t) est une fonction périodique de période fondamentale T , i.e. le plus petit réel positif
vérifiant x(t + T ) = x(t) ∀t, alors

+∞ Z T /2
X 1
x(t) = X(k) e j2πkt/T
avec X(k) = x(t) e−j2πkt/T dt
T −T /2
k=−∞

Propriétés similaires à celles de la TFTC.

2.3.2 Densité spectrale de puissance


Quelques commentaires relatifs à Parseval uniquement :

Z T /2 +∞
1 2
X
P = |x(t)| dt = |X(k)|2
T −T /2 k=−∞

On aimerait avoir une représentation de la distribution de la puissance suivant un axe fréquentiel


comme pour les signaux d’énergie finie. Utilisation des distributions nécessaire ici. Sinon les
distributions sont utiles dans le reste du cours uniquement pour obtenir quelques formules sous
une forme simple ⇒ utilisation très réduite et très formelle.
On appelle “densité spectrale de puissance” (ou spectre) du signal périodique x(t) de fréquence
fondamentale f0
+∞
X
SX (f ) = |X(k)|2 δ(f − kf0 )
k=−∞

Le spectre d’un signal périodique est un “spectre de raies”.

Nicolas Moreau 1er octobre 2007


page 7/78
Licence de droits d’usage
CHAPITRE 2. SIGNAUX À TEMPS CONTINU : ANALYSE ET FILTRAGE

Quelques transformées de Fourier “au sens des distributions” :

δ(t)
1(f )
j2πf0 t

δ(f − f0 )
e
a
a cos(2πf0 t)
[δ(f − f0 ) + δ(f + f0 )]
2
+∞ +∞ +∞
X 1 X j2π n t 1 X k
δ(t − nT ) = e T
δ(f − ).
n=−∞
T n=−∞ T T
k=−∞

2.3.3 Exemple
On désire avoir une représentation fréquentielle d’un signal d’horloge de la forme :
+∞
X
y(t) = x(t − pT2 )
p=−∞

avec
x(t) = rectT1 (t).

On a
+∞ Z T2 /2
X 1
y(t) = Y (k) e j2πkt/T2
avec Y (k) = y(t) e−j2πkt/T2 dt
T2 T2 /2
k=−∞
Z +∞
1 1 k T1 kT1
Y (k) = x(t) e−j2πkt/T2 dt = X(f = )= sinc( )
T2 −∞ T2 T2 T2 T2
ce qui donne comme transformée de Fourier “au sens des distributions”
+∞
T1 X kT1 k
Y (f ) = sinc( ) δ(f − ).
T2 T2 T2
k=−∞

1.2

0.8

0.6

0.4

0.2

−0.2

−0.4
−4 −3 −2 −1 0 1 2 3 4

Fig. 2.1 – Transformée de Fourier de la fonction rectangle x(t) et de la fonction rectangle


périodisée y(t).

Nicolas Moreau 1er octobre 2007


page 8/78
Licence de droits d’usage
2.4. FILTRAGE LINÉAIRE

2.4 Filtrage linéaire


2.4.1 Définitions
Un système est un organe physique (une boîte noire) qui transforme un signal d’entrée en un
signal de sortie. Existence d’une relation fonctionnelle : y(t) = F{x(t)}.
Exemples : système phonatoire (entrée : mouvement des cordes vocales, sortie : variation
d’une pression acoustique), chaîne Hi-Fi (ce qui va nous intéresser c’est surtout l’égaliseur et non
l’ampli qui est modélisable par un simple gain).
Représentation par des diagrammes fonctionnels. Vision très ingénieur (ce que cache la boîte
est très relatif). Systèmes en chaîne directe. Systèmes asservis (exemple : rebouclage HP-micro
avec risque d’effet Larsen).

2.4.2 Propriétés d’un système


Sans/avec mémoire : Exemple : résistance/capacité. Premier cas très peu intéressant.
Linéarité : Il vérifie le principe de superposition.
Si x1 (t) → y1 (t) et si x2 (t) → y2 (t) alors αx1 (t) + βx2 (t) → αy1 (t) + βy2 (t) ∀α, β.
Notion très relative, exemple : ampli. Les amplitudes des signaux d’entrée ne doivent pas
être trop élevées sinon on risque des saturations, ni trop faibles sinon on risque de manquer
de précision. La notion de modèle est très relative. Un “bon” modèle est le résultat d’un
compromis simplicité-performances.
Invariance : Les caractéristiques de la relation d’entrée-sortie ne changent pas au cours du
temps. Si x(t) → y(t) alors x(t + τ ) → y(t + τ ) ∀τ . Exemple d’un système non-invariant :
système phonatoire.
Causalité : A chaque instant, la sortie ne dépend que de l’entrée à des instants présents ou
passés. Systèmes non-anticipatifs. Propriété respectée pour tous les systèmes réels mais
pas obligatoire lors de simulation. Exemple : amélioration d’enregistrements dégradés.
Stabilité : Nombreuses définitions : stabilité locale, globale, asymptotique. Ici, simplement dé-
finition de la stabilité au sens “entrée bornée - sortie bornée” :
Si ∀M, t tels que |x(t)| < M ⇒ ∃N, t0 tels que |y(t)| < N ∀t > t0 , alors le système sera
stable “EB-SB”.
Dans ce cours, uniquement étude des systèmes linéaires, invariants et stables ⇒ filtre.

2.4.3 Caractérisation d’un filtre


On a cherché à représenter un signal comme une combinaison linéaire de quelques signaux de
base. Il suffit de connaître la réponse d’un filtre à ces signaux de base. Pour connaître la sortie
à une entrée quelconque, il suffira d’appliquer le principe de superposition.

Réponse impulsionnelle
On appelle réponse impulsionnelle d’un filtre, la réponse du filtre à l’entrée particulière x(t) =
δ(t). On la notera par la suite h(t). Cette réponse est-elle suffisante pour caractériser le filtre,
i.e. connaissant h(t) peut-on calculer y(t) quelque soit x(t) ?
Pseudo-démonstration via l’intégrale de Riemann. On suppose x(t) “suffisamment régulière”
pour pouvoir définir
+∞ +∞
X X 1
xT (t) = x(nT ) rectT (t − nT ) = x(nT ) rectT (t − nT ) T.
n=−∞ n=−∞
T

Nicolas Moreau 1er octobre 2007


page 9/78
Licence de droits d’usage
CHAPITRE 2. SIGNAUX À TEMPS CONTINU : ANALYSE ET FILTRAGE

Appelons hT (t) la réponse du filtre à l’entrée (1/T )rectT (t). En appliquant la propriété de linéa-
rité et d’invariance, on obtient
+∞
X
xT (t) → yT (t) = x(nT ) hT (t − nT ) T
n=−∞

Par passage à la limite et en supposant des propriétés de “continuité”, on obtient


Z +∞
x(t) → y(t) = x(τ ) h(t − τ ) dτ = h(t) ∗ x(t).
−∞

On dit qu’un filtre est un convolueur.

Remarques :
– Filtre causal. La réponse impulsionnelle est nulle pour t < 0. On a
Z t
y(t) = x(τ ) h(t − τ ) dτ.
−∞
R +∞
– Filtre stable. On montre qu’une CNS est que −∞ |h(τ )| dτ < ∞.
– Mémoire d’un filtre causal et stable. Le support de h(t) est l’intervalle [0, ∞[ à cause
de la causalité et h(t) tend vers zéro à cause de la stabilité. Dans la pratique, la réponse
impulsionnelle peut être considérée comme nulle à partir d’un certain instant t0 . On appelle
alors l’intervalle [0, t0 ] la mémoire du filtre (notion très relative). Dans la construction de
y(t) seules les valeurs de l’entrée dans l’intervalle [t − t0 , t] interviendront.

Réponse en fréquence, gain complexe


Entrée particulière : x(t) = x0 ej2πf0 t . Sortie correspondante :
Z +∞ Z +∞
y(t) = x0 h(τ )ej2πf0 (t−τ ) dτ = x0 ej2πf0 t h(τ )e−j2πf0 τ dτ
−∞ −∞
Z +∞
y(t) = x(t) H(f0 ) avec H(f ) = h(τ )e−j2πf τ dτ.
−∞

On appelle H(f ) la réponse en fréquence ou gain complexe du filtre. C’est la transformée de


Fourier de la réponse impulsionnelle.

Remarques :
– La relation y(t) = x(t) H(f0 ) si x(t) = x0 ej2πf0 t veut dire que les exponentielles complexes
sont les “fonctions propres” des filtres et que H(f0 ) est la “valeur propre” correspondante.
– Il est possible de caractériser un filtre en faisant une série de mesure. On fait une “analyse
harmonique” en mesurant l’amplitude y0 et le déphasage φ0 de la sinusoïde de sortie et en
écrivant |H(f0 )| = y0 /x0 et Arg H(f0 ) = φ0 . Lieu de transfert.
– La relation Y (f ) = H(f ) X(f ) traduit l’opération de filtrage. Suivant la forme du module
de H(f ), on parle de filtres passe-bas, passe-haut, passe-bande ...
– Dans la pratique, c’est l’application qui impose la forme du module de H(f ). Le problème
consiste alors à déterminer les caractéristiques du filtre dans le domaine temporel. On
parle de synthèse de filtre. Définition d’un gabarit nécessaire car un filtre passe-bas du
type H(f ) = rectT (f ) impossible.
La réponse impulsionnelle du filtre n’est pas la bonne réponse car on ne sait pas “implanter”
un produit de convolution. On montre qu’une équation différentielle à coefficients constants

Nicolas Moreau 1er octobre 2007


page 10/78
Licence de droits d’usage
2.4. FILTRAGE LINÉAIRE

est une autre façon de caractériser un filtre et que cette forme (ou plutôt) l’équation
intégrale correspondante s’implante à l’aide d’amplificateurs opérationnels, de résistances
et de capacités.
– Tous ces problèmes seront abordés lors de l’étude des systèmes à temps discret (étude plus
simple dans ce cas) et on dira que pour les systèmes à temps continu c’est la même chose
(plutôt que l’inverse).

2.4.4 Exos
Intégrateur
Système défini par la relation
Z t
1
y(t) = x(τ ) dτ.
T t−T

Lissage des variations brusques de x(t).


Système linéaire, invariant, causal, stable.
Réponse impulsionnelle : h(t) = 1/T si t ∈ [0, T ], h(t) = 0 sinon.
Réponse en fréquence H(f ) = sinc(f T ) e−jπf T . Interprétation ...

Filtre passe-bas idéal

Etude du filtre passe-bas idéal H(f ) = rect2B (f ) e−j2πf t0 soumis à l’entrée x(t) = rectT (t).
Réponse impulsionnelle : h(t) = 2B sinc(2B(t − t0 )). Par la suite t0 = 0 pour simplifier.
Réponse à l’entrée rectangulaire :
Z T /2 Z 2B(t+T /2)
y(t) = 2B sinc(2B(t − τ )) dτ = sinc(u) du
T /2 2B(t−T /2)
Z x
y(t) = g(2B(t + T /2)) − g(2B(t − T /2)) avec g(x) = sinc(u) du
0
Le tracé de g(x) puis de y(t) permet de constater que le temps de montée d’un filtre passe-bas
est de l’ordre de 1/B.

Nicolas Moreau 1er octobre 2007


page 11/78
Licence de droits d’usage
CHAPITRE 2. SIGNAUX À TEMPS CONTINU : ANALYSE ET FILTRAGE

Nicolas Moreau 1er octobre 2007


page 12/78
Licence de droits d’usage
Chapitre 3

Echantillonnage

3.1 Introduction
Numérisation = discrétisation de l’axe temporel (échantillonnage) + discrétisation des am-
plitudes (quantification, codage de source).
Intérêt d’un signal mis sous forme numérique : grande immunité au bruit, stockage (reproduc-
tibilité indéfinie), transmission possible avec a priori un taux d’erreur aussi faible que l’on veut,
possibilité de traitements très élaborés, arrivée du “tout numérique” pour des raisons économiques
(existence de circuits performants et peu coûteux).
Inconvénients : distorsions introduites lors de la numérisation, plus large occupation spectrale
lors de la transmission.
Exemple : CD. Bande “Hi-Fi” [20 Hz - 20 kHz]. Echantillonnage caractérisé par un para-
mètre : la fréquence d’échantillonnage fe = 1/T = 44.1 kHz. Quantification caractérisée par un
paramètre : la résolution b = 16 bits/ech. Débit = 44.1 × 16 = 705 kbits/s (par voie). Nécessité
de codes correcteurs d’erreur (codage de canal). Lecteur de CD ≡ récepteur d’une chaîne de
communication. Exemple d’un traitement : “MPEG-Audio”. Actuellement possibilité de réduire
le débit d’un signal par un facteur 10 (ordre de grandeur) sans perte de qualité grâce à un
traitement sophistiqué.

3.2 Théorème d’échantillonnage


3.2.1 Problème
On part d’un signal à temps continu x(t) dont on suppose connue la transformée de Fourier
X(f ). On construit une suite {x(nT )} en prélevant des valeurs à des instants régulièrement
espacés (multiples de T ). Quelles conditions à imposer à x(t) (un paramètre important : sa
largeur de bande B) et à T (ou fe ) pour qu’il n’y ait pas perte d’information, i.e. que l’on garde
la possibilité de reconstruire x(t) à partir des échantillons {x(nT )} ?

3.2.2 Résultat intermédiaire : formule de Poisson


Si X(f ) est la transformée de Fourier du signal à temps continu x(t), quel que soit T et quel
que soit le support de X(f ), on a
+∞ +∞
1 X k X
X(f − ) = x(nT )e−j2πf nT (3.1)
T T n=−∞
k=−∞

à la condition que ces sommations aient un sens (aucune difficulté dans la pratique). On appellera
Y (f ) cette expression.

Nicolas Moreau 1er octobre 2007


page 13/78
Licence de droits d’usage
CHAPITRE 3. ECHANTILLONNAGE

Démonstration : Y (f ) est une fonction périodique, de période 1/T , donc développable en


série de Fourier
+∞ +∞
1 X k X
Y (f ) = X(f − ) = cn ej2πf nT
T T n=−∞
k=−∞
avec
Z 1/2T +∞ +∞ Z 1/2T −k/T
1 X k X
cn = T X(f − )e−j2πf nT df = X(u)e−j2πunT du
−1/2T T T −1/2T −k/T
k=−∞ k=−∞

Z +∞
cn = X(u)ej2πu(−nT ) du = x(−nT )
−∞

ce qui entraîne (3.1).

X(f ) Y (f )
6 6

- -
f −fe +fe f

Fig. 3.1 – Transformée de Fourier périodisée.

3.2.3 Conditions suffisantes


Si on est dans le cas particulier montré figure 3.2, on constate qu’il suffit de multiplier Y (f )
par H(f ) = T × rectfe (f ) pour obtenir X(f ).

Y (f ) H(f )
6

-
−fe −B +B +fe f

Fig. 3.2 – Conditions suffisantes pour qu’il n’y ait pas perte d’information.

Nicolas Moreau 1er octobre 2007


page 14/78
Licence de droits d’usage
3.3. REMARQUES ET INTERPRÉTATIONS

On a donc
+∞
1 X k
X(f ) = H(f ) X(f − ).
T T
k=−∞
En appliquant la formule de Poisson, on obtient
+∞
X
X(f ) = x(nT )H(f )e−j2πf nT
n=−∞

ce qui entraîne
+∞
X
x(t) = x(nT )h(t − nT )
n=−∞
avec Z 1/2T
t
h(t) = T ej2πf t df = sinc( ).
−1/2T T
On en déduit le “théorème d’échantillonnage ”. Si les deux conditions suivantes sont respectées :
– x(t) est un signal à bande limitée [−B, +B], i.e. si X(f ) = 0 pour |f | > B
– et si fe ≥ 2B,
alors il n’y a pas de perte d’information par échantillonnage parce que l’on peut reconstruire
exactement le signal x(t) à partir de ses échantillons
+∞
X t − nT
x(t) = x(nT )sinc( ). (3.2)
n=−∞
T
C’est la “formule d’interpolation” dont le principe est visualisé figure 3.3 en se limitant à une
somme de 3 termes. Interpolation au sens “signal à bande la plus limitée”. La plus petite fréquence
0.3

0.25

0.2

0.15

0.1

0.05

−0.05

−0.1

−0.15

−0.2
0 20 40 60 80 100 120 140 160 180 200

Fig. 3.3 – Interpolation.

d’échantillonnage possible s’appelle la “fréquence de Nyquist”.

3.3 Remarques et interprétations


3.3.1 Recouvrement des spectres
Si le signal x(t) n’est pas à bande limitée, alors
fe
Y (f ) 6= X(f ) pour |f | < .
2

Nicolas Moreau 1er octobre 2007


page 15/78
Licence de droits d’usage
CHAPITRE 3. ECHANTILLONNAGE

Phénomène de “recouvrement des spectres” (aliasing). Il s’agit la plupart du temps d’une dégra-
dation (phénomène désagréable pour des signaux de parole ou de musique par exemple) excepté
de rares cas où ce phénomène est directement exploité (exemple de la stroboscopie).
En effet (sous forme d’exo)
Soit x1 (t) = cos(2πf1 t). On choisit fe = 4f1 . Expression de x1 (nT ), X1 (f ) et Y1 (f ) ?
Soit x2 (t) = cos(2πf2 t) avec f2 = f1 + fe . Expression de x2 (nT ), X2 (f ) et Y2 (f ) ?
Quel est le signal dont la transformée de Fourier est rectfe (f )Y2 (f ) ? Quel est le signal dont
la transformée de Fourier est [rectfe /2 (f − 5fe /4) + rectfe /2 (f + 5fe /4)]Y2 (f ) ?
Conclusion : après l’échantillonneur, on ne fait pas la différence. La fréquence f2 est vue
comme la fréquence f1 . L’échantillonneur est un opérateur modulo fe .

3.3.2 Fréquence numérique


Connaissant la suite {x(nT )}, il existe une infinité de fonctions prenant les valeurs x(nT ) aux
instants d’échantillonnage. Exemple précédent : toutes les fonctions xk (t) = cos(2π(f1 + kfe )t)
prennent les valeurs {1, 0, −1, 0, · · · } aux instants nT . On dira que la suite {1, 0, −1, 0, · · · } a
toute sa puissance concentrée aux deux “fréquences” +1/4 et -1/4. On introduit le paramètre

ftc
ftd = .
fe

La fréquence “analogique” (à temps continu) ftc est exprimable en Hertz. La fréquence “numéri-
que” (à temps discret) ftd est sans dimension.

3.3.3 Conditions non nécessaires (sous forme d’exo)


Exemple : On considère un signal x(t) à bande limitée [−B, +B]. On multiplie ce signal par
cos(10πBt). La fréquence d’échantillonnage fe ≥ 12B est-elle nécessaire ?
Réponse : non. Si fe = 4B, 6B ou 7B, la reconstruction reste possible.

3.4 Dans la pratique


3.4.1 Deux problèmes
– Un signal à bande limitée n’existe pas. Contre-exemple apparent : une sinusoïde mais elle
est de durée infinie. Dans la pratique, on ne peut observer une sinusoïde que pendant une
durée finie ⇒ sinusoïde×fenêtre, i.e., dans le domaine fréquentiel, deux sinus cardinaux
centrés en −f1 et +f1 .
Il faut toujours filtrer avant échantillonnage. On montre que le filtre G(f ) minimisant
l’erreur quadratique
Z +∞
|x(t) − x̂(t)|2 dt
−∞

est G(f ) = rectfe (f ) mais il est lui-même irréalisable !


– La formule d’interpolation est irréaliste. Pour reconstruire le signal à temps continu à
l’instant t à partir de ses échantillons, il faut connaître tous les échantillons de −∞ à
+∞. Il en résulte deux nouveaux problèmes. La limitation à un nombre fini d’échantillons
entraîne une approximation

+N
X t − nT
x̂(t) = x(nT )sinc( ) ≈ x(t).
T
n=−N

Nicolas Moreau 1er octobre 2007


page 16/78
Licence de droits d’usage
3.5. AUTRE EXO

Le “délai de reconstruction” est égal à N T . C’est grave ou pas suivant les applications. Com-
munications mono-directionnelles (diffusion) : pas grave. Communications bi-directionnelles
(conversation téléphonique) : grave si N T est trop important.

3.4.2 Bloqueurs
Dans la pratique, on emploie des bloqueurs
Bloqueur d’ordre 0 : x̂(t) = x(nT ) pour nT ≤ t < (n + 1)T .
Bloqueur d’ordre 1 : x̂(t) = ax(nT ) + bx((n − 1)T ) pour nT ≤ t < (n + 1)T .
Etc.
On peut chercher à caractériser le type de dégradation apportée lorsque l’on utilise un bloqueur
d’ordre 0 (convertisseur N/A). Il s’exprime sous la forme
+∞
X
x̂(t) = x(nT )h(t − nT )
n=−∞

avec h(t) = rectT (t − T /2). Donc


+∞ +∞
X 1 X k
X̂(f ) = x(nT )H(f )e−j2πnT = H(f ) X(f − )
n=−∞
T T
k=−∞
avec
1
H(f ) = sinc(f T ).
T
Le tracé de la figure 3.4 montre que la dégradation apportée par l’utilisation d’un bloqueur
d’ordre 0 est une distorsion dans la bande [−fe /2, fe /2], une création de puissance hors de cette
bande. A comparer avec la formule d’interpolation. Cas dual. Autres formules d’interpolation :

-
−fe +fe f

Fig. 3.4 – Bloqueur d’ordre 0.

fonctions splines, ...

3.5 Autre exo


Phénomène de stroboscopie. On considère un signal “haute fréquence” périodique de période
T0 très petite. On l’échantillonne à la fréquence fe = 1/T avec T = T0 + ∆. Montrer qu’il est
possible de “reconstruire” le signal à l’aide des échantillons x(nT ) de telle façon que l’on obtienne
un signal y(t) de même forme mais plus “lent” i.e. y(t) = x(αt) avec α < 1.
Exemple : x(t) = a + b cos(2πt/T0 ).

Nicolas Moreau 1er octobre 2007


page 17/78
Licence de droits d’usage
CHAPITRE 3. ECHANTILLONNAGE

Nicolas Moreau 1er octobre 2007


page 18/78
Licence de droits d’usage
Chapitre 4

Transformées

4.1 Transformée de Fourier à temps discret (TFTD)


4.1.1 Définition
Existence de signaux à temps discret qui ne sont pas forcément le résultat de l’échantillonnage
d’un signal à temps continu. Notation de tels signaux {x(n)} avec n ∈ Z. Equivalence avec
{x(nT )} : il suffit de poser T = 1. On cherche un outil permettant de les analyser : on utilise la
formule de Poisson.
Définition. On appelle transformée de Fourier à temps discret d’un signal à temps discret
{x(n)} l’expression (si elle existe)
+∞
X
X(f ) = x(n)e−j2πf n . (4.1)
n=−∞

X(f ) est une fonction périodique de période 1 (il suffit de connaître X(f ) pour f ∈ [−1/2, 1/2]).
La formule précédente est simplement le développement d’une fonction périodique en série de
Fourier. On a donc Z 1/2
x(n) = X(f )e+j2πf n df. (4.2)
−1/2

4.1.2 Justification
Considérons un signal x(t) et sa transformée de Fourier à temps continu Xtc (f ). On a la
relation Z +∞
x(t) = Xtc (f )e+j2πf t df.
−∞
Si x(t) est à bande limitée [−B, +B] et si 1/T = fe ≥ 2B, on peut écrire
Z +fe /2
x(nT ) = Xtc (f )e+j2πf n/fe df.
−fe /2

Si on fait la distinction entre “fréquences à temps continu” exprimées en Hz et “fréquences à


temps discret” (des fréquences normalisées) ftd = ftc /fe , on obtient la relation
Z +1/2
x(nT ) = fe Xtc (ftd fe )e+j2πftd n dftd
−1/2

qui est cohérente avec l’équation (4.2).


Si, à partir d’un signal à temps discret {x(n)}, on construit un signal à temps continu x(t)
qui prenne les valeurs x(n) aux instants nT avec T = 1 avec la contrainte : il est à bande limitée
[−fe /2, fe /2] avec fe = 1, alors

Nicolas Moreau 1er octobre 2007


page 19/78
Licence de droits d’usage
CHAPITRE 4. TRANSFORMÉES

– pour f ∈ [−1/2, 1/2], Xtd (f ) = Xtc (f ) (à un coefficient fe près).


– pour f 3 [−1/2, 1/2], Xtc (f ) = 0 et Xtd (f ) = Xtd (f mod 1).

4.1.3 Propriétés
On donne, table 4.1, les principales propriétés de la transformée de Fourier à temps discret
en les comparant à celles de la transformée de Fourier à temps continu (sans démonstration car
ces démonstrations sont très simples).

Transformée de Fourier Transformée de Fourier


à temps continu à temps discret
Définition x(t)
X(f ) x(n)
X(f )
R +∞ −j2πf t dt
P+∞ −j2πf n
−∞ x(t)e n=−∞ x(n)e
Modulation e j2πf t
0 x(t) ej2πf n
0 x(n)

X(f − f0 ) X(f − f0 )
Translation x(t − t0 ) x(n − n0 )
temporelle e−j2πf t0 X(f ) e −j2πf n0 X(f )
R +∞ P+∞
Convolution −∞ x(τ )y(t − τ )dτ m=−∞ x(m)y(n − m)
X(f )Y (f ) X(f )Y (f )
Produit x(t)y(t) x(n)y(n)
R +∞ R +1/2
−∞ X(λ)Y (f − λ)dλ −1/2 X(λ)Y (f − λ)dλ
Signal réel X(−f ) = X ∗ (f ) X(−f ) = X ∗ (f )
R +∞ 2
R +∞ 2
P+∞ 2
R +1/2 2
Parseval −∞ |x(t)| dt = −∞ |X(f )| dt n=−∞ |x(n)| = −1/2 |X(f )| df

Tab. 4.1 – Principales propriétés de la transformée de Fourier à temps continu et à temps discret.

4.1.4 Exemples
“Impulsion à temps discret”
x(n) = 1 si n = 0, x(n) = 0 sinon. On a X(f ) = 1 ∀f .
Signal rectangulaire
x(n) = 1 pour 0 ≤ n ≤ N − 1, sinon 0.
N −1
X 1 − e−j2πf N sin(πf N )
X(f ) = e−j2πf n = −j2πf
= e−jπf (N −1) .
1−e sin(πf )
n=0

A comparer avec la transformée de Fourier à temps continu de la fonction x(t) = rectT (t)
qui est égale à T sinc(f T ).
Remarque : x(n) peut être vu comme le résultat de l’échantillonnage de la fonction rectN T (t−
(N − 1)T /2) avec T = 1. On a donc
+∞
X
X(f ) = U (f − k) avec U (f ) = N sinc(f N ) e−jπf (N −1) .
k=−∞

Que se passe-t-il lorsque N → ∞ ? Convergence particulière : convergence non-uniforme,


phénomène de Gibbs.
Cosinus
x(n) = cos(2πf1 n) pour 0 ≤ n ≤ N − 1, sinon 0.
N −1 N −1
X
−j2πf n 1 X −j2π(f −f1 )n
X(f ) = cos(2πf1 n)e = [e + e−j2π(f +f1 )n ]
2
n=0 n=0

Nicolas Moreau 1er octobre 2007


page 20/78
Licence de droits d’usage
4.1. TRANSFORMÉE DE FOURIER À TEMPS DISCRET (TFTD)
10

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Fig. 4.1 – TFTD (module) d’un signal rectangulaire avec N = 10.

N −1
1 X −jπ(f −f1 )(N −1) sin(π(f − f1 )N ) sin(π(f + f1 )N )
X(f ) = [e + e−jπ(f +f1 )(N −1) ].
2 sin(π(f − f1 )) sin(π(f + f1 ))
n=0

Que se passe-t-il lorsque N → ∞ ? On admettra le résultat suivant

+∞ +∞
1 X −j2π(f −f1 )n 1 X
X(f ) = [e + e−j2π(f +f1 )n ] = [δ(f − f1 − k) + δ(f + f1 − k)]
2 n=−∞ 2
k=−∞

Problèmes mathématiques délicats. Convergence “au sens des distributions”.

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Fig. 4.2 – TFTD (module) du signal x(n) = cos(2πf1 n) pour 0 ≤ n ≤ N − 1, sinon 0 avec
f1 = 0.123 et N = 10.

4.1.5 Exos
Soit x(n) un signal à temps discret de bande [−1/9, 1/9]. Tracer symboliquement le spectre
du signal y(n) = x(n) × cos(4πn/3).

Nicolas Moreau 1er octobre 2007


page 21/78
Licence de droits d’usage
CHAPITRE 4. TRANSFORMÉES

4.2 Transformée de Fourier discrète (TFD)


4.2.1 Introduction
Plusieurs présentations possibles
– Comme
P+∞ une approximation de la TFTD. La transformée de Fourier à temps discret X(f ) =
−j2πf n
n=−∞ x(n)e n’est pas implantable directement dans une machine pour deux rai-
sons : la sommation est infinie et le paramètre f est à valeurs continues. Il faut se limiter
à un nombre fini d’échantillons et discrétiser l’axe fréquentiel.
– Comme une transformation unitaire avec une interprétation fréquentielle.

4.2.2 La TFD vue comme une approximation de la TFTD


Limitation à un nombre fini d’échantillons

Formalisation simple en employant une “fenêtre de pondération” de durée N .

y(n) = x(n) × v(n) avec v(n) = 0 pour n 6= 0 · · · N − 1

(choix d’une fenêtre non-centrée pour simplifier les notations par la suite). L’approximation
réalisée peut être étudiée en exploitant la relation
Z +1/2
Y (f ) = X(λ)V (f − λ)dλ.
−1/2

Apparition d’“ondulations”, cf figure 4.2. Existence de nombreuses fenêtres ayant des “bords plus
doux” (Hamming, etc, cf chapitre 5).

Discrétisation de l’axe fréquentiel

On évalue Y (f ) pour M fréquences uniformément réparties entre 0 et 1.

N −1
k X k
Y (f = )= x(n)e−j2π M n pour k = 0 · · · M − 1.
M
n=0

Quel est le type d’approximation réalisée ? Problème difficile.

Exemple

On part d’un cosinus de durée infinie x(n) = cos(2πf1 n). Sa transformée de Fourier à temps
discret est un double peigne de Dirac

+∞
1 X
X(f ) = [δ(f − f1 − k) + δ(f + f1 − k)].
2
k=−∞

La discrétisation de l’axe fréquentiel de la courbe de la figure 4.2 avec M = 20 donne le tracé


de la figure 4.3. On observe que les M valeurs Y (0) · · · Y ((M − 1)/M ) paraissent une mauvaise
approximation de X(f ) sauf si la fréquence f1 est un multiple de 1/M et si M = N .

Nicolas Moreau 1er octobre 2007


page 22/78
Licence de droits d’usage
4.2. TRANSFORMÉE DE FOURIER DISCRÈTE (TFD)
6

0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Fig. 4.3 – TFD (module) du signal x(n) = cos(2πf1 n) avec f1 = 0.123, N = 10 et M = 20.

Conclusion

On montre que la TFTC, la TFTD et la TFD donnent des résultats “cohérents” si le signal
à temps continu x(t) est périodique (période T0 ), à bande limitée (nombre fini de raies), si la
fréquence d’échantillonnage fe = 1/T est supérieure à la fréquence de Nyquist et si N est choisi
de façon que T0 = N T . Alors la TFTC est composée de N raies espacées de 1/N T et les N
valeurs de la TFTD sont précisément les coefficients du développement en série de Fourier.
Dans la pratique, ces conditions ne sont jamais vérifiées. On cherche par exemple à détecter
plusieurs sinusoïdes dans du bruit. Problème de la résolution spectrale. Condition suffisante
(pas forcément nécessaire) pour discriminer deux sinusoïdes : que leur fréquences (normalisées)
respectives f1 et f2 vérifient |f1 − f2 | > 2/N .

4.2.3 Définition directe de la TFD


On appelle transformée de Fourier discrète de la séquence x(0) · · · x(N −1) une autre séquence
X(0) · · · X(N − 1) avec
N
X −1
X(k) = x(n)e−j2πnk/N .
n=0

Ces N relations s’écrivent de façon matricielle


    
X(0) 1 1 ··· ··· 1 x(0)
 X(1)   1 W1 ··· ··· W (N −1)  x(1) 
= (4.3)
    
 .. .. .. .. .. ..  .. 
 .   . . . . .  . 
X(N − 1) 1 W (N −1) · · · ··· W (N −1)(N −1) x(N − 1)

où W = e−j2π/N est la N eme racine de l’unité.


La matrice précédente, notée A, caractérise une transformation linéaire. C’est une matrice
carrée de dimension
√ N ∗ N à valeurs complexes. Elle est composée de vecteurs orthogonaux 2 à
2 et de norme N . C’est une matrice unitaire, à un coefficient près,

1
A(A∗ )t = I.
N

Nicolas Moreau 1er octobre 2007


page 23/78
Licence de droits d’usage
CHAPITRE 4. TRANSFORMÉES

On peut donc calculer les N valeurs x(0) · · · x(N −1) à partir de X(0) · · · X(N −1). La transformée
de Fourier discrète inverse a pour expression
N −1
1 X
x(n) = X(k)e+j2πnk/N .
N
k=0

On remarque que x(n + mN ) = x(n). Toute suite périodique de période N peut s’écrire sous la
forme d’une somme de N exponentielles complexes.
Prenons l’exemple du “peigne de Dirac” à temps discret de période N ,
+∞
X
x(n) = λ(n − mN )
m=−∞

avec λ(n) = 1 si n = 0 et λ(n) = 0 sinon (symbole de Kronecker)1 . Ce signal s’écrit aussi


N −1
1 X k
x(n) = X(k)ej2π N n (4.4)
N
k=0
avec
N −1
X k
X(k) = λ(n)e−j2π N n = 1.
n=0
On obtient
+∞ N −1
X 1 X j2π k n
λ(n − mN ) = e N .
m=−∞
N
k=0
La transformée de Fourier à temps discret de x(n) a pour expression
+∞ N −1 N −1 +∞
X 1 X −j2π(f − k )n 1 X X k
X(f ) = e N = δ(f − − m)
n=−∞
N N m=−∞
N
k=0 k=0
+∞
1 X k
X(f ) = δ(f − ).
N N
k=−∞

4.2.4 Propriétés
On retrouve, table 4.2, toutes les propriétés habituelles à la condition de calculer tous les
indices modulo N (lorsque l’indice courant sort de l’intervalle [0 · · · N − 1]).
Remarque concernant le produit de convolution. Ce qui nous intéressera par la suite c’est le
produit de convolution
+∞
X
u(n) = x(m)y(n − m)
m=−∞
de deux signaux à temps discret de durée infinie. La propriété de “convolution” précédente est
relative à deux signaux de durée finie avec un indice courant calculé modulo N
N
X −1
v(n) = x(m)y[(n − m) mod N ].
m=0

Comme on sera amené par la suite à vouloir calculer u(n) à l’aide de v(n) en exploitant X(k)Y (k)
pour bénéficier de l’algorithme de FFT, on distinguera le produit de convolution à temps discret
u(n) = x(n) ? y(n) et le produit de convolution circulaire v(n) = x(n) ⊗ y(n).
1
Au niveau des notations, on distingue volontairement le symbole de Kronecker λ(n) qui ne présente aucun
problème mathématique de l’impulsion de Dirac δ(t) (ou δ(f )) qui réclame de grandes précautions.

Nicolas Moreau 1er octobre 2007


page 24/78
Licence de droits d’usage
4.2. TRANSFORMÉE DE FOURIER DISCRÈTE (TFD)

Définition x(0) · · · x(N − 1)


X(0) · · · X(N − 1)
PN −1 k
−j2π N n
n=0 x(n)e
l
Modulation ej2π N n x(n)
X[(k − l) mod N ]
Translation x[(n − n0 ) mod N ]
k
temporelle e−j2π N n0 X(k)
PN −1
Convolution m=0 x(m)y[(n − m) mod N ]
X(k)Y (k)
Produit x(n)y(n)
PN −1
[ l=0 X(l)Y [(k − l) mod N ]]/N
Signal réel X(N − k) = X ∗ (k)
PN −1 2
PN −1 2
Parseval n=0 |x(n)| = [ k=0 |X(k)| ]/N

Tab. 4.2 – Principales propriétés de la transformée de Fourier discrète.

4.2.5 Algorithme FFT


La matrice A a des propriétés très particulières ce qui entraîne l’existence d’un algorithme
rapide (FFT : Fast Fourier Transform). Montrons le principe de l’algorithme en choisissant N = 8.
Le système matriciel s’écrit
    
X(0) 1 1 1 1 1 1 1 1 x(0)
 X(1)   1 W 1 W 2 W 3 W 4 W 5 W 6 W 7   x(1) 
 X(2)   1 W 2 W 4 W 6 W 8 W 10 W 12 W 14   x(2) 
    
 X(3)   1 W 3 W 6 W 9 W 12 W 15 W 18 W 21   x(3) 
     
 X(4)   1 W 4 W 8 W 12 W 16 W 20 W 24 W 28   x(4)  .
 =  
    
 X(5)   1 W 5 W 10 W 15 W 20 W 25 W 30 W 35   x(5) 
    
 X(6)   1 W 6 W 12 W 18 W 24 W 30 W 36 W 42   x(6) 
X(7) 1 W 7 W 14 W 21 W 28 W 35 W 42 W 49 x(7)

On réalise une première décomposition suivant les indices pairs et les indices impairs de x(n)
       
X(0) 1 1 1 1 x(0) 1 1 1 1 x(1)
 X(1)   1 W 2 W 4 W 6   x(2)   W 1 W 3 W 5 W 7   x(3) 
 X(2)  =  1 W 4 W 8 W 12   x(4)  +  W 2 W 6 W 10 W 14   x(5) 
       

X(3) 1 W 6 W 12 W 18 x(6) W 3 W 9 W 15 W 21 x(7)

1 W 8 W 16 W 24 W 4 W 12 W 20 W 28
       
X(4) x(0) x(1)
 X(5)   1 W 10 W 20 W 30   x(2)   W 5 W 15 W 25 W 35   x(3) 
 X(6)  =  1 W 12 W 24 W 36   x(4)  +  W 6 W 18 W 30 W 42   x(5) 
       

X(7) 1 W 14 W 28 W 42 x(6) W 7 W 21 W 35 W 49 x(7)


En exploitant la propriété W 8n = W n , on obtient
       
X(0) x(0) 1 0 0 0 x(1)
 X(1)  0  x(2)
   0 W1 0 0  0  x(3) 
 X(2)  = A  x(4)
  + A  
  0 0 W2 0   x(5) 
X(3) x(6) 0 0 0 W3 x(7)
       
X(4) x(0) 1 0 0 0 x(1)
 X(5)  0  x(2)
   0 W1 0 0  0  x(3) 
−
 X(6)  = A  x(4)
  A  
  0 0 W2 0   x(5) 
X(7) x(6) 0 0 0 W3 x(7)

Nicolas Moreau 1er octobre 2007


page 25/78
Licence de droits d’usage
CHAPITRE 4. TRANSFORMÉES

Posons        
U (0) x(0) U (1) x(1)
 U (2)  x(2)   U (3)  x(3) 
 = A0   = A0 
 
  et  .
 U (4)   x(4)   U (5)   x(5) 
U (6) x(6) U (7) x(7)
On réalise une deuxième décomposition pour calculer [U (0), U (2), U (4), U (6)] (ou pour calculer
[U (1), U (3), U (5), U (7)]).
    
U (0) 1 1 1 1 x(0)
 U (2)   1 W 2 W 4 W 6   x(2) 
 U (4)  =  1 W 4 W 8 W 12   x(4) 
    

U (6) 1 W 6 W 12 W 18 x(6)
       
U (0) 1 1 x(0) 1 1 x(2)
= +
U (2) 1 W4 x(4) W2 W6 x(6)
1 W8 W 4 W 12
         
U (4) x(0) x(2)
= −
U (6) 1 W 12 x(4) W 6 W 18 x(6)
On obtient       
U (0) V (0) 1 0 V (2)
= +
U (2) V (4) 0 W2 V (6)
      
U (4) V (0) 1 0 V (2)
= −
U (6) V (4) 0 W2 V (6)
On réalise une troisième décomposition pour calculer [V (0), V (4)] etc.

V (0) = x(0) + x(4)

V (4) = x(0) − x(4).


On en déduit le diagramme de la figure 4.4 symbolisant l’ensemble des opérations à effectuer. Le
X(0) U(0) V(0) x(0)

0
X(1) U(2) V(4) W x(4)
-
0
X(2) U(4) W V(2) x(2)
-
2 0
X(3) U(6) W V(6) W x(6)
- -
0
X(4) W U(1) V(1) x(1)
-
1 0
X(5) W U(3) V(5) W x(5)
- -
2 0
X(6) W U(5) W V(3) x(3)
- -
3 2 0
X(7) W U(7) W V(7) W x(7)
- - -

Fig. 4.4 – Diagramme de calcul de l’algorithme FFT.

nombre de multiplications/accumulations (complexes) est égal à N ×γ à la condition que N = 2γ .


Comme le nombre de multiplications/accumulations pour effectuer le calcul brutal est égal à N 2 ,
on voit que cet algorithme est particulièrement performant. Par exemple lorsque N = 1024, le
rapport vaut
N2
≈ 100.
N log2 (N )
Nombreuses variantes : entrelacement temporel, fréquentiel, radix 2, 4, etc.

Nicolas Moreau 1er octobre 2007


page 26/78
Licence de droits d’usage
4.3. TRANSFORMÉE EN Z

4.3 Transformée en z
4.3.1 Introduction
Présentation “minimale”. Dans ce cours, la présentation de cette transformée n’est pas néces-
saire. Il est bon tout de même d’avoir un peu de vocabulaire car cette transformée est largement
utilisée par les “traiteurs de signaux”. Formalisme assez commode pour étudier des filtres “numé-
riques” (à temps discret).

4.3.2 Définition
On appelle transformée en z de la suite {x(n)} la fonction de la variable complexe z définie
par
+∞
X
X(z) = x(n)z −n
n=−∞
à la condition que cette somme ait un sens. Nécessité de définir lePdomaine de convergence. Si
on appelle ρ1 le rayon de convergence de la série entière X+ (z) =P +∞ n=0 x(n)z
−n i.e. que X (z)
+
existe pour |z −1 | < ρ1 et ρ2 le rayon de convergence de X− (z) = −∞ n=−1 x(n)z −n i.e. que X (z)

existe pour |z| < ρ2 , alors la transformée en z est définie sur la couronne 1/ρ1 < |z| < ρ2 .
Pour une suite à valeurs réelles ρ1 et ρ2 sont des réels positifs. Existence d’une correspondance
biunivoque entre une suite {x(n)} et une série X(z) uniquement à l’intérieur du domaine de
convergence (s’il est non vide) à cause de l’unicité du développement en série de Laurent pour
une fonction holomorphe dans une couronne.
Intérêt de la transformée en z : être une représentation “compacte” de l’ensemble des valeurs
numériques si on réussit à sommer la série ⇒ manipulation très aisée.
Si le cercle unité appartient au domaine de convergence, alors la transformée de Fourier à
temps discret existe et on a la relation

X(f ) = X(z)|z=ej2πf

avec l’abus d’écriture habituel. La transformée en z apparaît comme une généralisation de la


transformée de Fourier à temps discret.

4.3.3 Exemples
Exemple d’une suite causale : x(n) = 0 pour n < 0, x(n) = an pour n ≥ 0 avec a ∈ R.
+∞
X 1
X(z) = an z −n =
1 − az −1
n=0

à la condition que
|az −1 | < 1 ⇒ |z| > |a|.
Le domaine de convergence est donc dans ce cas l’extérieur du cercle de rayon |a|.
Exemple d’une suite anticausale : x(n) = 0 pour n ≥ 0, x(n) = −an pour n < 0 avec a ∈ R.
−1
X +∞
X
X(z) = − an z −n = 1 − a−n z n
n=−∞ n=0

1 1
X(z) = 1 − −1
=
1−a z 1 − az −1
à la condition que
|a−1 z| < 1 ⇒ |z| < |a|.

Nicolas Moreau 1er octobre 2007


page 27/78
Licence de droits d’usage
CHAPITRE 4. TRANSFORMÉES

Le domaine de convergence est l’intérieur du cercle de rayon |a|.


Conclusion : la condition de convergence impose “le sens du temps”. On parle de suites
“causales” (x(n) = 0 pour n < 0) et de suites “anti-causales” (x(n) = 0 pour n ≥ 0). Par la suite,
le domaine de convergence sera implicite. Si la suite est causale, le domaine de convergence est
l’extérieur d’un cercle. Si la suite est anti-causale, le domaine de convergence est l’intérieur d’un
cercle. Pour une suite quelconque, le domaine de convergence est une couronne. Tout le problème
sera de savoir si le cercle unité appartient ou non au domaine de convergence.

4.3.4 Inversion de la transformée en z


De façon générale : utilisation de l’intégrale de Cauchy pour une fonction holomorphe dans
une couronne Z
1
x(n) = X(z)z n−1 dz
2πj Γ
où Γ est un contour fermé appartenant au domaine d’holomorphie (le cercle unité pour les “bons”
signaux).
Habituellement X(z) se présente sous la forme d’une fraction rationnelle de deux polynômes
en z. Expression habituelle pour une suite causale

b0 + b1 z −1 + · · · + bQ z −Q
X(z) = .
1 + a1 z −1 + · · · + aP z −P

avec P et Q deux entiers positifs quelconques. La méthode la plus standard consiste à réaliser
une décomposition en éléments simples.
Exemple : quel est le signal causal ayant comme transformée en z
z−a
X(z) = .
(z − b)(z − c)

On a
b−a c−a
X(z) = +
(b − c)(z − b) (c − b)(z − c)
Comme on cherche un signal causal, il faut développer suivant les puissances négatives de z. On
écrit
z −1 b−a a−c
X(z) = [ −1
+ ]
b − c 1 − bz 1 − cz −1
z −1
X(z) = [(b − a)(1 + bz −1 + b2 z −2 + · · · ) + (a − c)(1 + cz −1 + c2 z −2 + · · · )].
b−c
Donc
1
x(n) = [(b − a)bn−1 + (a − c)cn−1 ].
b−c

4.3.5 Quelques propriétés


Translation temporelle
Si y(n) = x(n − k), alors
+∞
X
Y (z) = x(n − k)z −n = z −k X(z).
n=−∞

Remarque : z −1 apparaît comme l’opérateur de retard élémentaire.

Nicolas Moreau 1er octobre 2007


page 28/78
Licence de droits d’usage
4.3. TRANSFORMÉE EN Z

Produit de convolution
Si
+∞
X
y(n) = x(k)h(n − k)
k=−∞

alors
+∞
X +∞
X
Y (z) = x(k)h(n − k)z −n = X(z)H(z).
n=−∞ k=−∞

Généralement pas de problème dû au domaine de convergence car pour les signaux habituels
le cercle unité appartient au domaine de convergence de x(n) et de h(n) donc de y(n).
Autres propriétés
Translation temporelle, produit, parité, Parseval.

Nicolas Moreau 1er octobre 2007


page 29/78
Licence de droits d’usage
CHAPITRE 4. TRANSFORMÉES

Nicolas Moreau 1er octobre 2007


page 30/78
Licence de droits d’usage
Chapitre 5

Signaux à temps discret : Filtrage

5.1 Filtres discrets (numériques)


5.1.1 Définition
Système discret : à partir d’un signal à temps discret {x(n)} création d’un autre signal à
temps discret {y(n)} obéissant à une relation fonctionnelle
y(n) = F(· · · , x(n − 1), x(n), x(n + 1), · · · ; n).
Comme pour les systèmes à temps continu, si un système est linéaire (il vérifie le principe de
superposition) et invariant (la relation fonctionnelle est indépendante de n), alors le système est
appelé un filtre. Un filtre peut être causal (y(n) ne dépend pas de x(n + 1), · · · ) ou non et stable
(entrée bornée ⇒ sortie bornée) ou non.
On se limite dans le cadre de ce cours à la sous-classe des filtres décrits par une équation
récurrente à coefficients (réels) constants
y(n) + a1 y(n − 1) + · · · + aP y(n − P ) = b−Q1 x(n + Q1 ) + · · · + b0 x(n) + · · · + bQ2 x(n − Q2 )
avec P, Q1 , Q2 trois entiers positifs quelconques.
Remarques
– Cette forme est restrictive (en théorie mais pas en pratique) car elle traduit le fait que l’on
s’intéresse uniquement à une solution causale (n croissant) : on calcule y(n) en fonction de
l’entrée et de la sortie uniquement aux instants précédents.
– Pour simplifier les notations par la suite, on supposera en plus le filtre causal
Q
X P
X
y(n) = bi x(n − i) − ai y(n − i). (5.1)
i=0 i=1

Un filtre dans une machine c’est constamment calculer cette expression (souvent en temps
réel, c’est à dire au rythme d’arrivée des x(n)).
– Si l’entrée est elle-même causale (on démarre la récurrence à l’instant n = 0), les conditions
initiales associées à l’équation récurrente doivent être prises nulles pour respecter le principe
de superposition. L’introduction de conditions initiales non nulles posent un problème
(représentation des systèmes en “variables d’état” non abordée dans ce cours).
– Terminologie
P = 0, Q > 0 : Filtre RIF (réponse impulsionnelle finie), MA (moving average), tout zéro,
non-récursif.
P > 0, Q = 0 : Filtre AR (auto-régressif), tout pôle.
P > 0, Q > 0 : Filtre RII (réponse impulsionnelle infinie), ARMA.
– Ce n’est qu’une façon de caractériser un filtre. Il en existe trois autres.

Nicolas Moreau 1er octobre 2007


page 31/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE

5.1.2 Réponse impulsionnelle, produit de convolution


On appelle réponse impulsionnelle d’un filtre la réponse de ce filtre à l’entrée particulière
λ(n) = 1 si n = 0 sinon 0 (symbole de Kronecker). On la notera par la suite systématique-
ment h(n). La connaissance de la réponse impulsionnelle est suffisante pour caractériser un filtre
puisque, connaissant h(n), quelle que soit l’entrée x(n), on peut calculer la sortie
+∞
X +∞
X
y(n) = h(n) ∗ x(n) = h(k)x(n − k) = h(n − k)x(k).
k=−∞ k=−∞

En effet si (cas d’une entrée causale, tous les indices commencent à 0)

{x(n)} = {1, 0, · · · } ⇒ {y(n)} = {h(0), h(1), · · · }

alors, à cause de la propriété d’invariance,

{x(n)} = {0, · · · , 0, 1, 0, · · · } ⇒ {y(n)} = {0, · · · , 0, h(0), h(1), · · · }.

On utilisant la propriété de linéarité, on voit que la sortie y(n) à l’instant n est la somme de la
contribution de x(n) pondérée par h(0), plus la contribution de x(n − 1) pondérée par h(1), plus
etc.
Exemple d’un filtre défini par (5.1) avec P = 0

h(n) = bn si 0 ≥ n ≥ Q
= 0 sinon.

La réponse impulsionnelle est à support fini d’où la terminologie filtre RIF. Dans ce cas particu-
lier, le produit de convolution et la relation de récurrence ont exactement la même forme
Q
X
y(n) = bi x(n − i).
i=0

Dans le cas général, il existe deux façons de calculer le même résultat. L’équation récurrente est
beaucoup plus intéressante dans la pratique car elle ne met en jeu que des sommations finies.

5.1.3 Réponse en fréquence, gain complexe


On soumet le filtre à l’entrée non causale x(n) = aej2πf1 n . On suppose connue la réponse
impulsionnelle h(n). On obtient
+∞
X +∞
X
y(n) = h(k)x(n − k) = aej2πf1 n h(k)e−j2πf1 k .
k=−∞ k=−∞

On reconnait la transformée de Fourier à temps discret de la réponse impulsionnelle évaluée à la


fréquence f = f1 . Cette transformée H(f ) s’appelle la réponse en fréquence ou le gain complexe
du filtre. On a
y(n) = aej2πf1 n H(f1 ) = a|H(f1 )|ej[2πf1 n+Arg(H(f1 ))] .
Les exponentielles complexes sont les fonctions propres des filtres.
Ce développement se généralise sans difficulté au cas d’une entrée comportant un nombre fini
d’exponentielles. Si
XL
x(n) = al ej(2πfl n+φl )
l=1

Nicolas Moreau 1er octobre 2007


page 32/78
Licence de droits d’usage
5.1. FILTRES DISCRETS (NUMÉRIQUES)

alors
L
X L
X
y(n) = al ej(2πfl n+φl ) H(fl ) = al |H(fl )|ej(2πf1 n+φl +Arg(H(f1 )) .
l=1 l=1
Ce développement se généralise aussi au cas d’une entrée comportant un nombre infini d’expo-
nentielles. On trouve
Z +1/2 Z +1/2
j2πf n
y(n) = Y (f )e df = H(f )X(f )ej2πf n df.
−1/2 −1/2

On retrouve la relation standard Y (f ) = H(f )X(f ).


La réponse en fréquence caractérise l’opération de filtrage. Elle est accessible à la mesure
(analyse harmonique).

5.1.4 Fonction de transfert


Cherchons la transformée en z de y(n) donnée par (5.1). On obtient
+∞
X Q
X +∞
X P
X +∞
X
−n −n
Y (z) = y(n)z = bi x(n − i)z − ai y(n − i)z −n
n=−∞ i=0 n=−∞ i=1 n=−∞

Q
X P
X
Y (z) = bi z −i X(z) − ai z −i Y (z)
i=0 i=1
PQ −i
Y (z) i=0 bi z
= .
X(z) 1 + Pi=1 ai z −i
P

Le rapport H(z) = Y (z)/X(z) est indépendant de l’entrée. On l’appelle la fonction de transfert


du filtre. Il apparaît comme une fraction rationnelle de deux polynômes en z −1 (filtre causal). Il
se factorise sous la forme QQ
(1 − βi z −1 )
H(z) = b0 QPi=1 .
−1
i=1 (1 − αi z )
On appelle {β1 · · · βQ } les zéros de la fonction de transfert et {α1 · · · αP } les pôles. Les coefficients
des polynomes en z −1 étant supposés à valeurs réelles, les zéros et les pôles apparaîssent par paires
(βi , βi∗ ) et (αi , αi∗ ).
H(z) s’interprète comme étant la transformée en z de y(n) lorsque X(z) = 1 i.e. x(n) = λ(n).
H(z) est donc la transformée en z de la réponse impulsionnelle.

5.1.5 Relations
– La réponse en fréquence est la transformée de Fourier à temps discret de la réponse impul-
sionnelle.
– La fonction de transfert est la transformée en z de la réponse impulsionnelle.
– Interprétation géométrique du module de la réponse en fréquence. Exemple

(1 − βz −1 )
H(z) = .
(1 − αz −1 )(1 − α∗ z −1 )

|(1 − βe−j2πf )| MC
|H(f )| = −j2πf ∗ −j2πf
=
|(1 − αe )|.|(1 − α e )| M A.M B
où M, A, B et C sot respectivement les images de ej2πf , des poles et du zéro comme
le montre le tracé de la figure 5.1. Cette interprétation est intéressante car elle permet

Nicolas Moreau 1er octobre 2007


page 33/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE
Im
6

A
M
+1/2 C 0-
Re
-1/2

Fig. 5.1 – Interprétation géométrique du module de la réponse en fréquence.

un tracé approximatif de |H(f )|. Toutes les propriétés des filtres peuvent être appréciées
qualitativement par simple inspection de la position des pôles et des zéros. Les pôles proches
du cercle unité entraînent des “pics”. Les zéros proches du cercle unité entraînent des
“vallées”.

5.1.6 Notion de stabilité


Problème difficile si on veut le traiter de façon rigoureuse. Notion de stabilité locale, globale,
asymptotique, non-asymptotique, etc. On dira simplement qu’un filtre est stable si à toute entrée
bornée correspond une sortie bornée.
Deux CNS :
– h(n) est de module sommable : +∞
P
n=−∞ |h(n)| < +∞
P+∞
En effet, la condition est suffisante. Si n=−∞ h(n) < A et si |x(n)| < B ∀n, alors
+∞
X +∞
X
|y(n)| = | h(k)x(n − k)| ≤ |h(k)||x(n − k)| ≤ AB
k=−∞ k=−∞

La condition est également nécessaire car si elle n’est pas vérifiée, le signal x(n) = signe[h(n)]
qui est borné fait diverger le filtre pour n = 0.
– On montre que pour un système causal, une autre CNS est que les pôles de la fonction de
transfert soient strictement à l’intérieur du cercle unité.
La réponse en fréquence d’un filtre n’a un sens que si le filtre est stable. La fonction de transfert
d’un filtre instable existe mais le domaine de convergence de H(z) n’inclut pas le cercle unité.
Existence de critères algébriques.

5.1.7 Exemples sous forme d’exo


Filtre AR du 1er ordre

y(n) + a1 y(n − 1) = x(n)


Réponse impulsionnelle :
{h(n)} = {1, −a1 , a21 , · · · , (−a1 )n , · · · }
Réponse en fréquence :
+∞
X
H(f ) = (−a1 )n e−j2πf n = 1/(1 + a1 e−j2πf ) si |a1 | < 1
n=0

Nicolas Moreau 1er octobre 2007


page 34/78
Licence de droits d’usage
5.1. FILTRES DISCRETS (NUMÉRIQUES)

Module :
|H(f )|2 = 1/(1 + 2a1 cos(2πf ) + a21 )
Comme
|H(0)| = 1/|1 + a1 | et |H(1/2)| = 1/|1 − a1 |
on remarque que l’on obtient un filtre passe-bas si −1 < a1 < 0 et un filtre passe-haut si
0 < a1 < 1.
Fonction de transfert :
H(z) = 1/(1 + a1 z −1 )
Pôle = −a1 (si à l’intérieur du cercle unité ⇒ stabilité sinon instabilité (la réponse impulsionnelle
est une suite divergente)).

Filtre RIF

y(n) = x(n) + 2x(n − 1) + 3x(n − 2) + 4x(n − 3) + 3x(n − 4) + 2x(n − 5) + x(n − 6)

Réponse impulsionnelle :
{h(n)} = {1, 2, 3, 4, 3, 2, 1, 0, · · · }
Fonction de transfert :

H(z) = 1 + 2z −1 + 3z −2 + 4z −3 + 3z −4 + 2z −5 + z −6 = [(1 + z −1 )(1 + z −2 )]2

Zéros doubles en -1, +j, -j.


Réponse en fréquence :

H(f ) = [(1 + e−j2πf )(1 + e−j4πf )]2 = [e−j3πf (ejπf + e−jπf )(e2jπf + e−j2πf )]2

H(f ) = 16 cos2 (πf ) cos2 (2πf ) e−j6πf .


Le tracé de la figure 5.2 montre que, très grossièrement, ce filtre laisse passer les basses fréquences
et coupe les hautes fréquences. C’est un filtre passe-bas (pas très efficace). On peut chercher à

16
20

14
10

12
0
[echelle lineaire]

10
−10
[dB]

8 −20

6 −30

4 −40

2 −50

0 −60
0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5

Fig. 5.2 – Module de la réponse en fréquence du filtre RIF suivant une échelle linéaire (à gauche)
ou une échelle en dB (à droite).

calculer la chute en dB entre le lôbe principal et le lôbe secondaire. On trouve 22 dB. C’est
peu. Pour filtrer du signal de parole en bande téléphonique, il faut au moins 50 dB, pour de

Nicolas Moreau 1er octobre 2007


page 35/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE

la musique 100 dB (l’oreille est sensible à des variations de puissance de plus de 100 dB). Par
contre, il supprime totalement la fréquence 1/4.
Exemple de l’aspect manipulatoire de la transformée en z

Y (z) = [(1 + z −1 )(1 + z −2 )][(1 + z −1 )(1 + z −2 )]X(z)

s’écrit

U (z) = [(1 + z −1 )(1 + z −2 )]X(z)


Y (z) = [(1 + z −1 )(1 + z −2 )]U (z)

ce qui entraîne

u(n) = x(n) + x(n − 1) + x(n − 2) + x(n − 3)


y(n) = u(n) + u(n − 1) + u(n − 2) + u(n − 3).

La factorisation de H(z) permet d’en déduire de façon exhaustive tous les systèmes d’équations
récurrentes équivalents à la relation initiale.

Filtre réjecteur

La fonction de transfert comportant deux zéros sur le cercle unité s’écrit

H(z) = (1 − z0 z −1 )(1 − z0∗ z −1 ) avec z0 = ej2πf0

soit
H(z) = 1 − 2 cos(2πf0 )z −1 + z −2

ce qui donne dans le domaine temporel

y(n) = x(n) − 2 cos(2πf0 )x(n − 1) + x(n − 2).

L’interprétation géométrique du module de la réponse en fréquence indique qu’un signal compor-


tant une composante à la fréquence f0 verra cette composante totalement supprimée. On peut
effectivement remarquer que la sortie du filtre est nulle quel que soit n pour une entrée de la
forme x(n) = 2 cos(2πf0 n + φ).

Filtre AR du 2ème ordre

y(n) + a1 y(n − 1) + a2 y(n − 2) = b0 x(n)

1 1
H(z) = = .
(1 − ρejφ z −1 )(1 − ρe−jφ z −1 ) 1 − 2ρ cos(φ)z −1 + ρ2 z −2

sin[(n + 1)φ]
h(n) = ρn
sin(φ)
On donne, figure 5.3, la réponse impulsionnelle et le module de la réponse en fréquence lorsque
ρ = 0.9 et φ = π/4.

Nicolas Moreau 1er octobre 2007


page 36/78
Licence de droits d’usage
5.2. SYNTHÈSE DES FILTRES
1.5 20

15
1

10

0.5

−0.5
−5

−1 −10
0 5 10 15 20 25 30 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5

Fig. 5.3 – Réponse impulsionnelle et module de la réponse en fréquence du filtre AR du 2ème


ordre lorsque ρ = 0.9 et φ = π/4.

5.2 Synthèse des filtres


Connaissant la réponse en fréquence H(f ) imposée par une application, on veut en déduire
les coefficients de l’équation récurrente correspondante pour pouvoir “implanter” le filtre dans un
processeur. Il existe de nombreuses méthodes. Dans le cadre de ce document, on se limite à la
méthode la plus simple (mais pas la plus efficace) : la méthode dite de la fenêtre. En particulier,
on impose que le filtre soit un filtre RIF (nécessaire si on veut “une phase linéaire”).
Exemple. On désire synthétiser un filtre passe-bas de fréquence de coupure 1/2M (utilisation
assez fréquente en TS). On parle de filtre “demi-bande” si M = 2, “quart de bande” si M = 4,
etc.

5.2.1 Principe
On rappelle que H(f ) est nécessairement une fonction périodique. Comme H(f ) est la trans-
formée de Fourier à temps discret de la réponse impulsionnelle, on en déduit directement
Z +1/2
h(n) = H(f )ej2πf n df.
−1/2

Exemple du filtre passe-bas (idéal) de fréquence de coupure 1/2M

H(f ) = M si |f | < 1/2M


= 0 si 1/2M < |f | < 1/2.

On en déduit Z +1/2M
sin(πn/M ) n
h(n) = M ej2πf n df = = sinc( ).
−1/2M πn/M M

En théorie, la sortie du filtre y(n) est directement donnée par

+∞
X
y(n) = h(k)x(n − k). (5.2)
k=−∞

On montre, figure 5.4, la réponse impulsionnelle d’un filtre quart de bande.

Nicolas Moreau 1er octobre 2007


page 37/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE
1.2

0.8

0.6

0.4

0.2

−0.2

−0.4
−20 −15 −10 −5 0 5 10 15 20

Fig. 5.4 – Réponse impulsionnelle d’un filtre quart de bande (M = 4) et fenêtre de Hamming
sur 41 échantillons.

5.2.2 Mise en œuvre


Limitation à un nombre fini de termes
Le produit de convolution (5.2) n’est pas directement implantable dans une machine. Il faut
se limiter à un nombre fini d’opérations. Pour garder la propriété de symétrie de h(n), on choisit

+N
X
ŷ(n) = h(k)x(n − k) (5.3)
k=−N

et non pas
+2N
X
ŷ(n) = h(k)x(n − k).
k=0

Pour évaluer l’approximation réalisée, il faut comparer la nouvelle réponse en fréquence Ĥ(f ) à
celle du filtre idéal. Comme
+∞
X +∞
X
ŷ(n) = h(k)v(k)x(n − k) = ĥ(k)x(n − k)
k=−∞ k=−∞

en utilisant une “fenêtre de pondération”, on obtient dans le domaine fréquentiel


Z 1/2
Ĥ(f ) = H(λ)V (f − λ)dλ.
−1/2

Le tracé de la figure 5.5 montre l’introduction d’ondulations lorsque l’on utilise une fenêtre
rectangulaire particulièrement au voisinage des sauts brusques de H(f ). Ces ondulations sont
dues à celles présentes dans V (f ) (cf figure 4.1). Pour atténuer l’importance de ces ondulations, on
peut être tenté d’augmenter le paramètre N . Les ondulations ne disparaîssent pas. Elles sont juste
localisées dans une bande de fréquence de plus en plus étroite. C’est le “phénomène de Gibbs”.
On voudrait que V (f ) ait un lobe principal le plus étroit possible et des lobes secondaires les plus
petits possibles. C’est contradictoire. Un bon compromis est donné par la fenêtre de Hamming,
montrée figure 5.4, comme on peut le constater, figure 5.5, où on compare la réponse en fréquence
Ĥ(f ) lorsque l’on utilise une fenêtre rectangulaire et une fenêtre de Hamming.
On donne, table 5.1, les “performances” de quelques fenêtres.

Nicolas Moreau 1er octobre 2007


page 38/78
Licence de droits d’usage
5.2. SYNTHÈSE DES FILTRES
5 20

4.5 10

4 0

3.5 −10

3 −20
[echelle lineaire]

[dB]
2.5 −30

2 −40

1.5 −50

1 −60

0.5 −70

0 −80
0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5

Fig. 5.5 – Réponse en fréquence d’un filtre quart de bande avec N 0 = 2N + 1 = 41 après
pondération par une fenêtre rectangulaire (en bleu) ou par une fenêtre de Hamming (en rouge).
Les amplitudes sont données en échelle linéaire (à gauche) ou en décibel (à droite).

Fenêtre Largeur Taux


du lobe principal d’ondulation
rectangulaire 2/N 0 22 %
Hamming 4/N 0 2%
Blackman 6/N 0 0.1 %

Tab. 5.1 – Performances de quelques fenêtres avec N 0 = 2N + 1.

Filtre causal

Le produit de convolution (5.3) caractérise un filtre non-causal. Si on veut rendre causal le


filtre, il suffit de différer le résultat de N échantillons en posant

ỹ(n) = ŷ(n − N ).

On obtient
{h̃(0) · · · h̃(2N )} = {ĥ(−N ) · · · ĥ(N )}

H̃(f ) = e−j2πf N Ĥ(f ).

C’est très artificiel. On n’a rien changé exceptée la définition de l’origine des indices n !

Gabarit

Le développement précédent montre que la réponse en fréquence d’un filtre ne peut pas être
parfaitement plate ni en bande passante, ni en bande affaiblie et que la bande de transition ne peut
pas être infiniment étroite. Il faut donc introduire des degrés de liberté : on définit un “gabarit”
en précisant le taux d’ondulation en bande passante, le taux d’ondulation en bande affaiblie (ou
le taux de réjection) et la largeur de la bande de transition en fonction de l’application.
Dans l’exemple précédent, le taux de réjection est de l’ordre de 50 dB si l’ordre du filtre RIF
est égal à N 0 = 41 et si on utilise une fenêtre de Hamming. C’est suffisant si on désire filtrer
du signal de parole en bande téléphonique (dynamique en puissance de l’ordre de 50 dB). C’est
insuffisant pour du signal de musique (dynamique en puissance de l’ordre de 100 dB). Il faut
alors augmenter l’ordre.

Nicolas Moreau 1er octobre 2007


page 39/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE

Implantation temps réel dans un processeur

Peut-on prendre des ordres très élevés ? Cela dépend de la fréquence d’échantillonnage du
signal et de la puissance de calcul du processeur utilisé. Par exemple, pour du signal de musique
échantillonné à 44.1 kHz et filtré dans un microprocesseur capable d’effectuer 107 multiplica-
tions/accumulations par seconde1 , on trouve que N 0 doit vérifier 44100N 0 ≤ 107 soit N 0 ≤ 200.

5.2.3 Autres méthodes de synthèse


Algorithme de Remez. La réponse impulsionnelle est obtenue en utilisant des algorithmes
d’optimisation
min[ max W (f )|H(f ) − Ĥ(f )|]
h(n) f ∈[0,1/2]

où W (f ) est une fonction de pondération. On montre, figure 5.6, la réponse en fréquence d’un
filtre quart de bande après pondération par une fenêtre de Hamming (en rouge) ou après utilisa-
tion de l’algorithme de Remez (en noir). L’ordre du filtre et la bande de transition ont été choisis
identiques. Conclusion !
20

10

−10

−20
[dB]

−30

−40

−50

−60

−70

−80
0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5

Fig. 5.6 – Réponse en fréquence d’un filtre quart de bande après pondération par une fenêtre de
Hamming (trait plein) ou après utilisation de l’algorithme de Remez (trait pointillé).

5.3 Un exemple applicatif


On dispose d’un signal de musique échantillonné à 48 kHz. On veut en obtenir une version
échantillonnée à 32 kHz minimisant la distorsion.
Techniques de sous et sur-échantillonnage largement utilisées en TS : bancs de filtres, codage
de source, changement de fréquences d’échantillonnage, etc.

5.3.1 Sous-échantillonnage par un facteur M


A partir d’un signal à temps discret x(n), on construit un nouveau signal y(m) en prélevant
un échantillon sur M : y(m) = x(mM ). On recherche la relation existant entre la transformée de
Fourier à temps discret X(f ) du signal x(n) et la transformée de Fourier à temps discret Y (f )
du signal y(m). Pour formaliser l’opération de sous-échantillonnage, il est commode de créer un
1
Ordre de grandeur pour les processeurs actuels ?

Nicolas Moreau 1er octobre 2007


page 40/78
Licence de droits d’usage
5.3. UN EXEMPLE APPLICATIF

signal intermédiaire v(n) qui prend les valeurs x(n) si n = mM et qui est égal à zéro sinon. Il a
perdu de l’information contenue dans x(n) mais il reste à la même cadence. On a

+∞
X
v(n) = x(n) λ(n − mM )
m=−∞

En exploitant la formule (4.4), on obtient

M −1
1 X j2π k n
v(n) = x(n) e M .
M
k=0

Donc
M −1 +∞ M −1
1 X X 1 X k
V (f ) = x(n)e−j2π(f −k/M )n = X(f − ).
M n=−∞
M M
k=0 k=0

La transformée de Fourier du signal à temps discret y(m) qui est le signal v(n) débarrassé des
termes nuls, s’écrit
+∞ +∞
X X f f
Y (f ) = y(m)e−j2πf m = v(mM )e−j2π M mM = V ( ).
m=−∞ m=−∞
M

Donc
M −1
1 X f −k
Y (f ) = X( ).
M M
k=0

Conclusion : on observe deux phénomènes, un phénomène de recouvrement des spectres et


un phénomène de dilatation des fréquences comme le montre le dessin de la figure 5.7.

Fig. 5.7 – Phénomène de recouvrement dû au sous-échantillonnage.

Il est possible de supprimer le phénomène de recouvrement des spectres en filtrant la séquence


x(n) par un filtre passe-bas de fréquence de coupure 1/2M mais bien sûr on perd alors de
l’information. On appelle ce filtre, le filtre décimateur.
Interprétation du phénomène de dilatation des fréquences en prenant l’exemple du signal à
temps discret {· · · , 1, 0, −1, 0, 1, · · · }. Toute sa puissance est localisée à la fréquence 1/4. Si on
le sous-échantillonne par un facteur 2, on obtient {· · · , 1, −1, 1, · · · }. Toute sa puissance devient
localisée à la fréquence 1/2.

Remarque
On considère un signal à temps continu x(t) admettant comme transformée de Fourier à temps
continu Xtc (ftc ). On suppose que x(n) est obtenu par échantillonnage de x(t) à la fréquence fe .

Nicolas Moreau 1er octobre 2007


page 41/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE

La transformée de Fourier à temps discret du signal x(t) échantillonné à la fréquence fe = 1/T


en fonction de la fréquence ftc relative à des signaux à temps continu est donnée par
+∞ +∞
1 X X f
−j2π ftc n
Xtd (ftc ) = Xtc (ftc − kfe ) = x(nT )e e .
T n=−∞
k=−∞

Le signal à temps discret y(m) est obtenu par échantillonnage de x(t) à la fréquence fe0 = fe /2.
On a
+∞ +∞ f
1 X fe X −j2π ftc
0 m
Ytd (ftc ) = Xtc (ftc − k ) = y(mT 0 )e e .
2T 2 m=−∞
k=−∞

En décomposant la sommation sur k en indices pairs et indices impairs, on obtient


1 fe
Ytd (ftc ) = [Xtd (ftc ) + Xtd (ftc − )].
2 2
Il reste une ambiguïté au niveau du paramètre f . Dans les formules précédentes, ftc reste une
fréquence exprimée en Hertz. On veut utiliser la fréquence normalisée f = ftc /fe0 . Comme
+∞ ftc +∞ +∞
−j2π m ftc ftc −fe /2
−j2π n −j2π n
X X X
0 fe0
y(mT )e = x(nT )e fe + x(nT )e fe

m=−∞ n=−∞ n=−∞

+∞ ftc +∞ f +∞ ftc −2fe0 /2


X −j2π m X −j2π 2ftc0 n X −j2π n
y(mT 0 )e fe0 = x(nT )e e + x(nT )e 2fe0

m=−∞ n=−∞ n=−∞

On en déduit
1 f f −1
Y (f ) = [X( ) + X( )].
2 2 2

5.3.2 Sur-échantillonnage par un facteur M


A partir de la séquence x(n), on crée la séquence y(m) en intercalant des zéros

y(nM ) = x(n)
y(nM + l) = 0 pour l ∈ {1 · · · M − 1}.

Contrairement au cas précédent, le sur-échantillonnage n’entraîne pas de perte d’information (ni


de création d’ailleurs).

Relation entre les transformées de Fourier


La transformée de Fourier à temps discret de y(m) est égale à
+∞
X +∞
X
−j2πf m
Y (f ) = y(m)e = y(nM )e−j2πf nM = X(M f ).
m=−∞ n=−∞

Il y a création d’images et un phénomène de rétrécissement des fréquences comme le montre la


figure 5.8. Le dessin semble montrer que les spectres sont “identiques”. Ils n’ont pas la même
interprétation. Si x(n) est obtenu par échantillonnage à la fréquence fe d’un signal à temps
0
continu x(t) et si y(m) est obtenu par échantillonnage à la fréquence fe = M fe d’un signal à
temps continu y(t), on remarque que les deux signaux x(t) et y(t) ne sont pas “identiques” comme
le montre le tracé de la figure 5.9. Si l’on veut que x(n) et y(m) aient des spectres “identiques”, il
faut filtrer y(m) par un filtre passe-bas de fréquence de coupure 1/2M appelé filtre interpolateur.

Nicolas Moreau 1er octobre 2007


page 42/78
Licence de droits d’usage
5.3. UN EXEMPLE APPLICATIF

Fig. 5.8 – Spectres avant et après sur-échantillonnage lorsque M = 2.

Fig. 5.9 – Comparaison entre x(t) et y(t).

Filtrage interpolateur
La sortie du filtre interpolateur est donnée par
+∞
X
v(m) = h(m − k)y(k)
k=−∞

+∞
X m−k
v(m) = y(k)sinc( ).
M
k=−∞

En posant m = nM + l et k = pM + q, on obtient
+∞ M −1
X X nM + l − pM − q
v(nM + l) = y(pM + q)sinc( ).
p=−∞ q=0
M

Puisque

y(pM ) = x(p)
y(pM + q) = 0 pour q ∈ {1 · · · M − 1}

on obtient
+∞
X l
v(nM + l) = x(p)sinc(n − p + ). (5.4)
p=−∞
M

On remarque que v(nM ) = x(n) puisque

sinc(n − p) = 0 pour p 6= n.

La formule (5.4) est simplement la formule d’interpolation (3.2)


+∞
X t − pT
x(t) = x(pT )sinc( )
p=−∞
T

évaluée aux instants


l
t = (n + )T.
M

Nicolas Moreau 1er octobre 2007


page 43/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE

5.3.3 Conclusion
Pour passer d’un signal de musique échantillonné à 48 kHz à un signal échantillonné à 32
kHz, il faut d’abord sur-échantillonner par un facteur 2 puis sous-échantillonner par un facteur
3. Le filtre interpolateur qui suit le sur-échantillonneur doit être un filtre demi-bande et le filtre
décimateur qui précède le sous-échantillonneur doit être un filtre tiers de bande. Il suffit donc
d’intercaler entre le sur-échantillonneur et le sous-échantillonneur un unique filtre tiers de bande.

5.3.4 Un autre exemple : convertisseur A/N “sigma-delta 1 bit”


Comment un convertisseur analogique/numérique sur 1 bit peut-il devenir aussi précis qu’un
convertisseur sur 16 bits ? Réponse : en utilisant une fréquence d’échantillonnage bien supérieure à
la fréquence de Nyquist et en réalisant un traitement numérique. Intérêt ? Le coût (prix, fiabilité,
etc) d’un traitement analogique est supérieur au coût d’un traitement numérique. On a donc
intérêt à déplacer au maximum le traitement de l’analogique vers le numérique.
Le schéma de la figure 5.10 donne le principe d’un convertisseur A/N “sigma-delta”. L’entrée
x(n) est un nombre réel de précision infinie (compris entre -1 et +1), la sortie xa (m) est une
version arrondie, le signal intermédiaire (à temps discret et à valeurs discrètes) y(n) ∈ {−1, +1}
est une représentation binaire.

x(n) u(n)
- −m - +m - v(n) - y(n) - w(n)- xa (m)
z −1 H(f ) ?M
-
6 6

Fig. 5.10 – Schéma de principe d’un convertisseur A/N “sigma-delta 1 bit”.

On a

v(n) = u(n − 1) + v(n − 1)


y(n) = 1 si v(n) ≥ 0
= −1 si v(n) < 0
u(n) = x(n) − y(n) (5.5)

On montre, figure 5.11, la “représentation binaire” d’une sinusoïde à temps discret de fréquence
1/50. On montre, figure 5.12, le signal w(n) lorsque H(f ) est un filtre passe-bas de fréquence de
coupure 1/2M avec M = 4 ou 8 (pour rendre le tracé lisible, on représente w(n) comme s’il était
un signal à temps continu).
Pourquoi ? Prématuré. Il faut remplacer dans (5.5) y(n) par y(n) = v(n) + q(n) où q(n) est
une source de bruit.

Nicolas Moreau 1er octobre 2007


page 44/78
Licence de droits d’usage
5.3. UN EXEMPLE APPLICATIF

1 2

0.8
1.5

0.6
1
0.4

0.5
0.2

0 0

−0.2
−0.5

−0.4
−1
−0.6

−1.5
−0.8

−1 −2
10 20 30 40 50 60 10 20 30 40 50 60

Fig. 5.11 – Représentation binaire d’une sinusoïde de fréquence 1/50.

0.8

0.6

0.4

0.2

−0.2

−0.4

−0.6

−0.8

−1
10 20 30 40 50 60

Fig. 5.12 – Versions arrondies (traits pleins) d’une sinusoïde (trait pointillé) de réquence 1/50
lorsque M = 4 ou 8.

Nicolas Moreau 1er octobre 2007


page 45/78
Licence de droits d’usage
CHAPITRE 5. SIGNAUX À TEMPS DISCRET : FILTRAGE

Nicolas Moreau 1er octobre 2007


page 46/78
Licence de droits d’usage
Chapitre 6

Processus aléatoires : une introduction

6.1 Introduction
Nécessité d’un modèle probabiliste. Exemple d’une chaîne de communication. Une source
émet de l’information. Cette information passe dans un canal (sous la forme d’un signal déter-
ministe). Le récepteur observe la sortie du canal et cherche à récupérer l’information émise. Ce
signal est imprévisible pour (au moins) deux raisons :
– le message est inconnu pour le récepteur (sinon quelle est l’utilité de cette transmission),
– le message a subi des perturbations
– de nature déterministe (le canal peut être assimilé à une opération de filtrage),
– de nature aléatoire (addition de bruits).
Un modèle pour le signal reçu : un processus aléatoire à temps continu.
Autre exemple : cf figures 1.1 et 1.2. Problème : comment synthétiser de la parole (par un
filtre AR excité par un bruit blanc), de la musique (quelques sinusoïdes plus du bruit).
Dans ce développement (très sommaire), étude de quelques propriétés d’un processus aléatoire
à temps discret X(n) :
– Que veut dire processus aléatoire à temps discret stationnaire (au 2ème ordre au sens large),
2 , de fonction d’autocovariance R (k), de
centré, gaussien (laplacien,...), de puissance σX X
densité spectrale de puissance SX (f ), ergodique dont une réalisation (trajectoire) est le
signal (observé) x(n) ?
– Quelles sont les propriétés des processus aléatoires à temps discret après une opération de
filtrage ?
Rappel : Variable aléatoire = application mesurable (de Ω → R telle que l’image réciproque
X −1 (B) ∈ A ∀ borélien de R) définie sur un espace probabilisé(Ω, A, P ) où Ω est l’ensemble des
évènements, A un sous-ensemble de Ω possédant une structure particulière et P une mesure de
probabilité.

6.2 Processus aléatoire


6.2.1 Définition
Processus aléatoire à temps discret : famille de v.a. indexée par n ∈ Z. Notation : X(n, ω).
– Interprétation statistique. On fixe un (ou plusieurs) instant d’observation n = n1 , on ob-
tient une (ou plusieurs) variable aléatoire X(n1 , ω) qui se prête bien à un calcul (théorique).
– Interprétation temporelle. On fixe une épreuve particulière ω = ω1 , on obtient une obser-
vation (réalisation, trajectoire) particulière x(n, ω1 ) qui a une signification physique (on
interprète le signal que l’on cherche à traiter comme la réalisation d’un processus aléatoire).

Nicolas Moreau 1er octobre 2007


page 47/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

6.2.2 Interprétation statistique


Statistique du 1er ordre
On fixe un instant particulier n = n1 . On obtient une variable aléatoire X(n1 , ω). Si on
connaît sa fonction de répartition

FX (x; n1 ) = P {X(n1) ≤ x}

ou la densité de probabilité pX (x; n1 ) (si la fonction de répartition est différentiable par rapport
à x), on peut calculer tous les moments d’ordre M
Z +∞
E{X M (n1 )} = xM pX (x; n1 )dx.
−∞

Dans la pratique, on utilise les deux premiers moments, la moyenne

mX (n1 ) = E{X(n1 )}

et la variance (le moment d’ordre 2 centré)


2
σX (n1 ) = E{(X(n1 ) − mX (n1 ))2 }.

Statistique du 2ème ordre


On fixe deux instants particuliers n = n1 et n = n2 . On obtient deux variables aléatoires
X(n1 , ω) et X(n2 , ω). Si on connaît la fonction de répartition conjointe

FX (x1 , x2 ; n1 , n2 ) = P {X(n1 ) ≤ x1 et X(n2 ) ≤ x2 }

et la densité de probabilité conjointe pX (x1 , x2 ; n1 , n2 ), on peut en déduire tous les moments


d’ordre M et M 0
Z +∞ Z +∞
M M0 M0
E{X (n1 )X (n2 )} = xM
1 x2 pX (x1 , x2 ; n1 , n2 )dx1 dx2 .
−∞ −∞

Dans la pratique, on utilise le premier moment conjoint centré (fonction d’autocovariance)

RX (n1 , n2 ) = E{[X(n1 ) − mX (n1 )][X(n2 ) − mX (n2 )]}


Z +∞ Z +∞
RX (n1 , n2 ) = [x1 − mX (n1 )][x2 − mX (n2 )]pX (x1 , x2 ; n1 , n2 )dx1 dx2 .
−∞ −∞

Statistique d’ordre supérieur


Généralisation délicate pour un nombre infini (mais dénombrable pour les p.a. à temps dis-
cret) d’instants d’observation.

Processus stationnaire au sens large (SSL)


Si la moyenne et la variance ne dépendent pas de l’instant d’observation n, si la variance est
bornée et si l’autocovariance ne dépend que de l’écart entre les deux instants d’observation, alors
on dit que le processus X(n) est stationnaire au 2ème ordre au sens large. On supposera cette
propriété vérifiée dans toute la suite.
On notera la moyenne, la variance et la fonction d’autocovariance respectivement

mX = E{X(n)}

Nicolas Moreau 1er octobre 2007


page 48/78
Licence de droits d’usage
6.2. PROCESSUS ALÉATOIRE
2
σX = E{(X(n) − mX )2 }
RX (k) = E{(X(n) − mX )(X(n + k) − mX )}.
La puissance du processus est donnée par

E{X 2 (n)} = σX
2
+ m2X = RX (0) + m2X .

Processus stationnaire gaussien


Chaque variable aléatoire X(n1 ) suit une loi gaussienne, toute combinaison linéaire de X(n1 ) · · · X(nN )
avec n1 · · · nN et N quelconque suit une loi gaussienne. Un processus aléatoire stationnaire gaus-
sien est complètement caractérisé par sa moyenne mX et sa fonction d’autocovariance RX (k) ou
sa matrice de covariance
 
RX (0) RX (1) ··· RX (N − 1)
 RX (1) ··· ··· ··· 
ΓX (N ) =  
 ··· ··· ··· RX (1) 
RX (N − 1) ··· RX (1) RX (0)

pour N “suffisamment grand”. Si x = [x1 · · · xN ]t


1 1 t −1
pX (x) = p e− 2 (x−mX ) ΓX (x−mX )
(2π)N/2 det(Γ)

Pour N = 1, on retrouve la formule classique


1 2 /2σ 2
pX (x) = q e−(x−mX ) X
2
2πσX

6.2.3 Interprétation temporelle


On dispose d’une réalisation (trajectoire)
Moyenne temporelle
+N
1 X
mX (ω) = lim x(n)
N →+∞ 2N + 1
n=−N

Corrélation temporelle
+N
1 X
RX (k, ω) = lim x(n)x(n + k)
N →+∞ 2N + 1
n=−N

L’hypothèse d’ergodicité permet d’affirmer que mX (ω) et RX (k, ω) sont indépendants de ω et


que les moments statistiques sont presque sûrement égaux aux moyennes temporelles.

6.2.4 Densité spectrale de puissance


Propriétés de l’autocovariance
X(n, ω) : p.a. (SSL) à valeurs réelles et centré (pour simplifier)

RX (k) = E{X(n)X(n + k)}

– Expression déterministe
– symétrique RX (−k) = RX (k)

Nicolas Moreau 1er octobre 2007


page 49/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION
2
– bornée |RX (k)| ≤ RX (0) = σX
En effet
E{[X(n) + λX(n + k)]2 } ≥ 0 ∀λ
RX (0) + 2λRX (k) + λ2 RX (0) ≥ 0 ∀λ
2 2
RX (k) − RX (0) ≤ 0 ∀k
– limk→+∞ RX (k) = 0
– La matrice d’autocovariance
   
X(n) RX (0) ··· RX (N − 1)
ΓX (N ) = E{ ···  [X(n) · · · X(n+N −1)]} =  ··· ··· ··· 
X(n + N − 1) RX (N − 1) · · · RX (0)
PN −1
est semi-définie positive puisque E{[ k=0 λk X(n + k)]2 } ≥ 0 ∀λk

Définition de la densité spectrale de puissance


On appelle densité spectrale de puissance (ou spectre de puissance ou spectre tout court) la
transformée de Fourier à temps discret de la fonction d’autocovariance
+∞
X
SX (f ) = RX (k)e−j2πf k .
k=−∞

Existence si RX (k) est une suite de module sommable.

Propriétés de la densité spectrale de puissance


– Fonction périodique (période 1)
Z +1/2
RX (k) = SX (f )e+j2πf k df.
−1/2

– Fonction réelle et paire : SX (f ) = RX (0) + 2 +∞


P
k=1 RX (k) cos(2πf k).
– Fonction positive : théorème de Bochner SX (f ) ≥ 0 ∀f .
– La puissance (totale) est donnée par
Z +1/2
E{X 2 (n)} = RX (0) + m2X = SX (f )df + m2X .
−1/2

– Relation avec la transformée en z de la fonction d’autocovariance


+∞
X
SX (z) = RX (k)z −k
k=−∞

(domaine d’existence : une couronne comprenant le cercle unité). Evaluation de la trans-


formée en z sur le cercle unité.

6.2.5 Exemples de processus aléatoires


Bruit blanc
Processus aléatoire centré vérifiant
2
RX (0) = σX
RX (k) = 0 pour k 6= 0.
La densité spectrale de puissance est donc constante sur tout l’axe des fréquences
2
SX (f ) = σX .

Nicolas Moreau 1er octobre 2007


page 50/78
Licence de droits d’usage
6.2. PROCESSUS ALÉATOIRE

Suite i.i.d. (suite indépendante et identiquement distribuée)


– Chaque élément de la suite a même loi de probabilité à chaque instant

pX (x; n) = pX (x) ∀n.

– ∀n1 , n2 , les v.a. X(n1 ) et X(n2 ) sont mutuellement indépendantes, donc décorrélées (réci-
proque non vraie). C’est donc un bruit blanc.

Sinusoïde à phase aléatoire

X(n, ω) = a cos(2πf1 n + Φ(ω))


avec a et f1 constants et Φ(ω) équirépartie entre 0 et 2π.
Moyenne statistique
Z 2π

mX (n) = E{X(n)} = a cos(2πf1 n + φ) =0 ∀n
0 2π

Moyenne temporelle
+N
1 X
mX (ω) = lim a cos(2πf1 n + Φ(ω)) = 0 ∀ Φ(ω)
N →+∞ 2N + 1
n=−N

Corrélation statistique
Z 2π

RX (n, k) = E{X(n)X(n + k)} = a2 cos(2πf1 n + φ) cos(2πf1 (n + k) + φ)
0 2π

a2 2π
Z Z 2π
dφ dφ
RX (n, k) = [ cos(2πf1 (2n + k) + 2φ) + cos(2πf1 k) ]
2 0 2π 0 2π
a2
RX (n, k) = cos(2πf1 k).
2
Corrélation temporelle
+N
1 X
RX (ω) = lim a2 cos(2πf1 n + Φ(ω)) cos(2πf1 (n + k) + Φ(ω))
N →+∞ 2N + 1
n=−N

+N
a2 1 X
RX (ω) = lim [cos(2πf1 (2n + k) + 2Φ(ω)) + cos(2πf1 k)]
2 N →+∞ 2N + 1
n=−N

a2
RX (ω) = cos(2πf1 k).
2
Conclusion : processus stationnaire et ergodique au 2ème ordre (au sens large) de densité spectrale
de puissance
+∞
a2 X
SX (f ) = [δ(f − f1 − k) + δ(f + f1 − k)].
4
k=−∞

Deux cas extrêmes : une suite i.i.d. est totalement imprévisible (source “sans mémoire”), une
sinusoïde à phase aléatoire est presque totalement prévisible.

Nicolas Moreau 1er octobre 2007


page 51/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

6.2.6 Intérêt d’un modèle probabiliste : un exemple


Quantification = codage de source = discrétisation des amplitudes. Paramètre fondamental
= b = résolution = nombre de bits par échantillon (en moyenne).
Considérons un signal à temps discret x(n) prenant ses valeurs dans l’intervalle [−A, +A].
La démarche la plus naturelle pour définir un quantificateur consiste à
1. partitionner l’intervalle [−A, +A] en L = 2b intervalles distincts de même longueur ∆ =
2A/2b ,
2. numéroter chaque intervalle,
3. définir un représentant par intervalle, par exemple le milieu de l’intervalle.
La procédure d’encodage consiste à décider à quel intervalle appartient x(n) puis à lui associer
le numéro i(n) ∈ {1 · · · L = 2b } correspondant. C’est le numéro de l’intervalle choisi, le symbole
canal, qui sera transmis ou stocké. La procédure de décodage consiste à associer au numéro i(n)
le représentant correspondant x̂(n) = x̂i(n) choisi parmi l’ensemble des représentants {x̂1 · · · x̂L }.
On appelle l’ensemble des représentants un dictionnaire (codebook). Les procédures d’encodage
et de décodage de ce quantificateur sont schématisées figure 6.1.

6
+A
11 o
o
x(n)
o q(n)
10 o
-
n
01 o

00 o
o
-A

Fig. 6.1 – Quantificateur scalaire uniforme.

L’erreur de quantification a pour expression

q(n) = x(n) − x̂(n).

Pour caractériser la dégradation apportée par l’opération de quantification, il faut définir un


critère et proposer un modèle simple pour les signaux intervenant dans ce critère. On suppose
que x(n) est la réalisation d’un processus aléatoire X(n). Comme critère, on choisit le rapport
signal sur bruit
E{X 2 (n)}
Rsb = .
E{[X(n) − x̂(n)]2 }
A priori l’erreur de quantification q(n) n’est pas de nature probabiliste puisque q(n) est une fonc-
tion déterministe de x(n). Pour simplifier cette étude, on prend comme modèle pour représenter
cette erreur de quantification q(n), un processus aléatoire Q(n) avec les hypothèses suivantes.
– Il prend ses valeurs de façon équiprobable dans l’intervalle [−∆/2, +∆/2].
– Q(n), Q(n − 1), ... sont indépendants entre eux.
– Q(n) et X(n) sont indépendants.

Nicolas Moreau 1er octobre 2007


page 52/78
Licence de droits d’usage
6.3. FILTRAGE D’UN PROCESSUS ALÉATOIRE

Le processus aléatoire Q(n) est une suite de variables aléatoires indépendantes et identiquement
distribuées (une suite i.i.d.). La moyenne de l’erreur de quantification est nulle, sa variance est
donnée par
∆/2 2
∆2 A2 −2b
Z 
2 2 1 2 1 2A
σQ = E{Q (n)} = x dx = = = 2 .
−∆/2 ∆ 12 12 2b 3

Si on suppose X(n) uniformément réparti dans l’intervalle [−A, +A], hypothèse irréaliste pour
un signal quelconque mais cela entraîne un calcul simple, sa moyenne est nulle et sa variance a
pour expression Z A
2 2 1 A2
σX = E{X (n)} = x2 dx = .
−A 2A 3
On obtient la relation qui donne la puissance de l’erreur de quantification en fonction de la
puissance du signal et de la résolution b
2
σQ 2
= σX 2−2b .

Le rapport signal sur bruit a pour expression

E{X 2 (n)}
10 log10 = 10 log10 22b = 6, 02 b.
E{Q2 (n)}

Le fait de rajouter un bit revient donc à augmenter le rapport signal sur bruit de 6 dB. Exemple
du CD : entre le seuil d’audition absolu et le seuil de douleur, l’oreille a une dynamique en
puissance de l’ordre de 120 dB. Il aurait donc fallu que le CD réalise une quantification scalaire
sur 20 bits (au lieu de 16) !
Si on suppose X(n) gaussien, on obtiendrait
2
σQ 2
= c σX 2−2b avec 10 log10 (c) ≈ 4.3 dB

6.3 Filtrage d’un processus aléatoire


6.3.1 Problème
Connaissant les propriétés (au 2ème ordre) du processus stationnaire X(n), i.e. sa moyenne
mX , sa variance σX2 , sa fonction d’autocovariance R (k) et/ou sa densité spectrale de puissance
X
SX (f ), quelles sont les propriétés du processus filtré Y (n) par le filtre (stable) de réponse en
fréquence H(f ) ? Le processus Y (n) est-il stationnaire (au 2ème ordre au sens large) ? Si oui,
expression de mY , σY2 , RY (k) et SY (f ).
Problème mathématique délicat pour définir proprement le filtrage d’un processus aléatoire.
On se contentera d’admettre que l’opération de filtrage d’une réalisation
+∞
X
y(n) = h(k)x(n − k)
k=−∞

se généralise
+∞
X
Y (n) = h(k)X(n − k)
k=−∞

ou que l’on se limite à l’étude d’un filtrage d’un p.a. par un filtre RIF (combinaison linéaire
d’un nombre fini de v.a.). On rappelle qu’une condition de stabilité du filtre est que la réponse
impulsionnelle soit de module sommable.

Nicolas Moreau 1er octobre 2007


page 53/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

6.3.2 Formule de filtrage


Moyenne

+∞
X +∞
X
E{Y (n)} = E{ h(k)X(n − k)} = h(k)E{X(n − k)}
k=−∞ k=−∞

+∞
X
E{Y (n)} = mX h(k) = mX H(0).
k=−∞

Cette moyenne ne dépend pas de l’instant d’observation n. Si le processus d’entrée est centré, le
processus de sortie l’est aussi.

Intercovariance entrée-sortie
On supposera X(n) centré par la suite (pour simplifier). Appelons
+∞
X
RXY (n, k) = E{X(n)Y (n + k)} = E{X(n) h(l)X(n + k − l)}
l=−∞

+∞
X +∞
X
RXY (n, k) = h(l)E{X(n)X(n + k − l)} = h(l)RX (k − l)}.
l=−∞ l=−∞

Cette intercovariance ne dépend pas de n, elle ne dépend que de k. Comme

RXY (k) = h(k) ∗ RX (k)

on obtient
SXY (f ) = H(f )SX (f ).

Autocovariance de la sortie

+∞
X
RY (n, k) = E{Y (n)Y (n + k)} = E{Y (n) h(l)X(n + k − l)}
l=−∞

+∞
X +∞
X
RY (n, k) = h(l)E{Y (n)X(n + k − l)} = h(l)RY X (k − l).
l=−∞ l=−∞

L’autocovariance est indépendante de n, elle ne dépend que de k. La propriété de stationnarité


se conserve donc par filtrage. On a

RY (k) = h(k) ∗ RY X (k)

SY (f ) = H(f )SY X (f )
Comme
RY X (k) = E{Y (n)X(n + k)} = E{X(n)Y (n − k)} = RXY (−k)
et que
+∞
X +∞
X
SY X (f ) = RXY (−k)e−j2πf k = RXY (k)e−j2π(−f )k = SXY (−f )
k=−∞ k=−∞

Nicolas Moreau 1er octobre 2007


page 54/78
Licence de droits d’usage
6.3. FILTRAGE D’UN PROCESSUS ALÉATOIRE

on obtient
SY (f ) = H(f )SXY (−f ) = H(f )H(−f )SX (−f )

finalement
SY (f ) = H(f )H(−f )SX (f ) = |H(f )|2 SX (f )

puisque SX (f ) est une fonction paire et que H(f ) a la propriété de symétrie hermitienne.

6.3.3 Remarques et interprétation


Généralisations

Généralisation à un processus à temps continu

SY (f ) = T F T C{E{X(t)X(t + τ )}} = |H(f )|2 SX (f ).

Utilisation de la transformée en z

SY (z) = H(z)H(z −1 )SX (z).

Processus gaussien

Le caractère gaussien se maintient par filtrage

Densité spectrale de puissance

Considérons un filtre passe-bande (idéal) autour d’une fréquence f1

H(f ) = 1 si |f | ∈ [f1 − ∆/2, f1 + ∆/2]


= 0 sinon.

Comme

SY (f ) = SX (f ) si |f | ∈ [f1 − ∆/2, f1 + ∆/2]


= 0 sinon,

on obtient
Z +1/2
SY (f )df ≈ 2SX (f1 )∆
−1/2

pour ∆ suffisamment petit. Comme le premier membre est nécessairement positif (ou nul), on
remarque que SX (f ) ≥ 0 ∀f . SX (f ) s’interprète comme une densité spectrale de puissance.

Processus ARMA

Si le processus d’entrée X(n) est un bruit blanc centré de variance σX2 , la sortie Y (n) est un
2 2
processus de densité spectrale de puissance SY (f ) = |H(f )| σX . Par exemple, si le filtre est un
filtre AR(4) avec deux fois deux pôles complexes conjugués de module légèrement inférieur à 1
et d’argument +/ − π/8 et +/ − π/4, on obtiendra en sortie un processus aléatoire ayant de la
puissance surtout concentrée dans deux bandes de fréquence voisines de f = 1/8 et f = 1/4. On
voit qu’il est facile de créer un signal synthétique ayant des composantes spectrales particulières.

Nicolas Moreau 1er octobre 2007


page 55/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

6.3.4 Exemple sous forme d’exo


Filtre défini par la relation de récurrence

y(n) − ay(n − 1) = bx(n).

Processus d’entrée X(n) : processus stationnaire centré de fonction d’autocovariance


2 2
{RX (0), RX (1), RX (2), · · · } = {σX , cσX , 0, · · · }.

Caractéristiques du filtre

b
H(z) = si |z| > |a|
1 − az −1
b
H(f ) = si |a| < 1
1 − ae−j2πf
b2
|H(f )|2 =
1 − 2a cos(2πf ) + a2
{h(0), h(1), h(2), · · · } = {b, ba, ba2 , · · · }

Caractéristiques du processus d’entrée

2
SX (z) = σX (cz + 1 + cz −1 )
2
SX (f ) = σX (1 + 2c cos(2πf ))

Caractéristiques du processus de sortie

cz + 1 + cz −1
SY (z) = b2 σX
2
(1 − az)(1 − az −1 )
1 + 2c cos(2πf )
SY (f ) = b2 σX
2
1 − 2a cos(2πf ) + a2
Souvent, on désire connaître la puissance du processus de sortie
Z 1/2
2 2
σY = E{Y (n)} = RY (0) = SY (f )df
−1/2

Calcul plus facile en utilisant l’expression en z


Z
1 dz
σY2= RY (0) = SY (z)
2πj Γ z

cz 1 + 1 + cz −1 dz
Z
2 1
σY2 = b2 σX
2πj Γ (1 − az)(1 − az −1 ) z
cz 2 + z + c
Z
2 2 2 1
σY = b σX dz
2πj Γ z(1 − az)(z − a)
et en appliquant le théorème des résidus (deux pôles à l’intérieur du cercle unité z=0 et z=a)

c ca2 + a + c 2 1 + 2ac
σY2 = b2 σX
2
[ + ] = b2 σ X
(−a) a(1 − a2 ) 1 − a2

Nicolas Moreau 1er octobre 2007


page 56/78
Licence de droits d’usage
6.4. THÉORIE DE L’ESTIMATION : UNE INTRODUCTION

2ème méthode : décomposition en éléments simples puis développement en série de Laurent de

cz + 1 + cz −1
SY (z) = b2 σX
2
(1 − az)(1 − az −1 )

3ème méthode
X X
σY2 = E{ h(k)X(n − k) h(l)X(n − l)}
k l

X
σY2 = 2
h(k)[ch(k − 1) + h(k) + ch(k + 1)]σX
k

6.4 Théorie de l’estimation : une introduction


X(n) p.a. stationnaire et ergodique veut dire

+N
1 X
mX = E{X(n)} = lim x(n)
N →∞ 2N + 1
n=−N

+N
1 X
RX (k) = E{X(n)X(n + k)} = lim x(n)x(n + k).
N →∞ 2N + 1
n=−N

En pratique, on ne dispose qu’une seule réalisation comportant un nombre fini N d’échantillons


[x(0) · · · x(N − 1)]. Comment utiliser “au mieux” ces données ?

6.4.1 Moyenne “empirique”


On pose
N −1
1 X
m
bX = x(n).
N
n=0

m
b X est une v.a. puisqu’elle dépend d’un tirage aléatoire. Elle a une moyenne et une variance

N −1 N −1
1 X 1 X
b X } = E{
E{m X(n)} = mX = mX
N N
n=0 n=0

On dit que l’estimateur est “sans biais”.

N −1 N −1
1 X 2 1 X
b X } = E{[(
var{m X(n)) − mX ] } = 2 E{[ (X(n) − mX )]2 }
N N
n=0 n=0

N −1 N −1 N −1
1 X X 1 X |n|
b X} = 2
var{m RX (n − m) = [1 − ]RX (n).
N N N
n=0 m=0 n=0

Donc var{mb X} =
6 0 mais on montre que limN →∞ var{m
b X } → 0. On dit que cet estimateur tend
“en moyenne quadratique” vers mX .

Nicolas Moreau 1er octobre 2007


page 57/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

6.4.2 Covariance “empirique”


On supposera X(n) centré pour simplifier les notations. On pose

−1−k
NX
bX (k) = 1
R x(n)x(n + k) pour k = 0 · · · N − 1.
N
n=0

Comme
−1−k
NX
bX (k)} = 1 N −k
E{R E{X(n)X(n + k)} = RX (k)
N N
n=0

cet estimateur est biaisé mais “asymptotiquement sans biais”. Nécessité de calculer aussi var{R
bX (k)}
mais calcul compliqué.
Remarque : on aurait pu choisir

−1−k
NX
1
R
bX (k) = x(n)x(n + k)
N −k
n=0

ce qui aurait rendu l’estimateur sans biais mais on montre que la variance est alors beaucoup
plus importante. Globalement l’estimateur biaisé est bien préférable (on lui trouvera par la suite
de nouvelles qualités). Le fait de pondérer la sommation qui comporte en fonction de k de moins
en moins de termes par toujours 1/N permet de tenir de moins en moins compte de cette somme
qui est de moins en moins fiable !

6.4.3 Un estimateur spectral : le périodogramme


La détermination de la densité spectrale de puissance

+∞
X ∞
X
SX (f ) = RX (k) exp(−j2πf k) = RX (0) + 2 RX (k) cos(2πf k)
k=−∞ k=1

réclame la connaissance de RX (k) pour k = 0 · · · ∞. Ne disposant que de N données observées


[x(0) · · · x(N − 1)] on ne peut disposer qu’au plus N termes de la fonction d’autocovariance (plus
ou moins fiables).
En remplaçant dans l’expression précédente RX (k) par R bX (k), on obtient

N −1 −1−k
NX
X 1
SbX (f ) = [ x(n)x(n + |k|)] exp(−j2πf k)
N
k=−N +1 n=0

ce qui donne
N −1
1 X
SbX (f ) = | x(n) exp(−j2πf n)|2 .
N
n=0

En effet

N
X −1 N
X −1 N
X −1
2
| x(n) exp(−j2πf n)| = x(n) exp(−j2πf n) x(m) exp(+j2πf m)
n=0 n=0 m=0

Nicolas Moreau 1er octobre 2007


page 58/78
Licence de droits d’usage
6.4. THÉORIE DE L’ESTIMATION : UNE INTRODUCTION

ou en posant W = exp(−j2πf ) et en se limitant à N = 3


N
X −1
| x(n) exp(−j2πf n)|2 = x(0)[x(0)W 0 + x(1)W 1 + x(2)W 2 ]
n=0
+ x(1)[x(0)W −1 + x(1)W 0 + x(2)W 1 ]
+ x(2)[x(0)W −2 + x(1)W −1 + x(2)W 0 ]
= [x(2)x(0)]W −2 + [x(2)x(1) + x(1)x(0)]W −1
+ [x2 (0) + x1 (1) + x2 (2)]W 0
+ [x(0)x(1) + x(1)x(2)]W 1 + [x(0)x(2)]W 2
Si on évalue cette expression aux fréquences f = fk = k/L, on obtient la formule du “périodo-
gramme”
N −1
k 1 X k
SX ( ) = |
b x(n) exp(−j2π n)|2 pour k = 0 · · · L − 1.
L N L
n=0
Les représentations fréquentielles des signaux de parole et de musique des figures 1.1 et 1.2 ont
été calculées de cette façon.
Compléments
– Si L = N , le périodogramme est donné directement par le module au carré de la transformée
de Fourier discrète. Si L > N , on effectue également une TFD mais avec “zero padding”.
– Emploi d’une fenêtre de pondération dans le domaine temporel :
N −1
k 1 X k
SbX ( ) = | v(n)x(n) exp(−j2π n)|2
L N L
n=0
ou dans le domaine fréquentiel :
N −1
k X
bX (n) exp(−j2π k n).
SbX ( ) = W (n)R
L L
n=−N +1

– La densité spectrale de puissance SbX (f ) reste, pour toutes les valeurs de f possibles, une
variable aléatoire dont il peut être intéressant de déterminer la moyenne et la variance.
Dans le cas d’une pondération par une fenêtre dans le domaine temporel, la moyenne est
donnée par
+∞ +∞
X 1 X
E{SbX (f )} = [ v(n)v(n − k)E{X(n)X(n − k)}] exp(−j2πf k)
N n=−∞
k=−∞
+∞
X
E{SbX (f )} = q(k)RX (k) exp(−j2πf k)
k=−∞
avec
+∞
1 X
q(k) = v(n)v(n − k).
N n=−∞
On a donc Z 1/2
E{SbX (f )} = Q(f ) ∗ SX (f ) = Q(λ)SX (f − λ)dλ.
−1/2
Il existe donc un biais qui disparait si N → ∞. Le périodogramme est un estimateur
“asymptotiquement sans biais”.
Il faudrait donner l’expression de la variance E{[SbX (fk ) − E{SbX (fk )}]2 }. On peut montrer
que cette variance ne tend pas vers 0 lorsque N → ∞. On dit que cet estimateur n’est pas
“consistant”. Pour le rendre consistant, il faut pondérer le signal observé par une fenêtre ...

Nicolas Moreau 1er octobre 2007


page 59/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

6.5 Sinusoïde bruitée sous forme d’exo


6.5.1 Problème
On mesure un signal y(0) · · · y(N − 1) à l’extrémité d’un canal en sachant a priori que le
signal émis à l’entrée est une sinusoïde pure de la forme x(n) = a cos(2πf1 n + φ). On ignore la
valeurs numériques prises par a, f1 et φ. De plus, le signal a été perturbé par son passage dans le
canal. On désire calculer les valeurs numériques prises par a et f1 à partir des valeurs observées.
Problème d’“estimation statistique”.

6.5.2 Modélisation
On modélise le signal observé comme étant la réalisation d’un processus aléatoire

Y (n) = a cos(2πf1 n + Φ(ω)) + B(n)

où Φ(ω) est une variable aléatoire équirépartie entre 0 et 2π et où B(n) est un bruit blanc centré
de variance σ 2 non-corrélé avec X(n) = a cos(2πf1 n + Φ(ω)).
On a vu que X(n) est un processus aléatoire stationnaire au 2ème ordre (et ergodique) de
fonction d’autocovariance
a2
RX (k) = cos(2πf1 k).
2
On en déduit que Y (n) est un processus aléatoire stationnaire centré de fonction d’autocovariance

E{Y (n)Y (n + k)} = E{X(n)X(n + k)} + E{B(n)B(n + k)}

a2
cos(2πf1 k) + σ 2 δ(k).
RY (k) =
2
Etudions quelques propriétés de la matrice d’autocovariance
 
RY (0) · · · RY (N − 1)
ΓY =  ··· ··· ··· 
RY (N − 1) · · · RY (0)

On obtient
ΓY = ΓX + σ 2 I.
La matrice ΓX est de rang non complet. En effet, on a

cos(2πf1 k) + cos(2πf1 (k − 2)) = 2 cos(2πf1 (k − 1)) cos(2πf1 )

ce qui entraîne
RX (k) + α1 RX (k − 1) + RX (k − 2) = 0
donc     
RX (0) RX (1) RX (2) 1 0
 RX (1) RX (0) RX (1)   α1  =  0 
RX (2) RX (1) RX (0) 1 0
Seuls deux vecteurs colonne de ΓX sont linéairement indépendants ce qui montre que le rang de
ΓX est égal à = 2.
Appelons λ1 et λ2 les deux valeurs propres non-nuls (et positives puisque ΓX est semi-défini
positif) de ΓX et v 1 et v 2 les deux vecteurs propres correspondants. Appelons w3 · · · wN les
vecteurs engendrant l’espace nul. Comme

ΓY v = ΓX v + σ 2 v = (λ + σ 2 )v

Nicolas Moreau 1er octobre 2007


page 60/78
Licence de droits d’usage
6.5. SINUSOÏDE BRUITÉE SOUS FORME D’EXO

ΓY w = ΓX w + σ 2 w = σ 2 w
on en déduit que les valeurs propres de ΓY sont égales à

λ1 + σ 2 , λ 2 + σ 2 , σ 2 , · · · , σ 2

et que le vecteur propre associé à la plus petite valeur propre doit être de la forme [1, α1 , 1]t .

6.5.3 Algorithme
L’étude précédente donne la démarche. Il faut, dans la pratique, se limiter à la connaissance
d’un nombre fini de valeurs observées. On en déduit l’algorithme :
1. A partir de y(0) · · · y(N − 1), calculer
−1−k
NX
1
R̂X (k) = y(n)y(n + k) pour k = 0, 1, 2.
N
n=0

2. Construire la matrice Γ̂X de dimension 3.


3. Calculer le vecteur propre associé à la plus petite valeur propre. Il doit être de la forme
[1, α1 , 1]t .
4. En déduire
1 −α1
f1 = Arccos( ).
2π 2

6.5.4 Autre méthode


1. A partir de y(0) · · · y(N − 1), calculer
N
X −1
Y (k) = y(n)e−j2πnk/N pour k = 0, · · · , N − 1.
n=0

2. En déduire
arg max |Y (k)|
f1 = .
N
Comparaison avec la méthode précédente ?

Nicolas Moreau 1er octobre 2007


page 61/78
Licence de droits d’usage
CHAPITRE 6. PROCESSUS ALÉATOIRES : UNE INTRODUCTION

Nicolas Moreau 1er octobre 2007


page 62/78
Licence de droits d’usage
Annexe A

Codeur de parole “LPC10”

A.1 Introduction
– Parole en bande téléphonique : fe = 8 kHz ⇒ [300 Hz - 3.3 kH].
– Débit de référence ≈ 13 bits × 8 = 104 kbit/s. Nombreux codeurs existants. But : diminuer
le débit sans trop dégrader la qualité. Qualité acceptable de 64 kbit/s (RNIS) jusqu’à 6.3
kbit/s (partie son du visiophone). Pour un débit plus faible (débit visé = 2.4 kbit/s)
⇒ faible qualité. Codeur “LPC10” complètement démodé maintenant mais grand intérêt
pédagogique.
– Schéma très général d’un codeur de parole : schéma de la figure A.1. x(n) est le signal
de parole original, y(n) le signal en sortie du filtre d’“analyse”, ŷ(n) l’entrée du filtre de
“synthèse” et x̂(n) le signal de parole reconstruit avec

A(z) = 1 + a1 z −1 + · · · aP z −P .

x(n) y(n) Flux binaire y(n) 1 x(n)


A(z) −1 A(z)
Q Q
LPC a a

Emetteur Canal Récepteur

Fig. A.1 – Schéma très général d’un codeur de parole.

– Traitement en deux étapes : calculer les coefficients du filtre A(z) à partir du signal ori-
ginal puis déterminer l’entrée du filtre de synthèse de telle sorte que la qualité du signal
reconstruit soit la meilleure possible tout en respectant la contrainte de débit (2.4 kbit/s).
– Actualisation des coefficients du filtre et détermination de l’entrée du filtre de synthèse
toutes les vingtaine de ms en supposant le signal “localement stationnaire” dans chacune
de ces fenêtres. x(n) est de la forme x(n) = x(mN + l) mais on notera par la suite

x = [x(0) · · · x(N − 1)]t

les N = 160 échantillons (pour du signal de parole échantillonné à fe = 8 kHz) dans chacune
de ces fenêtres.

Nicolas Moreau 1er octobre 2007


page 63/78
Licence de droits d’usage
ANNEXE A. CODEUR DE PAROLE “LPC10”

A.2 Détermination des coefficients du filtre


PP
On appelle v(n) = − i=1 ai x(n − i) la valeur prédite à l’instant n et

P
X
y(n) = x(n) − v(n) = x(n) + ai x(n − i)
i=1

l’erreur de prédiction.

A.2.1 Minimisation au sens des moindres carrés


Problème : connaissant x = [x(0) · · · x(N − 1)]t déterminer a = [a1 · · · aP ]t minimisant la
puissance (empirique) de l’erreur de prédiction

aopt = arg min σ̂Y2


a

avec
N −1
1 X 2 1
σ̂Y2 = y (n) = y t y.
N N
n=0

Comme (en choisissant P = 2 pour simplifier l’écriture)


     
y(0) x(0) x(−1) x(−2)
y(1) x(1) x(0) x(−1)  
 a1
    
= +
    
.. .. .. ..
 a2
 
 .   .   . .
y(N − 1) x(N − 1) x(N − 2) x(N − 3)

soit
y = x + Γa
on peut écrire

1 t 1
σ̂Y2 = (x + (Γa)t )(x + Γa) = (xt x + 2(Γt x)t a + at Γt Γa).
N N

Le vecteur aopt est celui qui vérifie


∂ 2
σ̂ = 0
∂a Y
Γt x + Γt Γa = 0

aopt = −[Γt Γ]−1 Γt x


si Γt Γ est inversible (problème délicat).

A.2.2 Approche “propre”


On suppose que le signal x(n) peut être interprété (modélisé) comme la réalisation d’un p.a.
stationnaire de fonction d’autocovariance RX (k) = E{X(n)X(n − k)}. On recherche le vecteur
a minimisant la puissance de l’erreur de prédiction
P
X
σY2 2
= E{Y (n)} = E{[X(n) + ai X(n − i)]2 }
i=1

Nicolas Moreau 1er octobre 2007


page 64/78
Licence de droits d’usage
A.2. DÉTERMINATION DES COEFFICIENTS DU FILTRE
    
RX (1) Rx (0) ··· RX (P − 1) a1
2 2
σY = σX + 2[a1 · · · aP ]  .
. .
. .. .. .. 
 + [a1 · · · aP ] 
   
. . . .  . 
RX (P ) RX (P − 1) · · · RX (0) aP
σY2 = σX
2
+ 2rt a + at Ra
La minimisation de σY2 relativement à a entraîne

r + Raopt = 0

aopt = −R−1 r.
On a aussi
(σY2 )min = σX
2
+ 2at r − at r = σX
2
+ at r.

A.2.3 Remarques
– Les deux solutions sont très comparables ce qui n’est pas étonnant puisque N1 Γt x est une
estimée du vecteur r et N1 Γt Γ est une estimée de R. La différence (essentielle mais subtile)
est que la matrice de covariance exacte étant définie positive (sauf dans le cas limite où le
processus est harmonique), cette matrice est toujours inversible alors que la matrice

x2 (1) + · · · x2 (N − 2)
 
t x(0)x(1) + · · · + x(N − 3)x(N − 2)
ΓΓ =
x(0)x(1) + · · · + x(N − 3)x(N − 2) x2 (0) + · · · x2 (N − 3)
reste symétrique mais elle n’est plus forcément définie positive. On peut montrer aussi que
cette propriété de positivité entraîne que tous les zéros du polynôme A(z) sont à l’intérieur
du cercle unité ce qui assure la stabilité du filtre de synthèse (propriété très importante
dans la pratique).
– Dans la pratique où on ne dispose que de N données observées, on désire maintenir cette
propriété de positivité. On peut observer que remplacer une moyenne statistique RX (k) =
E{X(n)X(n − k)} par une moyenne temporelle
N −1
1 X
R̂X (k) = x(n)x(n − k)
N
n=k

en s’appuyant sur l’hypothèse d’ergodicité maintient l’inversion toujours réalisable et la


propriété de stabilité du filtre de synthèse.
– On peut montrer que l’erreur de prédiction Y (n) est blanche (plus exactement que le signal
X(n) a été blanchi). En effet, on remarque que

∂E{|Y (n)|2 }
=0 ⇒ E{Y (n)X(n − i)} = 0 ∀i = 1 · · · P.
∂ai
Supposons P grand. Comme Y (n) est non corrélé avec tous les X(n − i) précédents et que
Y (n−i) est une combinaison linéaire de ces X(n−i), on en déduit que Y (n) est non corrélé
avec Y (n − i). L’erreur de prédiction Y (n) est donc un bruit blanc mais cette propriété
n’est vérifiée a priori que si P → ∞ (comportement asymptotique). On appelle le filtre
donnant Y (n) à partir de X(n) le “filtre blanchissant”.
– Si Y (n) a été totalement blanchi, on peut écrire

SY (f ) = |A(f )|2 SX (f ) = σY2 .

On a donc
σY2
SX (f ) = .
|A(f )|2

Nicolas Moreau 1er octobre 2007


page 65/78
Licence de droits d’usage
ANNEXE A. CODEUR DE PAROLE “LPC10”

On rappelle que l’estimateur spectral le plus standard est le périodogramme : on applique


la TFD aux N données observées [x(0) · · · x(N − 1)] puis on élève au carré les modules.
La formule précédente suggère un deuxième estimateur spectral : à partir des N données
observées, on estime les P premiers coefficients de la fonction d’autocovariance, on résout
les équations normales puis on exploite la formule précédente.
On observera sur les tracés de droite de la figure A.2 le périodogramme d’un signal de parole
en bleu (un spectre de raies avec une fréquence fondamentale de 250 Hz correspondant à un
locuteur féminin) et le spectre LPC en noir qui ne donne que “l’enveloppe spectrale”. Ces
deux estimateurs ne sont pas équivalents. Le choix entre les deux dépend essentiellement
de l’application.

A.3 Détermination de l’entrée du filtre de synthèse


A.3.1 Cas des sons non voisés
Supposons que l’on choisisse comme entrée du filtre de synthèse une réalisation quelconque
d’un bruit blanc1 de puissance σŶ2 = σY2 . On voit alors que l’on a la propriété suivante au niveau
des densités spectrales de puissance

σY2
SX̂ (f ) = = SX (f ).
|A(f )|2

On reconstruit un signal qui a la même distribution de la puissance en fonction de la fréquence


mais les formes d’onde sont différentes. Comme l’oreille est relativement insensible à des modifi-
cations de la phase, on reconstruit par ce mécanisme un signal qui est perçu approximativement
identique au signal original. On appelle ce type de codeur un “vocodeur”.

A.3.2 Cas des sons voisés


Les tracés de la figure A.2 sont relatifs à un son voisé où l’on montre aussi bien dans le
domaine temporel (à gauche) que dans le domaine fréquentiel (à droite) le signal original x(n)
et l’erreur de prédiction y(n). Le filtre A(z) n’est manifestement pas totalement blanchissant.
Il reste dans le signal y(n) une périodicité assez marquée, visible aussi bien dans le domaine
temporel que dans le domaine fréquentiel. Pendant une durée de 32 ms, on a approximativement
7.5 périodes. La fréquence fondamentale est donc de l’ordre de f0 ≈ 7.5/0.032 ≈ 250 Hz. On
observe bien dans le domaine fréquentiel un spectre de raies avec une fréquence fondamentale de
250 Hz (correspondant à un locuteur féminin) et les différents harmoniques.
Le problème qui se pose maintenant est de trouver un modèle ŷ(n) pour y(n) qui permette
par filtrage d’obtenir SX̂ (f ) ≈ SX (f ) et qui soit très économe en débit. Un peigne de la forme

+∞
X
ŷ(n) = α λ(n − mT0 + ϕ)
m=−∞

est un bon candidat. Dans cette expression, λ(n) est le symbole de Kronecker qui vaut 1 si n
= 0, 0 sinon, T0 = fe /f0 est la période fondamentale exprimée en nombre d’échantillons et ϕ
une valeur appartenant à {0, · · · , T0 − 1} traduisant notre incertitude sur la phase. Le signal
ŷ(n) peut être alors interprété comme la réalisation d’un processus aléatoire Ŷ (n) dont il est
intéressant de déterminer les propriétés.
1
en exploitant la fonction matlab randn par exemple

Nicolas Moreau 1er octobre 2007


page 66/78
Licence de droits d’usage
A.3. DÉTERMINATION DE L’ENTRÉE DU FILTRE DE SYNTHÈSE
0.4 10

0.3 0

0.2
−10

0.1

Puissances en dB
−20

0
−30

−0.1

−40
−0.2

−50
−0.3

−60
−0.4 0 0.5 1 1.5 2 2.5 3 3.5 4
50 100 150 200 250 300 Fréquences en kHz

Fig. A.2 – Cas de sons voisés.

On remarque que la moyenne a pour expression


0 −1
TX +∞
1 X α
E{Ŷ (n)} = α λ(n − mT0 + ϕ) =
T0 m=−∞ T0
ϕ=0

et que la fonction d’autocorrélation2 qui est donnée par


+∞
X +∞
X
2
RŶ (k, n) = E{Ŷ (n)Ŷ (n − k)} = α E{ λ(n − mT0 + ϕ) λ(n − k − lT0 + ϕ)}
m=−∞ l=−∞

vaut
+∞
X α2
RŶ (k, n) = α2 E{ λ(n − mT0 + ϕ)} =
m=−∞
T0
si k est un multiple de T0 , 0 sinon. La moyenne et la fonction d’autocorrélation étant indépen-
dantes de l’instant d’observation n, le processus Ŷ (n) est stationnaire de moyenne 1/T0 et de
fonction d’autocorrélation
+∞
α2 X
RŶ (k) = λ(k − mT0 ).
T0 m=−∞

La densité spectrale de puissance SŶ (f ) est la transformée de Fourier à temps discret de RŶ (k).
On a donc
+∞ +∞ +∞
α2 X X α2 X
SŶ (f ) = λ(k − mT0 ) exp(−j2πf k) = exp(−j2πf mT0 )
T0 T0 m=−∞
k=−∞ m=−∞

ou de façon équivalente
+∞
X m
SŶ (f ) = α2 δ(f − )
m=−∞
T0
en se rappelant qu’un train d’impulsions de dirac admet une sorte de développement en série de
Fourier et que l’on a donc la formule générale
+∞ +∞
X 1 X j2π n t
δ(t − nT ) = e T .
n=−∞
T n=−∞
2
Son expression est plus simple que celle de la fonction d’autocovariance.

Nicolas Moreau 1er octobre 2007


page 67/78
Licence de droits d’usage
ANNEXE A. CODEUR DE PAROLE “LPC10”

Comme le tracé de droite de la figure A.2 a été obtenu à partir d’un nombre fini de données
observées, on “voit” sur ce tracé le lobe principal de la transformée de Fourier à temps discret
de la fenêtre de pondération. Dans cette simulation, on a pris une fenêtre de Hamming de durée
2N . Le lobe principal a pour largeur 4fe /2N ce qui correspond à 100 Hz, valeur cohérente avec
celle observée sur le tracé.
Le filtrage de ŷ(n) par le filtre dont le module au carré de la réponse en fréquence est
représenté en noir sur la figure A.2 donnera bien un signal dont la densité spectrale de puissance
se superposera avec le tracé en bleu de cette figure.

A.3.3 Détermination sons voisés/sons non voisés


Le codeur LPC10 doit déterminer s’il s’agit d’un son voisé ou non voisé dans la fenêtre
d’analyse courante. Cette détermination se fait en réalisant au préalable une estimation de la
fonction d’autocorrélation normalisée de l’erreur de prédiction y(n) pour de nombreuses valeurs
de k. Si cette fonction décroît vers 0 rapidement (comparaison à un seuil), on en déduit que le
son est non voisé. Dans le cas d’un signal voisé, cette fonction est presque périodique de période
T0 ce qui fournit un estimateur de la période fondamentale. Dans la pratique, il faut prendre
quelques précautions car un algorithme rudimentaire a tendance à déterminer comme fréquence
fondamentale, le 2ème ou même le 3ème harmonique ... Dans la pratique, il faudra assurer des
transitions pas trop brutales lorsque l’on passera d’une fenêtre d’analyse à la suivante. Il faudra,
en particulier, faire en sorte que le choix de ϕ (information non codée) entraîne une transition
douce entre les peignes successifs.

A.4 Conclusion
Toutes les 20 ms, on transmet les paramètres suivants.
– Les coefficients du filtre A(z). Dans le codeur LPC10, P = 10. En prenant comme ordre
de grandeur 3 ou 4 bits pour coder un coefficient, on obtient un débit de l’ordre de 1.5 à
1.8 kbit/s.
– La puissance σY2 dans le cas d’un son non voisé ou de façon équivalente le coefficient α
dans le cas d’un son voisé ce qui coûte 50 × 6 ≈ 300 bit/s (pour couvrir 50 dB par pas de
1 dB).
– La distinction voisé/non voisé soit 50 bits/s.
– Enfin la période fondamentale soit 50×log2 (T0max −T0min ) = 350 bits/s si on prend T0min =
20 (f0max = 400 Hz correspondant à une voix d’enfant très aigue) et T0max = 147 (f0min = 54
Hz correspondant à une voix d’homme très grave).
En sommant toutes ces contributions, on obtient un débit global de l’ordre de 2.4 kbit/s.

Nicolas Moreau 1er octobre 2007


page 68/78
Licence de droits d’usage
Annexe B

Transformée de Fourier à court terme

B.1 Introduction
Les signaux que l’on est amené à traiter dans la pratique (du signal de parole, des signaux
biomédicaux, des signaux radar ...) ne sont ni périodiques, ni stationnaires. Les caractéristiques
spectrales du signal à analyser évoluent avec le temps. On cherche une représentation temps-
fréquence adaptée.

B.2 Définitions
On appelle transformée de Fourier à temps discret à court terme la fonction de deux variables
+∞
X
X(n, f ) = h(n − l)x(l)e−j2πf l
l=−∞

où n est un entier relatif caractérisant le temps, f un réel représentant la fréquence et h(n)


une fenêtre dite d’analyse dont les valeurs seront supposées nulles à l’extérieur de l’intervalle
[0 · · · N − 1]. Comme la transformée de Fourier est appliquée à un signal de durée finie, il suffit
d’évaluer cette expression pour N fréquences multiples de 1/N . On obtient la transformée de
Fourier discrète à court terme
+∞
k X k
Xk (n) = X(n, f = ) = h(n − l)x(l)e−j2π N l . (B.1)
N
l=−∞

On ne considérera par la suite que ce cas.

B.3 Interprétation par bancs de filtres


Si dans l’expression précédente (B.1) on fixe la fréquence, c’est-à-dire l’indice k, Xk (n) s’in-
terprète comme un signal. Il est de la forme
+∞
X k
Xk (n) = h(n − l)yk (l) avec yk (l) = x(l)e−j2π N l .
l=−∞

Il s’agit d’une opération de modulation du signal x(n) par l’exponentielle complexe à la fréquence
(−k)/N suivie d’une opération de convolution. La transformée de Fourier discrète à court terme
s’interprète comme une translation vers la gauche du spectre du signal x(n) suivie par un filtrage
caractérisé par la réponse en fréquence H(f ), transformée de Fourier à temps discret de h(n).

Nicolas Moreau 1er octobre 2007


page 69/78
Licence de droits d’usage
ANNEXE B. TRANSFORMÉE DE FOURIER À COURT TERME

Lorsque l’on utilise une fenêtre rectangulaire, il s’agit d’un filtrage passe-bas comme le montre le
tracé de la figure 4.1. Il en est de même pour tout autre type de fenêtre (Hamming, Kaiser ...).
Ce filtrage passe-bas est d’autant plus sélectif que la fenêtre d’analyse est longue. La résolution
fréquentielle est donc directement fonction de la longueur de l’observation. Le choix de N est le
résultat d’un compromis car on désire généralement avoir aussi des fenêtres courtes pour pouvoir
analyser des phénomènes courts.
Si l’on met en parallèle N filtres de même réponse en fréquences H(f ), on obtient un banc
de filtres appelé banc de filtres d’analyse. Il extrait du signal sa contribution fréquentielle dans
toutes les sous-bandes centrées aux fréquences multiples de 1/N .
Il existe une deuxième interprétation en écrivant Xk (n) sous la forme
+∞
k X k
Xk (n) = e−j2π N n [h(n − l)ej2π N (n−l) ]x(l)
l=−∞

+∞
k X k
−j2π N n
Xk (n) = e [h(l)ej2π N l ]x(n − l).
l=−∞

La transformée de Fourier discrète à court terme s’interprète maintenant comme un filtrage passe-
bande de x(n) par un filtre de réponse en fréquence H(f − k/N ) suivi par une translation du
spectre du signal. Les deux interprétations sont équivalentes. Le banc de filtres est dit uniforme.

B.4 Fenêtre glissante et reconstruction


Connaissant Xk (n) pour k = 0 · · · N − 1, à l’instant n, il est possible de calculer h(n − l)x(l)
pour toutes valeurs de l ∈ [n − N + 1 · · · n] pourvu que h(n − l) soit différent de zéro. On obtient

N −1
1 1 X k
x(l) = [ Xk (n)ej2π N l ].
h(n − l) N
k=0

Cette expression montre qu’il n’est pas nécessaire de calculer Xk (n) à tous les instants n. La
connaissance de Xk (n) pour n multiple de N semble suffire a priori pour reconstruire le signal.
On donne dans ce paragraphe un premier élément de réponse à ce problème de reconstruction
[?].
On appelle
N −1
1 X k
x̃(n, p) = Xk (n)ej2π N p
N
k=0

l’échantillon à l’instant p provenant d’une transformée de Fourier discrète inverse de Xk (n)


obtenu à l’instant n. Si l’on pondère cet échantillon par la valeur f (p − n) où f (n) est une
nouvelle fenêtre de pondération et que l’on additionne les différentes contributions à l’instant p,
on obtient le signal
+∞
X
x̂(p) = f (p − n)x̃(n, p)
n=−∞

+∞ N −1
X 1 X k
x̂(p) = f (p − n) Xk (n)ej2π N p . (B.2)
n=−∞
N
k=0

Le signal reconstruit x̂(p) à l’instant p est obtenu par ”recouvrement et addition” (overlap-add).
Donnons à la fenêtre d’analyse la possibilité de glisser à chaque étape de M échantillons avec

Nicolas Moreau 1er octobre 2007


page 70/78
Licence de droits d’usage
B.4. FENÊTRE GLISSANTE ET RECONSTRUCTION

1 ≤ M ≤ N , c’est-à-dire que l’on évalue une transformée de Fourier discrète à court terme tous
les mM échantillons. Les relations (B.1) et (B.2) deviennent
+∞
X k
Xk (m) = h(mM − l)x(l)e−j2π N l
l=−∞

+∞ N −1
X 1 X k
x̂(p) = f (p − mM ) Xk (m)ej2π N p .
m=−∞
N
k=0

Cherchons les conditions que doivent remplir les deux fenêtres de pondération h(n) et f (n) de
façon que l’on puisse reconstruire exactement le signal. On a
+∞ N −1 +∞
X 1 X X k k
x̂(p) = f (p − mM ) h(mM − l)x(l)e−j2π N l ej2π N p
m=−∞
N
k=0 l=−∞

+∞ +∞ N −1
X X 1 X −j2π k (l−p)
x̂(p) = [ f (p − mM )h(mM − l) e N ]x(l).
N
l=−∞ m=−∞ k=0

La condition de reconstruction parfaite est donc


+∞ N −1
X 1 X −j2π k (l−p)
f (p − mM )h(mM − l) e N = λ(p − l).
m=−∞
N
k=0

Comme
N −1
1 X −j2π k (l−p)
e N = 1 si l = p + qN
N
k=0
= 0 sinon

cette condition se traduit par


+∞
X
f (p − mM )h(mM − p − qN ) = λ(q).
m=−∞

En pratique, comment remplir cette condition ? Supposons que les deux fenêtres sont de durée
finie et de même durée N . On peut vérifier que cette condition est remplie, par exemple,
– si les deux fenêtres sont des fenêtres rectangulaires avec N = M ,
– si M = N/2 et si
πn
f (n)h(−n) = sin2 ( ) pour 0 ≤ n ≤ N − 1
N
= 0 sinon.

La théorie dite des bancs de filtres à reconstruction parfaite permet de généraliser cette étude.

Nicolas Moreau 1er octobre 2007


page 71/78
Licence de droits d’usage
ANNEXE B. TRANSFORMÉE DE FOURIER À COURT TERME

Nicolas Moreau 1er octobre 2007


page 72/78
Licence de droits d’usage
Annexe C

Egalisation

C.1 Introduction
– Emetteur : d(n) ∈ {0, 1}, a(n) ∈ {−1, 1}, he (t) = filtre d’émission.
X
xe (t) = a(n)he (t − nT )
n

T = temps bit, 1/T = débit.


– Canal : Canal convolutif bruité
X
xr (t) = hc (t) ∗ xe (t) + b(t) = a(n)h(t − nT ) + b(t)
n

avec h(t) = hc (t) ∗ he (t).


– Récepteur : Structure “optimale” dans le cas d’un bruit gaussien : xr (t) filtré par hr (t) puis
échantillonné à la cadence T ⇒ x(nT ) puis une prise de décision ⇒ â(n) ∈ {−1, 1} ou
ˆ ∈ {0, 1}.
d(n)
Conclusion : “Canal numérique équivalent” = boite noire recevant en entrée a(n), caractérisé
par une réponse impulsionnelle globale (à temps discret) g(n), bruitée par b(n), fournissant en
sortie X
x(n) = g(k)a(n − k) + b(n).
k
Hypothèses :
– b(n) = bruit blanc centré gaussien de puissance σB2.

– g(n) = filtre causal de durée finie L prenant en compte l’ensembles filtre d’émission/filtre
caractérisant le canal/filtre de réception.

x(n) = g(0)a(n) + g(1)a(n − 1) + · · · + g(L − 1)a(n − L + 1) + b(n).

x(n) dépend non seulement de a(n) mais aussi des symboles précédents. On parle d’Inter-
férence Entre Symboles (IES).
Problème (cf figure C.1) : Détermination d’un filtre supplémentaire de réponse impulsionnelle
w(n) atténuant l’IES de façon à simplifier la prise de décision.
– 1ère idée : “égaliseur zero forcing” = “inverser” le filtre g(n). Si le bruit est de très faible
puissance, il suffit de filtrer x(n) par le filtre w(n) tel que W (z) = 1/G(z) pour obtenir
y(n) ≈ a(n).
– 2ème idée : “égaliseur de Wiener” = minimiser E{|A(n) − Y (n)|2 }.
Préalable : estimer g(n) à partir de x(n) connaissant a(n). En effet, dans la pratique, on
connaît les filtres d’émission et de réception mais pas le filtre caractérisant le canal (qui peut
varier fortement au cours du temps).

Nicolas Moreau 1er octobre 2007


page 73/78
Licence de droits d’usage
ANNEXE C. EGALISATION

b(n)
a(n) s(n)
- +mx(n)
? y(n) â(n)
- g(n) - w(n) - Seuil -

Fig. C.1 – Canal numérique équivalent.

C.2 Propriétés des signaux


– a(n) : réalisation d’une suite de v.a. A(n) ∈ {−1, 1}, équiprobables et indépendantes. On
a donc E{A(n)} = (−1) ∗ 1/2 + 1 ∗ (1/2) = 0 et E{A(n)A(n + k)} = δ(k).
– x(n) : réalisation d’un p.a. X(n).
X
E{X(n)} = g(n − k)E{A(k)} + E{B(n)} = 0
k

X X
E{X(n)X(n + k)} = E{[ g(n − l)A(l) + B(n)][ g(n + k − p)A(p) + B(n + k)]}
l p
XX
= g(n − l)g(n + k − p)E{A(l)A(p)} + E{B(n)B(n + k)}
l p
X
2
= g(l)g(l + k) + σB δ(k).
l

C.3 Identification du canal


Exploitation d’une “séquence d’apprentissage” (système GSM : 23 symboles tous les 126 ?).
Deux méthodes :

C.3.1 Expression de l’intercorrélation

X
RAX (k) = E{A(n)X(n + k)} = E{A(n) g(n + k − l)A(l) + A(n)B(n + k)}
l
X
RAX (k) = g(n + k − l)E{A(n)A(l)} = g(k).
l

L’hypothèse d’ergodicité veut dire


+N
1 X
g(k) = RAX (k) = lim a(n)x(n + k).
N →∞ 2N + 1
n=−N

Dans la pratique, le nombre d’échantillons observés est fini. On utilise l’“estimateur”


N −1
1 X
ĝ(k) = a(n)x(n + k) pour k = 0 · · · L − 1.
N
n=0

Pour garder un effet de moyenne, il faut que N  L par exemple N ≈ 10L.

Nicolas Moreau 1er octobre 2007


page 74/78
Licence de droits d’usage
C.4. EGALISEUR ZERO FORCING

C.3.2 Estimateur “des moindres carrés”


On connaît {a(0) · · · a(N − 1)}. On mesure au récepteur {x(0) · · · x(N − 1)}.
    
a(0) .. .. (0)

x(0)
...
 .. .. .. 
g(0)
  .. 


 
  . . .   . 
 y(L)  
= a(L) ··· a(0)
 .. +
 
(L)

.
.
     
 ..   .. ..

..

. g(L − 1)
  
   . ··· .   . 
y(N − 1) a(N − 1) · · · a(N − L + 1) (N − 1)
x = Ag + .
Système sur dimensionné N  L. Critère : minimiser
||||2 = xt x − 2xt Ag + g t At Ag.

At x + At Ag opt = 0
g opt = [At A]−1 At x.
Interprétation géométrique : opt est orthogonal au sous espace engendré par les vecteurs colonnes
de A.

C.3.3 Relation entre les deux résultats


g 1 = [RAX (0) · · · RAX (L − 1)]t et g 2 = [At A]−1 At x. Exemple pour L = 2 :
 
  a(0) a(−1)    
1 t 1 a(0) · · · a(N − 1)  . . RAA (0) RAA (1) 1 0
AA= .. .. →

=
N N a(−1) · · · a(N − 2) RAA (1) RAA (0) 0 1

a(N − 1) a(N − 2)
 
RAX (0)
1 t
A x →  RAX (1)  .
 
N ..
.
Conclusion : g 1 et g 2 tendent vers la même limite lorsque N tend vers l’infini.

C.4 Egaliseur zero forcing


En absence de bruit, il suffit de prendre w(n) tel que w(n) ∗ g(n) = δ(n) soit W (z) = 1/G(z).
La puissance du bruit en sortie est alors donnée par
Z +1/2 2
σB
σ2 = 2
df.
−1/2 |G(f )|

C.4.1 1er problème


Si |G(f )| ≈ 0 dans une certaine bande de fréquences, la puissance du bruit en sortie peut
devenir très importante. Exemple du “canal à évanouissement” (trajets multiples)
h(t) = 1 + a1 δ(t − t1 ) + a2 δ(t − t2 ) + · · ·
H(f ) = 1 + a1 exp(−j2πf t1 ) + · · ·
|H(f )|2 = 1 + 2a1 cos(2πf t1 ) + a21 + · · ·
Si le premier écho est presque aussi puissant que le trajet direct, a ≈ 1 et |H(f )|2 ≈ 0 pour
f t1 ≈ 1/2.

Nicolas Moreau 1er octobre 2007


page 75/78
Licence de droits d’usage
ANNEXE C. EGALISATION

C.4.2 2ème problème


“Propriété de phase” de G(z). Exemples :
– Si G(z) = 1 − (1/2)z −1 ⇒ W (z) = 1 + (1/2)z −1 + (1/2)2 z −2 + · · · si |z| > 1/2 ⇒

y(n) = x(n) + (1/2)x(n − 1) + · · · + (1/2)Q x(n − Q).

– Mais si G(z) = 1 − 2z −1 le pôle de W (z) est à l’extérieur du cercle unité. Il faut trouver un
développement (en série de Laurent) dont le domaine d’existence contienne le cercle unité.
Il suffit d’écrire W (z) sous la forme

1 z 2 2 z z2
W (z) = = [1 + (1/2)z + (1/2) z + · · · ] = − − + ···
−2z −1 (1 − (1/2)z) 2 2 4

La réponse implusionnelle du filtre est anticausale :


1 1 1
y(n) = − x(n + 1) − x(n + 2) − · · · − Q x(n + Q).
2 4 2
Pas de problème particulier excepté le fait qu’il faille attendre l’arrivée de l’échantillon à
l’instant n + Q pour prendre une décision en fonction du signe de y(n).

C.5 Egaliseur de Wiener


w(l)X(n − l)|2 } ⇒ 2E{[A(n) − l w(l)X(n − l)]X(n − k)} = 0 ∀k.
P P
Minimiser E{|A(n) − l
X
E{A(n)X(n − k)} = w(l)E{X(n − l)]X(n − k)}
l

Comme
X
E{A(n)X(n − k)} = E{A(n)S(n − k)} = g(p)E{A(n)]A(n − k − p)} = g(−k)
p

, on obtient X
w(l)RXX (k − l) = g(−k) ∀k
l
ce qui donne dans le domaine fréquentiel :

G(−f )
W (f )SXX (f ) = G(−f ) ⇒ W (f ) = 2 .
|G(f )|2 + σB

Concrètement on résout
..
 
   .
w(−Q))

 0 

RXX (0) ..

RXX (Q)  .
..  
g(L − 1)

  
.. .. ..    ..


.  w(0)  =  .

 . .     .
RXX (Q) · · · RXX (0) 
 .
..  
g(0)

 



w(Q)

 0 

..
.

où RXX (k) est estimé à partir des données x(n) et g(k) par la procédure précédente d’identifi-
cation du canal.

Nicolas Moreau 1er octobre 2007


page 76/78
Licence de droits d’usage
C.5. EGALISEUR DE WIENER

La puissance de l’erreur est égale à


Z +1/2 Z +1/2
2 2 2
E{|A(n) − Y (n)| } = |1 − G(f )W (f )| df + σB |W (f )|2 df
−1/2 −1/2

ou en remplaçant W (f ) par sa valeur


Z +1/2 2 Z +1/2 2
2 σB σB
σ = 2 df ≤
|G(f )|2 + σB |G(f )|2
df.
−1/2 −1/2

La puissance du bruit en sortie est plus faible que dans le cas zero forcing mais au détriment de
l’IES.

Nicolas Moreau 1er octobre 2007


page 77/78
Licence de droits d’usage
ANNEXE C. EGALISATION

Contexte public } avec modifications


Par le téléchargement ou la consultation de ce document, l’utilisateur accepte la licence
d’utilisation qui y est attachée, telle que détaillée dans les dispositions suivantes, et s’en-
gage à la respecter intégralement.

La licence confère à l’utilisateur un droit d’usage sur le document consulté ou téléchargé, totalement
ou en partie, dans les conditions définies ci-après et à l’exclusion expresse de toute utilisation commer-
ciale.
Le droit d’usage défini par la licence autorise un usage à destination de tout public qui comprend :
– le droit de reproduire tout ou partie du document sur support informatique ou papier,
– le droit de diffuser tout ou partie du document au public sur support papier ou informatique, y compris
par la mise à la disposition du public sur un réseau numérique,
– le droit de modifier la forme ou la présentation du document,
– le droit d’intégrer tout ou partie du document dans un document composite et de le diffuser dans ce
nouveau document, à condition que :
– L’auteur soit informé.

Aucune modification du document dans son contenu, sa forme ou sa présentation n’est autorisée.

Les mentions relatives à la source du document et/ou à son auteur doivent être conservées dans leur
intégralité.

Le droit d’usage défini par la licence est personnel et non exclusif.


Tout autre usage que ceux prévus par la licence est soumis à autorisation préalable et expresse de l’auteur :
sitepedago@enst.fr

Nicolas Moreau 1er octobre 2007


page 78/78
Licence de droits d’usage

Vous aimerez peut-être aussi