Brevet FR2849928A1 — publié le 20040716

i Procédé de reconnaissance d'hélicoptère par analyse en ondelettes du signal radar rétrodiffusé par les pales du rotor principal La présente invention concerne un procédé de reconnaissance d'hélicoptère par analyse en ondelettes du signal radar rétrodiffusé par les pales du rotor principal. Le signal radar rétrodiffusé par un hélicoptère est principalement constitué des échos de la cellule et des pièces en rotation à savoir, rotor 10 et turbines éventuelles, auxquels se superposent périodiquement les flashes des pales du rotor ainsi que ceux de l'anti-couple qui provoquent des signaux transitoires à chaque instant de passage d'une pale dans le plan orthogonal aux faisceaux du radar. L'onde modulée qui en résulte peut servir de base à des procédés 15 d'identification permettant de reconnaître le type de l'hélicoptère. Cependant, une analyse temporelle directe du signal s'avère insuffisante. De même, une analyse de Fourier simple de ce signal ne convient pas du fait du caractère transitoire du phénomène à analyser. Le but de l'invention est de pallier les inconvénients précités. A cet effet, l'invention a pour objet un procédé de reconnaissance d'hélicoptère par analyse en ondelettes du signal radar rétrodiffusé par les pales du rotor principal caractérisé en ce qu'il consiste à effectuer une analyse temps fréquence par transformée en ondelettes du signal rétrodiffusé, à détecter les flashes des pales de l'hélicoptère dans un 25 espace temps-fréquence, à analyser à long terme l'image temps- fréquence obtenue pour estimer le nombre de pales ainsi que la vitesse de rotation du rotor et à analyser de façon fine les flashes des pales repérés dans l'espace temps-fréquence pour affiner la classification obtenue par l'analyse à court terme. L'invention a également pour objet un dispositif pour la mise en oeuvre du procédé précité. D'autres caractéristiques et avantages de l'invention apparaîtront dans la description qui suit faite en regard des dessins annexés qui représentent: - la figure 1 un schéma général d'un système de reconnaissance d'hélicoptère selon l'invention; - la figure 2 un schéma synoptique général pour illustrer le procédé de reconnaissance selon l'invention; la figure 3 un exemple d'une configuration de rétrodiffusion d'une pale éclairée par un faisceau radar; - la figure 4 un exemple de localisation temps-fréquence d'une ondelette analysante; - les figures 5 et 6 des graphes représentatifs d'ondelettes SINC10 GAUSS mises en oeuvre dans le procédé selon l'invention; - la figure 7 un graphe temps-fréquence représentant schématiquement la variation du module d'une transformée d'ondelette mise en oeuvre dans le procédé selon l'invention; - la figure 8 un exemple de réalisation d'un dispositif pour la mise 1 5 en oeuvre du procédé selon l'invention; - la figure 9 un mode de localisation dans l'espace temps-fréquence de fenêtres d'énergie de flashes avançant et reculant; - la figure 10 une illustration d'une localisation grossière de flashes; - la figure 11 un diagramme schématique du module de l'image temps- fréquence autour du flash; - la figure 12 une coupe en fréquence d'une image d'analyse au niveau d'un flash; - la figure 13 un graphe pour illustrer le suivi de la courbe cosinu25 soldale dans l'espace temps-fréquence (composante à long terme du signal rétrodiffusé). Le dispositif selon l'invention est destiné, comme le montre la figure 1, a être relié à un dispositif connu de traitement distance 2 composant habituellement les récepteurs radars. Le procédé se déroule de la 30 manière représentée aux étapes référencées de 3 à 7 sur la figure 2. Il consiste selon l'étape 3 à effectuer un filtrage des échos fixes sur le signal temporel représentant le signal rétrodiffusé par l'hélicoptère dans la case distance correspondant à sa position pour effectuer ensuite à l'étape 4 un calcul de transformée en ondelettes sur le signal obtenu en sortie du filtrage réalisé à l'étape 3. Le signal obtenu à la fin du traitement de transformée en ondelettes de l'étape 4 est exploité à l'étape 5 pour détecter les flashes et aux étapes 6 et 7 pour effectuer respectivement une analyse à long terme et une analyse à court terme du signal transformé en ondelette et des flashes détectés. Cette analyse est basée sur le fait que le signal radar rétrodiffusé par un hélicoptère est principalement constitué par les échos de la cellule ou de la carlingue et des pièces en rotation à savoir, le rotor et les turbines, auxquels se superposent périodiquement les flashes des pales du 10 rotor ainsi que ceux produits par l'hélice d'anti-couple formés par les signaux transitoires correspondant aux instants de passage d'une pale dans le plan orthogonal au faisceau du radar. Le phénomène de rétrodiffusion d'une onde radar sur les pales du rotor principal peut être décrit suivant le modèle de la figure 3, qui assimile chaque pale du rotor princi15 pal d'un hélicoptère à un barreau rectiligne 8, de longueur L et de coefficient de rétrodiffusion uniforme, tournant autour d'une de ses extrémités o à une vitesse angulaire n constante, et en considérant que la pale est éclairée de façon uniforme par l'onde radar, le plan d'onde 9 formant sur la figure 3 un angle 0 avec le barreau rectiligne 8. En écrivant que le 20 signal émis par le radar est un signal de la forme E(t) = S(t). e 2,ifot (1) o S(t) est un signal réel formé d'une série de puises rectangulaires de fréquence Fe et fo est la fréquence porteuse de l'onde radar. Le signal réfléchi par chaque élément dx du barreau 8 est dans ce cas un signal de 25 la forme dux(t) = c.cos. Qt.e(t- x(t)).dx (2) dans laquelle Rx(t) représente la distance entre le radar et le point de coordonnées X correspondant à l'élément rayonnant dx et a le coefficient de rétrodiffusion du barreau 8. En désignant par R la distance radar30 hélicoptère la distance Rx(t) s'écrit encore Rx (t) = R + x. sin( Qt) (3) et la relation (2) devient (dans l'hypothèse bande étroite après démodulation par la porteuse fo) --du2 - (e -. 2ejfo i du x(t) =c. cos 0. s(t - -2 Jf e -2 c X Xsin e - 4. jfoe c dx (4) La relation (4) montre que le signal dux(t) est équivalent à un signal rétrodiffusé par une cible stationnaire de surface équivalente de rétrodiffusion égale à a(, située à une distance R, du radar modulé par une fonction dmx(t) définie par la relation _ 4l jx sin(n t) A dmx(t) = cos(nt)e U (5) o X0 = c/f0 est la longueur d'onde du radar. Dans ces conditions la contribution totale d'une pale s'obtient en effectuant la somme des contributions élémentaires conduisant au signal modulant global. Ce 10 signal est de la forme L L - 2njL sin(a t) o 2n L sin(út) m(t) = J dmx(t) = L.cos(2t).e.sinc X L) (s dans laquelle sinc représente l'expression de la fonction sinus cardinal. Le signal réfléchi global obtenu u(t) apparaît, en développant les relations 15 (4) et (6) précédentes, comme un signal égal au signal incident décalé du retard hélicoptère-radar et pondéré par la fonction modulante m(t) encore appelé ci-après fonction d'ouverture suivant la relation s( t2R)-2. fo2R u(t) = .s( t--e c. m(t) (7) La fonction d'ouverture m(t) de la relation (6) est composée du 20 produit de trois termes qui sont: - un terme en cos(Qt) faisant apparaître une fréquence pure réelle fR = Q/2x qui traduit la variation de la surface équivalente rétrodiffusée de la pale pendant une période de rotation - une exponentielle complexe de fréquence instantanée 25 f(t) = -Lcos(E2t) /X0 (8) qui suit une loi de variation cosinus. Cette fréquence est maximale lorsque la pale passe dans le faisceau radar pour 0 égal 0 ou X c'est-àdire au moment de l'éclair de pale. Elle est positive lorsque la pale est avançante et négative lorsque la pale est reculante. L'estimation de la 6) période de la fréquence f(t) (TR=1/fR) permet de remonter à la vitesse de rotation du rotor dans l'analyse à long terme des images temps-fréquence effectuée à l'étape 6. - un terme en sinus cardinal dont la largeur du lobe principal est 5 approximativement A = Xo/(Ln) ce terme contribue à la formation de l'éclair de pale et donne une réflection maximale au passage dans le faisceau radar avec une atténuation très rapide ensuite à -13 dB dès le deuxième lobe. A titre d'exemple pour une longueur d'onde Xo = 15cm et une 10 vitesse en bout de pale (v = QL) de l'ordre de 200 mètres par seconde la largeur du lobe principal est alors de 0,75 ms. Cette valeur lorsqu'elle est comparée à la période de rotation du rotor TR = 170 ms (pour une vitesse moyenne Q de 350 tours par minute) permet de souligner l'aspect transitoire et non stationnaire du phénomène. En posant dans la relation (8) i1 = 2QL/XO = 2V/10 l'expression du signal modulant m(t) peut s'écrire approximativement pendant l'intervalle de temps de passage de la pale dans le faisceau radar m(t) = L.e-7fint. sinc(7ntt) (9) Son spectre est constant dans la bande de fréquence (-rl,0) et nul 20 ailleurs. Les fréquences qui entrent en jeu sont négatives lorsque la pale considérée est reculante. Dans le cas d'une pale avançante le spectre est non nul dans la bande [0,il]. Dans ces conditions l'analyse temps-fréquence du signal de pale aux instants de passage des pales dans le faisceau radar donne une information dans le demi plan des fréquences posi25 tives pour une pale avançante et une information dans le demi plan des fréquences négatives pour une pale reculante. Naturellement les considérations précédentes ne valent que pour le cas o le radar n'intercepte qu'une seule pale à la fois, soit avançante soit reculante comme c'est généralement le cas pour des hélicoptères possédant un nombre impair 30 de pales. Dans le cas d'hélicoptères à nombre pair de pales, le faisceau radar intercepte à la fois une pale reculante et une pale avançante. Dans ce cas le signal rétrodiffusé global qui est pris en considération est la somme des signaux réfléchis par chacune des pales au moment du passage dans le faisceau radar et ceci donne dans le plan temps-fréquence de l'information sur les fréquences positives et négatives. Naturellement l'analyse qui vient d'être faite fait appel à une modélisation très idéalisée de la réflexion radar d'une pale en mouvement notamment en ce qui concerne la géométrie de la pale du fait que certaines déformations peu5 vent se produire à grande vitesse et la valeur constante du coefficient de rétrodiffusion le long de cette même pale. On a aussi considéré que la rétrodiffusion était identique sur les bords d'attaque et de fuite des pales ce qui généralement n'est pas vérifiée en pratique. Néanmoins les considérations précédentes permettent de dégager les composantes essen10 tielles du signal réfléchi qui permettent aussi de fixer des ordres de grandeur tant qu'à l'étendue fréquentielle et temporelle d'analyse des signaux ainsi que la cadence d'échantillonnage qu'il est nécessaire de mettre en oeuvre. Dans la pratique le signal total réfléchi par un hélicoptère supposé 15 en vol stationnaire se compose des échos de la cellule et des pièces en rotation auxquelles sont superposés les signaux engendrés par les pales du rotor telles que celles-ci viennent d'être décrites. A ce signal reçu par le radar il faut bien sr ajouter les bruits de sol qui correspondent à la case distance de l'hélicoptère ainsi que le bruit du récepteur. Après 20 démodulation complexe le signal reçu est un signal de la forme r(t) = s[t - 2-11e c (h(t) + f(t)) + n(t) (10) o s(t) représente le train d'impulsions émis par le radar, h(t) désigne le signal complexe qui caractérise la réflexion de l'hélicoptère, f(t) le signal caractérisant le fouillis du sol, n(t) désigne le signal qui caractérise le 25 bruit du récepteur et R est la distance radar-hélicoptère. La première opération effectuée au niveau du radar est le traitement distance c'est-à-dire le filtrage adapté à l'impulsion élémentaire dont le but est de maximiser le rapport signal à bruit dans chacune des cases distance balayées. Ce traitement est effectué à la cadence Fe de la fréquence de récurrence du 30 radar. Comme la durée de l'impulsion radar élémentaire est de quelques microsecondes l'effet Doppler de chaque pale est considéré comme constant sur cet intervalle. En conséquence le signal recueilli dans la case distance correspondant à un hélicoptère n'est autre en très bonne approximation qu'un échantillonnage à la fréquence Fe du signal h(t)+f(t) auquel il faut bien entendu ajouter un bruit dont la puissance a été rendue minimale. C'est ce signal qui sert de base à la reconnaissance. Comme il est intéressant de s'intéresser exclusivement au signal 5 de pale dans le cas d'un vol stationnaire il est nécessaire dans un premier temps pour des raisons de dynamique d'éliminer le fouillis de sol et les échos basses fréquences de la cellule. Ceci est obtenu de la façon figurée à la figure 2 par un filtrage de bande passe-haut des données. Ce filtre doit posséder une fréquence de coupure suffisamment faible de 10 l'ordre Fc = 0,01 en fréquence normalisée ce qui correspond à une vitesse de coupure de 15 n/s pour Fe = 20 KHz et X0 = 15cm. Ce filtrage est réalisable très simplement à l'aide d'un filtre récursif du premier ordre. L'analyse en ondelettes qui est effectuée à l'étape 4 sur la figure 2 consiste à effectuer une analyse temps-fréquence avec un rapport Af/f 1 5 constant adapté au signal filtré. Cela revient à étudier le signal au travers d'une fenêtre temporelle de taille variable en fonction de la fréquence analysée. La taille est réduite pour les hautes fréquences pour lesquelles c'est la résolution temporelle qui prime, la résolution fréquentielle étant alors réduite, la taille est importante pour les basses fréquences pour les20 quelles au contraire c'est la résolution fréquentielle qui est importante la résolution temporelle étant réduite. La transformation en ondelettes utilisée est une transformation de la forme T x (,ra) = 1\ J x(t). g (-). dt (Il1) o g(t) représente rondelette analysante, t le paramètre temporel et a, 25 un nombre réel appartenant à R*, le paramètre d'échelle. Cette transformée est inversible pourvu que g(t) soit suffisamment régulière et du type passe- bande. La relation précédente permet de calculer la projection du signal x(t) sur un ensemble de fonctions analysantes du type g(t, Y, a) = g(t * 9) (12) Ces fonctions analysantes sont localisées dans le plan temps- fréquence sur un pavé rectangulaire centré sur un point de coordonnées (r, FO/a) sur la figure 4 o F0 est la fréquence centrale de r'ondelette mère définie par la relation J Vj|g(V)l 2dv F0 = (13) I I(v)I 2dv o0 Les dimensions de ce pavé qui déterminent les résultats temporels et fréquentiels de l'analyse sont suivant l'axe temporel At = 21a1. F0 (14) et suivant l'axe fréquentiel Av = 2( F01al (15) a to et C F 0 représentent respectivement la dispersion temporelle 10 et fréquentielle de r'ondelette mère dans l'hypothèse o la valeur absolue de g(t) est une fonction centrée temporellement sur 0. La dimension globale de la case de résolution At.Av est constante et ne dépend que de r'ondelette mère soit: At.Av = 4.a t 0' a F 0 L'adaptation de la formulation continue précédente en signaux 15 échantillonnés a lieu suivant la relation T x (,a) = 1 x(nT, e)( a (16) ja1n a o Te = 1/Fe est la période d'échantillonnage du signal. La discrétisation de l'espace temps échelle a lieu de la manière suivante. La discrétisation naturelle de l'axe temporel consiste à prendre t = k.Te (k étant un 20 entier). Le paramètre d'échelle a est choisi de manière que a = 2am pour a > 0 a = -2am pour a 0 est alors déterminé en fonction de la résolution d'analyse le long de l'axe fréquentiel qui correspond à la dispersion fréquentielle de 25 I'ondelette mère. La transformée finale est - pour les échelles positives a Tx(k,m) = 2 2E, x(n).g ((n- k) 2-am.Te) (17) n - pour les échelles négatives Tx(k,m) = 2 2T x(n).g (k- n)2-a r.Te) (18) n Tx(km) s'obtient ainsi par convolution du signal à analyser avec un échantillonnage de r'ondelette mère à la période Tc = 2-am.Te. Le sup5 port temporel de r'ondelette qui fixe la fenêtre signal à prendre en compte est maximal pour les grandes valeurs de m c'est-à-dire pour l'analyse des basses fréquences. La fréquence nulle ne peut être analysée car elle correspond en effet à une valeur infinie de m. L'analyse s'effectue dans ces conditions sur 2 M points fréquentiels répartis sur des 10 intervalles (-Fmax, Fmini) et (Fmini, Fmax) en valeurs normalisées avec M points sur chaque intervalle. Fminj doit être choisi en pratique égal à la fréquence de coupure du filtre d'échos fixes. La fréquence Fmax quand à elle est déterminée en fonction de la fréquence maximale présente dans le flash qui dépend de la vitesse en bout de pale et de la fréquence de 15 récurrence du radar. Par exemple pour v = 200m/s et Fe = 20 KHz, Fmax peut être égal à 0,25. En notant par F0 la fréquence centrale de r'ondelette l'intervalle de variation de m est donnée par la relation (-1/ ().Log2(FMAx.Fe / F0) 1 + Log 2(FMAX / FMIN) / L0g 2((Fo + 6 FO) / Fo - Fo) (22) Dans la pratique M peut être choisi compris entre 128 et 256 pour que la condition précédente soit largement vérifiée. Cependant la fonction précédente ne peut être utilisée comme ondelette analysante du fait de sa mauvaise localisation dans le plan temps-fréquence et il est néces5 saire de la pondérer par une fenêtre. Une bonne localisation conjointe dans le plan temps-fréquence peut être obtenue au moyen d'une fenêtre gaussienne et rondelette résultante baptisée ondelette SINCGAUSS a alors la forme g(t) = e 2irjFot sin(ict) et2 (23) nt Cette ondelette est centrée temporellement sur 0 et sa fréquence centrale est Fo. Sa dispersion fréquentielle est de l'ordre de a FO = 0,3. Cette ondelette est pratiquement admissible dès que sa fréquence Fo est supérieure à 0,6. Cependant l'analyse à long terme et à court terme effectuée aux étapes 6 et 7 de la figure 2 nécessite de travailler sur deux 15 fréquences d'ondelettes différentes, en remarquant que le fait d'augmenter la valeur de Fo revient en fait à accroître la résolution fréquentielle de l'analyse et donc à diminuer sa résolution temporelle. Les parties réelles de cette ondelette sont représentées aux figures et 6 pour les fréquences F0 respectivement égales à 0,64 et 2,56. Pour la mise en oeuvre de l'invention d'autres ondelettes telle que par exemple rondelette de Morlet peuvent tout aussi bien être utilisées, bien que celles-ci conduisent à de moins bons résultats du point de vue de la classification notamment en ce qui concerne l'analyse à court terme. L'ondelette SINCGAUSS telle que définie précédemment à partir d'un modèle de signal rétrodiffusé par une pale seule, convient de la même manière par linéarité pour une analyse de signaux issus d'un hélicoptère à nombre pair de pales. Dans la suite du procédé seul le module de la fonction transformée 30 est utilisé pour la classification des hélicoptères. Ce module est quantifié sur 256 niveaux logarithmiquement pour l'analyse à long terme et linéairement après écrêtage pour l'analyse à court terme. Ceci conduit à une image de la forme de celle qui est représentée à la figure 7. 1 1 L'analyse à long terme qui est effectuée à l'étape 6 a pour but de mettre en évidence dans le plan temps-fréquence la composante exponentielle de fréquence variable et d'estimer ces paramètres. Cette fréquence pure suit une loi de variation en cosinus de la forme: 5 F(t) c cos n t et de période TR = 2n/0. Elle est maximale lorsque la pale passe dans le faisceau radar, c'est-à-dire au moment du flash, positive lorsque la pale est avançante et négative lorsqu'elle est reculante. Cette fonction laisse donc une trace en forme de cosinusode sur l'image du module de la transformée temps-fréquence lorsque les axes temporels et 10 fréquentiels sont gradués en échelle linéaire. Cette trace passe par le "sommet" S du flash, correspondant à la fréquence maximale du flash. Néanmoins l'amplitude de cette trace décroît très vite en 1/t au fur et à mesure que le point sur la courbe s'éloigne du flash c'est pourquoi, il est nécessaire pour bien faire apparaître le phénomène, de quantifier loga15 rithmiquement son module. Pour effectuer cette analyse il est aussi nécessaire de disposer d'images issues d'analyse temps-fréquence possédant une bonne résolution fréquentielle pour que le pic qui correspond à la courbe F(t) soit suffisamment contrastée par rapport au bruit ce qui signifie que la fréquence d'ondelette fo soit suffisamment grande, d'un 20 ordre supérieur à 3 par exemple. L'analyse à court terme qui est effectuée à l'étape 7 et qui se concentre sur le flash permet d'extraire dans le plan temps-fréquence certains paramètres de formes qui caractérisent ce signal transitoire. Cette analyse vient en complément de l'analyse à long terme pour obte25 nir une discrimination suffisamment fine des hélicoptères. Elle a lieu sur le module de la transformée temps-fréquence quantifiée linéairement après écrêtage bas et haut entre 30% et 60% de sa dynamique totale. Plus précisément, si XS(f,t) représente le module de la transformée tempsfréquence, bornée supérieurement par Xmax, l'opération de seuil30 lage préliminaire s'exprime sous la forme (en désignant par Xo le signal de sortie) XO(f,t) = Xmax si XS(f,t) > 0,6 Xmax XO(f,t) = si XS(f,t) Il résulte de ce qui précède que chacune des analyses à long terme et court terme effectuée aux étapes 6 et 7 nécessite l'utilisation d'une transformée en ondelettes particulière dont seuls les paramètres diffèrent. Elles nécessitent par ailleurs la détection préalable des flashes de pales dans le plan temps-fréquence. Celle-ci est effectuée sur les images 15 servant l'analyse à court terme du fait que les flashes y sont mieux mis en évidence. Les indices temporels des flashes détectés sont ensuite transmis à l'analyse à long terme comme schématisé à la figure 2. Un dispositif pour la mise en oeuvre du procédé qui vient d'être décrit est représenté à la figure 8. Celui-ci comporte, un dispositif de 20 transformation en ondelettes 10 couplé à un dispositif 11 d'écrêtage de la dynamique totale du module de la transformée en ondelettes obtenue fournie par le dispositif 10. Un dispositif de quantification linéaire 12 est couplé par son entrée à la sortie du dispositif d'écrêtage 11 et par sa sortie à l'entrée d'un dispositif de détection de flash 13. Un dispositif de 25 calcul de transformée en ondelettes à fréquence Fo élevée 14 est couplé par sa sortie à l'entrée d'un dispositif de quantification logarithmique 15 à grande dynamique de 60dB par exemple par rapport au maximum de la transformée. Un dispositif d'analyse à long terme 16 est couplé par une première entrée à la sortie du dispositif de quantification logarithmique 30 15 et par une deuxième entrée à la sortie du dispositif de détection de flashes 13. Un dispositif d'analyse à court terme 17 est couplé par son entrée à la sortie du dispositif de détection de flashes 13. Le dispositif d'analyse à long terme 16 délivre sur sa sortie des signaux indiquant le nombre de pales détectés de l'hélicoptère et la vitesse angulaire Ä des pales. Le dispositif d'analyse à court terme 17 fournit la position du centre de gravité du flash ainsi que les paramètres de courbure du flash. Avant d'effectuer les opérations de transformée en ondelettes, le signal radar est normalisé de façon que sa puissance soit égale à l'unité. Cette normalisation permet de traiter de la même façon les signaux issus d'hélicoptères à des distances variées et évite tout problème de mise à l'échelle dans la phase d'analyse. Le dispositif de détection de flashes 13 travaille sur le module de la transformée en ondelettes à fréquence Fo faible après écrêtage et 10 quantification linéaire effectués par les dispositifs 11 et 12. La détection a lieu suivant deux étapes composées, d'une détection et d'une localisation fine. La détection est basée sur l'utilisation d'un critère énergétique et la localisation fine est basée sur un filtrage adapté au modèle théorique de rétrodiffusion développé précédemment. 1 5 La première étape de détection permet de détecter et de localiser grossièrement dans le temps les flashes des pales reculantes ou avançantes. La précision de localisation est de l'ordre de 5 à 10 échantillons. Elle procède de manière représentée à la figure 9 par calcul de l'énergie moyenne sur une fenêtre glissante s'étendant sur 2L+1 points en tem20 porel et sur l'ensemble du demi-plan fréquentiel concerné positif pour les flashes avançant et négatif pour les flashes reculant, soit sur M = 128 points en pratique. Cette approche résulte des considérations du modèle théorique de la rétrodiffusion décrit précédemment qui considère que le flash est un phénomène qui est localisé dans le temps et dont le spectre 25 est théoriquement plat et maximal dans une bande de fréquence [0,,] ou [-i',0]. La valeur de l'énergie Ei est obtenue en effectuant la valeur moyenne de l'image à l'intérieur de la fenêtre d'analyse par la relation 1 L Tp E(i) = M + )E L x(i - k, j) (24) f 2IL + 1) k= -L j=:ri Si l'énergie ainsi calculée est supérieure à une valeur de seuil 30 déterminée, par exemple S1 pour les flashes positifs et S2 pour les flashes négatifs, le dispositif se place dans un état d'alarme et y demeure tant que l'énergie reste supérieure au seuil. Le dispositif quitte l'état d'alarme dès que l'énergie retombe en-dessous du seuil. Un flash est alors détecté et une première estimation de sa localisation temporelle est obtenue correspondant à l'indice i du maximum d'énergie à l'intérieur de la zone d'alarme. Cette procédure est résumée par le graphe de la figure 10. La largeur temporelle de la fenêtre est optimisée pour régler le compromis, probabilité de détection par rapport à la probabilité de fausse alarme pour une valeur fixée de seuil. Une fenêtre trop étroite a tendance à créer des paquets de fausses alarmes autour de la position réelle du flash et à l'opposé une fenêtre trop large a tendance à masquer certains 10 flashes. Pratiquement, la largeur de la fenêtre est choisie en fonction de la résolution temporelle de l'analyse par ondelettes, déterminée par la fréquence Fo de rondelette, afin de correspondre à la largeur attendue duflash. Cependant comme cette étendue temporelle n'est pas constante et diminue lorsque la fréquence augmente, un compromis est à trouver 1i5 entre les largeurs, minimale et maximale. En ce qui concerne le seuil de détection, il peut être fixé en supposant que lorsque la fenêtre d'énergie se trouve centrée sur le flash, tous les pixels correspondant à une valeur de fréquence inférieure à TI pour un flash positif doivent avoir une valeur très proche de 255, tous les autres ayant une valeur très proche de 0. 20 Cependant à cette étape de détection la localisation du flash n'est pas assez précise car le maximum d'énergie ne concide pas exactement avec le centre de phase du flash. Du fait que les traitements ultérieurs de classification ne sont pas tous invariants par translation temporelle, il est nécessaire de procéder à une localisation plus fine du flash par une 25 corrélation fine autour de la position détectée avec une image de référence obtenue par filtrage adapté. L'indice temporel donnant le maximum de corrélation définit alors la localisation finale du flash. En cas d'ambiguté se produisant par exemple pour deux positions de fenêtres donnant un même maximum de corrélation, c'est la position la plus proche de 30 celle fournie par l'étape de détection qui est privilégiée. L'image de référence est obtenue à partir du modèle théorique du signal rétrodiffusé par une pale déterminée précédemment. Seule la partie du module correspondant aux fréquences positives est conservé pour former l'image de référencee. Cette référence est utilisée telle qu'elle pour localiser les flashes positifs après retournement autour de l'axe fréquentiel pour le traitement des flashes négatifs. Cette étape de localisation fine améliore l'estimation de la position du flash dans 30 à 40% des cas. Le dispositif d'analyse à long terme 7 permet d'estimer la période 5 TR de la courbe en cosinus de la figure 7 qui représente la trace dans le plan temps-fréquence de la fréquence pure variable du signal de pale, avec pour objectif d'une part, d'estimer la vitesse de rotation Q du rotor et d'autre part, d'estimer son nombre de pales NP. Le nombre de pales du rotor est déterminé par la relation N = nint( R / TF) (25) o nint est la fonction "plus proche entier", tR la période estimée de la cosinusoide et tF l'estimation de l'intervalle de temps TF séparant deux flashes avançant ou reculant successifs. L'idéal pour estimer TR est de procéder au suivi de la cosinusode 15 dans l'image temps-fréquence à partir du sommet d'un flash et sur une période complète de rotation. Mais une telle approche n'est pas envisageable en pratique, car le suivi s'avère impossible à basse fréquence. En effet lorsque sa fréquence décroît en valeur absolue au fur et à mesure que l'on s'éloigne du flash, l'amplitude de sa trace décroît également si 20 bien qu'en-dessous d'une certaine fréquence celle-ci se trouve noyée dans le bruit. Par ailleurs, comme la transformée en ondelettes amène à utiliser un axe fréquentiel gradué logarithmiquement la trace théorique de F(t) présente une asymptote lorsque la pale passe dans l'axe du faisceau radar (F(t) = 0). Enfin, les cosinusodes correspondant à deux flashes 25 successifs se croisent dans l'image au milieu d'un intervalle les séparant ce qui risque de gêner considérablement le suivi de ces courbes. En fait l'estimation de la période TR n'est réalisable en pratique que sur la partie haute fréquence de la cosinusode dans le voisinage direct du flash, plus précisément sur une fenêtre temporelle (-t, +X) centrée sur le flash pour 30 t inférieur à TF/2. Pratiquement t est choisi aussi proche que possible de TF/2 afin d'obtenir une estimation suffisamment précise de TR. Il est nécessaire en effet de travailler suffisamment loin du flash pour que la courbure de la cosinusode puisse produire un effet. Néanmoins afin de limiter la charge de calcul, il est encore possible de travailler sur une analyse par ondelettes sous-échantillonnées temporellement. L'estimation de TR procède alors en deux étapes, une première étape qui consiste à rechercher le sommet du flash, c'est-à-dire le point 5 qui correspond à la fréquence maximale d'un flash avançant ou à la fréquence minimale d'un flash reculant et suivant une deuxième étape à estimer TR par identification de la cosinusode à partir de son sommet par un suivi de la courbe dans l'image associé à une identification par moindres carrés. Pour le traitement des flashes positifs, par exemple, la recherche du sommet du flash a pour but d'estimer l'ordonnée du point S correspondant à la fréquence maximale du flash comme décrit sur la figure 11. Un exemple tiré d'un enregistrement réel d'hélicoptère de la variation du module de l'image d'analyse en fonction de la fréquence au niveau du 1 5 centre de phase d'un flash (fréquence positive: ligne 0 = Fmax, ligne 128 = Fmini) est figuré à la figure 12. La courbe représentée apparaît peu bruitée, néanmoins il n'en est pas toujours ainsi en pratique et on peut voir apparaître des pics parasites en haute fréquence ou voir le niveau chuter relativement bas à l'intérieur du flash mettant ainsi en évi20 dence les difficultés de définition du sommet S. En fait, il est préférable de définir la position du sommet S comme le point à partir duquel l'amplitude se met à chuter rapidement sous un certain seuil S3 quand l'axe fréquentiel est parcouru de la fréquence Fmini vers la fréquence Fmax. Pour assurer une meilleure détection du sommet S il est aussi possible de 25 limiter le balayage en fréquence en fonction d'informations a priori sur la vitesse en bout de pale. En pratique le seuil S3 peut être pris égal à 1 80 ce qui représente un seuil à 3 dB par rapport à la valeur maximale 255 du module de la transformée qui peut être supposée atteinte dans le flash. La notion de rapidité de chute quant à elle, est évaluée par un cal30 cul de gradient le long de l'axe fréquentiel par filtrage de type CannyDériche. Ce filtre est un filtre à réponse impulsionnelle infinie définie par la relation h(n) = K.n.e-I|InI o K est une constante de normalisation choisie de telle sorte que l'énergie du filtre soit égale à l'unité soit Eh2(n) = 1 soit: K= (1- e2) (26) 2e-2a (1+ e2a) a peut être pris égal à 1. Ce filtrage peut être effectué sous une forme récursive. En notant Xn le signal d'entrée du filtre et Yn le signal à sa sortie, Yn peut être considéré comme la somme d'une partie causale yp(n) et d'une partie anticausale yM(n) définies comme suit: yp(n) = a.x(n-1)-b1 yp(n-1)-b2 yp(n-2) et yM(n) = -a.x(n + 1)-b1 yM(n + 1l)-b2 YM(n + 2) avec a = K e-a, b1 = 2ea et b2 = e-2e Après filtrage monodimensionnel de la colonne de l'image tempsfréquence correspondant au centre de phase du flash, obtenue par l'opé15 ration de détection, le sommet S correspond au point donnant la valeur Y(S) maximale sous la contrainte que X(S) est supérieur au seuil S3. A partir de ce sommet S, TR est estimé à partir d'une identification par moindres carrés de la cosinusode. Cette estimation a lieu en utilisant l'expression de la fréquence normalisée en fonction du numéro de ligne 20 dans l'image du module de la transformée en ondelettes. Cette image possède NFREQ lignes fréquentielles numérotées de 1 à NFREQ. La ligne 1 correspond à Fmax et la ligne Nfreq correspond à Fmini. La fréquence normalisée associée à la ligne m s'écrit alors ( ( I (FMAX.Fe) FO -aIm-1+l+int.I-.og2 F0 I f m- 1+ nint 109 2 (27) Fe o F0 est la fréquence de l'ondelette et a = 1o92(FMAX/FMiN)/(NFREQ-1). Néanmoins il est possible d'utiliser également l'expression approchée suivante (en négligeant la fonction nint) fm = FMAX.2-a(m-1) L'erreur commise est en effet de l'ordre de 0,5%. L'estimation de TR consiste alors à identifier la courbe F(t) = AcosQlt o l'instant t = 0 correspond au centre de phase du flash, à partir d'un ensemble de P points de mesures Mi = (Ti, Fi) relevé 5 dans l'image autour du flash. Ayant estimé A et úQ, TR s'obtient alors par la relation TR = 27r/Q. On effectue un calcul normalisé par rapport à la cadence d'échantillonnage Fe en remplaçant A par AN = A/Fe et O par ON = O/Fe. Les points de mesure Mi ont alors pour coordonnées: - Ti = K.ni, indice temporel normalisé o K est le facteur de sous10 échantillonnage de l'analyse et ni varie dans un intervalle [-L,L], o L est le nombre de points de chaque côté du centre du flash. et Fi désigne la fréquence normalisée du point Mi obtenue selon l'expression précédente à partir de l'indice de ligne mi du point Mi. Enfin pour rendre la résolution du problème linéaire il est néces15 saire de remplacer la fonction cosinus par son développement limité à l'ordre 2 par la relation F(t) = AN - ANn N.t 2 (28) Cette approximation même si elle ne se justifie pas toujours, étant donné la taille de la fenêtre temporelle prise en compte, donne toutefois 20 des résultats tout à fait acceptables. En posant BN = QN2 le problème se réduit à trouver AN et BN minimisant l'erreur C(ANBN) = ú(F - - AN AN.BN.Ti2 / 2)2 (29) Mij Pour un flash négatif, il suffit de remplacer AN par -AN. La minimisation de C est alors immédiate et conduit à 25 AN = -al/a2 avec a1 = .n fi.n - in (30) i i i i a2 P n P.nn2)et 2 Xn 2 Xf1.n1 -niZ Kz.-AN P = 2L + 1 représente le nombre de points de mesure pris en 30 compte. Les points de mesure Mi sont obtenus par un suivi de la cosinusolde dans l'image à partir du sommet estimé S du flash. L'algorithme utilisé consiste à faire ni = i entre -L à L; en désignant par mo la ligne dans l'image correspondant au sommet estimé, pour i variant de 1 à L 5 respectivement de -1 à -L et mi est déterminé à partir de mi-1 respectivement mi+1 comme étant l'indice correspondant au maximum de luminance dans la fenêtre fréquentielle d'indice temporel i est de taille [-G,+G] centrée sur la ligne mi-1 respectivement mi+1 de la façon représentée à la figure 13. La valeur G peut être prise égale à 1 à 3 pixels. Les valeurs mi inférieures à mi-1 sont privilégiées lors de la recherche du maximum en cas d'ambiguté du fait que la courbe théorique décroît avec i. Par ailleurs afin de prendre en compte l'étendue temporelle du flash à son sommet, et éviter que le suivi ne diverge à l'intérieur du flash, les Pa 15 premiers points de mesure de chaque côté de S sont forcés à la valeur de l'indice mo. Le procédé qui vient d'être décrit permet d'estimer Q2 et par la même TR qui est la période de la cosinusode dans le plan tempsfréquence. Cette estimation est cependant relativement grossière mais suffit pour déterminer le nombre de pales de l'hélicoptère. Il est en effet 20 commode de mesurer sur l'image temps-fréquence le nombre d'échantillons séparant deux flashes positifs ou négatifs successifs ( F). Il est alors possible d'estimer le nombre de pales de l'hélicoptère par la relation Np = nint(tR / IF) (31) En admettant une estimation suffisamment précise de TF c'est l'estimation de TR qui conditionne celle de NP. Plus précisément Np est correctement estimé pourvu que l'erreur commise sur TR soit inférieure à TF/2. Néanmoins en cas d'imprécision trop importante sur la mesure de TR le résultat peut être affiné par un test de parité du nombre de pales. 30 La parité du nombre de pales peut en effet être simplement obtenue en évaluant la distance entre deux flashes positifs et négatifs successifs. Si cette distance tombe en dessous d'un seuil (de l'ordre de la largeur de la fenêtre de calcul d'énergie de la procédure de détection) il peut être sup- posé raisonnablement que l'hélicoptère possède un nombre pair de pales. Ce test supplémentaire permet de corriger des erreurs d'estimation de TR allant jusqu'à TF. Le nombre de pixels Np est alors estimé comme étant le nombre pair ou impair le plus proche de TR / F* Une fois le nombre de pales estimé, une estimée plus précise de la vitesse de rotation peut être calculée par la relation Q COR = 2n (NP. F) (32) De la sorte, c'est l'erreur commise sur TF qui devient prépondérante sur la précision d'estimation de Q corrigé. Suivant ce procédé l'es10 timation de la vitesse de rotation peut être correctement estimée à deux tours/minute près environ. Bien que le nombre de pales et la vitesse de rotation du rotor soient considérés comme des paramètres discriminants très importants pour la classification des hélicoptères et que l'estimation précédente permette d'effectuer une première classification, elle peut 1 5 malgré tout s'avérer encore insuffisante dans certains cas pour des hélicoptères possédant le même nombre pales et des vitesses de rotation voisines, de sorte qu'il est nécessaire de procéder à une classification plus fine par analyse fine à court terme du flash de pale par le dispositif d'analyse à court terme 1 7. L'analyse effectuée par ce module travaille 20 sur le module de la transformée en ondelettes à haute résolution temporelle après écrêtage et quantification linéaire. L'application de cette transformation sur des enregistrements réels permet notamment de mettre en évidence une stabilité relativement bonne des flashes associés à un même hélicoptère. Cette condition de stabilité est bien évidemment 25 nécessaire pour le développement d'une technique de classification. D'autre part, il existe pour chaque hélicoptère une forme spécifique du flash ayant notamment une dissymétrie dans le temps. Cette courbure du flash traduit à priori une courbure de la pale lors de sa rotation. L'analyse à court terme travaille de préférence sur les pales avançantes car elles 30 présentent une rétrodiffusion bien meilleure et surtout bien plus stable dans le temps. Plus précisément, après détection des flashes, le dispositif d'analyse à court terme travaille sur une fenêtre centrée temporellement sur le centre de phase détecté et s'étendant sur le demiplan fré- quentiel positif analysé. Le module du flash est alors considéré comme une fonction de deux variables f(x,y) o x désigne la variable temporelle et y la variable fréquentielle x e [1, Nx] et y e [1, Ny] avec. Le dispositif d'analyse à court terme détermine quatre paramètres discrimi5 nants du flash qui sont ceux définissant son centre de gravité dans la fenêtre et ceux définissant la courbure. La détermination du centre de gravité a lieu en calculant ses coordonnées xG et YG dans la fenêtre par les relations Nx Ny xG = ' 1 x.f(x,y) / S (33) x= 1y=1 Nx=l yl Nx Ny YG = X T y.f(x,y) / S (34) x= y= 1 Nx Ny oS= L Ef(x,y) (35) x= 1 y= 1 L'identification de la courbure du flash s'effectue en deux phases, une première phase consistant à rechercher par chacune des fréquences analysées des centres de gravité temporels Gi et la seconde phase consistant à effectuer une régression parabolique de ces barycentres. Les barycentres Gi sont définis par leurs coordonnées xi et yi selon les 20 relations yi = i pour i variant de 1 à Ny et Nx Nx Xi = j.f(xj, Yi) / f(xj,Yi) (36) j=1 i=-1 Une fois les points G(i) déterminés un polynôme degré 2 x = g(y) 25 est recherché parmi ceux qui correspondent au mieux au sens de la norme euclidienne à cet ensemble. g est défini par la relation x = A(y- 1) 2 + B(y- 1) + C (37) dans laquelle C représente l'indice temporel du barycentre à fréquence maximale. Sa détermination revient à effectuer une localisation très fine 30 du flash et rend donc les deux autres paramètres A et B invariants par translation temporelle et de ce fait insensibles à toute erreur de localisation. Les paramètres A et B caractérisent la courbure du flash. En notant par 0 = (A, B, C)T la solution du problème des moindres carrés précédent et par H une matrice à Ny lignes et 3 colonnes définie par: H(i, 1) = (i-1)2 H(i,2) = (i-1) H(i,3) = 1 pour i variant de 1 à Ny 0 vérifie la relation matricielle 0 = (HT.H)-1.HT.X o X est le vecteur de taille Ny défini par X = (x1, x2. xNy)T...DTD: Naturellement, les lignes fréquentielles Yi pour lesquelles ú f(xj, y i) = 0 ne sont pas prises compte dans l'évaluation des para15 rhètres de courbure A et B. La classification proprement dite peut ainsi s'effectuer de façon hiérarchique, la détermination de Np et Q permettent dans un premier temps une classification grossière qui est ensuite affinée par les calculs de XG, YG, A et B. A cette fin, il pourra être utilisé des classifieurs tout 20 à fait classiques, tels que des classifieurs linéaires ou de réseaux de neurones par exemple. REVENDICATIONS 1. Procédé de reconnaissance d'hélicoptère par analyse en ondelettes du signal radar rétrodiffusé par les pales du rotor principal caracté5 risé en ce qu'il consiste à effectuer une analyse (4) temps/fréquence par transformée en ondelettes du signal rétrodiffusé, à détecter (5) les flashes des pales du rotor principal dans un espace temps-fréquence, à analyser à long terme (6, 16) l'image obtenue dans l'espace tempsfréquence pour estimer le nombre de pales ainsi que la vitesse de rotation 10 du rotor et à analyser (7, 17) de façon fine les flashes des pales repérées dans l'espace temps-fréquence pour affiner la classification obtenue par une analyse à court terme. 2. Procédé selon la revendication 1 caractérisé en ce qu'il utilise une transformée en ondelette du type sinus cardinal fenêtré par une 1 5 Gaussienne. 3. Procédé selon l'une quelconque des revendications 1 et 2 caractérisé en ce qu'il consiste pour détecter (5) les flashes, au niveau du module de la transformée en ondelettes du signal rétrodiffusé servant de base à l'analyse à court terme (7) pour transmettre des valeurs 20 d'indices de localisation aux analyses à court terme (7) et long terme (6) à calculer à partir du module de la transformée en ondelette l'énergie du signal transformé à l'intérieur d'une fenêtre glissante, à comparer l'énergie obtenue à un seuil d'énergie prédéterminée et à localiser grossièrement dans l'espace temps les flashes lorsque l'énergie du signal dans la 25 fenêtre est supérieure au seuil d'énergie prédéterminée et à localiser finement la position temporelle du flash par corrélation fine autour de la position détectée du signal avec une image de référence du flash obtenue à partir d'un modèle théorique (8, 9) du signal rétrodiffusé par une pale. 4. Procédé selon la revendication 3 caractérisé en ce que l'estima30 tion du nombre de pales (16) a lieu sur la transformée à haute fréquence d'ondelette (14) quantifiée logarithmiquement (15) en faisant le quotient de la période TR de la trace estimée dans le plan temps fréquence de la fréquence pure du signal de pale par une estimation de l'intervalle de temps Tf séparant deux flashes successifs dans le même domaine de fréquence positif ou négatif. 5. Procédé selon la revendication 4 caractérisé en ce qu'il consiste à effectuer un test de parité du nombre de pales pour choisir le nombre entier pair ou impair le plus proche du quotient tR / TF. 6. Procédé selon les revendications 4 et 5 caractérisé en ce qu'il consiste pour estimer la période tR, à rechercher le sommet d'un flash puis à identifier dans le plan temps-fréquence de la cosinusode représentant la fréquence pure du signal de pale, la cosinusode à partir de son 10 sommet. 7. Procédé selon la revendication 6 caractérisé en ce qu'il consiste à rechercher le sommet des flashes, à partir d'une valeur d'indice temporel obtenu lors de la détection (5) ainsi que par filtrage et seuillage couplés du signal fréquentiel correspondant à cet indice. 8 Procédé selon l'une des revendications 6 et 7 caractérisé en ce que l'identification de la cosinusode à partir de son sommet a lieu par application de la méthode des moindres carrés sur un ensemble de points obtenus par suivi de cette cosinusolide. 9. Procédé selon l'une quelconque des revendications 1 à 7 carac20 térisé en ce qu'il consiste dans l'analyse à court terme (7) à calculer les coordonnées de son centre de gravité à l'intérieur d'une fenêtre d'analyse centrée temporellement sur le centre de phase détecté, et à identifier la courbure du flash par une recherche pour chacune des fréquences des centres de gravités temporels et régression parabolique des bary25 centres.