SAT ENS ENSSAT

ENSSAT EII2
TP de Traitement Numrique du Signal
Il est demand de prparer ce TP
Premire partie
I. Analyse et synthse de ltres numriques
Nous allons tudier dans la premire partie de ce TP, lanalyse et la synthse de ltres numriques
rponses impulsionnelles nie et innie laide des logiciels Matlab et dune applet JAVA dveloppe
lENSSAT.
1 Analyse de ltres numriques passe-bas du second ordre
1.1 Analyse dun ltre numrique RII passe-bas du second ordre
Le ltre numrique est spci par le gabarit en dB de la gure 1.(b) avec comme paramtres :
Fp = 1kHz, Fs = 3kHz,
1
= 20 log(1 +
1
) = 1dB,
2
= 20 log(
2
) = 20dB, Fe = 10kHz
.

p
1
|H(e
j
)|
1-
1
1+
1

p
0 dB
|H(e
j
)| (dB)
20log(1-
1
)
20log(1+
1
)
20log
2

a
a) Gabarit frquentiel linaire b) Gabarit frquentiel en dB
Fig. 1 Gabarit des ltres passe-bas tudier
La synthse par la mthode bilinaire du ltre de Chebyshev H(z) aboutit la fonction :
H(z) =
0.079(z + 1)
2
z
2
1.2z + 0.516
Questions : ltre numrique RII en prcision innie
crire le programme Matlab permettant de visualiser la rponse frquentielle en module et en phase
du ltre numrique (voir lannexe A dcrivant ce programme Matlab). Faites la relation entre lchelle
des abscisses de la commande plot de Matlab et les frquences. Quelles sont donc les relations entre
frquence dchantillonnage et frquence de coupure ? Le ltre respecte-t-il le gabarit spci ?
Visualisez la rponse impulsionnelle.
Tracez dans le plan complexe les ples et les zros du ltre. Le ltre est-il stable ?
1.2 Analyse dun ltre numrique RIF passe-bas du second ordre
Le gabarit du ltre passe-bas tant dni comme la gure 1, la synthse du ltre donne deux fonctions
de transfert H(z) selon le nombre de points considr (N = 11 puis N = 7) :
H
11
(z) = 0, 0309396.(z
1
+z
9
)0, 0390182.(z
2
+z
8
)+0, 0766059.(z
3
+z
7
)+0, 288307.(z
4
+z
6
)+0, 4.z
5
H
7
(z) = 0, 0409365.(1 + z
6
) + 0, 078369.(z
1
+ z
5
) + 0, 289996.(z
2
+ z
4
) + 0, 4.z
3
Questions : ltre numrique RIF en prcision innie
Reprendre les questions du 1.1 (voir lannexe B dcrivant ce programme Matlab).
2 Synthse dun ltre passe-bas du second ordre
On dsire raliser un ltre numrique H(z) quivalent un ltre analogique de Butterworth passe-
bas H(j) du deuxime ordre qui prsente une frquence de coupure Fc de 1KHz. La frquence
dchantillonnage Fe sera de 8kHz. Son ltre de Butterworth normalis est le suivant :
H
Normalise
(j
N
) =
1
1 +

2j
N

2
N
H
Normalise
(p
N
) =
1
1 +

2p
N
+ p
2
N
2.1 Synthse dun ltre RII
Prparation
Le gabarit dun ltre passe-bas est dni comme indiqu sur la gure 1.(b). Calculer les paramtres
1 et 2 pour que le ltre passe bas de frquence de coupure Fc = 1KHz entre dans le gabarit,
lorsque Fp = 0.7Fc et Fs = 3Fc.
Faire la synthse par la mthode bilinaire du ltre H(j) an dobtenir H
b
(z). On tudiera lin-
uence de la distorsion en frquence implique par la mthode. On rappelle que la transformation
bilinaire est obtenue par p = f(z) 2Fe
1z
1
1+z
1
.
Faire la synthse par la mthode dinvariance impulsionnelle an dobtenir H
i
(z).
Tracez dans le plan complexe les ples et les zros des ltres. Les ltres sont-il stables ?
Manipulations sur Java_Filtre
Un logiciel de synthse de ltres numriques a t dvelopp lENSSAT. Il sagit dune applet JAVA
que lon peut lancer soit laide de la commande tns sous un serveur SUN, soit en vous connectant
laide dun navigateur lURL : http ://lasti.enssat.fr/GroupeArchi/enseignements/Tns/Tns.php.
Ce logiciel vous permet de choisir entre la synthse de ltres RII (Butterworth ou Chebychev) selon
la mthode bilinaire, et la synthse de ltres RIF par les mthodes de fentrage (Rectangle, Triangle,
Hamming, Hanning, Blackman, Kaiser), dchantillonnage frquentiel, ou de Parks MacClellan. Une
synthse directe par les ples et les zros est possible. Aprs avoir entr les paramtres 1, 2, Fp, Fs
(cf gabarit de la gure 1.(b)), la synthse seectue automatiquement et vous donne les valeurs des
coecients du ltre et la rponse impulsionnelle. Une aide est disponible la mme URL.
Dterminez laide de Java_Filtre les coecients des ltres prcdents.
Manipulations sur Matlab
crire le programme Matlab permettant de faire la synthse des ltres H
b
(z) et H
i
(z). Lannexe C
est complter pour rpondre cette question.
Visualiser la rponse frquentielle en module et en phase des ltres numriques. Les ltres respectent-
ils le gabarit spci ?
Observer leurs rponses impulsionnelles.
Comparer les deux mthodes dans les domaines frquentiel et temporel.
2.2 Synthse dun ltre RIF
Prparation
Soit le ltre passe bas idal de frquence de coupure c. Dterminez la rponse impulsionnelle h(n)
de ce ltre. On prendra comme application numrique une frquence de coupure Fc = 1kHz pour
une frquence dchantillonnage Fe = 8kHz et un f de 2kHz.
Quel est la longueur thorique du ltre respectant ce gabarit ?
Le gabarit du ltre passe-bas est dni comme indiqu sur la gure du 1.(b). On utilisera les options
suivantes de lapplet Java_Filtre qui correspondent 3 mthodes direntes de synthse :
Synthse selon la mthode de fentrage. Le gabarit idal H
i
() est chantillonn, et la rponse
impulsionnelle correspondante h
i
(n) est calcule puis tronque sur un nombre ni dchantillons M.
Cette rponse est ensuite multiplie par la fentre choisie. On obtient h(n) = h
i
(n).w(n). Le nombre
de coecients du ltre est donn par les estimations vues en cours.
Synthse par la mthode de lchantillonnage frquentiel.
Synthse base sur la mthode de Parks-MacClellan. Cest une procdure itrative qui ap-
proxime la fonction de transfert par une somme de cosinus jusqu convergence. Le nombre de points
de h(n) est donn par :
M
pm
=
10 log(12) 13
14.6fTe
Les deux ltres obtenus sont rponses impulsionnelles symtriques.
Manipulation sur Java_Filtre
Observez et comparez les coecients obtenus par les trois mthodes partir de la spcication dans
Java_Filtre correspondant au gabarit.
Observez galement le comportement de la synthse pour dautres types de fentres.
Manipulations sur Matlab
crire le programme Matlab permettant de faire la synthse du ltre RIF H(z) par la mthode
du fentrage pour des fentres rectangulaire et de Hamming. Le programme de lannexe D est
complter.
Visualisez la rponse frquentielle en module et en phase, observez la rponse impulsionnelle.
tudiez linuence de la fentre sur le ltre numrique. Peut on mettre en vidence le phnomne
de Gibbs.
crire le programme Matlab permettant de faire la synthse du ltre RIF H(z) par la mthode de
lchantillonnage frquentiel, en faisant attention ce que le ltre soit phase linaire.
Vriez que le ltre passe par les points frquentiels spcis en module et en phase.
Entrez les spcications pour synthtiser un ltre passe bande puis un ltre drivateur.
Deuxime partie
II. Algorithmes de traitement numrique du signal
3 Filtrage numrique en prcision nie
3.1 Etude dun ltre numrique RII passe-bas du second ordre en prcision nie
Ltude suivante reprend lexemple du ltre de la question 1.1 de la premire partie.
Les coecients sont arrondis sur 7 puis 3 bits aprs la virgule. laide du programme Matlab fourni
en annexeE, tudier les modications de la rponse frquentielle du ltre quanti, en particulier
par rapport au gabarit. Quelles sont les problmes apports par la quantication sur la stabilit du
ltre ? On placera un nombre susant de bits avant la virgule pour viter les dbordements.
On se basera sur le diagramme des ples et des zros du ltre en vriant linuence de la quanti-
cation sur les valeurs des coecients du ltre puis sur les valeurs des ples. Tracez le diagramme
des ples et des zros avant aprs quantications sur 8 et 4 bits.
tudiez les problmes de dbordements du ltre. A quelles conditions le ltre peut il diverger ?
3.2 tude dun ltre RII du quatrime ordre
Soit le ltre :
H(z) = 0.0599
1 2.1832z
1
+ 2.9976z
2
2.1832z
3
+ z
4
1 3.1914z
1
+ 4.17z
2
2.5854z
3
+ 0.6443z
4
Reprendre les questions prcdentes an de comparer limplmentation de ce ltre en prcision in-
nie et en prcision nie. Que remarquez vous ?
La fonction de transfert peut tre factorise comme suit :
H(z) =

0.0599
1 0.6511z
1
+ z
2
1 1.5684z
1
+ 0.6879z
2

0.0599
1 1.5321z
1
+ z
2
1 1.623z
1
+ 0.9366z
2
Cette factorisation aecte-t-elle la rponse frquentielle du ltre en prcision innie ?
Montrer que lon peut stabiliser le ltre en prcision nie en le mettant sous forme dune cascade de
cellules du second ordre. Quen concluez vous ?
3.3 tude dun ltre numrique RIF passe-bas en prcision nie
Reprendre les questions du 3.1 sur le ltre :
H
11
(z) = 0.0309396(z
1
+z
9
)0.0390182(z
2
+z
8
)+0.0766059(z
3
+z
7
)+0.288307(z
4
+z
6
)+0.4z
5
On vriera particulirement les points suivants : - caractristique en phase du ltre ; - problmes
de dbordement de laccumulateur.
4 Convolution Discrte
Soit le ltre h(n) dni la gure 2 pour M = 7, nous allons tudier le ltrage dun signal x(n) de
dure N limite, dni par une convolution discrte
y(n) =
n

i=0
h(i)x(n i), n = 0 . . . N + M 1,
au moyen de deux mthodes : une mthode temporelle directe puis une mthode frquentielle rapide.
h(n)
n
0.1
0.2
0.3
-1 0 1 2 3 4 5 6 7
M-1

Fig. 2 Signal h(n)
x(n) est une impulsion rectangulaire de longueur N = 11.
crire le programme Matlab ralisant la fonction de convolution dans le domaine temporelle laide
du programme de lannexe F. La fonction conv(x, h) prend 2 vecteurs de longueur N et M et
retourne un vecteur de longueur N + M 1.
crire le programme Matlab ralisant la fonction de convolution dans le domaine frquentiel par
lutilisation de la fonction FFT.
Vrier que le rsultat est identique dans les deux cas.
Calculer prcisment les complexits des deux approches. Comparer les complexits des mthodes
frquentielle et temporelle pour direntes valeurs de N et M. On pourra par exemple regarder les
cas N = 11, M = 7 ; N = M = 1024 ; N = 1024, M = 11.
5 Analyse spectrale
5.1 tude des fentres temporelles
tudiez les rponses temporelles et frquentielles des fentres Rectangulaire, Bartlett, Hamming,
Hanning, Blackman laide du programme de lannexe H.
Retrouvez les caractristiques du cours (largeur du lobe principal, attnuation du lobe secondaire)
et tudiez linuence du nombre de points N de la fentre.
5.2 Analyse spectrale dun signal sinusodal
Tracez la caractristique idale du spectre dune sinusode chantillonne. On prendra les valeurs
suivantes : frquence dchantillonnage Fe = 300 ; frquence de la sinusode F
0
= 50.
Lanalyse spectrale dun signal ne peut en pratique seectuer que sur un nombre limit de points. On
appelle donc T
0
= N.T lhorizon dobservation du signal. La multiplication temporelle du signal par
une fentre de longueur T
0
implique une convolution dans le domaine des frquences.
Donner les nesses en frquence et en amplitude pour les fentres rectangulaire et Hamming pour N
= 32 et N = 256.
crire le programme Matlab ralisant lanalyse spectrale de ce signal sinusodal dans le cas des
fentres rectangulaire et Hamming pour N = 32 et N = 256.
Que se passe-t-il lorsque la priode de la sinusode est un multiple entier de lhorizon dobservation
(T = K.T
0
= K.N.T).
crire le programme Matlab ralisant lanalyse spectrale dun signal form de la somme de deux
sinusodes de frquence 50Hz et 60Hz dans le cas des fentres rectangulaire et Hamming pour N =
32 et N = 256.
6 Dcimation et Interpolation
6.1 Dcimation
La dcimation dun signal consiste diminuer sa frquence dchantillonnage. Ce signal doit tre
bande troite pour respecter le thorme dchantillonnage. Ceci implique un ltrage passe bas avant
lopration de dcimation.
N = 192; n = 0:1:N-1;
x = sin(2*pi*0.02*n)+sin(2*pi*0.09*n);
tudiez les eets sur le spectre dun signal dune dcimation par 2 puis par 4. On exprimera la
dcimation par M par : y1 = x(1 :2 :N) ;.
Utilisez ensuite la fonction decimate(x,M) de Matlab qui eectue en plus un ltrage anti-repliement
avant dcimation. Comparez avec la version sans ltrage.
tudiez maintenant les eets sur le spectre dun signal dune dcimation par 6. Expliquez le rsultat
obtenu.
6.2 Interpolation
Linterpolation dun signal consiste augmenter sa frquence dchantillonnage. Les deux principales
mthodes sont linsertion de L1 zros dans la rponse temporelle ralisant une interpolation par un
facteur L (voir cours et gure 3.a), et lajout de Z zros dans le domaine frquentiel impliquant une
interpolation par un facteur L =
N+Z
N
o N est la taille de la FFT (voir gure 3.b).

L
x(nT) x
e
(nT')
Filtre Passe Bas
Gain = L

x
i
(nT')
L
x(nT) x
e
(nT')
Filtre Passe Bas
Gain = L

x
i
(nT')


x(n)
TFD
N pts
Fe
x(n)
ajout
Z zeros
TFD
-1
ZN
N
(ZN).Fe
x(n)
TFD
N pts
Fe
x(n)
ajout
Z zeros
TFD
-1
ZN
N
(ZN).Fe

(a) Interpolation par insertion de zros et ltrage (b) Interpolation par TFD et ajout de zros
Fig. 3 Interpolation
A partir du signal utilis dans la section 6.1 dcim dun facteur 4, ralisez linsertion de zros
(lvateur de frquence) et vriez son impact sur le spectre selon la mthode prsente en cours.
On prendra un facteur dinterpolation L = 4.
Raliser linterpolateur par insertion de zros et ltrage RIF passe bas.
Raliser linterpolateur par par TFD et ajout de zros.
Comparer qualitativement les deux mthodes dinterpolation.
6.3 Application au ltrage par dcimation-interpolation
La ralisation de ltres numriques peut tre simplie par une srie de traitements interpolation-
ltrage-dcimation. Lexemple le plus connu est le ltre de lissage du Compact Disc qui ralise un
ltrage trs slectif Fe/2. Il est en pratique ralis par le schma ci dessous.

x(n)
Fe/2
y(n) x(n)
Fe
y(n)
4 4
Fe/2 Fe
x(n)
Fe/2
y(n) x(n)
Fe
y(n)
4 4
Fe/2 Fe

Fig. 4 Filtrage par interpolation-dcimation
Comparer les fonctions de transfert et les complexits des deux types de ltrage par la mthode directe
et la mthode multicadence.
A Analyse dun ltre numrique RII passe-bas du second ordre
close all ;
% Premiere partie
% 1 Analyse de filtres numriques passe-bas du second ordre
% 1.1 Analyse dun filtre numrique RII passe-bas du second ordre
b = [0.079 2*0.079 0.079]; %Numrateur
a = [1 -1.2 0.516]; %Dnominateur
N = 32; n=0:N-1;
delta = [1; zeros(N-1,1)]; %Impulsion
step = ones(N,1); %Echelon unit
h = filter(b, a, delta); %Rponse du filtre (ou impz(b,a);)
figure(1); stem(n,h);
L = 256;
[h,w] = freqz(b,a,L);
m = abs(h); p = angle(h);
figure(2); plot(20*log10(m)); title(Log Magnitude Filtre RII 2nd ordre (echelle Matlab non modifiee));
figure(3); plot(w(1:L-1),20*log10(m(1:L-1))); title(Log Magnitude Filtre RII 2nd ordre);
axis([0 pi -60 2]); grid
figure(4); plot(w,m); title(Magnitude); figure(5); plot(w,p); title(Phase);
figure(6); zplane(b,a);
zero = roots(b); pole = roots(a);
Delta1 = 1; Delta2 = -20;
Fe = 10000;
Fp = 1000; NFp = round(L*Fp/(Fe/2));
Fs = 3000; NFs = round(L*Fs/(Fe/2));
gabh = [Delta1*ones(NFs,1); Delta2*ones(L-NFs,1)];
gabl = [-Delta1*ones(NFp,1); -5000*ones(L-NFp,1)];
figure(3); hold on; plot(w,gabh,r); plot(w,gabl,r);
B Analyse dun ltre numrique RIF passe-bas du second ordre
% 1.2 Analyse dun filtre numrique RIF passe-bas du second ordre
b11 = [0 -0.0309396 -0.0390182 0.0766059 0.288307 0.4 0.288307 0.0766059 -0.0390182 -0.0309396 0];
b7 = [-0.0409365 0.078369 0.289996 0.4 0.289996 0.078369 -0.0409365];
[h11,w] = freqz(b11,1,L);
[h7,w] = freqz(b7,1,L);
figure; plot(w,20*log10(abs(h11)),g); title(Log Magnitude Filtre RIF 2nd ordre);
axis([0 pi -60 2]); hold on;
plot(w,20*log10(abs(h7)),b); grid
plot(w,gabh,r); plot(w,gabl,r);
figure; plot(w,unwrap(angle(h11)),g); title(Phase Filtre RIF 2nd ordre);
grid hold on; plot(w,unwrap(angle(h7)),b);
h11 = filter(b11, 1, delta);
h7 = filter(b7, 1, delta);
figure; subplot(2,1,1); stem(n,h11,g); subplot(2,1,2); stem(n,h7,b)
figure; zplane(b11,1); figure; zplane(b7,1);
C Synthse dun ltre RII passe-bas
% 2 Synthse dun filtre passe-bas du second ordre
% 2.1 Synthse dun filtre RII
Fe = 8000; Fc = 1000; Fp = 700; Fs = 3000;
[N, Wn] = buttord(?,?,?,?);
disp(sprintf(Ordre du filtre de butterworth : %d \n,N));
[b a] = butter(?, ??,s); % Filtre analogique de Butterworth
[bb ab] = bilinear(b, a, Fe); % Transformation bilinaire sans prdistorsion
[bbp abp] = bilinear(b, a, Fe, Fc);% Transformation bilinaire avec prdistorsion en Fc
[bi ai] = impinvar(b, a, Fe); % Invariance impulsionnelle
L = 256; [hb,w] = freqz(bb,ab,L); figure;
plot(w,20*log10(abs(hb)),g); title(Log Magnitude Filtre RII Butterworth);
axis([0 pi -40 2]); grid
[hbp,w] = freqz(bbp,abp,L);
hold on; plot(w,20*log10(abs(hbp)),b);
[hi,w] = freqz(bi,ai,L);
plot(w,20*log10(abs(hi)),c);
D Synthse dun ltre RIF passe-bas
% 2.2 Synthse dun filtre RIF
% Rponse impulsionnelle thorique : h(n) = Wc/pi sinc(nWc/pi)
Fe = 8000; Fc = 1000; Wc = 2*pi*Fc/Fe;
% Mthode du fenetrage
N = 33; % Longueur du filtre RIF, doit etre impaire pour la suite
alpha = (N-1)/2; n=0:N-1;
b = ?????; %RI du filtre idal
bh = b.*hamming(N);
[h,w] = freqz(b,1,L); figure; plot(w,20*log10(abs(h)),b);
title(Log Magnitude Filtre RIF fenetrage); axis([0 pi -60 2]); grid
[hh,w] = freqz(bh,1,L);
hold on; plot(w,20*log10(abs(hh)),g);
figure; stem(n,b,b); hold on; stem(n,bh,g); grid
% Methode de lechantillonnage frequentiel avec une TFD sur N points
N = 33;
Hef = [ones(???,1); zeros(???,1); ones(???,1)] .* expression_de_la_phase;
%Mettre eventuellement un terme de phase
figure; subplot(2,1,1); stem(n,abs(Hef));
subplot(2,1,2); stem(n,angle(Hef));
bef = ifft(Hef,N);
[hef,w] = freqz(bef,1,L);
figure; subplot(2,1,1);
plot(w,20*log10(abs(hef)),b); title(Log Magnitude Filtre RIF echantillonage frequentiel);
axis([0 pi -40 2]); subplot(2,1,2); plot(w,unwrap(angle(hef)),b);
title(Phase Filtre RIF echantillonage frequentiel);
figure; stem(n,abs(bef),b); grid
E Filtrage numrique en prcision nie
b = [0.079 2*0.079 0.079]; %Numrateur
a = [1 -1.2 0.516]; %Dnominateur
L = 256;
[h,w] = freqz(b,a,L); figure; plot(w,20*log10(abs(h))); title(Log Magnitude Filtre RII 2nd ordre);
axis([0 pi -40 6]); grid
gabh = [Delta1*ones(NFs,1); Delta2*ones(L-NFs,1)];
gabl = [-Delta1*ones(NFp,1); -5000*ones(L-NFp,1)]; hold on;
% Arrondi sur b bits en Ca2 (1 bit de signe)
bit=8; %Nombre de bits
q=2^(1-bit);
bQ8= ???; %Quantification des coefficients
aQ8= ???; %Quantification des coefficients
[hQ8,w] = freqz(bQ8,aQ8,L); plot(w,20*log10(abs(hQ8)),g);
bit=4; %Nombre de bits
q=2^(1-bit); bQ4=???; aQ4= ???; [hQ4,w] = freqz(bQ4,aQ4,L);
plot(w,20*log10(abs(hQ4)),c);
plot(w,gabh,r); plot(w,gabl,r); legend(Prcision infinie,8 bits,4 bits);
figure; zplane(b,a); zero = roots(b); pole = roots(a); figure;
zplane(bQ8,aQ8); figure; zplane(bQ4,aQ4);
% 3.2 tude dun filtre IIR du quatrime ordre
% Forme directe
% a(1)*y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb)
% - a(2)*y(n-1) - ... - a(na+1)*y(n-na)
b = 0.0599*[1 -2.1832 2.9976 -2.1832 1];
a = [1 -3.1914 4.17 -2.5854 0.6443];
......... A COMPLETER
% Forme cascade
b1 = sqrt(0.0599)*[1 -0.6511 1]; b2 = sqrt(0.0599)*[1 -1.5321 1];
a1 = [1 -1.5684 0.6879]; a2 = [1 -1.623 0.9366 ];
......... A COMPLETER
% 3.3 Filtre RIF en prcision finie
b = [0 -0.0309396 -0.0390182 0.0766059 0.288307 0.4 0.288307 0.0766059 -0.0390182 -0.0309396 0];
......... A COMPLETER
F Convolution Discrete
close all; %Ferme les figures en cours
% II.4 Convolution
N=11; M=7;
L=N+M-1; n=0:L-1;
h = [0; 0.1; 0.2; 0.3; 0.2; 0.1; 0; zeros(L-M,1)];
x = [ones(N,1); zeros(L-N,1)];
y1 = conv(x, h);
figure; stem(n,h,r); hold on; stem(n,x,b); stem(n,y1(1:L),k);
......... A COMPLETER
G Analyse spectrale
% II.5 Analyse spectrale
N=64;
hr=boxcar(N);
figure; plot(hr,y.); hold on; grid on
ht=bartlett(N); plot(ht,r:)
hn=hanning(N); plot(hn,g--)
hm=hamming(N); plot(hm,m-.)
hb=blackman(N); plot(hb,b-)
xlabel(n);
ylabel(w(n)); title(Reponses temporelles de diffrentes fenetres);
legend(Rectangulaire,Bartlett,Hanning,Hamming,Blackman)
AXIS([0 N 0 1.1])
figure
subplot(2,3,1)
grid on
HR = fft(hr,1024);
plot((0:510)/1024,20*log10(abs(HR(1:511)/HR(1))))
AXIS([0 0.5 -80 0]);
ylabel(W dB)
Title(Fentre rectangulaire)
subplot(2,3,2)
grid on
HT = fft(ht,1024);
plot((0:510)/1024,20*log10(abs(HT(1:511)/HT(1))))
AXIS([0 0.5 -80 0]);
ylabel(W dB)
Title(Fentre de Bartlett)
subplot(2,3,3)
grid on
HN = fft(hn,1024);
plot((0:510)/1024,20*log10(abs(HN(1:511)/HN(1))))
AXIS([0 0.5 -80 0]);
ylabel(W dB)
Title(Fentre de Hanning)
subplot(2,3,4)
grid on
HM = fft(hm,1024);
plot((0:510)/1024,20*log10(abs(HM(1:511)/HM(1))))
AXIS([0 0.5 -80 0]);
ylabel(W dB)
Title(Fentre de Hamming)
subplot(2,3,5)
grid on
HB = fft(hb,1024);
plot((0:510)/1024,20*log10(abs(HB(1:511)/HB(1))))
AXIS([0 0.5 -80 0]);
ylabel(W dB)
Title(Fentre de Blackman)
subplot(2,3,6)
grid on
HR = fft(hr,1024);
plot((0:510)/1024,abs(HR(1:511)/HR(1)))
AXIS([0 0.25 0 1]);
ylabel(|W| lineaire)
Title(Fentre rectangulaire)
% Analyse dun sinus
N = 32; % nombre dechantillons
n = 0 : N-1; % vecteur des indices (nombres entiers)
a = 12; % amplitude de la sinusode
nu_o = 50; % frquence de la sinusode
nu_e = 300; % frquence dchantillonnage
f_o = nu_o/nu_e; % frquence normalise du signal
T_e = 1/nu_e; % priode dchantillonnage
nT_e = n*T_e; % vecteur des instants (en seconde)
y_n = cos(2*pi*f_o*n);
figure; subplot(3,1,1) ; stem(0:N-1+5,[y_n, zeros(1,5)]); axis([0 N-1+5 -1 1]);
title(Fonction cos(2pi f_o n) chantillonne tronque sur N=32 points); ylabel(x_N(n)); grid
fmin = -1/2; fmax = 1/2;
pas_f = (fmax-fmin)/1024 ; f=(fmin:pas_f:fmax-pas_f) ;
Y = fft(y_n,1024);
subplot(3,1,2) ; plot(f,abs(fftshift(Y)))
title(X_N(e^{j2pif})=TF\{x_N(n(\}) ; grid
Y = fft(y_n.*hamming(N),1024);
hold on; plot(f,abs(fftshift(Y)),r);
pas_f = (fmax-fmin)/N ; f=(fmin:pas_f:fmax-pas_f) ;
Y = fft(y_n);
subplot(3,1,3) ; stem(f,abs(fftshift(Y)))
title(X_N(k)=TFD\{x_N(n)\}) ; grid
Y = fft(y_n.*hamming(N));
hold on; stem(f,abs(fftshift(Y)),r)
......... A COMPLETER
H Dcimation et Interpolation
% II.6 Dcimation et interpolation
% 6.1 Dcimation
N = 192;
n = 0:1:N-1;
x = sin(2*pi*0.02*n)+sin(2*pi*0.09*n);
figure; subplot(4,1,1); stem(x); grid
subplot(4,1,2); plot(x)
y1 = x(1:2:N); subplot(4,1,2); stem(y1); grid
y2 = x(1:4:N); subplot(4,1,3); stem(y2); grid
y3 = x(1:6:N); subplot(4,1,4); stem(y3); grid
y1 = decimate(x,2); subplot(4,1,3); stem(y1)
y2 = decimate(x,4); subplot(4,1,4); stem(y2)
X = fft(x,300);
Y1 = fft(y1,300);
Y2 = fft(y2,300);
Y3 = fft(y3,300);
figure; subplot(4,1,1); plot(abs(X)); grid
subplot(4,1,2); plot(abs(Y1)); grid
subplot(4,1,3); plot(abs(Y2)); grid
subplot(4,1,4); plot(abs(Y3)); grid
% 6.2 Interpolation
L=4;
x = y2; %sinus dcim par 4
y = zeros(N, 1);
y(1:L:N) = x;
figure; stem(y); grid
X = fft(x,300);
Y = fft(y,300);
figure; subplot(2,1,1); plot(abs(X)); grid
subplot(2,1,2); plot(abs(Y)); grid
Nrif = 33; % Longueur du filtre RIF, doit etre impair pour la suite
...... A COMPLETER
y1 = filter(.....);
figure; plot(y1); grid
Y1 = fft(y1,300);
figure; plot(abs(Y1)); grid
% Interpolation par TFD
..... A COMPLETER
% 6.3 Filtrage
%filtre sans int/dec
Nrif = 21; % Longueur du filtre RIF, doit etre impair pour la suite
Wc = pi/2;
b = Wc/pi*sinc(Wc/pi*(-(Nrif-1)/2:(Nrif-1)/2));
bh = b.*hamming(Nrif);
N=256;
[h,t]=impz(bh,1,N);
H=fft(h,N);
figure; plot(abs(H)); grid
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%filtre quivalent avec int/dec
%%%% Interpolation
.......
%%%% Filtre
.......
%%%% Decimation
......