Brevet FR2960318A1 — publié le 20111125

Procédé d'estimation de l'état d'un système par résumé des mesures effectuées sur le système La présente invention concerne, de manière générale, un procédé d'estimation statistique de l'état d'un système. Elle s'applique notamment dans le domaine de la détection et de la localisation de cibles à partir de mesures fournies par un ou plusieurs capteurs. L'application de traitements d'estimation statistique à un nombre très élevé de données, qu'il s'agisse de mesures fournies par un capteur ou de tout autre type de données, permet d'obtenir une grande précision dans l'estimation des paramètres qui caractérisent la cible. On parle alors de traitement des données « par blocs » ou encore de traitement « batch ». Mais travailler sur toutes les mesures simultanément pose un double problème. Il s'agit d'une part d'un problème de mémorisation : l'espace mémoire alloué pour stocker les mesures étant forcément limité, comment prendre en compte de nouvelles mesures lorsque cet espace est déjà rempli ? Le cas échéant, quelles mesures supprimer ? Il s'agit d'autre part d'un problème de coût des calculs liés à la procédure d'estimation : plus l'algorithme d'estimation traite de mesures simultanément, plus sa durée d'exécution est longue, le coût de calcul augmentant en général de manière non linéaire. Une première solution connue consiste à opter pour un algorithme d'estimation récursif. L'inconvénient majeur d'une telle solution est qu'un algorithme récursif est réputé produire une estimation moins précise dans le cas général, ou être beaucoup plus lourd à mettre en oeuvre pour atteindre une performance similaire au traitement batch. Une autre solution connue consiste à suivre des règles de choix ad hoc pour supprimer des mesures, parmi lesquelles par exemple la suppression de la mesure d'écart-type le plus grand, de la mesure la plus rapprochée d'une autre temporellement ou la plus ancienne et d'appliquer l'algorithme d'estimation sur les mesures conservées. Un inconvénient majeur de ces règles de choix est qu'elles font travailler l'algorithme d'estimation sur un sous-ensemble des mesures. Elles entraînent une perte d'information et donc de précision, qui peut éventuellement s'amplifier au cours du temps. Un autre inconvénient majeur de ces règles de choix est qu'elles sont arbitraires et ne s'appuient pas sur des fondements statistiques. Elles font donc courir le risque d'éliminer des mesures de valeurs ajoutées importantes. L'invention a notamment pour but de limiter les pertes d'information et de précision dans les méthodes de traitement par blocs, tout en conservant les mêmes contraintes sur la charge calcul et la taille de mémorisation autorisées. Pour cela, l'invention propose de résumer N mesures réelles, où N entier strictement positif, par N' mesures virtuelles, où N' entier strictement positif et N' 3 Par exemple, l'état X du système peut être un vecteur d'état (X,,...,Xd) de dimension d, où d est un entier strictement supérieur à 1, l'étape de calcul pouvant calculer un vecteur d'état estimé (X,,...,Xd) de dimension d et sa matrice de covariance R de dimension d°x°d. Dans un mode de réalisation préférentiel, l'étape de calcul peut calculer le vecteur d'état estimé (kl,...,kd) par application de la méthode du Maximum de Vraisemblance, la matrice R pouvant être obtenue par inversion d'une matrice de Fisher. Par exemple, l'étape de compression peut inclure une étape de calcul, à partir des N mesures réelles, de la matrice de Fisher J1 associée au vecteur d'état estimé (X,,...,Xd), ainsi qu'une étape de détermination des N' mesures virtuelles satisfaisant à l'équation J1=J2, où J2 est la matrice de Fisher associée au vecteur d'état estimé (X,,...,Xd) exprimée en fonction des N' mesures virtuelles. Par exemple, le capteur peut être un capteur de détection embarqué sur un porteur et le système peut être une cible détectable par le capteur. Par exemple, si d=2, le vecteur d'état peut correspondre à la position (Xcible,Ycible) de la cible dans un plan et le vecteur d'état estimé correspondant à une estimation (Xcible, Ycible) de la position (Xcible, Ycible) - Les matrices J1 et J2 peuvent être données par : J, = J2 = N; cos z (Azi ) Lr 2 2 i=1 Dh .6i Nmx sin(Azi).cos(Azi ) 2 2 i=1 Dh .6i N; cos2(Azi) Z 2 2 i=1 Dh, .6i Ïsin(Azi).cos(AZi) D2 !Y z =1 h, i Nraa. sin(Az; ). cos(Azi ) L 2 2 i=1 Dh..Qi Ne sin 2 (Azi ) 2 2 i=1 Dh 6i N' sin(Azi ). cos(Azi ) Lr D2 i_1 h 6z i sin 2 (Azi ) 1 n2 2 i=1 Dh .6i et les azimuts Azi et les distances horizontales Dhi dans la matrice J1 pouvant être calculés en supposant que X eihie = xporteur (t) + Dh (t).sin(Az(t)) et Ycihle = yporteur (t) + Dh (t).cos(Az(t)) , où (xporteu, (t), yporteur (t)) peut être la position connue du porteur dans le plan à un instant t de mesure, où Dh(t) peut être la distance horizontale du porteur à la cible à l'instant t et où Az(t) peut être l'azimut entre le porteur et la cible à l'instant t, la variance o de la lème mesure pouvant être connue. La présente invention a encore pour principal avantage d'être aisément implémentée dans les systèmes actuels effectuant des traitements d'estimation statistique sur des blocs de données, ceci sans en modifier les interfaces, les mesures virtuelles étant construites sur le même format que les données initiales. De plus, un autre avantage de l'invention est qu'elle permet une mise en oeuvre optimale des traitements d'estimation statistique sur des blocs de données, en les exécutant sur un nombre de mesures ajustable en fonction des ressources disponibles (mémoire, unité de calcul). D'autres caractéristiques et avantages de l'invention apparaîtront à l'aide de la description qui suit, faite en regard de dessins annexés qui représentent : - la figure 1, par un diagramme fonctionnel, un exemple de système d'estimation statistique comportant un module de compression selon l'invention; - la figure 2, par un diagramme fonctionnel, un exemple de système d'estimation statistique distribué comportant plusieurs modules de 25 compression selon l'invention. La figure 1 illustre, par un diagramme fonctionnel, un exemple de système d'estimation statistique comportant un module 1 de compression de 30 données selon l'invention. Avantageusement, le module de compression 1 met en oeuvre l'invention de façon logicielle. Le module de compression 1 selon l'invention est inséré dans un système passif de détection et de localisation de cibles fixes comportant un module 4 d'estimation statistique de type traitement par blocs. Le module 4 35 effectue des traitements classiques d'estimation statistique sur des blocs de mesures d'azimut fournies par un capteur de détection passif 2, les mesures étant stockées dans un espace mémoire 3 pouvant contenir N_max mesures au maximum. Au sens de la présente invention, la notion de mesure recouvre non seulement la valeur mesurée, mais elle peut également recouvrir des éléments nécessaires à son interprétation ou encore des éléments nécessaires à son exploitation. Par exemple, dans l'exemple de réalisation de la figure 1, une mesure comporte une direction d'arrivée (azimut) mesurée par le capteur 2, la position du capteur 2 au moment de la mesure ainsi que la variance de la mesure. Le module d'estimation 4 peut par exemple implémenter la méthode du Maximum de Vraisemblance. Ainsi, pour une cible fixe de coordonnées (Xeible,Ycible) inconnues a priori, le module d'estimation 4 fournit des états estimés notés Xcible et Ycible , ainsi qu'une estimation de la qualité de l'estimation en (Xcible, Ycible) sous la forme d'une matrice 2x2 de covariance notée R, obtenue par inversion de la matrice de Fisher calculée au point supposé (5(Cib1e,Ycible) de localisation de la cible dans le plan XY. La matrice R fournit une estimée de la Borne de Cramer-Rao (BCR). Dans l'exemple de réalisation de la figure 1, (Xcible, Ycible) est un vecteur d'état estimé du système de dimension d=2 et R correspond à une ellipse autour 20 du point de coordonnées (kcjble,Ycible). Le module 4 pourrait implémenter toute autre méthode de localisation passive que le Maximum de Vraisemblance, dès lors que cette autre méthode fournit un point de localisation de la cible et une matrice de covariance. Ainsi, selon l'invention, tant que le nombre N de mesures reçues 25 du capteur 2 et stockées dans l'espace mémoire 3 est inférieur à N_max, le module 4 calcule notamment (Xeible,Ycible) et R, alors que le module 1 de compression demeure inactif. Lorsque N=N_max, alors non seulement le module 4 calcule (Xcible,Ycible) et R, mais le module 1 de compression détermine également N' mesures virtuelles, où N' Le principe général de l'invention est de déterminer les N' mesures virtuelles à partir de la matrice J1 = R-1, appelée matrice de Fisher. Cette matrice est une matrice carrée, symétrique et définie positive, dont la dimension est celle du vecteur d'état estimé. Ainsi, si le vecteur d'état est de dimension d, et si les N valeurs mesurées sont de dimension p chacune, il faut trouver N' tel que N_max > N' E(d)+1 mesures, où E désigne la partie P entière, qui donne, avec le même algorithme, le même résultat que les N mesures réelles. En supposant que le vecteur d'état est le même que celui calculé avec les N mesures réelles, la matrice de Fisher J2 associée aux N' mesures virtuelles inconnues peut facilement s'exprimer en fonction dudit vecteur d'état et des N' mesures virtuelles. Le système d'équations obtenu peut être résolu en faisant l'hypothèse que J2 = J1. On obtient alors les N' mesures virtuelles recherchées. Ainsi, dans l'exemple de réalisation de la figure 1, si la cible est observable par le capteur 2, la matrice de Fisher au point supposé (Xcible,tible) de localisation de la cible dans le plan XY est une matrice 2x2, symétrique, définie positive. En prenant N_max=3, il s'agit donc de déterminer deux mesures virtuelles, constituées chacune d'un azimut, d'une position du capteur 2 et d'une variance en azimut, telles que ces deux mesures virtuelles conduisent, en utilisant le même algorithme de localisation basé sur le Maximum de Vraisemblance, à la même localisation (Xcible,Ycible) et à la même matrice de Fisher. Les deux mesures virtuelles sont alors à 30 repositionner par rapport à la localisation (Xcible, Ycible) supposée de la cible. Dans le présent exemple de réalisation, l'invention permet une exécution rapide avec 3 mesures d'azimuts seulement : les deux mesures virtuelles qui résument le passé et la nouvelle mesure. En l'occurrence, ceci correspond à une mise en oeuvre récursive d'un traitement par blocs. Mais il faut bien comprendre que, en faisant l'hypothèse que le résultat trouvé avec ces trois azimuts est identique au résultat qui aurait été trouvé avec tous les azimuts, une approximation est faite. En effet, cette hypothèse n'est pas rigoureusement exacte du fait de la procédure non linéaire appliquée dans la méthode du Maximum de Vraisemblance. Analytiquement, la matrice de Fisher pour la localisation passive d'une cible fixe à partir d'azimuts, se calcule par la formule suivante pour N mesures : N cos 2 (Azi) N sin(Azi ). cos(Azi ) 2 2 2 i=1 Dh..6i2 i=1 Dh .6i sin(Azi ). cos(Azi) sin 2 (Azi ) 2 2 2 2 i=1 Dh .6i i=1 Dh .6i L'algorithme de localisation basé sur le maximum de 15 vraisemblance calcule cette matrice J1, en utilisant comme azimuts Azi et distances horizontales Dhi ceux et celles reconstitués à partir de la localisation (kçjble,tible) qu'il a déjà calculée, ainsi que les positions du porteur et les variances de mesures aux instants i de mesures, qui sont connues. 20 I I s'agit de trouver les azimuts Az1 et Az2, les distances horizontales Dh et Dhz (qui permettent de remonter aux positions du capteur associées aux azimuts Az1 et Az2), et les écarts-types 61 et ci-2, qui conduisent aux mêmes valeurs numériques que celles obtenues avec 25 l'ensemble des mesures : J1 J2 = cos2 (Azi ) L~ 2 2 i=1 Dh..6i sin(Azi ). cos(Azi ) 2 2 i=1 Dh, i 2 sin(Azi).cos(Azi ) 2 2 i=1 Dh,, .6i 2 1 sin 2 (Az ) D2 2 i=1 Comme la matrice J, est symétrique, elle ne contient que trois valeurs distinctes. Par conséquent, on peut limiter le nombre d'inconnues à 3 au maximum. Pour se faire, on fixe 61 = 62 = a et D,, = DhZ = D . Il reste donc à trouver Az1, Az2, et D. II faut noter qu'on pourrait aussi fixer directement le 5 produit D.6. De plus, la matrice J, est symétrique et définie positive. Elle est donc diagonalisable et ses valeurs propres sont réelles et positives. II existe donc une matrice P = /P' 12 et une matrice diagonale J; telles que \ 21 22) 10 J, = p * J id * PT , avec J; =(K- K \ . Les éléments diagonaux de Jd sont y/ les valeurs propres de J1 , les colonnes de P sont les vecteurs propres de J1 . Au passage, il faut noter que P est une matrice orthogonale, d'où l'équivalence J, = p* Jd * pT Jd = pT * J1 * P , car P * PT = pT *p = I donc PT = P-1. II faut également noter que P est une matrice de rotation 15 d'angle B , tel que 6= Atan P21 =-A tan P12 . P1 P1 On cherche alors à résoudre les équations (1), (2) et (3) qui suivent dans la base définie par les vecteurs propres. On note Az; , Az2 les angles et D la distance (indépendante du repère) correspondants : 2 cos2(Ai. ) cos2(Azi) cos2(Az2) =K 1 + ( ) 1 z - z z i=1 Dh .6i D2472 D .6 . sin2(Azi) = sin2(Azi) sin2(Az2) (2) =Ky E z z Z z+ z z i=1 Dh .6i D .6 D .6 cos(Azi ). sin(Azi) - cos(Az1). sin(Az1) + cos(Az2 ). sin(Az2) _ 0 (3) D2 62 D2.(72 D2.62 i=1 h; En posant X. = sin(Az;) et Y; = cos(Az;) , on a : -X, -Xz +2-Kx.62.D2 =0 (1) 20 25 X1 +X -Ky. r2.D2 =0 (2) d'où : D=i 2 6 Kx +Ky De plus, l'équation (3) implique : X1.Y1 + X2.Y2 =0 I I faut donc résoudre le système d'équations X1 + X 2 - Ky.62.D2 =cl X, .Y1 + X 2 .Y2 = 0 10 sachant les contraintes xi +Y1z -1=0 X2 +Y2?-1=0 La solution est : JKy.(Kx +Ky) X1 (Kx + Ky) Ky K x + Ky 15 _ K, K _ Y 1 J1-X1 11 Kx +Ky Kx +Ky X2= XI Y2 = -Y1 Donc : Kxj+yK ' Y1 KxI+K Y y X2 =X1,Y2 =-Y1 D=1 2 6 K + K X y 20 5 II faut ensuite repasser dans l'espace initial, par application inverse de la matrice de passage P. Or, en prenant xporteur (t) X cible - Dh (t).sin(Az(t)) et Y porteur (t) = Ycible - Dh (t).cos(Az(t)) , ce qui s'écrit sous forme matricielle comme suit : xporteur (t) - D/ sin(Az(t)) _ able Dh (t) yporteur (t) Ycible l cos(AZ(t)) sin(arctan Y1 - 0) 1 cos(arctan Y1 - 0) 1 10 On en déduit que : [X x pl (t) Yplorteur (t) cible - D. Ycible X sin(arctan - 0) Y2 cos(arctan Y2 2 - 0) - _ Xcible - D. _ Ycible 2 xporteur (t) _ Y porteur (t) et 15 La première mesure virtuelle est donc définie par : - 61 = 6 choisi arbitrairement en fonction d'un critère opérationnel ou technique; X - Azi = arctan y' - 0 ; X pl X 1 2 sin(Azl ) porteur(t) cible Yporteur (t) Ycible 6 Ky + Ky . cos(Azi ) Et la deuxième mesure virtuelle est donc définie par : - 62 = 6 identique à celui de la première mesure; - Az2 = arctan X 22 - 0 = - arctan X' - 0 = -Azl - 20 ; Y2 Yi 20 xporteur (t) _Y porteur (t)_ X cible "cible _ 6 "\/Kx + Ky fcos(Az2) 1 2 sin(Az2 ) On peut remarquer que la deuxième mesure peut être calculée à partir de la première mesure, sous réserve d'avoir mémorisé 0. Par conséquent, on peut se limiter à la mémorisation d'une seule mesure. Pour activer l'algorithme de localisation à partir d'azimuts, il faut néanmoins transmettre les deux mesures, et donc reconstituer la deuxième mesure. Il faut noter que la valeur minimale de N' est liée au problème d'estimation à résoudre, notamment à la dimension du vecteur d'état recherché ainsi qu'à la dimension d'une mesure. Dans le présent exemple de réalisation, l'espace mémoire nécessaire est très faible: c'est l'espace nécessaire au stockage de deux mesures. De plus, la perte d'information dépend de l'application, mais elle est quasiment nulle dans le cas du présent exemple de localisation passive mono-porteur. Bien-sûr, le temps de calcul nécessaire à l'estimateur est d'autant plus faible que le nombre de mesures sur lequel l'estimateur travaille est faible. La figure 2 illustre, par un diagramme fonctionnel, un autre exemple de système d'estimation statistique distribué comportant cette fois deux modules 11 et 21 de compression de données selon l'invention. Avantageusement là encore, les modules de compression 11 et 21 mettent en oeuvre l'invention de façon logicielle. Les modules 11 et 21 sont insérés dans des systèmes passifs de détection et de localisation de cibles comportant des modules 14 et 24 respectivement d'estimation statistique. Les modules 14 et 24 effectuent des traitements classiques d'estimation statistique sur des blocs de mesures d'azimut fournies par des capteurs de détection passifs 12 et 22 respectivement embarqués sur des porteurs mobiles 10 et 20 respectivement, les mesures étant stockées dans des espaces mémoires 13 et 23 respectivement. Par exemple, les modules 14 et 24 peuvent implémenter des algorithmes de type Maximum de Vraisemblance. De manière similaire au module 1 de la figure 1, les modules 11 et 21 résument l'information contenue dans N mesures réelles et la synthétisent en N' mesures virtuelles, où N' REVENDICATIONS1. Procédé d'estimation d'un état X d'un système, le procédé incluant une étape de calcul, à partir de N mesures réelles effectuées sur le système par un capteur (2), où N est un entier supérieur ou égal à 1, d'un état estimé du système ainsi que d'une qualité estimée R de l'estimation, le procédé étant caractérisé en ce qu'il comporte une étape de compression permettant, si N=Nmax où Nmax est un entier prédéfini supérieur ou égal à 1, de déterminer N' mesures virtuelles, où N' est un entier supérieur ou égal à 1 et strictement inférieur à N=Nmax, de telle sorte que l'étape de calcul appliquée aux N' mesures virtuelles fournit le même état estimé ainsi que la même qualité estimée R que si elle est appliquée aux N=Nmax mesures réelles effectuées par le capteur. 2. Procédé selon la revendication 1, caractérisé en ce que les mesures sont stockées dans un espace mémoire (3) pouvant stocker au maximum N_max mesures, les mesures réelles étant remplacées, dans l'espace mémoire, par les N' mesures virtuelles dès lors que N=N_max. 3. Procédé selon la revendication 1, caractérisé en ce que les mesures sont envoyées via un réseau de données pouvant véhiculer au maximum N_max mesures, les mesures virtuelles étant envoyées via le réseau de données à la place des mesures réelles dès lors que N=N_max. 4. Procédé selon la revendication 1, caractérisé en ce que l'état X du système est un vecteur d'état (X,,...,Xd) de dimension d, où d est un entier strictement supérieur à 1, l'étape de calcul calculant un vecteur d'état estimé (X,,...,Xd) de dimension d et sa matrice de covariance R de dimension d°x°d. 5. Procédé selon la revendication 4, caractérisé en ce que l'étape de calcul calcule le vecteur d'état estimé (X,,...,Xd) par application de la méthode du Maximum de Vraisemblance, la matrice R étant obtenue par inversion d'une matrice de Fisher. 6. Procédé selon la revendication 5, caractérisé en ce que l'étape de compression inclut les étapes suivantes: une étape de calcul, à partir des N mesures réelles, de la matrice de Fisher J1 associée au vecteur d'état estimé (Xl,...,Xd) ; une étape de détermination des N' mesures virtuelles satisfaisant à l'équation J1=J2, où J2 est la matrice de Fisher associée au vecteur d'état estimé (Xl,...,Xd) exprimée en fonction des N' mesures virtuelles. 7. Procédé selon la revendication 4, caractérisé en ce que le capteur est un capteur de détection (2) embarqué sur un porteur et le système est une cible détectable par le capteur. 8. Procédé selon les revendications 5 et 7, caractérisé en ce que d=2, le vecteur d'état correspondant à la position (Xeible, Ycible) de la cible dans un plan et le vecteur d'état estimé correspondant à une estimation (X cible ,Ycible) de la position (Xcible,Ycible) 9. Procédé selon les revendications 6 et 8, caractérisé en ce que les matrices J1 et J2 sont données par: Jl = J2 = Ne cos2 (Azi ) Lr i=l Dz h .6i Nrna. sin(Azi ). cos(Azi ) D2 2 i=l hi ' i Ncos2(Azi) 2 z i=l Dh X i Ïsin(Azi).cos(Azi) Dz6-z i=1 hi .. Nrna. sin(Azi).cos(Azi ) Ld Dz i=l hi .6-2 i N; sin 2 (Azi ) 2 i=l Dht' i N' sin(Azi).cos(Azi ) i=l Dh 2.6iz N sin 2 (Azi ) D2 i=l h: '6z et les azimuts Azi et les distances horizontales Dhi dans la matrice J1 étant calculés en supposant que X,ihie =xporteur(t)+Dh(t).sin(Az(t)) et Ycihle = yporteur (t) + Dh (t).cos(Az(t)) , où (x porteur (t), yporteur (t)) est la positionconnue du porteur dans le plan à un instant t de mesure, où Dh(t) est la distance horizontale du porteur à la cible à l'instant t et où Az(t) est l'azimut entre le porteur et la cible à l'instant t, la variance 6`2 de la ième mesure étant connue.5