FR 2665316 A1 19920131 FR 9009614 A 19900727 i PROCEDE DE TRAITEMENT D'IMAGES POUR LA DETECTION DE L'HORIZON ET DISPOSITIF DE MISE EN OEUVRE L'invention se rapporte au domaine du traitement d'images en vue de détecter automatiquement l'horizon qui apparaît sur ces images. La détection automatique de l'horizon est déterminante pour certaines applications: systèmes de surveillance, aide au pilotage, aide à la détection de cibles, etc. Parmi les procédés de détection connus on peut distinguer les procédés d'élaboration d'horizon artificiel et les procédés de détection d'horizon dans une image à partir d'un signal vidéo délivré par une caméra Les premiers procédés sont basés sur un repérage de l'aéronef en position et en direction grâce à une centrale inertielle et les seconds sont basés sur des traitements d'images exploitant différents attributs ou paramètres caractéristiques de l'horizon (par exemple étude des gradients de luminance). Les procédés connus manquent de précision tant au plan de la localisation spatiale qu'au plan de la détermination temporelle Pour la première catégorie de procédés, la détermination de l'horizon est absolue, sans tenir compte du relief réel, car elle ne sert qu'à définir l'assiette de l'aéronef Le manque de précision se traduit, pour la seconde catégorie de procédés, par un bruit important pouvant aller jusqu'à fournir une information aberrante. L'invention vise à remédier à ces inconvénients en proposant un procédé utilisant des traitements d'images permettant la détection automatique de l'horizon avec une grande précision, en temps réel, et pour toutes les positions possibles de l'horizon L'invention concerne également un dispositif de mise en oeuvre de ce procédé. Pour atteindre ces objectifs, le procédé de traitement d'images selon l'invention, pour la détection de l'horizon à partir d'un signal vidéo délivré par une caméra thermique, est caractérisé en ce qu'il comporte les étapes suivantes: traitement numérique du signal vidéo pour obtenir les valeurs de luminance échantillonnées de l'image analysée rangée par rangée; détermination, pour chaque position i d'une rangée k donnée, d'une fonction de modélisation en échelon dont le saut est situé en i et dont l'amplitude est calculée de façon à donner une erreur minimale de l'ensemble des valeurs prises par la fonction par rapport à l'ensemble des valeurs du signal de luminance de la rangée considérée; sélection de la valeur ihk de i, correspondant à la position de l'horizon sur la rangée considérée, pour laquelle la fonction d'échelon déterminée à l'étape précédente a une amplitude calculée de façon à donner la plus faible de toutes les erreurs minimales; reproduction des étapes précédentes pour toutes les rangées de l'image analysée et visualisation d'une courbe de lissage, représentant l'horizon, à partir de l'ensemble des valeurs ihk sélectionnées. La présente invention sera mieux comprise et d'autres caractéristiques et avantages apparaîtront à la lecture de la description qui va suivre, en référence aux figures annexées qui représentent, respectivement: la figure 1, un schéma illustrant la structure d'une caméra thermique panoramique; les figures 2 a à 2 c, des modélisations du signal de luminance d'une colonne donnée de l'image analysée sous forme de fonctions Heaviside; les figures 3 a et 3 b, deux schémas synoptiques illustrant la formation de l'erreur commise en modélisant le signal de luminance d'une colonne donnée, la figure 4, le dispositif de détection d'horizon selon l'invention. Le champ d'observation est analysé rangée par rangée par une caméra thermique, c'est-à-dire colonne par colonne ou ligne par ligne Dans l'exemple de réalisation illustrant une première mise en oeuvre de l'invention, l'image est délivrée par une caméra thermique panoramique Une telle caméra, en référence à la figure 1, est munie habituellement d'un capteur linéaire 1 constitué d'une barrette de N cellules photosensibles, disposée sur un support cylindrique 2, et sur laquelle est focalisé un secteur d'observation grâce à un système optique de focalisation 3 Les charges intégrées dans chaque cellule sont lues dans un registre 4 et délivrées sous forme d'un signal de luminance qui, après traitement électronique dans un circuit 5, fournit le signal vidéo Les variations du signal vidéo traduisent donc les variations en température du champ d'observation analysé, colonne par colonne, sur p colonnes de N pixels Dans ces conditions, la détection de l'horizon est opérée à partir des valeurs de luminance d'une image du champ d'observation analysée colonne par colonne. L'horizon terrestre peut se définir comme la frontière entre le ciel et la terre En imagerie thermique, le ciel et la terre ont des températures moyennes différentes, le ciel étant toujours beaucoup plus froid que la terre L'horizon est donc un lieu continu o une séparation brutale se produit en température Cette variation brutale est traduite, dans le cadre de la présente invention, par un saut. L'invention repose, en effet, sur la modélisation de la courbe représentative des variations de luminance de chacune des rangées de valeurs par une fonction d'échelon ou fonction d'Heaviside La démarche suivie par l'invention est donc une démarche globale qui tend, par un traitement d'optimisation tenant compte de l'ensemble des valeurs de luminance d'une rangée donnée, à sélectionner une fonction d'échelon dont le saut est situé à une position indicative de la position de l'horizon sur la colonne considérée Le procédé selon l'invention peut être appliqué indifféremment à une caméra d'analyse d'image en colonnes ou en lignes Dans le développement qui suit, ce procédé est d'abord appliqué à une analyse d'images en colonnes effectuée par la caméra panoramique, puis appliqué à une analyse d'image en lignes avec des variantes de mise en oeuvre. Plus précisément, en référence aux figures 2 a à 2 c, qui représentent par trois fois la courbe Sk des valeurs de luminance 1 jk le long d'une colonne de rang k en fonction des positions I à N des points de la colonne, l'invention propose la modélisation de la courbe Sk selon les étapes successives suivantes: échantillonnage du signal de variation de luminance en n valeurs réparties le long de la colonne de rang k; modélisation de la courbe Sk à chaque position i dans la colonne k, i variant de I à n, par une fonction d'échelon i fi définie à partir: d'une fonction d'échelon normée f dont le saut se situe à la position i, et, de l'amplitude de saut i O i qui permet de s'ajuster au mieux aux variations du signal de luminance de la colonne: l'amplitude Q kième colonne. Pour accéder à la fonction d'échelon vérifie que l'erreur ei, élevée au carré, est minimale. La fonction d'échelon normée f utilisée doit être "centrée" sur une colonne, c'est-à-dire que la somme des valeurs prises par cette fonction sur cette colonne doit être nulle Il est en effet impératif de conserver le niveau moyen de luminance de l'image, colonne par colonne, du fait de la puissance limitée du signal émis par la caméra Un exemple de fonction d'échelon "centrée" est donné par les expressions suivantes: fi (j) = -(n-i) pour j appartenant à l 1,il fi (j) = i pour j appartenant à li+ 1,nl qui vérifient bien la condition imposée par le centrage: -(n-i)l + i(n-i) = O Symboliquement, la figure 3 a matérialise la mise en forme de l'erreur ei par différence entre les valeurs de luminance 1 jk de la courbe Sk et les valeurs respectives de la fonction d'échelon normée f à laquelle est appliquée le coefficient variable o En notation matricielle, soit L la matrice des valeurs de luminance correspondant à la courbe Sk pour l'ensemble des niveaux de la klème colonne, T le symbole de transposition de matrice, F la matrice des valeurs de la fonction d'échelon a normée f (i) L'erreur e a alors pour expression 2 T e = (L F(E) (L F) (I) La valeur de i qui rend e minimale vérifie que ? i la dérivée de ei est nulle La valeur de l'amplitude i de la fonction d'échelon qui lisse au mieux la courbe Sk est déterminée en annulant cette dérivée La valeur de ( est alors donnée par l'expression iL F/ FTF (II) En remplaçant i( par sa valeur o( déterminée par la formule (II) dans l'expression de l'erreur e (formule (I)), l'erreur minimale e m pour le niveau i correspondant est alors donné par: a T T 2 T e L L -(L F) /F E (III) Pour chaque niveau i, la formule (II) fournit donc la valeur de l'amplitude de la fonction d'échelon lissant au mieux, au sens des moindres carrés, la courbe Skl et la formule (III) fournit la valeur de l'erreur correspondant au lissage ainsi déterminé D'après la structure de la formule (III), cette erreur mesure la différence d'énergies entre l'énergie lumineuse LTL de la colonne k et une énergie lumineuse "modélisée", notée Q, obtenue par intervention de la fonction d'échelon normée. La plus petite de toutes les erreurs ejm, appelée em, correspond à la fonction d'échelon fhk dont la courbe représentative (figure 2 c) s'ajuste au plus près à la courbe Sk pour l'ensemble des i de la kme colonne La position ihk du saut de cette fonction fhk situe la position de l'horizon sur la colonne considérée La valeur minimale e de lim l'erreur eim correspond, d'après la formule (III), à la valeur maximale de la quantité Q: Ainsi, en notant respectivement ai la somme des valeurs de luminance 1 jk des points j situés entre les positions 1 et i, et bi celle des valeurs 1 jk des points j situées entre les positions 1 + 1 et N de la colonne considérée, la valeur ihk de i, qui indique la position de l'horizon sur la colonne considérée, rend maximale l'expression analytique qik de l'énergie lumineuse Q, à savoir: (- (n-i)a + ib)2 / i(n-i) (IV) Pour accéder plus explicitement au type de répartition d'énergie lumineuse qui opère dans l'expression (IV), une variante du procédé selon l'invention consiste à modéliser la courbe représentative des valeurs de luminance de chaque colonne de l'image analysée par une combinaison de deux fonctions d'échelon unité f 1 et f 2, ce qui fournit une plus grande souplesse de modélisation. Les deux fonctions d'échelon sont définies de manière à être décorellées l'une de l'autre Un exemple de définition est fourni par les expressions suivantes: f 1 (= 1 et f 2 (j) = O pour j appartenant à l 1,1 l f 1 (j) = O et f 2 (j) = 1 pour j appartenant à l 1 + 1, nl La suite des étapes du procédé selon l'invention se déroule alors suivant le même type de schéma que dans le cas de l'approximation par une seule fonction d'échelon: détermination de deux coefficients multiplicateurs (ou amplitudes) ( et 13 i qui minimisent l'erreur e.i= 1 k (o( f 1 + 313 f) ei ljk (i f 3 i f 2 et sélection de la valeur ihk de i correspondant à la plus petite des erreurs ei. La figure 3 b symbolise par un schéma synoptique la formation de l'erreur ei obtenue par différence entre les valeurs de luminance ijk de la courbe Sk et les valeurs respectives des deux fonctions d'échelon f et f 2 auxquelles sont appliquées les coefficients i et 13 i. Dans ces conditions, par calcul direct ou par analogie avec l'étude à partir d'une seule fonction d'échelon, l'erreur minimale eim, pour un i donné, s'exprime à l'aide des notations matricielles par l'expression: eim 2 =LTL-l (LT I 1) 2/,T F 1 + (LT F 2)2/F 2 T F 2 l La recherche de la valeur ihk qui correspond à la position de l'horizon sur la colonne considérée, est identifiée, dans le cadre de l'invention, à la recherche de la valeur de i pour laquelle l'erreur ei est la plus petite de toutes sur cette même colonne Or, cette plus petite erreur, appelée em, correspond d'après l'expression précédente au maximum de: T T T T (L F 1) /Fl F 1 + (L F 2) /F 2 F 2 Ainsi, en développant ces expressions matricielles à l'aide des notations ai et b précédemment définies il vient, en notant 1 k l'énergie lumineuse de la colonne k eim lk qik avec 1 k (ljk) et qk = a 2/i+b 2/(n-i) (V) k= ( 4) La valeur ih de i est donc celle qui correspond à la valeur maximale q Mk des quantités q 1 k Cette expression est à rapprocher de celle (expression (IV)) obtenue dans le cadre d'une modélisation à l'aide d'une seule fonction d'échelon : la différence entre les expressions (IV) et (V) est en effet une constante égale à l'énergie lumineuse moyenne de la colonne considérée: sk /n, avec sk = a + b = somme des valeurs de luminance de la colonne k. Ainsi, chacune des valeurs ihk déterminée est la position pour laquelle la somme qik des énergies lumineuses moyennes en-deçà et au-delà de cette position est maximale dans chacune des p colonnes; comme chaque valeur ihk est la position de l'horizon sur la colonne correspondante, l'ensemble des p valeurs ihk donne la position de l'horizon modélisé point par point sur l'ensemble de l'image. Afin de visualiser l'horizon ainsi modélisé par une courbe continue, il est connu de lisser un ensemble de points, par une courbe ou une droite (en utilisant, par exemple, la transformée de Hough ou une régression linéaire) Cependant, le procédé précédemment décrit fournit systématiquement une position d' horizon pour chaque colonne, même lorsqu'aucun horizon n'est présent sur cette colonne (par exemple lorsque l'horizon fait un angle de plus de 450 avec l'horizontale) Un lissage à partir de l'ensemble des positions déterminées par le procédé selon l'invention serait perturbé par les positions erronées qui en résultent, en particulier pour une approximation par une droite. Pour neutraliser l'effet perturbateur de ces positions aberrantes, il est utile d'introduire un coefficient de pondération définissant le degré de certitude des différentes positions déterminées par le procédé de l'invention Cette pondération peut être avantageusement représentée par un rapport signal sur bruit (S/B), fonction de la valeur de l'erreur minimale em commise en modélisant la courbe Sk par la courbe représentative de la fonction d'échelon fhk sélectionnée pour la colonne considérée Cette approche est justifiée par le fait qu'une fonction d'échelon présentant une erreur relativement importante a peu de chance de correspondre à une transition réelle entre le ciel et la terre Le signal S/B peut, par exemple, être défini pour la kème colonne comme l'écart relatif entre l'énergie lumineuse 1 k de cette colonne et le carré de l'erreur minimale e sélectionnée pour cette m colonne, soit 2 2 S/B (k) = ( 'k em) / em ou encore, d'après les expressions (V) S/B (k) = q Mk / (lk -q Mk) En définitive, il convient de rechercher la droite (D) dont la position et la pente minimisent l'expression suivante Z k S/B (k) x distance lniveau ihk, Dl (VI), ihk étant la position sélectionnée comme représentant la position de l'horizon selon le procédé de l'invention pour la kième colonne, et S/B (k) le degré de certitude de ce niveau. Le procédé selon l'invention peut également s'appliquer à une caméra effectuant une analyse en ligne de l'image du champ d'observation Pour appliquer le procédé de l'invention à une telle caméra, plusieurs variantes sont possibles: l'application directe: la recherche de la position de l'horizon s'effectue directement en ligne, de la même manière qu'elle s'effectue en colonne; le résultat obtenu est de qualité inférieure à celui obtenu à partir de la détermination en colonne pour les horizons proches de l'horizontale, du fait de l'imprécision qui apparaît alors dans cette zone; mise en oeuvre d'une mémoire d'image afin de disposer de toutes les informations colonne par colonne; cette solution est coûteuse et peut introduire un temps de retard peu souhaitable dans la détermination de l'horizon; l'utilisation de l'énergie lumineuse de la trame précédente: l'expression (V), dont la valeur maximale détermine la position de l'horizon sur une colonne k donnée, peut encore s'écrire en faisant intervenir la somme sk des valeurs de luminence sur cette colonne: qi = ai/ i + (sk ai (n-i) (V) bis L'image étant relativement stable d'une trame à l'autre, il suffit de prendre pour S k la valeur obtenue sur la trame précédente Cette approximation est justifiée par le fait que la variation de l'énergie lumineuse d'une trame à la suivante peut être considérée comme négligeable. Un mode de réalisation du procédé selon l'invention, adapté à une analyse d'images en mode colonne ou en mode ligne, est décrit ci-après, en référence à la figure 4. En mode colonne, les variations de luminance des niveaux d'une même colonne k, après échantillonnage numérique dans l'échantillonneur 6, sont accumulées successivement dans un accumulateur 7, puis ces sommations partielles (alk, a 2 k, aik, ank Sk) sont mémorisées dans une mémoire double flipflop 8 Un processeur 9 calcule alors, pour chaque niveau i de la colonne k considérée, la somme des énergies lumineuses moyennes au-delà et en deçà du niveau 1, soit la quantité q 1 déjà rencontrée: q.= a 2/i + (a -a) l = alk ? (ank ik / n-i (V) ter Les quantités ainsi calculées sont appliquées à un comparateur 10 qui permet, par comparaison successive de celles-ci deux à deux, de déterminer le rang ih correspondant au maximum de ces quantités et donc, d'après le procédé de l'invention, à la position de l'horizon sur la colonne considérée Les valeurs ihk, pour l'ensemble des colonnes de l'image analysée, sont mémorisées dans une mémoire 11 Il est avantageux, pour conserver la valeur de pondération S/B des niveaux ainsi sélectionnés, de mémoriser également les quantités maximales q, correspondantes dans une autre mémoire 12 Le processeur 9 calcule alors, à partir des valeurs fournies par les mémoires 11 et 12, les paramètres de la droite de régression linéaire qui lisse les niveaux déterminés par le procédé de l'invention à partir de la formule (VI) Le processeur 9 délivre alors un signal de forte luminance correspondant au niveau du point de la droite de régression dans la colonne considérée Ce signal est appliqué, à travers un soustracteur 13, au niveau correspondant du signal de luminance numérisé de telle sorte que la valeur de luminance qui en résulte au niveau considéré soit inversée en valeur par rapport à la valeur initiale Ainsi la droite de synthèse, matérialisant l'horizon déterminée par le procédé de l'invention, apparaît en contraste sur le paysage tel que visualisé sur un écran de contrôle à partir du signal vidéo fourni par le caméra. En mode ligne, la double mémoire flip flop 8 est utilisée de la manière suivante: les valeurs de luminance cumulées, colonne par colonne jusqu'à une ilème ligne (a 1 I a 12, aik,,a i) sont stockées dans une des mémoires flip flop, et dans l'autre sont stockées les valeurs de luminance (a nî, an 2,, ank, '*, an) des colonnes de la trame précédente Le processeur 9 calcule alors, de façon similaire au mode colonne, les quantités ayant pour expression la formule (V) adaptée selon (V) ter, et envoie les signaux de luminance maximale correspondant à la droite de régression, calculée à partir des valeurs mémorisées dans les mémoires 11 et 12 au soustracteur 13. L'invention n'est pas limitée au mode de réalisation décrit et représenté En particulier, le fonctionnement en temps réel du dispositif de mise en oeuvre selon l'invention est rendu possible par l'utilisation d'un processeur de signal de capacité adapté, par exemple un VLSI (circuit intégré à très grande échelle), ou par un traitement statistique des valeurs de luminance fournies par la caméra; par exemple l'évaluation des fonctions d'échelon peut se faire à partir des valeurs de luminance prélevées tous les x pixels, par réduction de la zone de recherche de l'horizon, ou par recherche de l'horizon toutes les y colonnes D'autres part, l'implantation de la carte de détection de l'horizon peut être aménagée en différents sites: directement au niveau de la délivrance du signal de luminance d'une caméra panoramique, ou au niveau du boîtier électronique contenant les différentes cartes de traitement Le signal vidéo est alors préalablement échantillonné de manière à pouvoir être exploité à partir de valeurs numériques. REVENDICATIONS 1 Procédé de traitement d'images pour la détection de l'horizon, à partir du signal vidéo délivré par une caméra formant l'image du champ d'observation analysée rangée par rangée, caractérisé en ce qu'il comporte les étapes suivantes: traitement numérique du signal vidéo pour obtenir les valeurs de luminance échantillonnées de l'image analysée rangée par rangée; détermination, pour chaque position i d'une rangée k donnée, d'une fonction de modélisation en échelon dont le saut est situé en i et dont l'amplitude est calculée de façon à donner une erreur minimale par rapport à l'ensemble des valeurs du signal de luminance de la rangée considérée; sélection de la valeur ihk de i, correspondant à la position de l'horizon sur la rangée considérée, pour laquelle la position d'échelon déterminée à l'étape précédente a une amplitude calculée de façon à donner la plus faible de toutes les erreurs minimales; reproduction des étapes précédentes pour toutes les rangées de l'image analysée et visualisation d'une courbe de lissage, représentant l'horizon, à partir de l'ensemble des valeurs lhk sélectionnées. 2 Procédé de traitement d'images pour la détection de l'horizon selon la revendication 1 caractérisé en ce que la fonction d'échelon de modélisation est la somme de deux fonctions d'échelon décorrélées dont les amplitudes ont des valeurs telles que, pour chaque position i d'une rangée donnée, l'erreur entre les valeurs de la fonction de modélisation, et l'ensemble des valeurs échantillonnées prises respectivement au niveau correspondant par le signal de luminance, soit minimisée. 3 Procédé de traitement d'images selon les revendications 1 ou 2 caractérisé en ce que le niveau ihk de la position de l'horizon sur la colonne de rang k est celui pour lequel la somme des énergies lumineuses moyennes en-deçà et au-delà de ce niveau est maximale dans la rangée considérée. 4 Procédé de traitement d'images selon les revendications 1 ou 2 appliqué à une analyse d'image par colonnes, caractérisé en ce que l'image est découpée en p colonnes, en ce que i varie en lignes de 1 à N et en ce que les valeurs ihk sélectionnées indiquent directement le traitement de l'horizon sur l'image. 5 Procédé de traitement d'images selon les revendications 1 ou 2 appliqué à une analyse d'images par lignes, caractérisé en ce que l'image est découpée en N lignes, en ce que i varie en colonnes de 1 à p, et en ce que les valeurs ihk sélectionnées indiquent soit directement, soit indirectement après traitement des valeurs de luminance par colonne, la hauteur de l'horizon sur l'image. 6 Dispositif de traitement d'images pour la mise en oeuvre du procédé selon la revendication 4 caractérisé en ce qu'il comporte un échantillonneur de signal (b) dont les valeurs numérisées sont successivement adressées à un accumulateur ( 7) pour sommet ces valeurs numérisées au fur et à mesure de leur apparition, une mémoire flip flop ( 8) couplée à l'accumulateur ( 7) et qui mémorise les sommations cumulées, un processeur de signal ( 9) qui détermine, à partir des valeurs transmises par la mémoire ( 8) et pour chaque position i de chaque colonne k de valeurs de luminance, la quantité q 1 égale à la somme des énergies lumineuses moyennes au-delà et en deçà de la position i, un comparateur ( 10) couplé au processeur ( 9) pour sélectionner la quantité q Mk ayant la valeur maximale pour la colonne de rang k, des mémoires ( 11) et ( 12) couplées au comparateur ( 10) pour mémoriser respectivement les quantités q Mk et les valeurs ihk de i correspondantes, et un soustracteur ( 13) monté sur le circuit délivrant le signal vidéo échantillonné et sur le processeur ( 9) pour soustraire su signal provenant de l'échantillonneur ( 6) un signal de luminance maximale provenant du processeur ( 9) pour le niveau de chaque colonne o est positionné le point de la droite représentant la position de l'horizon sur la colonne considérée, cette droite étant calculée dans le processeur ( 9) par régression linéaire par application des valeurs de niveaux mémorisés dans la mémoire ( 12) et des valeurs de pondération correspondantes mémorisées dans la mémoire ( 11). 7 Dispositif selon les revendications 4 et 5 caractérisé en ce que l'une des mémoires flip flop ( 8) mémorise successivement les sommations partielles des valeurs de luminance en ce que l'une des deux mémoires flip flop ( 8) cumule et stocke les valeurs de luminance colonne par colonne et que l'autre mémoire stocke les valeurs de luminance des colonnes de la trame d'image précédente, en ce que les quantités q 1 sont calculées dans le processeur ( 9), à partir des sommations partielles de luminance provenant de la première des mémoires flip flop ( 8) et des sommations des valeurs de luminance des colonnes provenant de la seconde des mémoires flip flop ( 8).