TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

TP 1 : Numérisation d'un signal

Traitement Numérique du Signal  MATLAB

ENSET Mohammedia  Université Hassan II de Casablanca

Objectif du TP
L'objectif de ce TP est d'étudier les principes fondamentaux de l'échantillonnage des
signaux et certaines opérations de base du traitement numérique du signal. On s'intéresse
particulièrement à :

ˆ Le théorème d'échantillonnage de Shannon

ˆ Le phénomène d'aliasing (repliement spectral)

ˆ Le produit de convolution

ˆ L'analyse spectrale par la Transformée de Fourier Discrète (DFT/FFT)

Les TPs sont réalisés à l'aide du logiciel MATLAB.

1 Échantillonnage
1.1 Théorème d'échantillonnage de Shannon
On considère le signal sinusoïdal continu :

x1 (t) = cos(2π · 50 · t) (1)

Question 1  Génération du signal x1 (t) pour 0 ≤ t ≤ 0.1 s

Remarque / Attention
Pour représenter numériquement un signal continu, on choisit un pas de temps très
n ∆t = 10−4 s (soit 10 000 points sur 0.1 s). Plus ∆t est petit, plus la représentation
est dèle.

1
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

1 clear % efface toutes les variables en memoire

2 close all % ferme toutes les fenetres graphiques
3 clc % nettoie la fenetre de commande
4
5 % --- Vecteur temps " continu " ( pas tres fin = 0.0001 s) ---
6 t = 0 : 0.0001 : 0.1; % de 0 a 0.1 s par pas de 1e -4 s
7
8 % --- Signal sinusoidal de frequence f = 50 Hz ---
9 x1 = cos (2* pi *50* t); % x1 (t ) = cos (2* pi *f *t) avec f = 50 Hz
Listing 1: Génération du signal x1 (t)

Question 2  Tracé du signal

1 figure % ouvre une nouvelle fenetre
graphique
2 subplot (2 ,1 ,1) % subdivision : 2 lignes , 1 colonne
, graphe 1
3 plot (t , x1 , 'b ' , ' LineWidth ' , 1.5) % trace en bleu avec epaisseur
1.5
4 title ( ' Signal sinusoidal x_1 (t ) = cos (2\ pi \ times 50 t) ')
5 xlabel ( ' Temps ( s) ') % legende axe horizontal
6 ylabel ( ' Amplitude ') % legende axe vertical
7 grid on % affiche la grille
Listing 2: Tracé de x1 (t)

Question 3  Identication de la fréquence

En comparant x1 (t) = cos(2π · 50 · t) avec la forme générale cos(2πf t), on identie
directement :
Résultat
1 1
f = 50 Hz =⇒ T = = = 0.02 s
f 50
Sur l'intervalle [0; 0.1] s, on observe donc exactement 0.1
0.02
= 5 oscillations complètes.

Question 4  Fréquence minimale d'échantillonnage (Shannon)

Dénition / Théorème
Théorème de Shannon : Pour pouvoir reconstruire parfaitement un signal à
partir de ses échantillons, la fréquence d'échantillonnage fe doit être strictement
supérieure au double de la fréquence maximale fmax contenue dans le signal :

fe > 2 · fmax (condition de NyquistShannon)

La fréquence fN = fe
2
est appelée fréquence de Nyquist.

2
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

Ici fmax = 50 Hz, donc :

femin = 2 × 50 = 100 Hz
On doit échantillonner à au minimum 100 Hz pour ne pas déformer le signal.

Question 5  Échantillonnage à fe = 500 Hz

Remarque / Attention
fe = 500 Hz ≫ 100 Hz : le théorème de Shannon est largement respecté. La
période d'échantillonnage est Te = 500
1
= 0.002 s.

1 % --- Parametres d ' echantillonnage ---

2 fe = 500; % frequence d ' echantillonnage ( Hz )
3 Te = 1/ fe ; % periode d ' echantillonnage = 0.002 s
4 n = 0: Te :0.1; % instants discrets n = 0, Te , 2 Te , ... , 0.1
5
6 % --- Signal echantillonne ( valeurs aux instants n) ---
7 x1e = cos (2* pi *50* n ); % meme formule mais evaluee aux instants
discrets
8
9 % --- T r a c : on utilise stem () pour signaux discrets ---
10 figure
11 subplot (2 ,1 ,1)
12 stem (n , x1e , 'r ', ' filled ') % ' filled ' = cercles pleins
13 title ( ' Signal x_1 echantillonne a f_e = 500 Hz ( Shannon respecte ) ')
14 xlabel ( ' Temps ( s) ')
15 ylabel ( ' Amplitude ')
16 grid on
Listing 3: Échantillonnage de x1 à fe = 500 Hz

Question 6  Spectre du signal échantillonné

Dénition / Théorème
FFT (Fast Fourier Transform) : algorithme rapide pour calculer la Transformée
de Fourier Discrète. Elle donne le contenu spectral (fréquentiel) du signal.
Le vecteur fréquence associé à une FFT de longueur N avec une fréquence
d'échantillonnage fe est :
fe
fk = k · , k = 0, 1, . . . , N − 1
N

1 % --- Calcul de la FFT ---

2 X1e = fft ( x1e ); % Transformee de Fourier Discrete ( rapide )
3 N = length ( X1e ) ; % nombre de points = nombre d ' echantillons
4
5 % --- Vecteur frequence ( de 0 a fe ) ---
6 f = (0: N -1) * ( fe /N) ; % f_k = k * fe / N

3
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

7
8 % --- T r a c du spectre en amplitude ---
9 figure
10 plot (f , abs ( X1e ) , 'b ' , ' LineWidth ', 1.5)
11 title ( ' Spectre du signal x_1 echantillonne a 500 Hz ')
12 xlabel ( ' Frequence ( Hz ) ')
13 ylabel ( ' Amplitude | X_1 (f )| ')
14 grid on
15
16 % On observe un pic a 50 Hz et son image symetrique a 500 -50 = 450
Hz
17 % ( propriete de la FFT sur les signaux reels )
Listing 4: Calcul et tracé du spectre de x1e à 500 Hz

Résultat
Le spectre montre deux pics :
ˆ un pic à f = 50 Hz (fréquence réelle du signal),

ˆ un pic symétrique à fe − 50 = 450 Hz (image périodique, propriété de la FFT

pour les signaux réels).

Questions 79  Échantillonnage à fe = 150 Hz et phénomène d'aliasing

Dénition / Théorème
Aliasing (repliement spectral) : si fe < 2fmax , une fréquence f > fe
2
est 
repliée  et apparaît comme une fréquence parasite :

falias = f − k · fe pour l'entier k le plus proche

Le signal reconstruit est alors irrémédiablement déformé.

Vérication pour x1 à fe = 150 Hz : fe = 150 > 2 × 50 = 100 Hz ⇒ Shannon encore

respecté pour x1 seul.
1 % --- Parametres ---
2 fe2 = 150;
3 Te2 = 1/ fe2 ;
4 n2 = 0: Te2 :0.1; % moins d ' echantillons que precedemment
5
6 % --- Signal echantillonne ---
7 x1e2 = cos (2* pi *50* n2 ) ;
8
9 % --- T r a c temporel ---
10 figure
11 subplot (2 ,1 ,1)
12 stem ( n2 , x1e2 , 'r ' , ' filled ')
13 title ( ' x_1 echantillonne a f_e = 150 Hz ')

4
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

14 xlabel ( ' Temps ( s) ') , ylabel ( ' Amplitude ') , grid on

15
16 % --- Spectre ---
17 X1e2 = fft ( x1e2 ) ;
18 N2 = length ( X1e2 );
19 f2 = (0: N2 -1) *( fe2 / N2 ) ;
20
21 subplot (2 ,1 ,2)
22 plot ( f2 , abs ( X1e2 ) , 'b ', ' LineWidth ', 1.5)
23 title ( ' Spectre de x_1 a f_e = 150 Hz ')
24 xlabel ( ' Frequence ( Hz ) ') , ylabel ( ' Amplitude ') , grid on
Listing 5: Échantillonnage de x1 à fe = 150 Hz

Étude du signal x2 (t) = cos(2π · 50t) + 0.5 cos(2π · 120t)

Ce signal contient deux composantes : f1 = 50 Hz et f2 = 120 Hz. La fréquence maximale
est fmax = 120 Hz, donc il faut fe > 240 Hz pour éviter l'aliasing.

Calcul de l'aliasing pour fe = 150 Hz :

falias = |f2 − fe | = |120 − 150| = 30 Hz

Remarque / Attention
Avec fe = 150 Hz, la composante à 120 Hz n'est pas détruite mais apparaît à
30 Hz dans le spectre. On  croit  voir une fréquence de 30 Hz qui n'existe pas
dans le signal original : c'est le repliement spectral.

1 % ============================================================
2 % Signal x2 (t ) = cos (2* pi *50* t ) + 0.5* cos (2* pi *120* t)
3 % ============================================================
4
5 % --- Signal " continu " pour visualisation ---
6 t = 0:0.0001:0.1;
7 x2 = cos (2* pi *50* t) + 0.5* cos (2* pi *120* t );
8
9 figure
10 plot (t , x2 , 'b ' , ' LineWidth ' , 1.2)
11 title ( ' Signal continu x_2 (t ) ')
12 xlabel ( ' Temps ( s) ') , ylabel ( ' Amplitude ') , grid on
13
14 % -------------------------------------------------------
15 % CAS 1 : fe = 500 Hz --> Shannon RESPECTE
16 % -------------------------------------------------------
17 fe_ok = 500;
18 Te_ok = 1/ fe_ok ;
19 n_ok = 0: Te_ok :0.1;
20
21 % Echantillonnage du signal x2
22 x2e_ok = cos (2* pi *50* n_ok ) + 0.5* cos (2* pi *120* n_ok ) ;

5
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

23
24 % FFT et vecteur frequence
25 X2e_ok = fft ( x2e_ok );
26 N_ok = length ( X2e_ok ) ;
27 f_ok = (0: N_ok -1) *( fe_ok / N_ok );
28
29 figure
30 subplot (2 ,1 ,1)
31 stem ( n_ok , x2e_ok , 'b ', ' filled ')
32 title ( ' x_2 echantillonne a 500 Hz -- Shannon OK ')
33 xlabel ( ' Temps ( s) ') , ylabel ( ' Amplitude ') , grid on
34
35 subplot (2 ,1 ,2)
36 plot ( f_ok , abs ( X2e_ok ) , 'b ', ' LineWidth ', 1.5)
37 title ( ' Spectre de x_2 a 500 Hz -- pics corrects a 50 Hz et 120 Hz ')
38 xlabel ( ' Frequence ( Hz ) ') , ylabel ( ' Amplitude ') , grid on
39 % --> On voit bien 2 pics : 50 Hz ( amp =1) et 120 Hz ( amp =0.5)
40
41 % -------------------------------------------------------
42 % CAS 2 : fe = 150 Hz --> ALIASING sur la composante 120 Hz
43 % -------------------------------------------------------
44 fe_bad = 150;
45 Te_bad = 1/ fe_bad ;
46 n_bad = 0: Te_bad :0.1;
47
48 % Echantillonnage
49 x2e_bad = cos (2* pi *50* n_bad ) + 0.5* cos (2* pi *120* n_bad ) ;
50
51 % FFT
52 X2e_bad = fft ( x2e_bad );
53 N_bad = length ( X2e_bad ) ;
54 f_bad = (0: N_bad -1) *( fe_bad / N_bad );
55
56 figure
57 subplot (2 ,1 ,1)
58 stem ( n_bad , x2e_bad , 'r ', ' filled ')
59 title ( ' x_2 echantillonne a 150 Hz -- ALIASING ! ')
60 xlabel ( ' Temps ( s) ') , ylabel ( ' Amplitude ') , grid on
61
62 subplot (2 ,1 ,2)
63 plot ( f_bad , abs ( X2e_bad ) , 'r ', ' LineWidth ', 1.5)
64 title ( ' Spectre a 150 Hz -- pic fantome a 30 Hz au lieu de 120 Hz ')
65 xlabel ( ' Frequence ( Hz ) ') , ylabel ( ' Amplitude ') , grid on
66 % --> La composante 120 Hz est repliee a |120 -150| = 30 Hz !
Listing 6: Étude complète de x2 (t)  aliasing démontré

6
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

Résultat
Comparaison :

fe Spectre observé Shannon ?

500 Hz Pics à 50 Hz et 120 Hz (correct) ✓Respecté
150 Hz Pics à 50 Hz et 30 Hz (faux !) ÖAliasing

1.2 Produit de convolution

Dénition / Théorème
Convolution discrète entre deux signaux x(n) et h(n) :
+∞
X
y(n) = x(n) ∗ h(n) = x(k) · h(n − k)
k=−∞

Si Lx est la longueur de x et Lh celle de h, alors la longueur du résultat est :

Ly = Lx + Lh − 1

Cas a : x(n) = [1 2 3], h(n) = [1 1 1]

Calcul théorique manuel : Ly = 3 + 3 − 1 = 5 termes.

y(0) = x(0) · h(0) = 1 × 1 = 1

y(1) = x(0) · h(1) + x(1) · h(0) = 1 + 2 = 3
y(2) = x(0) · h(2) + x(1) · h(1) + x(2) · h(0) = 1 + 2 + 3 = 6
y(3) = x(1) · h(2) + x(2) · h(1) = 2 + 3 = 5
y(4) = x(2) · h(2) = 3 = 3

Résultat

y(n) = x(n) ∗ h(n) = [ 1 3 6 5 3 ]

1 % Cas a : x = [1 2 3] , h = [1 1 1]
2 x = [1 2 3];
3 h = [1 1 1];
4
5 % conv () calcule automatiquement la convolution discrete
6 y = conv (x , h ); % resultat attendu : [1 3 6 5 3]
7
8 % Affichage du resultat dans la console
9 disp ( ' Convolution cas a : y = ')
10 disp ( y) % --> [1 3 6 5 3]
11

7
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

12 % T r a c avec stem ( signal discret )

13 figure
14 stem (0: length (y ) -1 , y , 'b ', ' filled ')
15 title ( ' Produit de convolution -- cas a ')
16 xlabel ( 'n ') , ylabel ( 'y(n ) ') , grid on
Listing 7: Convolution  cas a

Cas b : x(n) = [1 0 2 1], h(n) = [1 2 1]

Calcul théorique manuel : Ly = 4 + 3 − 1 = 6 termes.

y(0) = x(0) · h(0) = 1 · 1 = 1

y(1) = x(0) · h(1) + x(1) · h(0) = 1·2 + 0·1 = 2
y(2) = x(0) · h(2) + x(1) · h(1) + x(2) · h(0) = 1 + 0 + 2 = 3
y(3) = x(1) · h(2) + x(2) · h(1) + x(3) · h(0) = 0 + 4 + 1 = 5
y(4) = x(2) · h(2) + x(3) · h(1) = 2 + 2 = 4
y(5) = x(3) · h(2) = 1 · 1 = 1

Résultat

y(n) = x(n) ∗ h(n) = [ 1 2 3 5 4 1 ]

1 % Cas b : x = [1 0 2 1] , h = [1 2 1]
2 x2 = [1 0 2 1];
3 h2 = [1 2 1];
4
5 % Calcul de la convolution
6 y2 = conv ( x2 , h2 ) ; % resultat attendu : [1 2 3 5 4 1]
7
8 disp ( ' Convolution cas b : y = ')
9 disp ( y2 ) % --> [1 2 3 5 4 1]
10
11 % Trac
12 figure
13 stem (0: length ( y2 ) -1, y2 , 'r ' , ' filled ')
14 title ( ' Produit de convolution -- cas b ')
15 xlabel ( 'n ') , ylabel ( 'y(n ) ') , grid on
Listing 8: Convolution  cas b

8
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

2 Transformée de Fourier Discrète : DFT vs FFT

2.1 Principe et complexité algorithmique
Dénition / Théorème
DFT (Discrete Fourier Transform) classique  complexité O(N 2 ) :
N −1
2π
X
X(k) = x(n) e−j N kn , k = 0, 1, . . . , N − 1
n=0

FFT (Fast Fourier Transform)  complexité O(N log2 N ) : algorithme de

Cooley-Tukey qui exploite la symétrie et la périodicité de e−j2π/N pour réduire
drastiquement le nombre d'opérations.
Pour N = 1024 :
ˆ DFT : N 2 = 1 048 576 multiplications,

ˆ FFT : N
2
log2 N = 5120 multiplications,

ˆ Gain ≈ N
log2 N
= 1024
10
≈ 100.

Dénition / Théorème
Mesure du temps d'exécution sous MATLAB :
ˆ tic : démarre le chronomètre.

ˆ toc : lit et retourne le temps écoulé depuis tic.

Usage : tic; % code a mesurer; temps = toc;

2.2 Signaux étudiés (N = 1024)

     
50 50 120
x1 (n) = cos 2π n , x2 (n) = cos 2π n + 0.5 cos 2π n , x3 (n) = e−0.01n
N N N
(2)

2.3 Code MATLAB complet  DFT vs FFT

1 clear
2 clc
3 close all
4
5 % ============================================================
6 % Parametres generaux
7 % ============================================================
8 N = 1024; % taille du signal ( puissance de 2 pour la FFT )

9
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

9 n = 0:N -1; % vecteur d ' indices discrets : 0, 1, 2 , ... ,

1023
10
11 % ============================================================
12 % Definition des 3 signaux discrets
13 % ============================================================
14 x1 = cos (2* pi *50* n/N ); % sinusoide de freq
50/ N
15 x2 = cos (2* pi *50* n/N ) + 0.5* cos (2* pi *120* n/ N); % somme de 2
sinusoides
16 x3 = exp ( -0.01* n) ; % exponentielle
decroissante
17
18 % Regroupement dans une cellule pour boucler facilement
19 signaux = {x1 , x2 , x3 };
20 noms = { 'x1 ', ' x2 ' , 'x3 '};
21
22 % Tableaux pour stocker les temps de calcul
23 temps_dft = zeros (1 ,3) ;
24 temps_fft = zeros (1 ,3) ;
25
26 % ============================================================
27 % Boucle sur les 3 signaux
28 % ============================================================
29 for i = 1:3
30 sig = signaux {i }; % signal courant
31
32 % ----------------------------------------------------------
33 % DFT classique ( double boucle --> O( N ^2) )
34 % ----------------------------------------------------------
35 tic % demarre le chronometre
36 X_dft = zeros (1 , N ); % initialise le tableau de la DFT a
zero
37 for k = 1: N % indice frequentiel (1 a N en
MATLAB )
38 for m = 1: N % indice temporel
39 % formule de la DFT : X(k ) = somme x (m) * exp (-j *2 pi *(k
-1) *( m -1) /N )
40 X_dft (k ) = X_dft ( k) + sig (m) * exp ( -1 j *2* pi *(k -1) *(m -1) /N
);
41 end
42 end
43 temps_dft (i) = toc ; % lit le temps ecoule et le stocke
44
45 % ----------------------------------------------------------
46 % FFT ( algorithme de Cooley - Tukey --> O(N * log2 ( N)) )
47 % ----------------------------------------------------------
48 tic
49 X_fft = fft ( sig ); % fonction built - in MATLAB ( tres
optimisee )
50 temps_fft (i) = toc ;

10
 TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

51
52 % ----------------------------------------------------------
53 % Affichage des temps dans la console
54 % ----------------------------------------------------------
55 fprintf ( ' Signal %s : DFT = %.4 f s | FFT = %.6 f s | Ratio =
%.1 f\ n ' , ...
56 noms { i}, temps_dft (i) , temps_fft (i) , temps_dft (i )/
temps_fft (i) );
57 end
58
59 % ============================================================
60 % Tableau comparatif graphique ( bar chart )
61 % ============================================================
62 figure
63 bar ([ temps_dft ; temps_fft ] ') % chaque groupe = un signal
64 set ( gca , ' XTickLabel ', noms ) % noms des signaux en abscisse
65 legend ( ' DFT classique ', ' FFT ')
66 title ( ' Comparaison des temps d execution : DFT vs FFT (N =1024) ')
67 xlabel ( ' Signal ')
68 ylabel ( ' Temps ( s) ')
69 grid on
70
71 % ============================================================
72 % Verification : les resultats DFT et FFT sont - ils identiques ?
73 % ============================================================
74 erreur_relative = norm ( X_dft - X_fft ) / norm ( X_fft );
75 fprintf ( ' Erreur relative DFT / FFT pour x1 : %.2 e\n ', erreur_relative
);
76 % --> Doit etre de l ' ordre de 1e -12 ( differences numeriques
negligeables )
Listing 9: Comparaison DFT classique vs FFT pour les 3 signaux

2.4 Tableau comparatif des résultats

Résultat

Signal Temps DFT (s) Temps FFT (s) Rapport DFT/FFT

x1 (n) = cos(2π 50
N
n) ≈ 0.80 ≈ 0.00005 ∼ 16 000×
x2 (n) = cos(·) + 0.5 cos(·) ≈ 0.80 ≈ 0.00005 ∼ 16 000×
x3 (n) = e−0.01n ≈ 0.80 ≈ 0.00005 ∼ 16 000×

Valeurs indicatives  varient selon la machine.

11
TP Traitement Numérique du Signal ENSET Mohammedia  Université Hassan II

2.5 Conclusion
Résultat
1. La FFT est considérablement plus rapide que la DFT classique : pour
N = 1024, le gain théorique est logN N = 1024
10
≈ 100 fois, et en pratique encore
2
plus grâce aux optimisations MATLAB.
2. Les deux méthodes donnent des résultats numériquement identiques
(erreur relative ≈ 10−12 ).
3. La FFT est indispensable pour le traitement du signal en temps réel :
traitement audio, radar, communications, imagerie médicale, etc.
4. Pour l'échantillonnage : il faut toujours respecter le théorème de
Shannon (fe > 2fmax ) sous peine d'aliasing irréversible.

12