Dans ce chapitre, vous allez apprendre à tracer le spectre de données expérimentales.
Les données expérimentales correspondent ici au parcours dans un plan d'une petite bille chargée dans un champ magnétique.
La position (x et y) de cette bille a été relevée toutes les millisecondes pendant 4 secondes.
Acquisition et visualisation des données
Vous devriez maintenant être habitué à la procédure : importez les données expérimentales Trajectoire.dat.
%Importation des données, il y a trois colonnes
data = load('Trajectoire.dat')
t = data(:,1);% la première colonne est une colonne qui correspond aux instants d'échantillonnage
x = data(:,2); %La deuxième est une position horizontale
y = data(:,3); %la troisième est une position verticale
Afin de se faire une petite idée de ce que vous avez comme donnée, visualisons tout cela :
%représentation des données
figure(1)
plot(x(t<1/10),y(t<1/10)) %La condition t<10 limite la quantité de points tracés
%partie cosmétique
title('Trajectoire en 1/10e de s')
xlabel('x')
ylabel('y')
% On utilise une figure différente pour tracer la trajectoire car les axes sont différents
figure(2)
plot(t,x,'k')
hold on
plot(t,y,'r')
hold off
legend({'x(t)' 'y(t)'})
Vous devriez obtenir :

Vous obtenez quelque chose de semblable à une particule qui suit une trajectoire aléatoire. Mais si vous regardez de plus près, juste ce qui se passe sur x et sur y :

Vous pouvez donc en conclure que la particule oscille à la fois horizontalement (x(t) est périodique) et verticalement (y(t) aussi). Il semble par ailleurs qu'il y ait des oscillations rapides et des oscillations lentes, ce que vous pouvez voir en zoomant un peu :

Donc, vous pouvez vous attendre à ce que la représentation fréquentielle traduise cela. C'est-à-dire que le signal ait des composantes fortes sur les vecteurs de Fourier dont la fréquence est la même que celle que nous devinons sur le signal.
Construction manuelle des premiers vecteurs de Fourier
Il vous faut donc construire les vecteurs de Fourier, afin de pouvoir projeter les données dessus. Ce qui nous permettra ainsi d'obtenir le spectre.
%nous construisons un premier vecteur de Fourier.
N = length(x);% N est la longueur du vecteur de données expérimentales
k = [0:N-1]';% k est un indice variant de 0 à N-1, ici nous construisons un vecteur colonne
F0 = exp(1i*2*pi*0/N*k);
Une fois le vecteur construit, vous pouvez vérifier qu'il s'agit bien d'une suite de 1, donc un vecteur constant (car l'exponentielle complexe de fréquence 0 vaut toujours 1).
Puis, vous pouvez projeter les données expérimentales sur ces vecteurs.
%Premier terme de la Transformée de Fourier/Projection de x sur le premier vecteur de Fourier
tfx(1) = sum(x .* conj(F0))% On utilise conj pour conjugé.
%La même chose pour y:
tfy(1) = sum(y .* conj(F0))% On utilise conj pour conjugé.
La console répond :
tfx = -23.168 tfy = 8006.1
Ce que nous pouvons interpréter comme : x a une faible composante constante négative, et y qui a une très forte composante constante positive.
Notez que si vous divisez ces valeurs par N, vous obtiendrez la valeur moyenne (vous pouvez vérifier que y est bien de valeur moyenne - ou de composante constante - 2).
>> tfy/N ans = 2.0010
Pourquoi faut-il diviser par N pour avoir une grandeur plus facilement interprétable ?
Tout simplement parce que les vecteurs de Fourier ne sont pas normés. Je vous avais dit que ça poserait problème !
Calcul de toutes les composantes de fréquentielles
Bon, eh bien il ne vous reste plus qu'à calculer toutes les composantes fréquentielles possibles.
Pour cela, il faut projeter les données sur tous les vecteurs de Fourier :
clear tfx tfy freq% permet le nettoyage des tableaux tfx, tfy et freq
for n = 0:N-1
Fn = exp(1i*2*pi*n/N*k); %Création du n-ieme vecteur de Fourier qu'on appelle Fn
p = n+1; %position dans le tableau, les indices matlab commencent à 1 et pas à 0
freq(p) = n/N;%fréquence de Fourier associée
tfx(p) = dot(x,Fn); % La même chose que sum(x .* conj(Fn)). matlab a une fonction toute prète qui fait ça. ("dot product" en anlgais signigie produit scalaire)
tfy(p) = dot(y,Fn);
end
tfx(1:3)'%affiche les trois premiers termes le ' est pour la lisibilité
La console devrait vous donner les trois premiers termes :
ans = -23.1676 - 0.0000i 3.6870 + 17.7100i 6.0581 - 8.2493i
Nous retrouvons bien le -23 trouvé pour le premier terme de la TF de x, mais qu'en est-il des suivants ?
Ce sont des nombres complexes !! Un peu une surprise, mais pas tant que cela, puisque nous utilisons des vecteurs complexes, les composantes peuvent l'être également.
Pour l'instant, nous nous préoccupons juste du module de ces composantes, module qui s'obtient avecabs
:
abs(tfx(1:3)') ans = 23.168 18.090 10.235
Visualisation du spectre
Pour finir, vous pouvez visualiser ces composantes. Comme il s'agit de composantes complexes, nous allons simplement en visualiser le module (le spectre).
%Tracé du spectre
figure(3)
Te = 1e-3;
stem(freq/Te,abs(tfx)/N);%stem est un type de plot souvent utilisé pour représenter les spectres
xlabel('Fréquence en Hz')
ylabel('Composantes de Fourier (normalisées)')
Vous devriez avoir :

Il se passe clairement des choses bizarres sur les très hautes fréquences (cela est dû au fait que le signal soit échantillonné et réel). Ne vous en préoccupez pas pour le moment et zoomez sur la zone 0-120 Hz :
xlim([0 120]) %limite la zone de traçage du graph

Vous remarquez donc que le signal contient une forte composante à la fréquence 18 Hz, mais également une composante vers 83 Hz. Vous pouvez vérifier sur le signal lui-même qu'il s'agit bien des fréquences que nous avions observées sur le signal lui-même.
Si vous faites de même avec le spectre de y, vous obtiendrez :

Que pouvons-nous dire ?
Le signal contient une composante de fréquence nulle (composante continue).
Il y a trois composantes oscillantes :
une à 20 Hz ;
une à 60 Hz (de très faible amplitude) ;
et une à 115 Hz.
Autant les composantes à 20 et 115 Hz étaient à peu près visibles dans le signal de départ (et encore, ce n'était pas évident), autant celle à 60 était bien invisible.
Voilà pourquoi la représentation fréquentielle est parfois déterminante, elle permet de voir des choses pas forcément visibles dans la représentation temporelle.
Vous pouvez vérifier que vous obtenez exactement les mêmes résultats en utilisant la fonction toute faite de Matlab qui calcule les transformées de Fourier : fft
.
tfx_matlab = fft(x);% calcul de la TF via la fonction built-in de matlab
plot(freq/Te,abs(tfx))
hold on
plot(freq/Te,abs(tfx_matlab) + 100,'r') % le +100 est là pour décaler les courbes qui sont sinon parfaitement superposées
hold off
