Estimation spectrale de quoi ?

Il y a deux situattions totalement différentes dans laquelle on pourrait vouloir obtenir un spectre de signal: la premières est une estimation de la transformation de Fourier d'un morceau de signal continu dans le temps, comme une fonction déterministe du temps.  La deuxième nous concerne principalement: c'est l'estimation de la densité spectrale d'un ensemble de fonctions continues du temps.  Une erreur basique est de confondre ces deux problèmes.  Cela dit, avant de pouvoir comprendre l'estimation de la densité spectrale d'un ensemble, il faut déjà comprendre les problèmes qui se posent quand on veut estimer la transformation de Fourier d'un signal donné.  Cette contribution s'occupe donc de l'estimation d'une transformée de Fourier d'un signal déterministe.

Dans tous les cas, l'estimation spectrale est faite à partir de "morceaux de signal" limités dans le temps, en prenant un nombre fini d'échantillons à une fréquence finie d'échantillonnage.   Cela introduira des effets dans le spectre estimé: repliement (aliasing) et fenêtrage (windowing).

Repliement

Le repliement est une conséquence de la fréquence finie d’échantillonnage, comme décrit par le théorème de Nyquist: une fréquence d'échantillionnage de 2 fs divises le spectre en bandes successives de largeur fs et projette ces bandes l'une sur l'autre dans le sens direct ou inversé.  Le plus souvent, le théorème de Nyquist est formulé en disant que la fréquence maximale dans un signal échantillonné est fs.  En fait, ce théorème rend équlvalent toutes les fréquences qui diffèrent de 2 fs ou qui diffèrent de 2 fs si une des deux est rendu négative.

Si la fréquence d'échantillonnage est 1 KHz, alors un signal qui est un sinus de 10 Hz sera équivalent à un signal de 1010 Hz et à un signal de 990 Hz.  Un sinus de 490 Hz sera équivalent à un sinus de 510 Hz.  La bande spectrale de 500 Hz à 1000 Hz est projetée sur la bande de 0 à 500 Hz, mais dans le sens inverse: ce qui est proche de 500 Hz, restera proche de 500 Hz, et ce qui est près de 1 KHz sera proche de 0 Hz.  La bande de 1 KHz à 1.5 KHz est projetée sur la bande 0 - 500 Hz, dans le même sens: ce qui est proche de 1000 Hz sera proche de 0 Hz, et ce qui est proche de 1500 Hz sera proche de 500 Hz.  Et ainsi de suite. De façon plus formelle: les fonctions du temps suivants:

sin(2 π 10 t), - sin(2π 990 t), sin(2 π 1010 t) , ...  

ont les mêmes échantillons si on les prend tous les 1 ms.

Il faut ce rendre compte que cette équivalence des trains d'échantillons est indépendant du traitement ultérieur de ces échantillons, comme une analyse spectrale, du filtrage ou autre: c'est la nature de l'acte d'échantillonnage qui induit cette équivalence.

Si une fréquence de 10 Hz et une de 990 Hz ne peuvent pas être distinguées, cela implique que si nous échantillons un signal qui contient ces deux composantes, il ne sera plus possible, à partir des échantillons seul, de dire quelle partie vient des 10 Hz et quelle partie vient des 990 Hz: les échantillons des deux composantes seront additionnés.  Il ne sera donc pas possible de reconstituer le signal continu à partir de ces échantillons.  Par contre, le théorème de Nyquist nous apprend aussi que si nous savons que le signal continu était garanti de n'avoir des composantes dans la bande 0 - 500 Hz, qu'il n'y a que un seul signal continu qui correspond aux échantillons: ainsi, le signal continu est donc parfaitement déterminé et unique dans ce cas, et nous pouvons le reconstruire à partir des échantillons.  Dans la plupart des cas, nous nous limitons à la bande de base, ici 0 - 500 Hz.  Mais cela n'est pas nécessaire.  Si nous avons la garantie que le signal continu est limité à une seule bande de n. fs to (n+1) fs et nous connaissons cette bande, alors le signal peut être reconstruit à partir des échantillons, et est unique.

Par exemple, si nous avons un signal avec un contenu fréquentiel autour de 400 MHz, avec une largeur de bande de 100 KHz, de 400 MHz à 400.1 MHz, en principe nous pouvons échantillonner ce signal à 200 KHz, car le signal continu tombe dans la bande de 4000-ième harmonique.  C'est une bonne nouvelle: un flux de données à 200 KHz est facilement maîtrisable par une électronique modeste, tandis que un flux à 800.2 MHz (ce qui serait nécessaire si nous voulions garder le signal dans la bande de base, poserait une difficulté sauf pour des systèmes de traitement très puissants.

Mais il y a bien sûr une difficulté pratique: si nous voulons utiliser un échantillonneur de 200 KHz pour notre signal de 400 MHz, il faut que chaque échantillon soit pris avec suffisamment de précision (en dessous de la nanoseconde dans ce cas).  En fait, il faut que la prise d'échantillons soit aussi précise que si nous échantillons à 800 MHz, sauf que nous n'en prenons pas autant (un échantillon sur 4000).  Certains échantillonneurs sont effectivement capables d'échantillonner avec une précision bien plus grande que la fréquence d'échantillonnage, mais d'autres, pas.  Beaucoup d'échantillonneurs font l'hypothèse qu'ils seront utilisés dans leur bande de base, ou dans une bande harmonique de faible ordre.

Filtre anti-repliement

Le repliement est un phénomène intrinsèque à l'opération d'échantillonnage, indépendamment de ce que nous allons faire de ces échantillons.  Le repliement n'est pas une conséquence d'une estimation spectrale ou autre. C'est un phénomène universel de confusion entre différentes composantes fréquentielles dans les signaux continus produisant des échantillons identiques.  Si nous voulons éviter cette confusion, un filtre analogique anti-repliement est absolument nécessaire, sauf dans le cas où nous pouvons garantir que seulement les composantes dans la bande permise seront présentes - y compris le bruit.  Comme, surtout pour le bruit, cela est quasiment jamais le cas, il faut toujours un filtre de repliement.

La plupart du temps, quand on parle de filtre anti-repliement, on pense à un filtre passe-bas qui coupe à une fréquence fs, mais comme nous avons vu, toute bande harmonique de n.fs à (n+1).fs est possible.  Dans ce cas, le filtre anti-repliement est un filtre passe-bande avec ces fréquences comme fréquences de coupure.

Comme dans la pratique, il est difficile de construire un filtre analogique hautement sélectif, il vaut mieux utiliser une fréquence d' échantillonnage un peu plus élevée que nécessaire et construire un filtre passe-bas ou passe-bande auquel on peut donner un peu plus d'espace de transition entre la bande passante utile, et la bande passante numérique.

Par exemple, si nous avons un signal utile dans la bande 0 - 30 KHz, alors, strictement parlé, une fréquence d'échantillonnage de 60 KHz est suffisante.  Mais cela implique un filtre passe-bas parfait, qui laisse tout passer jusqu'à 30 KHz, et coupe parfaitement à 30.1 KHz, car la bande 30 KHz - 60 KHz sera projetée dans la bande 0 - 30 KHz.  Du bruit à 31 KHz apparaîtra comme un bruit à 29 KHz, dans notre partie utile, si le filtre ne coupe pas sévèrement à 31 KHz.  Il est donc bien plus judicieux d' échantillonner à 100 KHz, ce qui met la bande passante numérique à 50 KHz.  Maintenant, nous avons besoin d'un filtre qui laisse passer tout jusqu'à 30 KHz, et qui coupe fortement à 50 KHz.  Notre bruit à 31 KHz sera vu comme un signal à 31 KHz et peut donc être distingué numériquement du signal utile à 29 KHz.  Un bruit à 51 KHz se manifestera comme un signal à 49 KHz par repliement.  Il faudra un bruit à 71 KHz pour qu'il sera confondu avec un signal de 29 KHz.  Le filtre doit donc bien supprimer tout au-dessus de 70 KHz.  Cela donne de l'espace de coupure de 30 KHz à 70 KHz au filtre, ce qui est plus confortable que de devoir construire un filtre quasi-parfait.  Le prix à payer: un flux et un échantillonnage un peu plus rapide, est largement compensé par la facilité de la construction du filtre.

Fenêtrage et fuites spectrales.

Si nous avons un morceau de signal limité dans le temps, disons de 0 à T, et zéro en dehors de cette fenêtre, alors:

  1. La transformée de Fourier ce ce signal est continue
  2. La bande spectrale est infinie

Cela rend, strictement parlé, des morceaux finis de signal incompatibles avec l'exigence de bande limitée nécessaire pour l'anti-repliement.  Strictement parlé, il faudrait s'arrêter ici si on veut calculer le représentation spectrale d'un morceau de signal échantillonné, car c'est voué à l'échec comme nous voulons faire deux choses contradictoires.

Cela est bien sûr en partie le fait de la définition même de la transformation de Fourier, qui s'applique à des signaux continus de durée infinie.  Il est parfaitement possible de définir d'autres transformations de signal, compatibles avec un jeu fini d'échantillons.  Ces transformations existent et sont utilisés.  Mais notre but était d'avoir une estimation de la vraie transformation de Fourier d'un signal qui est réellement continu et de durée infinie.   Ce que nous apprenons ici, c'est que notre estimation ne sera jamais parfaite: il n'est simplement pas possible de retrouver un estimateur parfait qui rend la transformée de Fourier à base d'un nombre fini d'échantillons sur un morceau fini de signal.  Mais nous pourrons s'en approcher.  Il faut juste être conscient des erreurs que l'on va inévitablement commettre.  Il n'y a bien sûr pas moyen de se défaire de l'hypothèse d'un échantillonnage à fréquence finie et dans une fenêtre de temps finie, car c'est la seule façon d'obtenir un résultat pendant notre vie ! 

Comme toute partie du signal peut influencer toute partie de la transformée de Fourier, il est normal que se limiter à un morceau fini ne pourra jamais nous donner la vraie transformée de Fourier, car nous ignorons la partie du signal en dehors du morceau qu'on étudie.  Mais nous allons bien sûr prendre comme hypothèse que la partie que nous observons est représentatif pour le reste du signal.  Une façon de faire cela, est de prendre l'hypothèse que la partie finie qu'on décide d'observer, sera répétée infiniment en dehors de la fenêtre d'observation.  C'est mieux que de faire l'hypothèse que le signal sera zéro en dehors de cette fenêtre.

Nous illustrons cela avec un exemple simple.  Supposons que le signal d'origine est un sinus (mais nous l'ignorons, bien sûr).  Nous prenons un morceau fini de ce signal, de 0 à 2 secondes.  Nous mettons le reste à zéro.  En suite, nous répétons ce morceau infiniment (bien sûr, pour l'illustration, nous devons nous limiter aussi dans le temps...).  Il faut s'imaginer que le signal que nous allons étudier est l'étendu infinie de l'illustration à droite, qui est notre meilleure hypothèse du signal d'origine à gauche à partir du morceau fini que nous en avons pris au milieu.

Ainsi, le nouveau signal hypothétique de durée infinie, reconstruit à partir du morceau fini, que nous avons reconstruit, est un signal périodique.  Un signal périodique n'a de contenu spectral que sur une grille, à des multiples de la fréquence Fw = 1/T, où T est la période de répétition.  Il y a une relation entre la transformée de Fourier continue du morceau fini de signal et zéro en dehors (illustration du milieu) et les valeurs spectrales discrètes sur la grille du signal périodique qui est la répétition de ce morceau:

Les amplitudes du spectre discret au fréquences qui sont des multiples de Fw du signal périodique sont proportionnels aux échantillons du spectre continu du morceau de signal fini et zéro ailleurs, à exactement ces fréquences.  Ces amplitudes correspondent aux coefficients de la série de Fourier du signal périodique.

Nous pouvons illustrer cela dansle cas d'un signal d'origine qui est un sinus de 0.8 Hz.  La transformée de Fourier de ce signal d'origine, infini de durée, est la ligne verte: un seul pic à 0.8 Hz.  Le spectre de Fourier continu du morceau de ce sinus sur 2 secondes, et mis à zéro ailleurs, correspond à la ligne bleue.  Nous avons un spectre continu, d'extension spectrale infinie.  Si nous construisons un signal périodique où ce morceau se répète infiniment, nous obtenons une série de Fourier avec des coefficients correspondant aux multiples de l'inverse de la période: Fw = 0.5 Hz (car ce sont des morceaux de longueur de 2 s qui se répètent).  Ces coefficients, en rouge, sont proportionnels aux échantillons à ces fréquences (multiples de 0.5 Hz) de la courbe bleue.

D'une certaine façon, en répétant le seul morceau de signal que nous possédons, de façon infinie, nous avons la meilleure hypothèse de ce qu'était le signal original.  Un signal périodique peut avoir un contenu spectral fini.  Ce n'est pas nécessaire, mais c'est possible, contrairement à un bout fini de signal.  Comment est-ce possible ?  Comment s'élimine le contenu spectral infini du simple morceau, quand on le répète, quand l'un est l'échantillonnage de l'autre comme nous venons de l'annoncer ?  La seule possibilité pour que cela arrive, est que, bien que l'étendue du spectre continue est infinie, sa valeur est zéro aux multiples de Fw à partir d'une certaine fréquence.  Cela revient à dire que le signal périodique est composé d'une somme finie de sinus et de cosinus ayant des fréquences qui sont quelques multiples de Fw.  Si le plus haut multiple de Fw présent ainsi, est plus petit que fs, alors ce signal sera fidèlement représenté si on l'échantillonne à une fréquence 2 fs.   Dans ce cas, et seulement dans ce cas, un échantillonnage à une fréquence finie (2fs) pendant un temps fini (T) nous permettra d'estimer correctement la transformée de Fourier qui se limite alors à un nombre fini de coefficients de Fourier dans une série de Fourier.

Si nous limitons un signal de durée infinie à un morceau entre 0 et T, et nous mettons le signal à zéro en dehors de cette plage, nous pouvons considérer cela comme un produit de deux fonctions: le signal original, et une fonction qui est zéro pour t < 0 et t > T, et qui prend la valeur 1 entre 0 et T.  Nous appelons cette fonction: une fonction fenêtre.  La multiplication de deux fonctions dans le domaine temporaire correspond à une convolution des transformées de Fourier.  La transformée de Fourier de la fonction fenêtre est une fonction sinc avec une largeur de Fw.  Si notre signal d'origine était un sinus unique ayant un spectre de Fourier avec un seul pic, la convolution en fait une "copie" de la fonction sinc, mais centrée sur ce pic.  Si ce pic ne se trouve pas lui-même sur un multiple de Fw, les lobes auront des valeurs non-nulles aux multiples de Fw.  Par contre, si le pic de départ est lui-même un multiple de Fw, la fonction "sinc" déplacée et centrée sur ce pic, aura des zéros exactement à toutes les autres multiples de Fw.  Ainsi, les échantillons prises sur la grille au multiples de Fw seront zéro, sauf pour la valeur du pic d'origine.  Cela correspond à dire que si nous prenons un sinus, et nous en retenons un morceau fini qui est un multiple de sa période, la reconstruction périodique de ce morceau va reconstituer le sinus d'origine, et seulement dans ce cas-là.  Ce n'est que quand on coupe un morceau de sinus qui n'est pas sa période, que "les deux bouts ne collent pas" quand on va répéter ce morceau, et que nous allons générer, dans le nouveau signal périodique, pleins de composantes qui ne sont pas présentes dans le signal d'origine.

En réalité, nous sommes rarement dans le cas du signal répété égal au signal d'origine, qui était le seul et unique cas dans lequel nous pouvions estimer la transformée de Fourier (devenu une série de Fourier limitée) correctement, car le signal qui nous intéresse peut très bien ne pas être périodique et donc contenir des composantes fréquentielles qui ne sont pas un multiple d'une fréquence Fw.  Dans ce cas, la fenêtre finie va "disperser" le contenu spectral de chaque composante du signal d'origine selon la transformée de la fenêtre, à savoir la fonction sinc.   La répétition d'un bout fini de signal pour former un signal périodique artificiel sera alors représenté par une série de Fourier dont les coefficients sont l'échantillonnage de cette transformée convoluée avec la fonction sinc.  Comme les composantes dans le signal d'origine ne sont pas alignées à la grille des multiples de Fw les fonctions "sinc" centrées sur ces composantes n'ont pas leurs zéros sur la grille, et vont donc contribuer à la valeur des coefficients proches et même relativement loin de la composante d'origine.   Cela implique donc:

  1. les composantes de fréquences proches des composantes réelles dans le signal d'origine vont recevoir des contributions importantes de celles-ci, de la lobe principale, et des lobes auxiliaires du sinc aligné avec les composantes réelles.
  2. des composantes de fréquences élevés, au-delà des fréquences présentes dans le signal d'origine, auront des contributions non-nulles via les lobes auxiliaires des basses fréquences.  Ces composantes peuvent en suite souffrir d'un repliement qui n'aurait pas été présent vu le contenu fréquentiel original.

Ces effets sont appelés "fuite spectrale".

 Il est important de ce rendre compte que tant que le repliement proprement dit est indépendant de toute estimation spectrale ou non, la fuite spectrale est le résultat du fait que la transformation de Fourier est défini sur un signal continu de durée infinie, et que, sauf dans le cas particulier d'un signal déjà périodique et avec un nombre de termes finie dans la série de Fourier, il est impossible de déduire cette transformée d'un morceau fini de signal.  La fuite spectrale est donc un artefact de l'exercice impossible qu'on veut s'imposer: estimer une transformée de Fourier à partir d'un bout de signal seulement. 

Illustration avec un signal sinusoidal.

Dans l'exemple qui suit, nous allons utiliser une fréquence d'échantillonnage de 1 KHz.  Quand nous échantillons un signal sinus avec une fréquence de 1010 Hz, alors le graphique des échantillons successives se présente comme:

Nous constatons l'effet du repliement: les échantillons sont les mêmes que ceux d'un sinus à fréquence 10 Hz.

Quand nous échantillons un sinus de 300 Hz et nous prenons une fenêtre de temps très très longue, la transformée de Fourier calculée ressemble très fortement à la Transformée théorique: un pic à 300 Hz (nous ne regardons que les fréquences positives, et nous prenons la valeur absolue de la transformée pour les graphiques):

Par contre, si nous faisons la même chose, mais avec une fenêtre de 131 ms, ce qui impose une grille de fréquences qui sont des multiples de  7.63358... Hz, 300 Hz ne se trouve pas sur la grille, et nous observons la fuite spectrale:

Par contre, si notre fenêtre fait 150 ms, 300 Hz se trouve bien sur la grille des multiples de 6.6667 Hz, et le spectre est estimé correctement:

Nous constatons effectivement qu'il n'y ait que un seul composant non-nul, comme il se doit pour un sinus pur.

Nous constatons de nouveau que la seule façon de calculer correctement une transformation de Fourier à partir d'un morceau fini de signal, est quand le signal original est périodique, et que nous avons pris un multiple précis de la période de ce signal comme fenêtre.  Dans ce cas, nous calculons en fait les coefficients de la série de Fourier du signal périodique reconstruit.  Ce qui est fortement problématique, c'est que selon une composante, dans un vrai signal qui nous intéresse et qui n'est pas périodique, se trouve par hasard sur la grille des multiples de Fw ou non, le résultat dans le spectre calculé sera fortement différent.  Une composante sur la grille sera correctement représentée sans fuite spectrale, et une composante très légèrement différente provoquera un effet de fuite spectrale considérable.

Nous illustrons cet effet avec deux signaux qui seront tous les deux échantillonnés à 1 KHz, avec une fenêtre de temps de 200 ms, ce qui implique un Fw de 5 Hz.  Un signal composé d'un sinus de 300 Hz et un sinus de 310 Hz aura deux composantes sur la grille.  Un deuxième signal, avec une composante à 302 Hz et une autre composante à 312 Hz, sera analysé de la même façon.  Remarquez que dans les deux cas, les deux composantes diffèrent de 10 Hz.

On remarque que les spectres sont très différents, pour des signaux très semblables.  Le spectre du signal contenant le 300 Hz et le 310 Hz est parfaitement représenté, tandis que l'effet de fuite du deuxième signal est particulièrement sévère.  C'est très embêtant que des signaux si similaires ont des spectres si radicalement différents.  Le fait d'utiliser une fenêtre de 200 ms indique que nous sommes intéressés par une résolution spectrale de l'ordre de 5 Hz.  Nous ne devrions pas avoir une si grande différence entre deux signaux que nous ne voudrions même pas distinguer !

Fonctions de fenêtrage

La différence importante entre le signal de 300/310 Hz et le signal 302/312 Hz dans l'exemple précédent vient du fait qu'une fenêtre de 200 ms implique une grille de fréquences au pas de 5 Hz, et que tous les sinus sur cette grille seront représentés correctement, et tous les autres auront des fuites spectrales.   On pourrait être intéressé par une technique qui rend un peu plus "toutes les fréquences égal pour la loi".  Nous ne voulons pas une représentation parfaite de ces quelques fréquences sur la grille, si nous pouvons obtenir moins de fuites pour tous les autres.  Cela peut être obtenu en rendent les bords de la fenêtre "flous", en appliquant une fonction fenêtre non-triviale.

 La fonction fenêtre considéré jusqu'ici est la fenêtre triviale rectangulaire: une copie parfaite du signal entre 0 et T, et zéro ailleurs.  Le "zéro ailleurs" restera de mise pour toutes les fonctions fenêtre et n'est pas négociable.  Mais on pourrait rendre la transition vers zéro plus graduelle, en diminuant l'amplitude du signal près de 0 et près de T.  Il y a beaucoup de fonctions de fenêtrage.  Ces fonctions sont un peu comme les approximations des filtres: ils essaient de résoudre un problème impossible à résoudre parfaitement selon un ou autre critère, de façon optimale.

Ce que nous voulons, d'une fonction de fenêtrage, c'est que sa transformation de Fourier qui donne l'élargissement de toute composante suite à la restriction du signal dans le temps, aura un lobe centrale mince, et aura de très petits lobes secondaires.  Idéalement, la transformée serait un pulse de Dirac.  Mais alors la fenêtre est la fonction constante "1", et ne peut pas être 0 en dehors de l'intervalle 0 - T.

Une liste non-exhaustive de fonctions de fenêtrage:

  • la fenêtre rectangulaire (fenêtre triviale)
  • fenêtre triangulaire
  • fenêtre de Parzen (ou de la Vallée Poussin)
  • Welch
  • Hamming
  • Hann
  • Blackman
  • Nuttal
  • Blackman-Harris
  • Haut plat
  • Cosinus
  • Gaussien
  • Tukey
  • Slepian
  • Kaiser
  • Dolph-Chebyshev
  • Lanczos
  • ...

Le choix de la bonne fenêtre dépend de l'application et son choix nécessite une certaine expertise.   Dans beaucoup de cas, le choix importe finalement peu, mais il faut être conscient des effets et vérifier dans quelle mesure ces effets influent le résultat important pour l'application en question.  Dans beaucoup de cas, il y a un compromis entre la perte de résolution spectrale (la largeur du lobe central) et l'intensité et diminution des lobes secondaires.

Illustrons cela avec deux cas extrêmes: la fenêtre rectangulaire et la fenêtre Blackman.  Cette dernière a l'aspect suivant:

Nous appliquons exactement les mêmes opérations sur les signaux contenant 300 Hz et 310 Hz, et 302 Hz et 312 Hz, mais nous appliquons une fenêtre Blackman à la place de la fenêtre rectangulaire.  Voici les deux spectres:

L'effet est frappant !  Il n'y a presque pas de différence entre les spectres des deux signaux, sauf que le signal contenant 302 Hz et 312 Hz donne un spectre très légèrement décalé vers la droite - comme il se doit.  Il n'y a plus de fuite spectrale importante à grande distance.  Mais sans doute l'effet le plus frappant est que nous avons un seul pic par spectre.  Nous ne distinguons plus le 300 Hz du 310 Hz.  C'est le compromis habituel: la diminution des lobes secondaires à grande distance est compensée par un élargissement du lobe central (et donc une perte de résolution).  La fenêtre rectangulaire a une très bonne résolution spectrale, mais une fuite spectrale terrible.  La fenêtre Blackman est fantastique en ce qui concerne les fuites, mais a une résolution terrible.  D'où le besoin, parfois, d'expertise pour choisir la fonction de fenêtrage de façon judicieuse.  Il n'y a pas de fenêtre universellement meilleur que une autre.  Cela dépend de l'application.

DFT, FFT et transformée de Fourier

La transformée de Fourier rapide (FFT) n'est rien d'autre qu'un algorithme efficace pour implémenter la transformée de Fourier Discrète (DFT).  La DFT s'applique à un jeu fini d' échantillons, c'est à dire, un signal limité par une fenêtre dans le temps (entre 0 et T), et échantillonné à une fréquence 2 fs, de telle façon que nous avons N échantillons.  La DFT est simplement le calcul des lignes de la transformée de Fourier (réduit à une série de Fourier) du morceau de signal répété dans le temps pour fabriquer un signal périodique.  Les lignes du spectre sont sur la grille Fw, jusqu'à la fréquence fs.   La DFT sont N nombres complexes (y compris les fréquences négatives).  Si la transformée de Fourier du morceau de signal contient des fréquences supérieures à fs sur la grille Fw, un effet de repliement aura lieu, même si le signal d'origine ne contenait pas ces fréquences.  C'est à dire, la composante dans la DFT, disons, à 7 Fw sera l'échantillon de la valeur de la transformée du morceau unique à 7 Fw, diminuée de l'échantillon de cette transformée à 2 fs -  7 Fw, augmenté avec la valeur à 2 fs + 7 Fw etc...  C'est la meilleure estimation de la transformée à 7 Fw du signal original que nous pouvons obtenir, si nous n'avons pas d'autres informations que les échantillons, mais nous savons que cette estimation aura soufferte de repliement et de fuite spectrale.

Il y a quelques facteurs de normalisation qui doivent êtres prises en compte pour passet d'une DFT à N points vers un estimateur d' une véritable transformée de Fourier continue.

La DFT comme série de Fourier

Le seul cas où la DFT donnera un résultat exactement correct est le cas du signal périodique de période T, à série de Fourier finie, donc avec seulement un nombre de termes finis, sans terme au delà de fs.  S'il y a des petits coefficients non-nuls au-delà de fs ils apparaîtront comme des contributions aux plus basses fréquences: l'effet de repliement.  Dans la mesure où on peut accepter cela, la DFT calcule donc bien une série de Fourier.

Considérons un signal carré qui est 1 la moitié du temps, et 0 l'autre moitié, avec une période de 0.4 secondes.  Cette fonction périodique a une série de Fourier qu'on peut calculer théoriquement:

f(t) = 0.5 + 2/π Σ 1/n sin(2 π t / T)  

où la somme s’effectue sur les n impair, en commençant à 1.

Pour notre série spécifique, cela donne des coefficients suivants:

pour n = 0: 0.5
pour n = 1: 2/π = 0.6366...
pour n = 3: 2/3π = 0.2122...
pour n = 5: 2/π 1/5 = 0.1273...
pour n = 7: 2/π 1/7 = 0.0909...

etc...

Si nous prenons les 5 premiers termes de la série de Fourier, cela donne une fonction qui est une approximation de notre fonction carré.  Nous pouvons nous imaginer que si nous ajoutons de plus en plus de termes de la série, le résultat convergera bien vers le signal rectangulaire d'origine:

Nous allons maintenant échantillonner le vrai signal rectangulaire, de 3 façons différentes: avec 400 échantillons, avec 40 échantillons, et avec 8 échantillons.  Cela correspond à des fréquences d'échantillonnage différents:

Nous observons déjà un effet: bien que la fonction d'origine monte à 1 a t = 0 exactement, et descend à 0 à t = 0.2 exactement, les échantillons sont tels que la montée et la descente apparaissent un peu plus tôt.  Cela donnera un décalage de phase dans les coefficients de Fourier, comme nous allons le constater.

La DFT de la suite de 400 échantillons produira aussi 400 coefficients.  La deuxième moitié de cette suite consistera des conjugués complexes de la première moitié.  Pour l'illustration, nous allons représenter les valeurs absolus:

Les valeurs sont:

200
1 - 127.32 i
0
1 - 42.43 i
0
1 - 25.45 i
0
1 - 18.17 i

...

Les valeurs de DFT doivent être divisés par N/2 pour obtenir l'estimation des coefficients de la série de Fourier à la façon dont ceux-ci sont habituellement défini en mathématique.  La partie réelle est le coefficient du cosinus, et moins la partie imaginaire est le coefficient du sinus.  Le tout premier terme est la constante, et doit être divisé par N, et non par N/2.

En divisant les valeurs de notre DFT par N/2 (ici donc 200, sauf le premier terme qu'on divise par 400) et en prenant les valeurs absolus, nous trouvons:

0.5
0.6366
0
0.2122
0
0.1274
0
0.0910

...

Comme nous voyons, ces valeurs sont très proches des coefficients théoriques de la série de Fourier.  Seulement, dans la série théorique il n'y a que des termes sinus et tous les cosinus sont absents.  Dans notre DFT complexe, nous observons qu'il y a une petite partie réelle (égale à 1 avant division par 200), correspondant à des petits coefficients de termes cosinus.  Ceci est l'effet du décalage de phase qu'on a observé dans le graphique avec les échantillons: le fait que la descente se fait un peu avant 0.2 secondes.

Si nous passons sur notre suite de 40 échantillons au lieu de 400, et nous montrons les valeurs absolues de la DFT:

nous allons observer un comportement similaire:  Les valeurs DFT sont:

20
1 - 12.7062 i
0
1 - 4.1653 i
0
1 - 2.4142 i
0
1 - 1.6319 i

...

Nous constatons déjà que le décalage de phase sera plus important (le terme réel est le même, mais la partie imaginaire est plus petite).  Le décalage de phase est 10 fois supérieur au cas avec 400 échantillons.

En regardant les valeurs absolues, normalisés en divisant par N/2 = 20 (sauf le premier terme), nous obtenons:

0.5
0.6373
0
0.2142
0
0.1307
0
0.0957

...

Ces valeurs sont toujours proches des valeurs théoriques, mais nous constatons que pour les termes des fréquences plus élevés, les valeurs seront un peu supérieurs aux valeurs théoriques.  Par exemple, le coefficient de la 7ième harmonique devrait être 0.0909.  Dans la DFT à 400 points, c'est 0.0910.  Ici, avec 40 points, c'est devenu 0.0957.  Nous commençons à voir l'effet du repliement !

Cet effet sera très prononcé si nous ne prenons que 8 échantillons:

La DFT complète est:

4
1 - 2.414 i
0
1 - 0.4142 i
0
1 + 0.4142 i
0
1 + 2.414 i

Les 4 dernières points sont les conjugués complexes et ils ne nous apprendront rien de nouveau qui ne sera pas présent dans les 4 premiers.  Le décalage de phase est énorme: les termes en cosinus sont aussi importants que les termes sinus.

En normalisant, nous trouvons pour les valeurs absolues:

0.5
0.6533
0
0.2706

L'effet de repliement est très fort et se voit même dans le premier terme qui est 0.6533 au lieu de 0.6366.

Nous pouvons conclure que la DFT nous permet de calculer les coefficients de la série de Fourier d'une fonction réellement périodique, mais que l'échantillonnage fini peut introduire un décalage de phase, et introduit bien sûr un effet de repliement si la série d'origine ne s'arrête pas à la fréquence de Nyquist.  Nous voyons aussi que la normalisation ne dépend que du nombre de points N dans la DFT, et que l'échelle de temps absolue ne joue aucun rôle.

Finalement, nous allons illustrer ce qui se passe si nous prenons plusieurs périodes du signal dans la fenêtre.  Nous allons prendre notre même fonction carré, mais notre fenêtre sera 5 périodes, à savoir 2 secondes.  Nous allons considérer N = 2000, N = 200, et N = 20.

Finally, we illustrate what happens if we take a periodic signal of which we take a window that is an integer multiple of the signal period.  We take our same square function, but this time we take a window that is 5 times the period (that is to say, 2 seconds, instead of 0.4 seconds which is the signal period).  We will sample with N = 2000, N = 200 and N = 20.

Dans ce graphique répété, nous voyons encore mieux l'effet du décalage dans le temps (ce qui explique le décalage de phase): le signal vert a visiblement des pulses "plus tôt" que les signaux bleus.  Si nous regardons bien, nous pouvons même voir que la courbe rouge est un peu plus vers la gauche que la courbe bleue.

Les valeurs absolues de la DFT de 2000 points est montrée en suite:

Les valeurs sont 5 fois plus grands que dans le cas d'une simple période avec 400 échantillons à la même fréquence d'échantillonnage.  C'est normal, car il faudra normaliser par N/2, et N est 5 fois plus grand aussi.

Pour le cas N = 200:

Nous constatons qu'entre chaque valeur non-zéro il y a 9 valeurs zéro (avant, nous n'avions que 1 valeur zéro entre les non-zéro).   C'est normal, car Fw est un cinquième maintenant de la valeur de Fw avec une fenêtre d'une période.  Mais les composantes fréquentielles de la série de Fourier n'ont pas changé parce qu'on prend 5 périodes.  Il y a donc 4 termes sinus (et cosinus) avec un coefficient 0 avant de retomber sur une fréquence présente potentiellement dans la série de Fourier d'origine.  Effectivement, avec une période de 0.4 s, les fréquences de la série de Fourier sont des multiples de 2.5 Hz.  Pour notre signal carré, il y a une composante à 2.5 Hz, mais celui de 5 Hz est zéro dans notre cas.  Puis il y a un terme à 7.5 Hz.  Et ainsi de suite.

A 2 secondes de fenêtre, Fw est 0.5 Hz.  La DFT va donc calculer les coefficients correspondant à des multiples de 0.5 Hz.  Pour 0.5 Hz, 1 Hz, 1.5 Hz et 2 Hz, il n'y a pas de terme potentiel dans la série de Fourier bien sûr, donc ces 4 coefficients ne peuvent être que zéro.  Mais comme dans notre cas précis, un terme sur 2 de la série de Fourier est zéro (les harmoniques paires: 5 Hz, 10 Hz, 15 Hz...), nous avons des suites de 9 zéros successives (4 impossibles, une harmonique paire, 4 impossibles).

Les premières valeurs de notre DFT à 200 points est:

100
0
0
0
0
5 - 63.5310 i
0
0
0
0
0
0
0
0
0
5 - 20.8265 i
0
0

...

En normalisant par N/2, en éliminant les 4 zéros "impossibles" et en prenant les valeurs absolues:

0.5
0.6373
0
0.2142
0

...

Ces valeurs ne sont pas identiques aux valeurs théoriques, mais elles sont strictement égales aux valeurs trouvés avec une période et N = 40.  C'est très important de se rendre compte: une période à N = 40 est parfaitement équivalent à 5 périodes et N = 200, sauf pour le fait que nous ajoutons 4 zéros inutiles entre chaque valeur utile et qu'il faut diviser par 100 au lieu de diviser par 20 pour retrouver les coefficients de Fourier.

La DFT et la Transformée de Fourier.

Maintenant que nous savons que la DFT est une méthode de calcul de coefficients d'une série de Fourier, que reste-t-il de notre espoir de calculer des transformées de Fourier ?

D'abord, dans la mesure où il n'y a aucun doute sur la normalisation standard d'une série de Fourier, les transformées de Fourier ont plusieurs normalisations conventionnées - essentiellement selon qu'on est un ingénieur, un physicien, ou un mathématicien.  Nous allons adopter la convention standard pour les ingénieurs:

F(ω) = ∫ f(t) e-jωt dt

où ω est la fréquence angulaire.

Nous savons déjà que les valeurs de la DFT sont proportionnels aux échantillons F(2 π n Fw) où F est la transformée de Fourier du morceau fini du signal entre 0 et T, et 0 en dehors de cette fenêtre.  Quelle est la constante de proportionnalité ?

Il est facile de démontrer que le n-ième coefficient de la Serie de Fourier complexe est:  cn = 1/T X(2 π n Fw),où T est la longueur de la fenêtre en temps que nous avons pris pour calculer la transformée de Fourier X, et que nous avons répété infiniment pour avoir un signal périodique dont nous considérons la série de Fourier.

Le coefficient du cosinus et le coefficient du sinus sont respectivement deux fois la partie réelle et moins deux fois la partie imaginaire de cn qui sont eux-mêmes 2/N fois la DFT.

En d'autres termes: cn = 1/N * DFT(n + 1) pour n < N/2.  Il en suit:

X(2 π n Fw) = T/N * DFT(n + 1) for n < N/2

L'estimation de la valeur complexe de la transformée de Fourier d'un morceau fini de signal entre 0 et T aux échantillons, multiples de Fw et donc fréquence angulaire  2 π n Fest donnée par la (n+1)ième valeur de la DFT, multipliée par T/N, ce qui n'est rien d'autre que la période d'échantillonnage.

Vérifions cela avec notre morceau de signal carré, qui est 1 entre t = 0, et t = 0.2 s, mais zéro ailleurs.  Nous allons prendre une période longue: T = 5 secondes.  Ainsi, nous échantillonnerons la transformée de Fourier au pas de 0.2 Hz.

La transformée de Fourier d'un morceau de signal qui est 1 entre t = 0 et t = 0.2 est théoriquement:

F(ω) = i / ω ( cos(0.2ω) - i sin(0.2ω) - 1)

ce qui peut être vérifié facilement par simple intégration.

En échantillonnant notre fonction à une fréquence de 500 Hz, nous aurons N = 2500, dont 100 seront 1, et 2400 seront 0.  Il faut se souvenir que cette fois, nous ne considérons pas un signal périodique, mais un signal fini, qui est zéro en dehors de la fenêtre de temps.  La fenêtre étant 5 secondes, la grille sur laquelle nous allons échantillonner la transformée sera 0.2 Hz, ou 1.256 rad/s.  Nous allons comparer les valeurs absolues de la DFT normalisés avec un facteur 0.002 s (la fréquence d'échantillonnage) sur une partie de l'axe omega, avec la valeur absolue de la transformée théorique:

Nous pouvons comparer les valeurs complexes, cela marche bien aussi.  On pourrait croire que nous exagérons avec une fréquence d'échantillonnage de 500 Hz.  Nous avons, finalement, regardé le graphique des fréquences jusqu'à 10 Hz (62.5 rad/s).   Une fréquence d'échantillonnage de plus que 20 Hz serait donc OK...

Si nous répétons l'exercice avec une fréquence d'échantillonnage de 25 Hz et une période de 5 s, nous voyons l'effet de repliement: les valeurs de la DFT sont plus larges que les valeurs théoriques, surtout aux fréquences un peu plus élevées.

Conclusion

Si nous voulons calculer numériquement la transformée de Fourier ou la transformée de Fourier un outil est la DFT, qui va transformer un jeu fini d'échantillons d'un signal à fréquence d'échantillonnage limité sur une durée de signal limitée en un jeu de nombres complexes.  Ces nombres complexes peuvent, si nous considérons notre morceau fini de signal échantillonné comme une période d'un signal périodique, être considérés comme des estimations des coefficients de la série de Fourier de ce signal périodique.  Ils peuvent aussi, si nous considérons notre morceau de signal comme non-zéro seulement pendant la fenêtre observée, être considérés comme des échantillons de la transformée de Fourier de ce signal fini dans le temps.  Dans les deux cas, il faut appliquer une normalisation différente.

Dans ces deux cas, le repliement peut se manifester, un effet inévitable d'utiliser une fréquence d' échantillonnage fini.Dans la mesure où notre morceau de signal est représentatif pour un signal de durée infinie, de la fuite spectrale sera inévitable sauf dans des cas assez particuliers.  Cette fuite spectrale est dépendante du choix de la fonction de fenêtrage utilisée pour pondérer le morceau fini du signal.Dans cette contribution, nous avons parlé seulement de l'estimation de la transformée de Fourier d'un seul signal déterminé, et pas du tout d'un ensemble de signaux.Le choix de la fréquence d'échantillonnage, de la période et de la fonction fenêtre nécessitent une certaine expertise pour mettre en accord une estimation spectrale avec l'application en question.  Il n'y a pas de solution universelle qui donne la transformée de Fourier parfaite, car nous ne disposons que d'un nombre fini d'échantillons sur une durée limitée dans le domaine numérique.