Transformée de Fourier

1. Transformée de Fourier

Soit un signal u(t). Sa transformée de Fourier(TF) est :

S (f) = \int_{- \infty}^{\infty} u (t) exp (- j 2 π f t) d t

Si u(t) est réel, sa transformée de Fourier possède la parité suivante :

S (- f) = S {(f)}^{*}

Le signal s'exprime avec sa TF par la transformée de Fourier inverse :

u (t) = \int_{- \infty}^{\infty} S (f) exp (j 2 π f t) d f

Lors du traitement numérique d'un signal, on dispose de u(t) sur une durée T, par exemple sur l'intervalle [-T/2,T/2]. Une approximation de la TF est calculée sous la forme :

S_{a} (f) = \int_{- T / 2}^{T / 2} u (t) exp (- j 2 π f t) d t

Soit un échantillonnage de N points, obtenu pour $0 \leq k \leq N - 1$ :

\begin{matrix} t_{k} = - \frac{T}{2} + k \frac{T}{N} \\ u_{k} = u (t_{k}) \end{matrix}

Une approximation est obtenue par la méthode des rectangles :

S_{a} (f) ≃ \frac{T}{N} exp (j π f T) \sum_{k = 0}^{N - 1} u_{k} exp (- j 2 π k f \frac{T}{N})

On recherche la TF pour les fréquences suivantes, avec $0 \leq n \leq N - 1$ :

f_{n} = \frac{n}{T}

c'est-à-dire :

S_{a} (f_{n}) = T exp (j π n) \frac{1}{N} \sum_{k = 0}^{N - 1} u_{k} exp (- j 2 π k \frac{n}{N})

En notant S_n la transformée de Fourier discrète (TFD) de u_k, on a donc :

S_{a} (f_{n}) ≃ T exp (j π n) S_{n}

Dans une analyse spectrale, on s'intéresse généralement au module de S(f), ce qui permet d'ignorer le terme exp(jπ n)

Le spectre obtenu est par nature discret, avec des raies espacées de 1/T. C'est donc le spectre d'un signal périodique de période T. Pour simuler un spectre continu, T devra être choisi assez grand.

Le spectre obtenu est périodique, de périodicité fe=N/T, la fréquence d'échantillonnage.

2. Signal à support borné

2.a. Exemple : gaussienne

On choisit T tel que u(t)=0 pour |t|>T/2. Considérons par exemple une gaussienne centrée en t=0 :

u (t) = exp (- \frac{t^{2}}{a^{2}})

dont la transformée de Fourier est

S (f) = a \sqrt{π} exp (- π^{2} a^{2} f^{2})

En choisissant par exemple T=10a, on a $| u (t) | < 1 0^{- 10}$ pour t>T/2

a=1;
signal[t_]:=Exp[-t^2/a^2];

Plot[signal[t],{t,-5,5},PlotRange->{{-5,5},{0,1}},AxesLabel->{"t","s"}]

plotA.pdf

La fonction suivante calcule la TFD pour une durée T et un nombre points N

spectre[T_,npoints_]:=Module[{te,temps,echantillons,transFourier,frequences,spectre},
	te=T/npoints;
	temps[k_]:=-T/2+(k-1)*te;
	echantillons:=Table[signal[temps[k]],{k,npoints}];
	transFourier=Fourier[echantillons,FourierParameters->{-1,-1}]*T;
	frequences[k_]:=(k-1)/T;
	spectre=Table[{frequences[k],Abs[transFourier[[k]]]},{k,npoints}];
	Return[spectre];
]

Pour comparaison, on tracera aussi la TF théorique de la gaussienne :

tfSignal[f_]:=a*Sqrt[Pi]*Exp[-(Pi*a*f)^2];

Exemple :

T=20; npoints=100;
sp=spectre[T,npoints];
p1=ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{-npoints/T,npoints/T},{0,2}},
									PlotStyle->{RGBColor[1,0,0]}];
p2=Plot[tfSignal[f],{f,-npoints/T,npoints/T},PlotRange->{{-npoints/T,npoints/T},{0,2}}];

Show[p1,p2]

plotB.pdf

Le spectre obtenu doit être prolongé par périodicité de période N/T (5 sur l'exemple). Pour augmenter la résolution, il faut augmenter T. Il est intéressant de maintenir constante la fréquence d'échantillonnage N/T :

T=100; npoints=500;
sp=spectre[T,npoints];
p1=ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{-npoints/T,npoints/T},{0,2}},
									PlotStyle->{RGBColor[1,0,0]}];
p2=Plot[tfSignal[f],{f,-npoints/T,npoints/T},PlotRange->{{-npoints/T,npoints/T},{0,2}}];

Show[p1,p2]

plotC.pdf

La fréquence d'échantillonnage doit vérifier le critère de Shannon : N/T doit être supérieur à 2 fois la plus grande fréquence présente dans le spectre. Voici un exemple où ce critère n'est pas vérifié :

T=100; npoints=100;
sp=spectre[T,npoints];
p1=ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{-npoints/T,npoints/T},{0,2}},
									PlotStyle->{RGBColor[1,0,0]}];
p2=Plot[tfSignal[f],{f,-npoints/T,npoints/T},PlotRange->{{-npoints/T,npoints/T},{0,2}}];

Show[p1,p2]

plotD.pdf

2.b. Exemple : sinusoïde modulée par une gaussienne

On considère le signal suivant :

u (t) = exp (- t^{2} / a^{2}) cos (2 π \frac{t}{b})

avec $b ≪ a$ .

a=1;
b=0.1; 
signal[t_]:=Exp[-t^2/a^2]*Cos[2*Pi*t/b];

Plot[signal[t],{t,-5,5},PlotRange->{{-2,2},{-1,1}},
                        AxesLabel->{"t","s"}]

plotE.pdf

Dans ce cas, il faut choisir une fréquence d'échantillonnage supérieure à 2 fois la fréquence de la sinusoïde, c.a.d. N/T>2/b.

T=100; npoints=3000;
sp=spectre[T,npoints];

ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{0,npoints/T},{0,2}},
                        PlotStyle->{RGBColor[1,0,0]}]

plotF.pdf

2.c. Fenêtre rectangulaire

Soit une fenêtre rectangulaire de largeur a :

a=1;
signal[t_]:=If[Abs[t]<(a/2),1,0];

Voyons son spectre avec T=100 et N=1000 :

T=100; npoints=3000;
sp=spectre[T,npoints];

ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{0,npoints/T},{0,2}},
                        PlotStyle->{RGBColor[1,0,0]}]

plotG.pdf

3. Signal à support non borné

Dans ce cas, la fenêtre [-T/2,T/2] est arbitrairement imposée par le système de mesure. Par exemple sur un oscilloscope numérique, T peut être ajusté par le réglage de la base de temps. Considérons par exemple un signal périodique comportant 3 harmoniques :

b=1; (*periode*)
w0=2*Pi/b;
signal[t_]:=Cos[w0*t]+0.5*Cos[2*w0*t]+0.1*Cos[3*w0*t];

La fréquence d'échantillonnage doit vérifier N/T>6/b.

T=100; npoints=1000;
sp=spectre[T,npoints];

ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{0,npoints/T},{0,T}},
                    PlotStyle->{RGBColor[1,0,0]}]

plotH.pdf

On reconnait dans le spectre les raies correspondant aux 3 harmoniques. Les coefficients de Fourier se retrouvent en divisant par T. Modifions la fréquence du signal :

b=0.9;
w0=2*Pi/b;
signal[t_]:=Cos[w0*t]+0.5*Cos[2*w0*t]+0.1*Cos[3*w0*t];
T=100; npoints=1000;
sp=spectre[T,npoints];

ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{0,npoints/T},{0,T}},PlotStyle->{RGBColor[1,0,0]}]

plotI.pdf

La base des raies est notablement plus large. Le signal échantillonné est en fait le produit du signal périodique défini ci-dessus par une fenêtre h(t) rectangulaire de largeur T. La TF est donc le produit de convolution de S avec la TF de h :

H (f) = T \frac{sin (π T f)}{π T f}

qui présente des oscillations lentement décroissantes qui peuvent conduire à des interférences entre les différentes harmoniques du signal. Pour remedier à ce problème, on remplace la fenêtre rectangulaire par une fenêtre dont le spectre présente des lobes secondaires plus faibles, par exemple la fenêtre de Hamming :

b=0.9;
w0=2*Pi/b;
hamming[t_]:=0.54+0.46*Cos[2*Pi*t/T];
signal[t_]:=(Cos[w0*t]+0.5*Cos[2*w0*t]+0.1*Cos[3*w0*t])*hamming[t];
T=100; npoints=1000;
sp=spectre[T,npoints];

ListPlot[sp,AxesLabel->{"f","S"},PlotRange->{{0,npoints/T},{0,T}},
                    PlotStyle->{RGBColor[1,0,0]}]

plotJ.pdf

Transformée de Fourier discrète : transformée de Fourier