Transformées de Fourier

La transformée de Fourier est une formule mathématique qui transforme un signal échantillonné en temps ou en espace, en un même signal échantillonné en fréquence temporelle ou spatiale. En traitement du signal, la transformée de Fourier peut révéler d’importantes caractéristiques d’un signal, notamment ses composants de fréquence.

La transformée de Fourier est définie pour un vecteur $$x$$ avec $$n$$ points échantillonnés uniformément par

$$y_{k+1}=\sum _{j=0}^{n-1}{\omega }^{jk}{x}_{j+1}.$$

$$\omega =e^{-2\pi i/n}$$ est l’une des $$n$$ racines complexes de l’unité où $$i$$ est l’unité imaginaire. Pour $$x$$ et $$y$$, les indices $$j$$ et $$k$$ vont de $$0$$ à $$n-1$$.

La fonction fft de MATLAB® utilise un algorithme de transformée de Fourier rapide qui calcule la transformée de Fourier des données. Prenez un signal sinusoïdal x qui est une fonction du temps t avec des composants de fréquences de 15 et 20 Hz. Utilisez un vecteur temporel échantillonné en incréments de 1/50 seconde sur une période de 10 secondes.

Ts = 1/50;
t = 0:Ts:10-Ts;
x = sin(2*pi*15*t) + sin(2*pi*20*t);
plot(t,x)
xlabel('Time (seconds)')
ylabel('Amplitude')

Calculez la transformée de Fourier du signal et créez le vecteur f correspondant à l’échantillonnage du signal dans l’espace fréquentiel.

y = fft(x);
fs = 1/Ts;
f = (0:length(y)-1)*fs/length(y);

Tracez l’amplitude du signal transformé en tant que fonction de la fréquence.

plot(f,abs(y))
xlabel('Frequency (Hz)')
ylabel('Magnitude')
title('Magnitude')

Le tracé montre quatre pics fréquentiels, alors que le signal devrait avoir deux pics de fréquence, respectivement à 15 et 20 Hz. Ici, la deuxième moitié du tracé est le reflet miroir de la première moitié. La transformée de Fourier discrète d’un signal dans le domaine temporel est de nature périodique, où la première moitié de son spectre se trouve dans les fréquences positives et la deuxième moitié dans les fréquences négatives. Les composants de fréquence 30 et 35 Hz du tracé correspondent aux composants de fréquence –20 et –15 Hz. Pour mieux visualiser cette périodicité, vous pouvez utiliser la fonction fftshift, qui effectue un décalage circulaire centré sur zéro sur la transformée.

n = length(x);
fshift = (-n/2:n/2-1)*(fs/n);
yshift = fftshift(y);
plot(fshift,abs(yshift))
xlabel('Frequency (Hz)')
ylabel('Magnitude')

Signaux bruités

Dans les applications scientifiques, les signaux sont souvent corrompus par du bruit aléatoire, ce qui masque leurs composants fréquentiels. La transformée de Fourier peut retirer ce bruit aléatoire et révéler les fréquences. Par exemple, créez un nouveau signal, xnoise, en injectant du bruit gaussien dans le signal original, x.

rng('default')
xnoise = x + 2.5*randn(size(t));

La puissance du signal comme fonction de la fréquence est une métrique courante utilisée en traitement du signal. La puissance est l’amplitude au carré de la transformée de Fourier d’un signal, normalisée par le nombre d’échantillons de fréquence. Calculez et tracez le spectre de puissance du signal bruité centré à la fréquence zéro. Malgré le bruit, vous pouvez quand même distinguer les fréquences du signal grâce aux pics de puissance.

ynoise = fft(xnoise);
ynoiseshift = fftshift(ynoise);
power = abs(ynoiseshift).^2/n; 
plot(fshift,power)
title('Power')
xlabel('Frequency (Hz)')
ylabel('Power')

Efficacité computationnelle

Utiliser la formule de la transformée de Fourier directement pour calculer chacun des $$n$$ éléments de $$y$$ nécessite de l’ordre de $$n^2$$ opérations à virgule flottante. L’algorithme de transformée de Fourier rapide ne nécessite que de l’ordre de $\mathit{n}\log \mathit{n}$ opérations pour calculer. Cette efficacité computationnelle est un grand avantage pour le traitement des données comportant des millions de points de données. De nombreuses implémentations spécialisées de l’algorithme de transformée de Fourier rapide sont encore plus efficaces quand $$n$$ comporte de faibles facteurs premiers, par exemple lorsque $$n$$ est une puissance de 2.

Prenez des données audio collectées depuis des microphones sous-marins au large de la Californie. On peut trouver ces données dans une bibliothèque gérée par le programme de recherche en bioacoustique de l’Université de Cornell. Chargez et formatez un sous-ensemble des données de bluewhale.au, contenant un chant de baleine bleue du Pacifique. Comme les cris des baleines bleues sont des sons de basses fréquences, les humains peuvent à peine les entendre. L’échelle temporelle des données est compressée d’un facteur 10 pour augmenter la hauteur et rendre le cri plus facile à entendre. Vous pouvez utiliser la commande sound(x,fs) pour écouter le fichier audio en entier.

whaleFile = 'bluewhale.au';
[x,fs] = audioread(whaleFile);
whaleMoan = x(2.45e4:3.10e4);
t = 10*(0:1/fs:(length(whaleMoan)-1)/fs);

plot(t,whaleMoan)
xlabel('Time (seconds)')
ylabel('Amplitude')
xlim([0 t(end)])

Spécifiez une nouvelle longueur de signal étant la prochaine puissance de 2 supérieure à la longueur originale. Ensuite, utilisez fft pour calculer la transformée de Fourier avec la nouvelle longueur de signal. fft ajoute automatiquement des zéros aux données pour augmenter la taille d’échantillon. Ce remplissage peut rendre le calcul de la transformée beaucoup plus rapide, notamment pour les tailles d’échantillon ayant de grands facteurs premiers.

m = length(whaleMoan);
n = pow2(nextpow2(m));
y = fft(whaleMoan,n);

Tracez le spectre de puissance du signal. Le tracé indique que le gémissement est composé d’une fréquence fondamentale autour de 17 Hz et d’une séquence d’harmoniques, où le deuxième harmonique est accentué.

f = (0:n-1)*(fs/n)/10; % frequency vector
power = abs(y).^2/n;   % power spectrum
plot(f(1:floor(n/2)),power(1:floor(n/2)))
xlabel('Frequency (Hz)')
ylabel('Power')

Phase de sinusoïdes

Avec la transformée de Fourier, vous pouvez également extraire le spectre de phase du signal original. Par exemple, créez un signal constitué de deux sinusoïdes dont les fréquences sont de 15 et 40 Hz. La première sinusoïde est une onde cosinusoïdale de phase $-\pi /4$, et la deuxième est une onde cosinusoïdale de phase $\pi /2$. Échantillonnez le signal à 100 Hz pendant 1 seconde.

fs = 100;
t = 0:1/fs:1-1/fs;
x = cos(2*pi*15*t - pi/4) - sin(2*pi*40*t);

Calculez la transformée de Fourier du signal. Tracez l’amplitude de la transformée en tant que fonction de la fréquence.

y = fft(x);
z = fftshift(y);

ly = length(y);
f = (-ly/2:ly/2-1)/ly*fs;

stem(f,abs(z))
xlabel("Frequency (Hz)")
ylabel("|y|")
grid

Calculez la phase de la transformée en supprimant les valeurs de faible amplitude. Tracez la phase comme une fonction de la fréquence.

tol = 1e-6;
z(abs(z) < tol) = 0;

theta = angle(z);

stem(f,theta/pi)
xlabel("Frequency (Hz)")
ylabel("Phase / \pi")
grid

Voir aussi

fft | fftshift | nextpow2 | ifft | fft2 | fftn | fftw

Rubriques

• 2-D Fourier Transforms

Sites web externes

• Analyse de Fourier (Ressources pour l’enseignement MathWorks)