Maîtrise d'informatique - Module 18 Traitement du signal - Cours Chapitre 3 : 1. Introduction. 2. La discrétisation fait-elle perdre de l'information ?. 3. La transformation de Fourier discrète (TFD). 4. Algorithme de la transformation de Fourier rapide (TFR).

L'objet essentiel de ce chapitre est d'introduire le moyen de calculer des transformées de Fourier à l'aide de l'outil informatique.

1. Introduction.

1.1. Calcul d'intégrales grâce à l'outil informatique.

Les calculs de transformées de Fourier "à la main" doivent être abandonnés pour deux raisons :

• De tels calculs sont difficiles, comme le prouvent les quelques exemples effectués dans le chapitre 2.

• Ils sont impossibles si l'on n'a pas l'expression analytique du signal, donc en particulier dans le cas d'un signal réel !

Notons quand même qu'il existe, en "optique de Fourier", un système permettant de calculer la transformée de Fourier d'un signal bidimensionnel quelconque. Ce moyen constitue un "calculateur analogique".

1.1.a. L'outil informatique peut-il être un "calculateur analogique" ?

Contrairement aux idées reçues, l'outil informatique peut constituer un calculateur analogique. Il s'agit alors de ce qu'on appelle le calcul formel. A l'aide de certains logiciels, qui autorisent le calcul de l'intégrale d'une fonction entre deux bornes, on peut calculer la transformée de Fourier d'une fonction d'expression analytique connue.

Exemple :

On désire calculer, à l'aide d'un logiciel de calcul formel, la transformée de Fourier de la fonction "fenêtre" déjà vue au chapitre 2. Comme un tel logiciel sait calculer l'intégrale d'une fonction entre deux bornes, et que la fonction "fenêtre" est nulle en dehors de l'intervalle [-1,1], il suffit d'intégrer la fonction  entre les bornes -1 et 1, et le logiciel fournira comme résultat :



qui est bien la transformée de Fourier de la fonction "fenêtre".

Malheureusement, le calcul formel ne permet pas, lui non plus, de calculer la transformée de Fourier d'un signal réel.

1.1.b. L'outil informatique en tant que "calculateur discret".

Il s'agit dans ce cas de ce qu'on appelle le calcul numérique. L'intérêt du calcul numérique par rapport au calcul formel est qu'il permet de traiter tous les signaux, y compris les signaux réels. En revanche, le calcul numérique présente le défaut de faire perdre de l'information, et ceci à cause de deux phénomènes distincts et qu'il ne faut pas confondre :

• La "discrétisation" du signal.

• La "quantification" du signal discret.

Nous allons détailler successivement chacun de ces deux phénomènes.

1.2. Discrétisation et quantification du signal.

1.2.a. Discrétisation.

Soit par exemple le signal analogique f(t) suivant :

La discrétisation de ce signal revient à n'en garder qu'un certain nombre de valeurs (..., f_-1, f₀, f₁, ...) correspondant aux valeurs (..., t_-1, t₀, t₁, ...) de la variable t :

On ne dispose donc plus, après discrétisation, que d'une suite (f_n)_Z, où Z désigne l'ensemble des entiers relatifs. On appelle cette suite un signal discret associé au signal analogique f(t).

Dans la quasi-totalité des cas, l'échantillonnage est périodique, c'est-à-dire que ... = t₀ - t_-1 = t₁ - t₀ = ... = T_e, valeur qu'on appelle la période d'échantillonnage. On définit aussi la pulsation d'échantillonnage par :

Un problème majeur en traitement du signal est le choix de cette pulsation d'échantillonnage w_e :

• Si w_e est trop grand, on aura trop de données à traiter et à stocker.

• Si w_e est trop petit, la suite (f_n)_Z ne sera pas fidèle au signal initial f(t).

Nous verrons un peu plus loin qu'il existe un critère optimal pour choisir cette pulsation, qui s'appelle le "critère de Shannon".

1.2.b. Quantification.

La quantification est une opération qui modifie légèrement les valeurs du signal discret précédent. Ceci vient du fait qu'en informatique, les supports d'information ne sont pas analogiques mais binaires, et que donc toute valeur réelle ne peut pas être mémorisée. On peut par exemple représenter ce phénomène par le schéma suivant :

Sur ce schéma, les traits horizontaux représentent l'ensemble des valeurs possibles pour le signal.

On obtient donc, après quantification, une suite (f'_n)_Z telle que, pour tout entier n, .

Remarques :

• On dit parfois que la discrétisation est une "discrétisation en abscisse", alors que la quantification est une "discrétisation en ordonnée".

• La quantification introduit une erreur dans les valeurs du signal, ce que ne fait pas la discrétisation.

Dans toute la suite de ce cours, nous ne prendrons plus en compte que la discrétisation, et nous ne parlerons donc pas du moyen de tenir compte de l'erreur introduite par la quantification.

2. La discrétisation fait-elle perdre de l'information ?

On a vu qu'il était très simple de passer d'un signal analogique f(t) au signal discret correspondant (f_n)_Z, pourvu qu'une pulsation d'échantillonnage w_e ait été choisie. L'opération inverse est-elle réalisable ?

2.1. Notations.

Soit f(t) le signal analogique de départ et T_e la période d'échantillonnage. Le signal discret est la suite (f_n)_Z, où f_n = f(nT_e), ce qui signifie que l'origine des t correspond à n = 0.

On cherche une fonction f*(t) qui soit "proche" de ce signal discret. Une idée consiste à prendre une fonction nulle en dehors des valeurs t_n de l'échantillonnage, et valant f_n pour t = t_n, c'est-à-dire la fonction suivante :

Cette fonction n'est malheureusement pas satisfaisante, car sa transformée de Fourier est la fonction nulle, ce qui est facile à vérifier. En revanche, si on définit f*(t) de la manière suivante :

où  désigne l'impulsion de Dirac pour t = t_n, on n'a plus le problème précédent. Néanmoins, il ne faudra pas oublier que f*(t) n'est pas une vraie fonction, mais une pseudo-fonction (ou "distribution"), et qu'elle n'admet pas de transformée de Fourier, mais seulement une "transformée de Fourier généralisée" (on omettra, par la suite, l'adjectif "généralisée", et on parlera de "transformée de Fourier").

En utilisant la formule (2.11) du chapitre 2, il est évident que :

Par ailleurs, la même intégrale appliquée à f(t) peut être approximée ainsi :

ce qui veut dire qu'on a l'approximation suivante :

On peut se demander si, par analogie, l'approximation suivante ne serait pas également vraie :

c'est-à-dire si la transformée de Fourier de f(t) ne serait pas approximable par la "transformée de Fourier" de f*(t) multipliée par T_e.

2.2. Calcul de la "transformée de Fourier" de f*(t).

Supposons que la fonction f(t) admette une transformée de Fourier . La distribution f*(t) qui vient d'être introduite n'admet pas de transformée de Fourier, mais elle admet une "transformée de Fourier", que nous allons calculer.

2.2.a. "Transformée de Fourier" de f*(t).

L'implusion de Dirac  vérifie, par définition (voir chapitre 2) :

En remplaçant t par -t dans cette égalité, on obtient :

donc la "transformée de Fourier" de  est , fonction de la famille de Fourier, qui est périodique.

En utilisant la linéarité de la transformation de Fourier, on a donc :

Or, d'après ce qu'on vient de voir, en inversant les variables muettes t et w :

donc :

L'expression du second membre est ce qu'on appelle une série de Fourier.

Remarque :

La fonction  vérifie :



Or, comme , quel que soit l'entier relatif n, il en résulte que :



ce qui signifie que  est une fonction périodique, de période . Cela suffit d'ailleurs à prouver que  n'est pas contenue dans T.

2.2.b. "Transformée de Fourier inverse" de la "transformée de Fourier" de f*(t).

La "transformée de Fourier inverse" de la série de Fourier  est très différente de la transformée de Fourier inverse d'une fonction, telle qu'elle a été définie au chapitre 2. Au lieu d'associer une fonction f(t) à une fonction , elle associe à la série de Fourier  la suite des coefficients f_n de cette série de Fourier. On demande ici au lecteur d'admettre la relation suivante, vraie pour un entier relatif n quelconque :

On va maintenant revenir à notre intuition initiale, selon laquelle il pourrait y avoir une relation entre les transformées de Fourier de f(t) et f*(t).

2.3. Relation entre les transformées de Fourier de f(t) et f*(t).

Appliquons d'abord la transformation de Fourier inverse à  :

En particulier, si l'on remplace dans cette égalité t par t_n = nT_e :

A ce niveau du raisonnement, on utilise une astuce consistant à découper R en une infinité d'intervalles du type , qui sont disjoints deux à deux si k est un entier relatif. L'égalité précédente devient alors :

En faisant les changements de variables suivants :

il vient :

En remarquant que les facteurs  sont tous égaux à 1, et en remplaçant les variables muettes w_k par w, cette égalité devient :

Si on intervertit la somme et l'intégrale du second membre, et qu'on égale les expressions de f_n fournies par (3.3) et (3.4), on obtient l'égalité suivante :

En identifiant les deux intégrandes, on établit l'égalité fondamentale suivante :

On a bien trouvé une relation entre les transformées de Fourier de f(t) et f*(t), légèrement plus compliquée que celle qu'on attendait. L'égalité (3.5), qui traduit cette relation, signifie que  est égale à la somme de  et de toutes ses "translatées" de , où k désigne un entier relatif quelconque.

Remarque :

On retrouve bien, dans (3.5), le fait que  est périodique, de période .

2.4. Théorème de Shannon.

Soit une fonction f(t) admettant une transformée de Fourier, et telle que la fonction  soit nulle en dehors d'un intervalle [w_min,w_max] de largeur L, où L = w_max-w_min.

Plutôt que de rester dans le cas général, et pour mieux illustrer notre propos, prenons l'exemple de la fonction "sinus cardinal", dont on a vu au chapitre 2 que sa transformée de Fourier était égale, à un facteur  près, à la fonction "fenêtre", sauf aux deux points de discontinuité de cette fonction :

Pour la fonction , deux cas très différents se présentent :

• Si , alors , qui est la somme de   et de toutes ses "translatées" de  = kw_e, a l'allure suivante :



On peut retrouver la fonction  à partir de la fonction  en la "découpant" sur l'intervalle [w_min,w_max], et en complétant ce "tronçon" par la valeur 0.

• Si w_e < L, alors  a l'allure suivante :



Dans ce cas, à cause du chevauchement des motifs, on ne peut plus remonter de la fonction  à la fonction f(t).

On peut généraliser ce résultat en énonçant le :

Théorème de Shannon :

La discrétisation d'un signal analogique f(t), dont le spectre  est nul en dehors d'un intervalle [w_min,w_max], fait perdre de l'information si et seulement si la pulsation d'échantillonnage w_e est strictement inférieure à la largeur du spectre L = w_max-w_min.

2.5. Application du théorème de Shannon : enregistrement des disques compacts.

Le théorème de Shannon a une application très importante dans le domaine du son, et plus précisément pour les enregistrements de disques compacts, sur lesquels le signal sonore est stocké de manière discrète. Pour qu'un lecteur de disques compacts puisse restituer le son analogique initial de manière satisfaisante, il faut éviter que la discrétisation de ce signal fasse perdre de l'information. Or, d'après le théorème de Shannon, la discrétisation d'un signal ne fait perdre d'information que si la pulsation d'échantillonnage w_e est strictement inférieure à la largeur du spectre L = w_max-w_min. Le signal qui provient d'un enregistrement sonore a-t-il un spectre de largeur finie ? En général non, mais comme l'oreille humaine n'est sensible qu'aux fréquences comprises dans l'intervalle [20 Hz,20000 Hz], on peut très bien se passer de tenir compte de la partie du spectre qui se trouve en dehors de cet intervalle audible. Voici donc les étapes permettant l'enregistrement des disques compacts sans perte d'information :

• Le signal analogique f(t) recueilli à l'aide du micro est "filtré" (nous reparlerons du filtrage au chapitre 4), c'est-à-dire qu'on élimine de son spectre les fréquences supérieures à la plus grande fréquence audible 20000 Hz, et on obtient un nouveau signal analogique f'(t).

• La largeur du spectre de f'(t) n'est pas 20000 Hz, comme on pourrait le penser naïvement, car il faut également tenir compte des "fréquences négatives". A quoi correspondent ces "fréquences négatives" ? A rien d'un point de vue physique, mais dans toutes les définitions que nous avons données jusqu'à présent, la pulsation w est réelle, et donc il en va de même pour la fréquence n qui lui est liée par la relation . La largeur du spectre de f'(t) est donc 40000 Hz.

• On discrétise le signal f'(t) à la "fréquence d'échantillonnage" n_e = 44100 Hz, qui est supérieure à la largeur de son spectre. C'est cette fréquence d'échantillonnage de 44100 Hz qui est utilisée pour les disques compacts commerciaux !

3. La transformation de Fourier discrète (TFD).

Désormais, nous savons qu'un signal contient la même information sous forme discrète que sous forme analogique, pourvu que le critère de Shannon  soit vérifié, avec L = w_max-w_min. La transformation de Fourier discrète est censée être l'équivalent, pour un signal discret, de la transformation de Fourier, pour un signal analogique. C'est cette nouvelle transformation que nous allons définir dans ce paragraphe.

3.1. Introduction.

3.1.a. Récapitulatif des transformations déjà étudiées.

Une même information est censée pouvoir être représentée par les quatre "objets" suivants  :

• Objet 1 : un signal analogique, qui est en général le signal physique original.
• Objet 2 : un spectre analogique, qui est obtenu par transformation de Fourier du signal précédent.
• Objet 3 : un signal discret, qui est obtenu par discrétisation du signal analogique.
• Objet 4 : un spectre discret, qui est obtenu par discrétisation du spectre analogique.

Si les discrétisations sont faites de manière à respecter le critère de Shannon, ces quatre "objets" contiennent la même information. La question que nous nous posons maintenant est la suivante : sait-on passer de l'un quelconque de ces objets à un autre ? Remarquons d'abord qu'il existe douze transformations de ce genre, que l'on peut représenter par les relations suivantes :

Parmi ces douze transformations, nous en avons rencontré quatre jusqu'à présent, et ce sont ces quatre transformations qui apparaissent en rose. Faisons-en le détail :

• 1 -> 2 : c'est la transformation de Fourier du signal, vue au chapitre 2.
• 2 -> 1 : c'est la transformation de Fourier inverse du spectre, vue également au chapitre 2.
• 1 -> 3 : c'est la discrétisation du signal, dont on a parlé dans le paragraphe I.
• 3 -> 2 : c'est la "restauration" du spectre analogique à partir du signal discret, transformation qui a fait l'objet de l'étude du paragraphe II, et qui n'est possible que si le critère de Shannon est vérifié.

3.1.b. Nouvelles transformations.

Comme un spectre est une fonction complexe d'une variable réelle, on peut considérer un spectre comme un signal particulier et, de la même façon, on peut considérer un signal comme un spectre particulier, donc parmi les douze transformations citées précédemment, par cette simple symétrie, on sait en réaliser deux autres :

• 2 -> 4 : c'est la discrétisation du spectre, "symétrique" de la discrétisation du signal 1 -> 3.
• 4 -> 1 : c'est la restauration du signal analogique à partir du spectre discret, "symétrique" de la restauration du spectre analogique à partir du signal discret 3 -> 2.

Enfin, en appliquant successivement deux des transformations déjà connues (il y en a six maintenant), on peut en réaliser quatre de plus :

• 1 -> 4 peut être réalisée en "enchaînant" 1 -> 2 et 2 -> 4, transformations déjà connues : cette transformation ne porte pas de nom.
• 2 -> 3 peut être réalisée en "enchaînant" 2 -> 1 et 1 -> 3, transformations déjà connues : cette transformation non plus ne porte pas de nom.
• 3 -> 1 peut être réalisée en "enchaînant" 3 -> 2 et 2 -> 1, transformations déjà connues : on l'appelle "conversion digital/analogique" du signal.
• 4 -> 2 peut être réalisée en "enchaînant" 4 -> 1 et 1 -> 2, transformations déjà connues : on l'appelle "conversion digital/analogique" du spectre.

Il ne reste donc plus que les deux transformations symétriques suivantes à réaliser :

Ces deux transformations sont les seules qui associent à un objet discret un autre objet discret, donc ce sont les seules qui pourront être utilisables par un calculateur discret. A l'aide du raisonnement précédent, qui consiste à "enchaîner" des transformations déjà connues, il semble que l'on puisse les réaliser toutes les deux :

• 3 -> 4 semble être réalisable en "enchaînant" 3 -> 2 et 2 -> 4.
• 4 -> 3 semble être réalisable en "enchaînant" 4 -> 1 et 1 -> 3.

Or, on commettrait une grave erreur dans le raisonnement en procédant de la sorte, erreur qui ne saute pas aux yeux, et que nous allons exposer dans le prochain sous-paragraphe.

3.1.c. Signal discret et spectre discret ne peuvent pas contenir simultanément toute l'information.

Soit un signal analogique f(t), et  son spectre analogique. Soit w_e la pulsation d'échantillonnage du signal, et w_f la pulsation d'échantillonnage du spectre. Si l'on veut que le signal discret et le spectre discret contiennent tous deux la même information que le signal analogique (ou que le spectre analogique), alors, d'après le théorème de Shannon :

• La pulsation d'échantillonnage du signal w_e doit être supérieure ou égale à la largeur du spectre .

• La pulsation d'échantillonnage du spectre w_f doit être supérieure ou égale à la largeur du spectre, c'est-à-dire à la largeur de  f(-t), donc w_f doit être supérieure ou égale à la largeur du signal multipliée par .

Ceci implique en particulier que signal et spectre doivent tous deux avoir une largeur finie. Or, on peut montrer que ceci est impossible : il est impossible qu'un signal et son spectre soient simultanément de largeur finie. Par conséquent, il est impossible de discrétiser simultanément un signal et son spectre en respectant pour chacun d'eux le critère de Shannon. Ceci a la conséquence très importante suivante :

• Si le signal discret contient toute l'information (ce qui veut dire que w_e respecte le critère de Shannon), alors on pourra bien trouver une transformation du type 3 -> 4, mais il n'existera pas de transformation du type 4 -> 3, puisque l'objet 4 contiendra moins d'information que l'objet 3.

• Si le spectre discret contient toute l'information (ce qui veut dire que w_f respecte le critère de Shannon), alors on pourra bien trouver une transformation du type 4 -> 3, mais il n'existera pas de transformation du type 3 -> 4, puisque l'objet 3 contiendra moins d'information que l'objet 4.

Le seul moyen pour qu'il existe simultanément deux transformations 3 -> 4 et 4 -> 3 inverses l'une de l'autre, est que ni l'objet 3 ni l'objet 4 ne contiennent toute l'information. Les deux transformations 3 -> 4 et 4 -> 3 inverses l'une de l'autre seront alors respectivement la transformation de Fourier discrète et la transformation de Fourier discrète inverse. Nous allons les définir dans le prochain paragraphe.

3.2. La transformation de Fourier discrète.

3.2.a. Signal discret de longueur finie et spectre discret de longueur finie.

On appelle signal discret de longueur finie d'un signal analogique f(t), et on note (f_n)_[0,N-1], tout ensemble de N valeurs successives (N fini) provenant de l'échantillonnage de f(t) à une certaine pulsation w_e. Un signal discret de longueur finie est composé de valeurs f_n qui forment un sous-ensemble des valeurs f_n du signal discret (f_n)_Z. Alors que le signal discret (f_n)_Z peut contenir toute l'information du signal f(t), si w_e vérifie le critère de Shannon, un signal discret de longueur finie (f_n)_[0,N-1] ne peut bien sûr pas contenir toute cette information.

De la même façon, on appelle sous-spectre discret de longueur finie, et on note _[0,M-1], tout ensemble de M valeurs successives (M fini) provenant de l'échantillonnage de  à une certaine pulsation w_f. Un spectre discret de longueur finie est composé de valeurs  qui forment un sous-ensemble des valeurs du spectre discret _Z. Un spectre discret de longueur finie ne peut pas contenir, lui non plus, toute l'information de départ.

C'est entre ces deux nouveaux objets que l'on va essayer de définir deux transformations inverses l'une de l'autre, que l'on appellera "transformation de Fourier discrète" et "transformation de Fourier discrète inverse". La raison en est double :

• D'après ce qui a été dit plus haut, il ne peut exister de transformation inversible entre un signal discret et un spectre discret que si aucun des deux ne contient l'information totale. Or, c'est le cas d'un signal discret de longueur finie et d'un spectre discret de longueur finie.

• En informatique, on ne sait traiter que des ensembles de données de cardinal fini, ce qui est le cas d'un signal discret de longueur finie et d'un spectre discret de longueur finie.

Pour obtenir un signal (ou un spectre) discret de longueur finie, deux cas se présentent :

• Cas 1 : si le signal (ou le spectre) analogique de départ est de largeur finie L (on dit aussi "à support compact"), alors une discrétisation à la pulsation d'échantillonnage  fournira un nombre fini  de valeurs non nulles, qui constitueront un signal (ou un spectre) discret de longueur finie.
• Cas 2 : sinon, après discrétisation à une pulsation w, il faudra conserver seulement N valeurs successives, et éliminer l'infinité de valeurs restantes.

Or on sait que, si un spectre analogique  est de largeur L finie, alors le signal analogique f(t) qui lui est associé est forcément de largeur infinie, ce qui signifie que si un spectre analogique rentre dans le cas 1, alors le signal analogique associé rentre forcément dans le cas 2. C'est une situation de ce type que nous allons détailler dans le paragraphe qui suit.

3.2.b. Choix des pulsations d'échantillonnage w_e et w_f.

Soit f(t) un signal analogique, tel que le spectre  soit nul en dehors de l'intervalle [w_min,w_max], de largeur L = w_max - w_min finie. On sait que l'on peut discrétiser le signal f(t) de manière à vérifier le critère de Shannon. Il suffit pour cela de choisir la pulsation d'échantillonnage w_e supérieur ou égal à L. Supposons par exemple que l'on choisisse w_e = L. En revanche, on ne peut pas également vérifier le critère de Shannon lors de la discrétisation du spectre. En conséquence, comment choisir la pulsation d'échantillonnage du spectre w_f ? C'est à cette question que nous allons tenter de répondre dans ce paragraphe.

Puisque le spectre rentre dans le cas 1 mentionné précédemment, sa discrétisation à la pulsation d'échantillonnage w_f fournira  valeurs non nulles, qui constitueront un spectre discret de longueur  finie. En revanche, comme le signal rentre dans le cas 2, il ne faudra conserver que N valeurs successives parmi l'infinité de valeurs non nulles fournies par la discrétisation à la pulsation d'échantillonnage w_e = L, pour obtenir un signal discret de longueur finie N. Comme on souhaite trouver une transformation inversible entre ce spectre de longueur finie et un signal discret de longueur finie, il faut nécessairement que l'un et l'autre contiennent la même information, c'est-à-dire le même nombre de valeurs, donc que :

En remplaçant M par , puis L par w_e, l'égalité (3.6) devient :

soit encore :

Maintenant que les deux pulsations d'échantillonnage sont connues, nous allons chercher une transformation permettant de relier les N valeurs non nulles du spectre discret à N valeurs successives du signal discret.

3.2.c. Peut-on exprimer les N valeurs non nulles du spectre discret en fonction de N valeurs successives du signal discret ?

Après discrétisation du spectre analogique à la pulsation d'échantillonnage w_f calculée précédemment, on obtient un spectre discret ne comportant que N valeurs non nulles, que l'on note _[0,N-1]. Soit m un entier appartenant à l'intervalle [0,N-1]. On souhaiterait pouvoir calculer  en fonction de N valeurs successives du signal discret (f_n)_Z. Est-ce possible ? C'est ce que nous allons voir dans ce paragraphe.

En partant de la définition de  :

Or, comme w_e vérifie par hypothèse le critère de Shannon, on sait que sur l'intervalle [w_min,w_max], d'après l'égalité (3.5) :

où f*(t) désigne la distribution introduite au paragraphe II. Comme  est bien dans l'intervalle [w_min,w_max], ceci permet de réécrire (3.8) sous la forme suivante :

qui devient, d'après l'expression (3.2) de la "transformée de Fourier" de f*(t) :

Cette expression se simplifie encore en remplaçant w_f par son expression (3.7) :

Il est clair qu'on ne peut pas exprimer  en fonction de N valeurs successives du signal discret (f_n)_Z, puisque toutes ses valeurs (il y en a une infinité) interviennent dans l'expression (3.9). En restreignant l'intervalle de variation de l'indice n à [0,N-1] dans cette somme infinie, on n'a au mieux que l'approximation suivante :

Nous pouvons maintenant donner la définition mathématique de la transformée de Fourier discrète, qui s'inspire fortement de (3.10).

3.2.d. Transformation de Fourier discrète et transformation de Fourier discrète inverse.

On appelle transformation de Fourier discrète l'application qui, à un signal discret (f_n)_[0,N-1] de longueur N finie, associe la suite (F_m)_[0,N-1] de N valeurs complexes, appelée transformée de Fourier discrète de (f_n)_[0,N-1], et définie par :

En comparant (3.11) à (3.10), on voit que F_m est une approximation de , ce que l'on peut exprimer en disant que "la transformée de Fourier discrète du signal discrétisé est une approximation de la discrétisation de la transformée de Fourier du signal" (au facteur T_e près).

Comme on pouvait s'y attendre, cette transformation est inversible, et on demande au lecteur d'admettre que la transformation de Fourier discrète inverse associe à une suite (F_m)_[0,N-1] de N valeurs complexes, la suite (f_n)_[0,N-1] de même longueur N, appelée transformée de Fourier discrète inverse de (F_m)_[0,N-1], et définie par la relation suivante, inverse de (3.11) :

Remarque :

La transformation de Fourier discrète inverse ne diffère de la transformation de Fourier discrète que par le signe de l'exponentielle et le facteur .

Il nous reste maintenant à présenter un des algorithmes permettant le calcul de transformées de Fourier discrètes par ordinateur, ce qui était notre objectif initial. Cet algorithme est "l'algorithme de la transformation de Fourier rapide".

4. Algorithme de la transformation de Fourier rapide (TFR).

Pour calculer la transformée de Fourier discrète d'un signal discret (f_n)_[0,N-1] à l'aide de la formule (3.11), il semble que l'on ait à effectuer N multiplications pour chaque valeur F_m, c'est-à-dire N² multiplications en tout. Pour une image numérique usuelle, de taille 256 x 256, cela fait 256⁴ multiplications, soit à peu près 4 milliards de multiplications ! Cette seule constatation démontre la nécessité de trouver un algorithme "rusé".

Il existe de nombreux algorithmes "rusés" permettant le calcul de la transformée de Fourier discrète d'un signal discret de longueur finie. Nous n'en présenterons qu'un, appelé "algorithme de la transformation de Fourier rapide" (TFR, ou FFT en anglais), dont la programmation en langage C sera effectuée lors des deux premières séances de travaux pratiques.

4.1. Notations.

Dans l'égalité (3.11), le facteur suivant apparaît :

On note :

L'égalité (3.11) se réécrit donc :

Remarques :

• Le facteur W_N est appelé "facteur de rotation".

• Ce facteur vérifie :

Exemple :

Si N = 8, alors :



Il existe de nombreuses redondances dans le calcul des différentes puissances de W₈. En effet, on vérifie facilement que :





...

Il suffit de calculer les 8 premières puissances de W₈ pour connaître toutes les puissances de W₈ apparaissant dans (3.13).

4.2. Description de "l'algorithme de la transformation de Fourier rapide".

4.2.a. Définition de la matrice T_N.

Dans ce paragraphe, nous allons mener les calculs sous forme matricielle, ce qui nous permettra de leur donner une expression beaucoup plus compacte. Une telle écriture est possible si on remarque que, dans (3.13), les valeurs de la suite (F_m)_[0,N-1] sont des combinaisons linéaires des valeurs de la suite (f_n)_[0,N-1], et que les coefficients de ces combinaisons linéaires sont des puissances du facteur de rotation W_N. On peut donc réécrire (3.13) sous la forme matricielle suivante :

où T_N désigne la matrice carrée, de taille N x N, dont l'élément courant vaut :

si m désigne le numéro de ligne et n le numéro de colonne.

La matrice T_N est facile à calculer grâce aux redondances citées précédemment, mais la multiplication matricielle (3.14) nécessite bien N² multiplications, comme cela a déjà été mentionné.

4.2.b. Comment choisir le nombre de valeurs du signal discret ?

L'algorithme de la transformation de Fourier rapide étant récursif, comme nous allons le voir, on a intérêt à choisir un nombre N qui soit une puissance de 2 :

C'est la raison pour laquelle les images numériques sont en général de taille 128 x 128 (128 = 2⁷), 256 x 256 (256 = 2⁸), 512 x 512 (512 = 2⁹), ...

Dorénavant, nous fixons N à la valeur suivante :

ce qui conduira à une profondeur de 3 dans la récursion. Par ailleurs, nous notons W le facteur de rotation W₈.

Le fait de raisonner sur cette valeur particulière de N n'enlevera rien au caractère général de ce qui suit.

4.2.c. Premier étage de la récursion.

En découpant la matrice T₈ en quatre quarts de taille 4 x 4, on peut calculer les valeurs de la suite (F_m)_[0,7] en deux moitiés. Les quatre premières valeurs de cette suite sont données par :

où W est une notation "allégée" pour W₈.

A ce niveau du raisonnement, deux remarques très importantes (mais qui ne sautent pas forcément aux yeux !) doivent être faites :

• En posant Y = W², la première matrice 4 x 4 qui apparaît dans le membre droit de (3.16) se réécrit :

Comme, par ailleurs :



on constate, d'après (3.15), que cette matrice n'est autre que T₄. Ainsi donc, on voit poindre le caractère récursif de ce calcul.

• Une deuxième remarque consiste à noter que la deuxième matrice 4 x 4 qui apparaît dans le membre droit de (3.16) peut s'exprimer très facilement en fonction de la première. En effet :

c'est-à-dire, d'après ce qu'on vient de montrer :

On peut donc réécrire (3.16) sous la forme suivante :

On demande au lecteur d'admettre que, par un calcul strictement similaire, on obtiendrait pour les quatre autres valeurs de la suite (F_m)_[0,7] :

Le nombre de multiplications nécessaires pour calculer les expressions (3.17) et (3.18) vaut :

• 16 pour le produit :

                (3.19)

qui apparaît deux fois ;

• 16 pour le produit :

               (3.20)

qui apparaît deux fois ;

• 4 enfin pour le produit :



qui apparaît également deux fois.

On constate donc qu'au lieu des 64 multiplications annoncées, on n'en est déjà plus qu'à 36.

4.2.d. Suite et terminaison de la récursion.

Les produits matriciels (3.19) et (3.20) peuvent être interprétés comme les calculs des transformées de Fourier discrètes des deux signaux discrets (f₀,f₂,f₄,f₆) et (f₁,f₃,f₅,f₇) de longueur 4 chacun, lesquelles transformées de Fourier discrètes peuvent être notées (G_m)_[0,3] et (H_m)_[0,3]. On peut donc appliquer à nouveau le raisonnement précédent, consistant à découper la matrice T₄ en quatre quarts de taille 2 x 2. On obtient alors les écritures matricielles suivantes pour le calcul de la transformée de Fourier discrète (G_m)_[0,3], par exemple :

Les produits matriciels correspondant à des transformées de Fourier de signaux discrets de longueur 2 peuvent être aisément développés. Par exemple :

où on doit effectuer une addition et une soustraction, mais plus aucune multiplication. Or on sait que le temps de calcul d'une addition ou d'une soustraction est négligeable devant celui d'une multiplication.

Combien de multiplications faut-il faire au final ? Dénombrons-les :

• Au premier étage de la récursion, il ne reste que les 4 multiplications correspondant au produit :

• Au deuxième étage, il reste 2 multiplications correspondant au produit :



et 2 autres correspondant au produit :



qui intervient dans le calcul de (H_m)_[0,3].

On n'a donc plus que 8 multiplications à effectuer, au lieu des 64 prévues initialement !

4.2.e. Complexité de l'algorithme.

De manière générale, quel est le nombre de multiplications nécessaires à l'algorithme que nous venons de décrire dans le cas particulier où N = 8, algorithme dit de la transformée de Fourier rapide ?

Supposons, comme nous l'avons préconisé plus haut, que N = 2^P. Le nombre d'étages de la récursion, qui était 2 avec N = 8, vaut  de manière générale. A chaque étage, on peut montrer que le nombre de multiplications est le même, et vaut . Il y a donc en tout  multiplications, au lieu de N².

Exemple :

Pour une image 256 x 256, au lieu des 4 milliards de multiplications déjà mentionnés plus haut, il n'y en a plus que 400000 environ, c'est-à-dire 10000 fois moins ! Comprenez-vous pourquoi cet algorithme est appelé "algorithme de la transformation de Fourier rapide" ?

4.3. "Algorithme papillon".

L'algorithme de la transformation de Fourier rapide peut être programmé de multiples façons, et nous allons en détailler une version particulière, appelée "algorithme papillon".

4.3.a. Ordre dans lequel les valeurs du signal discret sont traitées.

Reprenons le cas particulier d'un signal discret (f_n)_[0,7] de longueur N = 8. Dans l'algorithme de la transformation de Fourier rapide, il y aura trois niveaux de traitement des 8 valeurs de ce signal :

• Au "niveau 0" de la récursion, correspondant à l'égalité (3.14), ces valeurs apparaissent dans leur ordre naturel, c'est-à-dire :

f₀, f₁, f₂, f₃, f₄, f₅, f₆, f₇

• Au premier étage de la récursion, correspondant aux égalités (3.19) et (3.20), elles apparaissent dans l'ordre suivant :

f₀, f₂, f₄, f₆, f₁, f₃, f₅, f₇

• Enfin, au deuxième et dernier étage de la récursion, c'est dans l'ordre suivant qu'elles apparaissent :

f₀, f₄, f₂, f₆, f₁, f₅, f₃, f₇

Comme on commence toujours une récursion par son dernier étage, c'est-à-dire par sa terminaison (se rappeler, par exemple, l'algorithme récursif du calcul de la factorielle), il faudra avant toute chose réordonner les valeurs du signal discret. Il existe un algorithme permettant de faire cela, et qui procède en trois étapes :

• On traduit en binaire (sur P bits, avec N = 2^P) les indices entiers n des termes f_n de la suite.

• On inverse l'ordre des bits, ce qui signifie qu'on lit de droite à gauche l'écriture binaire réalisée précédemment.

• On traduit en décimal les suites de bits obtenues.

Pour convaincre le lecteur de la justesse de cet algorithme, appliquons-le sur un exemple :

Exemple :

Avec N = 8, donc P = 3, les trois étapes qui viennent d'être citées procèdent comme suit :



ce qui fournit bien l'ordre souhaité pour les indices.

4.3.b. Suite du calcul.

On va à nouveau décrire ce qui se passe dans le cas particulier où N = 8, en supposant que les valeurs du signal discret ont été préalablement réordonnées grâce à l'algorithme qui vient d'être décrit. La suite du calcul peut en fait être résumée par un dessin, en convenant d'adopter les conventions suivantes :

• Lorsque deux flèches convergent, cela signifie une addition.

• Lorsque, parmi deux flèches convergentes, l'une des deux est accompagnée d'un caractère -, cela signifie une soustraction.

• Lorsqu'enfin une flèche est surmontée d'une inscription du type x W⁰, cela signifie une multiplication par la valeur indiquée.

Voici le dessin correspondant au cas N = 8 :

Remarques :

• Le nom "d'algorithme papillon" vient d'une certaine ressemblance de ce dessin avec un papillon.

• Le gros avantage de cet algorithme est qu'il permet d'effectuer ce qu'on appelle un "calcul en place" : en effet, à aucun moment, on n'a besoin d'une place mémoire supérieure à celle qui permet de stocker le signal discret de départ. De tels calculs sont bien sûr très prisés en informatique !

• Cet algorithme a enfin l'avantage d'être assez facile à programmer, chose que nous proposerons de faire durant les deux premières séances de travaux pratiques.

Il reste maintenant à convaincre le lecteur de l'intérêt de calculer des transformées de Fourier, discrètes ou non. Cela sera l'objet du quatrième et dernier chapitre de ce cours.