[Débutant] Démodulation de fréquence

Bonjour à tous et à toutes,
Après plusieurs tentatives, j'ai pas pu trouver la solution. Voici mon problème : j'ai un signal yb (bruité) modulé en amplitude et en fréquence que j'ai construit pour lui appliquer la démodulation en amplitude et en fréquence. La démod en amplitude marche bien avec l'enveloppe mais mon problème et dans la démod de fréquence que j'arrive pas à faire (fonction de modulation de fréqunce). Dans le signal il y trois fréquence : f1 =2000 hz (la porteuse), f2=100 (fréq mod d'amplitude) , f3=50 Hz (fréquence de mod de freq). il faut que je visualise cette dernière "f3" en analyse spectrale de la démod de fréquence c-à-d un seul pic de 50 hz !

voici mon programme

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
% Modulation  et démodulation d'amplitude et de fréquence-------------------------
clear all;close all;
A=1;  % amplitude maxi du signal (de la porteuse)
Ba=0.5  % amplitude de la modulante
Bf=0.5; %indice de modulation de fréquence
f1=2000; % fréquence de la porteuse
f2=100;  % fréquence de la modulante d'amplitude
f3=40;  % fréquence de la modulante de fréquence
Fs=44100;dt=1/Fs;  % Fréquence d'échantillonnage
tf=.5; % temps d'acquisition
cbf=2;  % amplitude maxi du bruit de fond (aléatoire)
t=0:dt:tf;
Nt=length(t);
y=A*(Ba*sin(2*pi*f2*t)+1).*sin(2*pi*f1*t+Bf*sin(2*pi*f3*t));
figure;plot(t(1:Nt/3),y(1:Nt/3),'k');
axis([0 t(round(Nt/3+5)) min(y)-0.2 max(y)+0.2]);
xlabel('temps en secondes');title('signal d''origine y(t)');
%---------signal bruité yb(t) et son spectre------------------------
bf=randn(Nt,1);sM2=max(bf);sm2=min(bf);
maxi=max(sM2,abs(sm2));bf=cbf*bf'/maxi;
yb=y+bf;
figure;
plot(t(1:Nt/5),yb(1:Nt/5),'k');
xlabel('temps en secondes');title('signal d''origine y(t) bruité');
[asp, freqs] = autosp(Fs, 'win', 'han', 'duree', tf,'format', 'PSD',yb);
asp=sqrt(asp*(freqs(2)-freqs(1)))*sqrt(2);
df=1/tf;
fmax=max(f1,f2)*2;Nmax=round(fmax/df);
figure;
plot(freqs(1:Nmax),asp(1:Nmax),'k')
xlabel('fréquence en Hz');title('Spectre du signal bruité');
axis([0 freqs(Nmax+5) 0 1.2*max(asp)]);
%-------------calcul des enveloppes------------------------------
%--------------modulation de fréquence---------------------------
z=hilbert(yb);
% fmf=f1+Bf*f3*cos(2*pi*f3*t);


 fmf = demod(yb,f1,Fs,'fm');
 

 
figure;plot(t(1:Nt/2),fmf(1:Nt/2),'k');
xlabel('temps en secondes');title('Fonction de modulation de fréquence');
grid on;
[asp, freqs] = autosp(Fs, 'win', 'han', 'duree', tf,'format', 'PSD',fmf);
asp=sqrt(asp*(freqs(2)-freqs(1)))*sqrt(2);
df=1/tf;
fmax=max(f1,f2)*2;Nmax=round(fmax/df);
figure;
plot(freqs(3:Nmax),asp(3:Nmax),'k')
xlabel('Fréquence en Hz');title('Spectre de la fonction de modulation de fréquence');
%--------------modulation d'amplitude---------------------------
Env=abs(z);
figure;plot(t(1:Nt/5),Env(1:Nt/5),'k');
xlabel('temps en secondes');title('Enveloppe du signal bruité');
[asp, freqs] = autosp(Fs, 'win', 'han', 'duree', tf,'format', 'PSD',Env);
asp=sqrt(asp*(freqs(2)-freqs(1)))*sqrt(2);
df=1/tf;
fmax=max(f1,f2)*2;Nmax=round(fmax/df);
figure;
plot(freqs(3:Nmax),asp(3:Nmax),'k')
xlabel('Fréquence en Hz');title('Spectre de l''enveloppe');

je vous remercie !