FR 3049354 A1 20170929 FR 1600503 A 20160325 1 PROCEDE DE GEOPOSITIONNEMENT PAR SATELLITES ET TERMINAL ASSOCIE La présente invention concerne le domaine du positionnement par satellites. L'invention concerne plus particulièrement un procédé de géopositionnement par satellites et un terminal associé. Un système de positionnement par satellites (ou GNSS pour « Global Navigation Satellite System" selon la terminologie anglosaxonne) utilise une constellation de satellites qui tournent autour de la Terre sur des orbites très précisément déterminées, c'est à dire pour lesquels il est possible de connaître à tout instant la position. Les orbites des satellites sont choisies de sorte qu'à toute heure, six à douze satellites soient visibles en tout point de la Terre. Chaque satellite, dit bi-fréquences, émet des signaux électromagnétiques de géopositionnement sur deux fréquences différentes (par exemple, L1=1 575,42 MHz et L2=1 227,6 MHz pour le système GPS et L1=1 575,42 MHz et E5=1 227,6 MHz pour le système GALILEO). Un récepteur GNSS reçoit les signaux émis par des satellites visibles et mesure la durée de propagation requise pour qu'une marque horaire transmise par un satellite lui parvienne. Les marques horaires sont codées sur des ondes porteuses par la technique de modulation de phase. Chaque satellite transmet ainsi un code pseudo-aléatoire qui lui est propre. Une réplique de la séquence du code est générée par le récepteur et le décalage que doit subir la réplique afin de coïncider avec le code reçu correspond à la durée de propagation du signal pour parcourir la distance satellite-récepteur. Cette durée multipliée par la vitesse de la lumière dans le milieu traversé donne une mesure de distance appelée pseudo-distance. A partir des mesures des pseudo-distances le séparant de chaque satellite visible, et de la connaissance de la position des satellites, le récepteur déduit sa position précise en latitude, longitude, et en altitude dans un repère terrestre par une résolution numérique voisine de la triangulation. Il peut également en déduire la date et l'heure précise dans le repère temporel du système GNSS. La référence de temps du récepteur, fournie par son horloge, ne 35 coïncide pas parfaitement avec la référence de temps des satellites de la constellation, ce qui induit un biais dans les mesures de temps de 3049354 2 propagation, donc de distance, égal au retard de la référence de temps du récepteur par rapport à la référence de temps des satellites. On emploie pour cela le terme de pseudo-distance ". Le biais de temps, commun à toutes les mesures, constitue une quatrième inconnue, outre les trois inconnues de 5 position, ce qui oblige à disposer d'au moins quatre mesures pour calculer la position. En outre, la position du terminal de réception est estimée en faisant un certain nombre d'approximations. La mesure de la pseudo-distance ne peut pas par exemple, s'affranchir des erreurs liées au système telles que le 10 manque de précision des éphémérides ou des horloges embarquées dans les satellites. La mesure de la pseudo-distance est également entachée d'erreurs liées aux interactions entre les signaux et les couches d'atmosphère (troposphère et ionosphère) qu'ils traversent. Le retard de propagation du signal dans la troposphère et l'ionosphère dépend de 15 l'inclinaison du trajet et de l'heure à laquelle il a lieu. Typiquement, les erreurs de positionnement du GNSS liées à l'atmosphère sont plus marquées de jour que de nuit et plus sensibles lorsqu'un satellite est proche de l'horizon qu'au zénith. Dans certaines applications, comme l'approche de précision en aéronautique, la précision de positionnement obtenue par une 20 mesure directe ou absolue de la pseudo-distance n'est pas suffisante. L'utilisation d'une mesure différentielle permet d'améliorer sensiblement la précision du positionnement. Elle consiste à transmettre par un canal dédié (VHF, UHF ou téléphonie cellulaire) des corrections des mesures de pseudo-distance élaborées à partir de mesures de pseudo- 25 distances provenant de récepteurs disposés dans des stations au sol et dont les positions sont connues très précisément et proches du récepteur embarqué. La mesure de la pseudo-distance séparant un récepteur au sol et un satellite est comparée à la distance théorique séparant ces deux dispositifs. La distance théorique est calculée à partir des coordonnées 30 spatiales respectives du récepteur au sol et du satellite qui sont connues à tout instant. La différence entre la mesure de distance et la distance théorique représente l'erreur de mesure, elle est calculée pour chaque satellite à chaque époque d'observation. Ces différences de distances constituent des termes correctifs (aussi appelés corrections différentielles) qui sont retranchés aux mesures de pseudo-distance réalisées par le 3049354 3 récepteur mobile. Ces corrections ont pour effet de supprimer pratiquement totalement les erreurs qui présentent une corrélation spatiale importante quelle que soit leur origine, système ou atmosphérique. Les corrections sont d'autant plus efficaces que les deux récepteurs sont proches. Toutefois, la 5 mesure différentielle n'élimine pas les erreurs liées aux réflexions du signal sur des objets proches de l'antenne du récepteur, ni les erreurs propres au récepteur (bruit thermique). Ces erreurs sont présentes sur le récepteur de référence ainsi que sur le récepteur embarqué, elles dégradent la mesure de positionnement lors de la correction différentielle ; la précision obtenue est de 10 l'ordre de quelques mètres. Pour améliorer la précision du positionnement, les récepteurs au sol et les récepteurs mobiles, embarqués sur des porteurs, peuvent aussi tirer profit d'une seconde information élaborée par le récepteur qui est la mesure de la phase de la porteuse, pour chaque signal satellite reçu. La mesure de 15 la phase instantanée de la porteuse reçue permet en effet de calculer une pseudo-distance, dite pseudo-distance de porteuse, entre le récepteur et le satellite, au même titre que la mesure de la phase instantanée du code pseudo-aléatoire. Cette pseudo-distance de porteuse subit les mêmes variations que la pseudo-distance de code, lorsque la distance entre le 20 récepteur et le satellite ou le biais de temps dû à l'horloge du récepteur varient. Cette pseudo-distance mesurée par la phase est a priori ambigüe puisque la phase est connue modulo 2u mais elle est beaucoup moins bruitée que les mesures de pseudo-distance de code. Une solution connue pour améliorer les mesures de pseudo-distance 25 consiste à lisser la mesure de pseudo-distance bruitée réalisée sur le code par les mesures de phase peu bruitées. Pour cela, le récepteur applique un filtre passe-bas, généralement du premier ordre, à la différence entre les mesures de pseudo-distance de code et de pseudo-distance de porteuse, puis ajoute cette différence filtrée à la mesure de pseudo-distance de 30 porteuse pour reconstituer la mesure de phase de code. Ce traitement est réalisé axe satellite par axe satellite. Si la mesure est différentielle, un lissage identique est appliqué sur les récepteurs de la station sol afin que l'erreur de traînage du filtre passe bas, due à la divergence entre le code et la porteuse liée aux fluctuations du retard ionosphérique, soit identique au sol et sur le 3049354 4 récepteur mobile, et ne perturbe pas la mesure de positionnement après de l'application de la correction. Le principal intérêt de la méthode de lissage axe par axe tient dans sa simplicité et dans l'absence d'effet de couplage entre les mesures des 5 pseudo-distances des différents satellites ou canaux, néanmoins elle n'est pas complètement satisfaisante. En effet, le gain sur la précision de la mesure est significatif uniquement lorsque le lissage est effectué avec une constante de temps longue ; et dans ce cas, la durée de réinitialisation pour récupérer la précision après une modification brusque de l'ensemble des 10 mesures disponibles, est également longue, par exemple lors de la disparition d'un satellite par masquage, d'une panne d'un satellite ou encore d'une panne d'un récepteur au sol dans le cas du GNSS différentiel. En présence de scintillation ionosphérique, notamment près de 15 l'équateur, les mesures de phase de porteuse sur les signaux émis par les satellites présentent des discontinuités. En effet, dans l'ionosphère, le rayonnement solaire provoque une ionisation partielle des atomes et molécules de gaz. On se trouve donc en présence de particules chargées (ions et électrons libres) qui sont susceptibles d'interagir avec une onde 20 électromagnétique incidente et d'influer sur sa propagation. La scintillation ionosphérique se manifeste sous forme d'évènements plus ou moins localisés, d'une durée typique de quelques heures et affectant les liaisons Terre-satellites. Durant un évènement de scintillation, en raison de la présence de particules ionisées dans l'ionosphère, le signal reçu sur Terre, 25 présente des fluctuations d'amplitude et de phase susceptibles de perturber fortement la liaison et de compromettre son exploitation. Du fait des fluctuations sur l'amplitude des signaux reçus par le terminal de géopositionnement, les boucles de poursuite qui assurent la démodulation du signal décrochent lorsque l'amplitude est trop faible. Ainsi, le terminal ne peut 30 plus se baser sur la continuité de la phase pour effectuer le lissage code porteuse. A chaque fois que la boucle décroche le filtre de lissage code porteuse doit être réinitialisé et on doit attendre que ce filtre converge de nouveau. A titre d'exemple, il faut attendre au moins une seconde pour ré-accrocher le signal et retrouver une mesure de pseudo-distance valide. Ces 3049354 5 scintillations ionosphériques entrainent donc une dégradation régulière de la précision. Dans le cas d'un récepteur standard ce phénomène peut ne pas être 5 gênant car les satellites de géopositionnement ne décrochent pas tous en même temps, et si les durées de décrochage sont courtes, le récepteur dispose toujours d'au moins quatre mesures disponibles. Dans le cas d'un récepteur aéronautique chaque décrochage est problématique car il faut réinitialiser le filtre de lissage code porteuse et 10 attendre que le filtre converge de nouveau. Pendant ce temps de convergence, qui est typiquement d'environ cent secondes, la mesure satellite n'est pas disponible pour calculer le point, ou très dégradée. Sachant que la durée moyenne entre deux décrochages est de dix secondes, on se retrouve avec des mesures satellites toujours indisponibles. 15 Avec les récepteurs aéronautiques actuels le positionnement devient indisponible en présence de fortes scintillations ionosphériques. Il est connu dans l'art antérieur des récepteurs conçus pour résister à la scintillation ionosphérique grâce à des bandes de boucle de poursuite 20 étroites. Cela fonctionne pour des récepteurs fixes, avec des horloges stables, mais pas pour des récepteurs aéronautiques. Un but de l'invention est notamment de corriger tout ou partie des inconvénients de l'art antérieur en proposant une solution permettant à un 25 récepteur de géopositionnement d'assurer la continuité de son positionnement même en présence de scintillation ionosphérique. A cet effet, l'invention a pour objet un procédé de géopositionnement par satellites mis en oeuvre par un terminal de géopositionnement à partir 30 d'au moins quatre satellites émettant chacun des signaux électromagnétiques de géopositionnement, dit bi-fréquences, sur deux fréquences différentes, chacun desdits signaux étant formés d'une fréquence porteuse modulée par un code d'étalement, ledit terminal comprenant au moins un module de réception configuré pour recevoir les signaux 35 électromagnétiques provenant de chaque satellite sur deux fréquences 3049354 6 différentes et au moins un module de calcul configuré pour traiter lesdits signaux, ledit procédé comprenant, pour chaque satellite de géopositionnement, une étape de calcul, en régime nominal, de quatre pseudo- 5 distances à partir des deux codes et des deux porteuses des signaux de géopositionnement bi-fréquences reçus, une étape de correction des retards ionosphériques sur chaque pseudo-distance calculée en appliquant un modèle de propagation de l'erreur ionosphérique, 10 une étape de lissage code porteuse à l'aide d'un filtre de Kalman afin de fournir une mesure de pseudo-distance sans bruit de mesure et de corriger le résidu d'erreur ionosphérique, la position du terminal de géopositionnement étant estimée en utilisant les pseudo-distances corrigées calculées pour chaque satellite. 15 Selon un mode de mise en oeuvre, le procédé comprend une étape de correction des retards troposphériques. Selon un mode de mise en oeuvre, le procédé comprend une étape de correction des erreurs système. Selon un mode de mise en oeuvre, le procédé comprend une étape de 20 compensation d'un biais inter-fréquences entre les signaux bi-fréquences. Selon un mode de mise en oeuvre, le filtre de Kalman est configuré pour propager un vecteur d'état comprenant la valeur de la pseudo-distance sans bruit de mesure, la valeur de l'erreur ionosphérique sur une fréquence après correction par le modèle de propagation de l'erreur ionosphérique et 25 les valeurs des l'ambiguïtés flottantes des mesures de phase de porteuse sur chacune des deux fréquences. Selon un mode de mise en oeuvre, le filtre de Kalman est configuré pour propager un vecteur d'état comprenant en outre la valeur du biais inter-fréquences entre les signaux bi-fréquences. 3049354 7 Selon un mode de mise en oeuvre, le modèle de propagation de l'erreur ionosphérique est un modèle de Kobuchar. Selon un mode de mise en oeuvre, le modèle de propagation de l'erreur ionosphérique est un modèle de Nequick. 5 Selon un mode de mise en oeuvre, le procédé comprend une étape de mémorisation de l'erreur ionosphérique. L'invention a également pour objet un terminal de géopositionnement par satellites comprenant au moins un module de réception configuré pour 10 recevoir des signaux électromagnétiques de géopositionnement, dit bifréquences, émis par au moins quatre satellites sur deux fréquences différentes et au moins un module de calcul configuré pour traiter lesdits signaux de géopositionnement et mettre en oeuvre le procédé de géopositionnement précédemment décrit. 15 D'autres particularités et avantages de la présente invention apparaîtront plus clairement à la lecture de la description ci-après, donnée à titre illustratif et non limitatif, et faite en référence aux dessins annexés, sur lesquels : 20 - La figure 1 représente un exemple de schéma synoptique d'un procédé de géopositionnement par satellites selon l'invention ; - La figure 2 représente un exemple d'organigramme de fonctionnement d'un filtre de Kalman récepteur selon l'invention. 25 Par la suite on supposera que chaque satellite de géopositionnement émet des signaux électromagnétiques de géopositionnement sur deux fréquences différentes. On parlera de signaux bi-fréquences et on notera ces deux fréquences, Fa et Fb. 30 On supposera également qu'il n'y a jamais de décrochage des deux porteuses en même temps, (ou pendant un laps de temps très court), et que le terminal de géopositionnement reçoit au moins un des deux signaux de géopositionnement. 3049354 8 La figure 1 représente un exemple de schéma synoptique, pour chaque axe satellite, d'un procédé de géopositionnement par satellites selon 5 un mode de mise en oeuvre préféré de l'invention. Ce procédé est mis en oeuvre par un récepteur ou terminal de géopositionnement à partir de signaux électromagnétiques de géopositionnement émis par au moins quatre satellites de géopositionnement. Comme énoncé précédemment, chacun de ces signaux 10 bi-fréquences est formé d'une fréquence porteuse modulée par un code d'étalement. Le procédé comprend, pour chaque satellite visible, une première étape Etp1 de mesure de pseudo-distances. Pour cela, le terminal de 15 géopositionnement comprend au moins un module de réception configuré pour recevoir ces signaux électromagnétiques provenant de chaque satellite sur deux fréquences différentes. De façon connu, ce module de réception peut comprendre au moins une antenne, un circuit analogique effectuant des amplifications et 20 changements de fréquence, un convertisseur analogique-numérique et au moins N canaux de traitement numérique. Chaque canal étant affecté à un satellite, l'entier N sera choisi supérieur au nombre de satellites a partir desquels on souhaite recevoir des signaux de géopositionnement. Chaque canal numérique reçoit un signal numérisé contenant 25 l'ensemble des signaux satellites qu'il soumet à une double boucle d'asservissement permettant d'une part de synchroniser en phase une porteuse localement générée avec la porteuse issue du satellite considéré et d'autre part de synchroniser un code pseudo-aléatoire généré localement avec un code identique présent dans le signal satellite et propre à ce 30 satellite. Les deux boucles d'asservissement peuvent comprendre chacune deux oscillateurs dont le contenu représente pour le premier oscillateur, la phase instantanée du code pseudo-aléatoire local (aligné sur le code présent dans le signal), qui constitue la mesure de la phase instantanée du code reçu 35 et pour le deuxième oscillateur, la phase instantanée de la porteuse locale 3049354 9 (alignée sur la phase de la porteuse présente dans le signal reçu du satellite, au déphasage près introduit par les circuits du récepteur), qui constitue la mesure de la phase instantanée de la porteuse reçue. La mesure de la phase instantanée du code dans chaque canal est 5 utilisée pour calculer une première donnée numérique appelée pseudodistance de code PDcode représentant une première mesure de la pseudodistance entre le récepteur et le satellite considéré. Cette mesure est non ambiguë mais assez bruitée. La mesure de la phase instantanée de la porteuse dans le canal 10 considéré est utilisée pour calculer une deuxième donnée numérique, appelée pseudo-distance de porteuse PDporteuse, représentant une deuxième mesure de la pseudo-distance entre le récepteur et le satellite considéré. Cette mesure est peu bruitée mais ambiguë. Les signaux étant bi-fréquences, en régime nominal, deux pseudo- 15 distances de code et deux pseudo-distances de porteuse peuvent être calculées pour chaque satellite visible du récepteur. On dispose donc, dans chaque canal numérique, de quatre mesures indépendantes (PDcode a , PDcode b , PDporteuse a , PDporteuse b) d'une même pseudo-distance séparant le récepteur de chaque satellite visible de ce dernier. Pour chaque satellite, les 20 deux types de mesures de pseudo-distances (PDcode et PDP°rteus obtenus sous forme de durées de propagation de signal entrelese) sont considéré et le récepteur le long de l'axe (axe satellite) joignant le satellite considéré et le récepteur. A l'instant n, les pseudo-distances PDc.ode a(n) , Ppoode b(n) et 25 PDporteuse a(n) , PDporteuse b(n), mesurées respectivement sur le code et la porteuse à la fréquence Fa et la fréquence Fb, sont données par les formules : PDcode a (n) = ( Tréception(n) - Tsat a(11) ) PDcode b (n) = ( Tréception(n) - Tsat b(n) ) 30 PDporteuse a (n) = ( Tréception(n) - (Pporteuse a(n) / Fa ) PDporteuse b (n) = ( Tréception(n) - (Pporteuse b(n) / Fb ) Dans lesquelles : - Treception(n) représente la date, à l'instant n de prise de mesure, donnée par l'horloge du récepteur, 3049354 10 - Tsat a(n) représente la date d'émission, par le satellite, du signal reçu à l'instant n par le récepteur, donnée par la phase de code locale, - (Pporteuse a(n) représente la phase de la porteuse locale, ramenée en fréquence de porteuse (cycles), 5 Pour chaque satellite visible du récepteur de géopositionnement, le procédé de géopositionnement comprend une étape Etp2 de correction des retards ionosphérique sur chaque pseudo distance calculée en appliquant un modèle de propagation de l'erreur ionosphérique. Pour cela, le terminal 10 comprend également au moins un module de calcul configuré pour traiter les signaux reçus. Le modèle de propagation de l'erreur ionosphérique peut être un modèle de Kobuchar, un modèle de Nequick ou tout autre modèle équivalent connu de l'homme du métier. 15 Le modèle de propagation de l'erreur ionosphérique fournit une estimation de l'erreur ionosphérique sur la mesure de pseudo-distance de code d'une fréquence, par exemple la fréquence Fa, en fonction de la position sur terre et de l'heure dans la journée. Cette erreur ionosphérique 20 représente le retard de groupe induit sur la propagation du signal lors de la traversée de l'ionosphère. Ce retard est proportionnel au contenu électronique total de la colonne d'atmosphère traversée par le signal et inversement proportionnel au carré de la fréquence de porteuse. Le retard se retrouve sur la mesure de code. L'effet sur la porteuse est une avance de 25 phase de la même valeur absolue. Si on considère l'erreur ionosphérique sur la fréquence Fa : sur la mesure de porteuse Fa, on a une avance de phase de - Biono sur la mesure de porteuse Fb, on a une avance de phase de -y.Biono, 30 sur la mesure de code Fa, on a un retard de Bionc, = CC . CET / Fa2 et sur la mesure de code Fb, on a un retard de y.Biono = a . CET / Fb2 avec y = Fa2/ Fb2 Donc pour corriger les mesures de code et de phase de porteuse on 35 applique les corrections suivantes : 3049354 11 Sur PDcode a (n) : - Biono modèle Sur PDporteuse a (n) : + Biono modèle Sur PDcode b (n) : - 'y Biono modèle Sur PDporteuse b (n) : + y Biono modèle 5 Bien entendu, un résonnement similaire peut être effectué en considérant un modèle fournissant une estimation de l'erreur ionosphérique sur la fréquence Fb. 10 Afin de traiter les signaux émis par les satellites de géopositionnement aux deux fréquences Fa et Fb, le module de réception comprend deux voies analogiques distinctes. Lorsque ce module de réception reçoit deux signaux de fréquences différentes, les temps de propagation dans chacune des voies analogiques du module peuvent être différents. Ceci a pour effet d'introduire 15 des différences de retard de groupe entre les deux fréquences Fa et Fb, appelés biais RF. Lorsque les mesures sont homogènes c'est-à-dire toutes mono-fréquences sur la même fréquence ou toutes bi-fréquences, les biais RF, communs à tous les satellites, n'ont pas d'effet sur la position résolue. Par 20 contre lorsque l'on mélange des mesures mono-fréquences sur des fréquences différentes ou mono-fréquence et bi-fréquences on introduit des distorsions entre les mesures satellite donc une erreur sur la position résolue. Ce cas de figure peut ce produire par exemple lorsqu'un des deux signaux bi-fréquences n'est pas reçus par le terminal de géopositionnement 25 par exemple suite à des problèmes de scintillation ionosphérique. Afin d'éviter cela, pour chaque satellite visible, le procédé peut comprendre une étape Etp3 de compensation d'un biais inter-fréquences entre les signaux bi-fréquences. Au cours de cette étape, le module de calcul identifie l'écart de biais 30 RF entre les deux fréquences sur les mesures de signaux bi-fréquences et corrige cet écart sur tous les satellites en compensant l'écart entre les deux fréquences uniquement sur les mesures de code et seulement sur une des deux fréquences. L'identification du biais RF et sa correction peut être mis en oeuvre en utilisant les techniques connues de l'homme du métier notamment 35 celle développé dans la demande de brevet FR 2 943 868. L'écart de biais 3049354 12 RF peut être estimé en faisant la moyenne de toutes les mesures de signaux bi-fréquences disponibles depuis le début (moyenne dans le temps et sur les satellites) et en considérant que les bais RF sont constants dans le temps. 5 Les mesures de code bi-fréquences de chaque satellite sont ensuite lissées par les mesures de porteuse bi-fréquences au cours une étape Etp4 de lissage des pseudo-distances compensées à l'aide d'un filtre de Kalman. Au cours de cette étape, le filtre de Kalman effectue simultanément le lissage code porteuse et la correction du résidu d'erreur ionosphérique. 10 En référence à la figure 2, on rappelle le principe de fonctionnement du filtrage de Kalman. Le filtrage de Kalman utilise un modèle d'état, établi sur la base de la 15 connaissance (avérée ou supposée) du comportement des grandeurs physiques inconnues que l'on cherche à déterminer et aux mesures disponibles. Ce modèle d'état est constitué : d'un vecteur d'état Xn représentant les grandeurs physiques 20 modélisées à l'instant tn, comportant un nombre de composantes Nétat d'un modèle de propagation, de la forme : Xn+1 = Fn. Xn + Vn Dans lequel : 25 Xn+1 est le vecteur d'état à l'instant tn+1 Fr, est la matrice de propagation sur l'intervalle [tu, tu+1], de dimension Nétat X Nétat (dans laquelle l'opérateur " x " représente le signe multiplié). , Vn estle vecteur bruit de propagation sur l'intervalle r je, 1.1+1 11 30 blanc, gaussien, de moyenne nulle, de matrice de covariance Qn = E[VnT . Vn] (où VuT est le vecteur transposé de Vu) et de dimension Nétat. d'un modèle d'observation, de la forme : 3049354 13 Zn = Hn . Xn Wn Dans lequel : Zn est le vecteur d'observation à l'instant tn, de dimension Nobs Hn est la matrice d'observation à l'instant tn, de dimension Nes x 5 Nétat Wn est le vecteur bruit de mesure à l'instant tn, blanc, gaussien, de moyenne nulle, de matrice de covariance Rn = E[WnT Wn] (où WnT est le vecteur transposé de Wn) et de dimension Nobs. 10 Dans ce modèle d'état, le vecteur d'état Xr, est de valeur a priori inconnue. Il n'est pas directement accessible par la mesure, contrairement au vecteur d'observations Zn, mais seulement par l'intermédiaire du modèle d'observation. Le filtre de Kalman réalise l'estimation du vecteur d'état par un calcul 15 de propagation, à partir du modèle de propagation, et par un calcul de recalage, à partir des observation et du modèle d'observation. Pour cela le filtre utilise deux variables : - le vecteur d'état estimé, noté Xnin après recalage à l'instant tn, noté Xn+11n après propagation entre les instants tn et tn+1, et noté Xn+11n+1 après 20 recalage à l'instant tn+1, de dimension Nétat- - la matrice de covariance de l'état estimé, notée Pn/n après le recalage à l'instant tn , notée P - n+1/n après la propagation entre les instants tn et tn+1, et noté Pn+11n+1 après recalage à l'instant tn+1 , de dimension Nétat X Nétat 25 Pour réaliser le calcul de propagation, le filtre utilise les formules : - Pour le vecteur d'état estimé : Xn+11n = Fn- Xn/n La matrice de propagation Fr, servant à établir une relation linéaire entre le vecteur d'état avant propagation et après propagation. 30 - Pour la matrice de covariance : Pn+1/n = Fn - l'en/n - FnT Qn (dans laquelle Fn T représente la matrice transposée de Fn) 3049354 14 Dans cette formule, les coefficients de la matrice de covariance Pn représentent la variance de chacune des composantes du vecteur d'état estimé (termes diagonaux) et la covariance des différentes paires de composantes de ce vecteur (termes non diagonaux). Cette matrice Pn 5 représente le degré de confiance qu'on attribue au vecteur d'état estimé. La matrice de covariance Qn du bruit de propagation permet de quantifier la part aléatoire et les approximations faites dans le modèle de propagation pour chaque composante du vecteur d'état Xn/n. 10 Pour réaliser le calcul de recalage, le filtre utilise les formules : - Pour le gain de recalage : 1 = Pn+1/n - Hn+1T - ( H n+1 - Pn+1/n - Hn+1T Rn+1 (où Hni-iT est la matrice transposée de Hn+1) La matrice d'observations Hn sert à établir une relation linéaire entre 15 le vecteur d'état et les mesures La matrice de covariance Rn caractérise le bruit propre des mesures. Rn est une matrice carrée de dimension Nobs x Nobs- Pour le vecteur d'état estimé : 20 Xn+1/n+1 = Xn+1/n Kn+1 - ( Zn+1 Hn+1 - Xn+11n ) Pour la matrice de covariance : Pn+1/n+1 = IdNétat Kn+1 . Hn+1 ) - Pn+1/n (dans laquelle Id Né. représente la matrice identité de dimension N état. Nétat) 25 Au départ, le vecteur d'état Xo est initialisé à zéro et la matrice Po est initialisée avec les variances et covariances représentatives de l'incertitude sur les grandeurs physique modélisées dans le vecteur d'état. Le calcul de propagation fait intervenir les matrices Fn et Qn pour 30 déterminer le vecteur d'état estimé propagé Xn+vr, à partir du vecteur d'état estimé recalé Xn1n et la matrice de covariance propagée P - n+11n à partir de la matrice de covariance recalée Pn/n. 3049354 15 Le calcul de recalage fait intervenir les observations Zn-F1 issues des mesures et les matrices Hn+1 et Rn+1 , pour déterminer le vecteur d'état estimé recalé Xn+11n+1 à partir du vecteur d'état estimé propagé Xn.vn et la matrice de covariance recalée P - n+1/n+1 à partir de la matrice de covariance 5 propagé Pnin. L'indice n est ensuite incrémenté de 1 (soit n = n+1) et la propagation, puis le recalage sont réitérés. On enchaine ainsi une alternance de propogation, recalage et incrémentation de l'indice. Sur la figure 2, on a représenté une mémoire recueillant Xn+1/n+1 et 10 fournissant Xn/n pour indiquer que la valeur de vecteur d'état introduite dans le modèle de propagation pour le calcul de Xn+11n daté de l'instant tn+1 est la valeur de vecteur d'état Xnin qui a été calculée à l'instant précédent tn. De même, on a représenté une mémoire pour la matrice de covariance Pn/n et pour la position estimée Gn/n- 15 Dans le filtre de lissage code porteuse selon le mode de mise en oeuvre préféré de l'invention le vecteur d'état X et le modèle d'évolution, représenté par les matrices F et Q, sont les suivants : X= PD F= 1 0 0 0 0 Q= fi PD 'I PD 0 0 0 0 B iono Sa Sb 0 /lion° 0 0 0 0 qIono 0 0 0 B RF 0 0 1 0 0 0 0 a n 0 0 0 0 0 1 0 0 0 o ab 0 0 0 0 0 1 0 0 0 . 0 0 20 Dans lesquels : PD représente la pseudo-distance sans bruit de mesure récepteur (bruit thermique, interférence, multi-trajets) et sans erreur ionosphérique ; PD = distance satellite récepteur + biais horloge récepteur x c + erreur troposphérique + biais satellite 25 - Biono représente l'erreur ionosphérique sur une fréquence (par exemple sur la fréquence Fa) après correction par le modèle (résidu) ; Sa représente l'ambiguïté flottante de la mesure de phase de porteuse sur la fréquence Fa ; 3049354 16 Sb représente l'ambiguïté flottante de la mesure de phase de porteuse sur la fréquence Fb ; BRF représente le biais inter-fréquences de la voie analogique de la fréquence Fb par apport à la voie analogique de la fréquence Fa ; 5 Xiono représente le facteur d'atténuation du modèle de Markov de l'erreur ionosphérique. = 1 - Art AT représentant la période de recalage du filtre et 'I la constante de temps du modèle de Markov (du 1ier ordre) de l'erreur ionosphérique 10 après correction par le modèle. On prendra AT = ls et i = 2000 s CIPD représente le bruit d'état sur la pseudo-distance vrai. On choisira une valeur très grande pour cette composante car on ne possède pas de modèle d'évolution fiable. On prendra par exemple qpip = (1000 m)2; qicino représente le bruit d'état du modèle de Markov (du 1 ier ordre) de 15 l'erreur ionosphérique après correction par le modèle ; Cliono = Giono2 - A-1-/ . i it Giono représentant l'écart-type de l'erreur ionosphérique après correction par le modèle. On prendra Giono = (20 m)2 qa représente le bruit d'état de l'ambiguïté de la mesure de phase de 20 porteuse sur la fréquence Fa. On prendra qa = (10-2 m) 2 qb représente le bruit d'état de l'ambiguïté de la mesure de phase de porteuse sur la fréquence Fb. On prendra qb = (10-2 m) 2 Suivant un mode mise en oeuvre alternatif, le vecteur d'état peut 25 comprendre quatre états. Dans ce cas, les matrices F et Q du modèle d'évolution sont également de dimension quatre. Le vecteur d'état X et le modèle d'évolution deviennent : PD 1 0 0 0 0 q PD 0 0 0 0 Biono O 2,',i0 O O O qt0b0 O O O X 8a F= 0 0 1 0 0 Q 0 qa 0 0 86 0 0 0 1 0 0 0 qb 0 3049354 17 En l'absence de saut de cycle de la boucle de phase de porteuse le bruit d'état sur les ambiguïtés est nul. Pour tenir compte du risque de saut de cycle (modèle non linéaire et non gaussien) on limite la constante de temps du filtre de lissage, qui tend naturellement vers l'infini. En effet, si la 5 constante de temps est trop grande l'effet du saut de cycle sur la pseudo- distance filtrée dure longtemps, le temps que le filtre converge et se recale sur la mesure de code (sans saut). En cas de saut de cycle à répétition, à faible rapport signal à bruit, les sauts vont se cumuler avant que le filtre n'ait le temps de converger sur la mesure de code, ce qui peut introduire une 10 erreur de mesure inacceptable. A titre d'exemple, cinq sauts de cycle donnent 1 mètre d'erreur. Pour limiter la constante de temps on utilise une valeur de bruit d'état non nulle. 15 Les matrices représentatives des mesures et du modèle d'observation s'écrivent : - PDcode a (n) PD porteuse a (n) - 2 PDcode ( n) codea 0 0 0 b 2 PDporteuse b (n) porteusea 0 0 Zu = R= 2 0- code b 0 0 0 0Cr 2 porteuse b +1 0 0 0 -1 1 0 0 +y 0 0 1 -y 0 1 1_ Dans lesquelles : 20 - PDcode a(n) représente la pseudo-distance mesurée sur le code de la fréquence Fa , à l'instant t(n) = n.AT - P Dporteuse a(n) représente la pseudo -distance mesurée sur le code de la fréquence Fa , à l'instant t(n) = n.AT - PDcode b(n) représente la pseudo-distance mesurée sur le code de la 25 fréquence Fb , à l'instant t(n) = n.AT H= 3049354 18 - n P . .... porteuse ID(n) représente la pseudo-distance mesurée sur le code de la fréquence Fb , à l'instant t(n) = n.AT code a et code b représentent respectivement l'écart-type de l'erreur de mesure de la phase de code propre au récepteur pour les fréquences a 5 et b. Gporteuse a et aporteuse b représentent respectivement l'écart-type de l'erreur de mesure de la phase de porteuse propre au récepteur pour les fréquences a et b. 10 Lors de l'initialisation, le vecteur d'état X et la matrice de covariance P ont pour valeur : 0 10002 0 0 0 0 0 0 Giono2 0 0 0 X0/0 P0/0 = 0 0 10002 0 0 0 0 0 0 10002 0 0 _ 0 0 0 0 CT 2 RF _ 15 Giono représente l'écart-type de l'erreur ionosphérique sur l'axe à vue du satellite considére. Sa valeur est donnée par des modèles de l'atmosphère, en fonction de l'inclinaison de l'axe satellite, de la latitude du récepteur, et de l'heure du jour. On prend Giono = 30 m. Giono représente l'écart-type de l'erreur sur le biais interfréquence sur 20 l'axe à vue du satellite considéré. On prend GRF = 3 m. Lorsqu'un satellite décroche sur une fréquence, c'est-à-dire qu'un des deux signaux bi-fréquences émis par un satellite de géopositionnement n'est pas reçu par le terminal ou est reçu avec une amplitude trop faible, le modèle 25 d'observation est modifié. Les deux lignes correspondant à ces mesures dans les matrices Z, R et H sont supprimées. Si par exemple, sur un axe satellite, les mesures sur la fréquence Fb deviennent indisponibles et seules les mesures sur la fréquence Fa sont disponible, les matrices représentatives du modèle d'observation 30 deviennent : 3049354 19 Z [PD code a (n) R- (3- 2code a 0 n PD porteuse a (n) [ 0 [1 +1 0 0 01 H= 1 -1 1 0 0 2 porteuse a De façon analogue, lorsque les mesures sur la fréquence Fa deviennent indisponibles le modèle d'observation s'écrit : 5 Z PD code b (n) R= n [PDporteuse b (n) 0 [62 code b 0 2 6 porteuse b [1 +y 0 0 11 H= 1 -y 1 0 1 Lorsque les mesures redeviennent disponibles, les deux lignes 10 associées à ces mesures dans les matrices Z, R et H sont rétablies. De plus, pour la fréquence réapparue, l'estimation de l'ambigüité de la mesure de phase du vecteur d'état X est réinitialisée à 0 avant le recalage sur les 4 mesures. Le coefficient diagonal (variance) associé à l'ambiguïté de phase de porteuse dans la matrice P est réinitialisé à une valeur très grande devant 15 les ordres de grandeur couramment utilisés. Cette valeur est augmentée de façon artificielle pour indiquer au filtre que l'estimation de l'ambigüité n'est plus bonne, et qu'il faut donc la réinitialiser. Cette variance peut être fixée à une valeur de 10002. 20 Lorsque les mesures sur la fréquence Fa redeviennent disponibles, le vecteur d'état X et la matrice de covariance P s'écrivent : PD n+11 n 0 0 0 0 0 Bionon+11 n 0 0 0 0 0 X n+11 n = 0 Pn+11 n = Pn+1In 0 0 10002 0 0 8 b n+11n 0 0 0 0 0 BRF 0 0 0 0 0 25 De façon analogue, lorsque les mesures sur la fréquence Fb redeviennent disponibles on obtient : 3049354 20 P D n+1 I n 0 0 0 0 0 B iono n+11 n 0 0 0 0 0 Xn+1In 8an+lln Pn+11 n Pn+11 n 0 0 0 0 0 0 0 0 0 10002 0 BRF 0 0 0 0 0_ Lorsque l'on observe une discontinuité sur la mesure de phase de 5 porteuse de l'une des deux fréquences, le vecteur d'état et la matrice de covariance, sont réinitialisés de la même manière avant le recalage sur la nouvelle mesure de phase de porteuse. Le Filtre de Kalman permet d'effectuer à la fois le lissage code porteuse et la 10 combinaison " iono-free ". En sortie du lissage code porteuse, la donné utile est la pseudo-distance PD filtrée et corrigée de l'erreur ionosphérique, correspondant à la première coordonnée du vecteur d'état estimé après recalage Xn-vun+1. 15 De façon avantageuse, du fait que l'on applique les corrections du modèle en amont du filtre de Kalman, celui-ci n'a plus qu'à identifier le résidu de l'erreur ionosphérique. L'erreur calculée grâce au modèle obéit au même modèle de Markov que l'erreur ionosphérique, donc à fortiori la différence aussi, avec une amplitude plus petite. La correction étant elle même 20 proportionnelle à l'inverse du carré de la fréquence, l'erreur résiduelle l'est aussi, donc le modèle d'observation liant l'état du filtre aux mesures bifréquence reste inchangé. Les modèles de Kobuchar ou de Nequick sont des modèles de prévision de l'erreur ionosphérique qui normalement sont utilisés pour 25 corriger des mesures mono fréquence. Dans le cas de mesures bifréquences, le fait de faire la combinaison " iono-free " de signaux sur deux fréquences différentes, permet normalement d'éliminer l'erreur ionosphérique. Dans le procédé selon l'invention, ce modèle est utile lorsque l'une des deux fréquences n'est plus disponible et que la mesure se retrouve 30 mono-fréquence.. 3049354 21 Lorsqu'une des deux fréquences n'est plus disponible sur une longue période, par exemple suite à un brouillage ou une amplitude trop faible, on passe d'un fonctionnement bi-fréquence à un fonctionnement mono-fréquence. Dans un premier temps le filtre va fonctionner naturellement sur 5 l'acquis du fonctionnement bi-fréquences avec le biais ionosphérique estimé lors de la présence des signaux bi-fréquences et petit à petit cette estimation du biais va se périmer parce que ce biais a changé. Le fait d'appliquer une correction de l'erreur ionosphérique en amont du filtre de Kalman, fait que le procédé va passer progressivement d'une combinaison " iono free » sur deux 10 fréquences à un fonctionnement équivalent à un cas mono-fréquence mais corrigé grâce au modèle de prévision de l'erreur ionosphérique. Une étape Etp5 de correction des retards troposphériques peut ensuite être appliquée sur la mesure de pseudo-distance en sortie du filtre de 15 Kalman pour chaque satellite. Cette compensation est obtenue en appliquant un modèle classique dépendant entre autre de l'heure de la journée, de la position géographique du satellite considéré. Une étape Etp6 de correction des erreurs " système " peut ensuite être 20 appliquée pour chaque axe satellite. Ces erreurs sont liées au principe du GNSS. Il s'agit par exemple de l'Effet Sagnac du au décalage temporel de la réception de deux signaux tournant en sens inverse. Parmi ces erreurs, on trouve également l'imprécision des horloges atomiques. Lors de cette étape l'effet relativiste est également corrigé. Ces corrections à effectuer sont 25 données par le satellite considéré à travers un message « navigation " contenant les termes correctifs. Ces termes correctifs sont communs aux quatre mesures issues d'un même satellite. Une fois ces corrections effectuées, le module de calcul estime la 30 position du terminal de géopositionnement en combinant les pseudodistances filtrées et corrigées calculées pour chaque satellite visible du terminal à l'aide de l'algorithme de résolution, PVT pour "Position, Velocity and Time" selon la terminologie anglosaxonne. Comme énoncé précédemment, au moins quatre pseudo-distances sont nécessaires. 35 3049354 22 Suivant un mode de mise en oeuvre alternatif, les étapes de correction des retards troposphériques Etp5 et des erreurs " système " Etp6 peuvent être réalisées en amont du filtrage de Kalman à condition d'appliquer ces corrections aux quatre mesures de pseudo distance. Le module de calcul peut comprendre un ou plusieurs microprocesseurs, processeurs, ordinateurs ou tous autres moyens équivalents programmés de façon opportune. REVENDICATIONS1. Procédé de géopositionnement par satellites mis en oeuvre par un terminal de géopositionnement à partir d'au moins quatre satellites émettant chacun des signaux électromagnétiques de géopositionnement, dit bi-fréquences, sur deux fréquences différentes, chacun desdits signaux étant formés d'une fréquence porteuse modulée par un code d'étalement, ledit terminal comprenant au moins un module de réception configuré pour recevoir les signaux électromagnétiques provenant de chaque satellite sur deux fréquences différentes et au moins un module de calcul configuré pour traiter lesdits signaux, ledit procédé étant caractérisé en ce qu'il comprend, pour chaque satellite de géopositionnement, une étape de calcul (Etp1), en régime nominal, de quatre pseudo-distances à partir des deux codes et des deux porteuses des signaux de géopositionnement bi-fréquences reçus, une étape de correction (Etp2) des retards ionosphériques sur chaque pseudo-distance calculée en appliquant un modèle de propagation de l'erreur ionosphérique, une étape (Etp4) de lissage code porteuse à l'aide d'un filtre de Kalman afin de fournir une mesure de pseudo-distance sans bruit de mesure et de corriger le résidu d'erreur ionosphérique, la position du terminal de géopositionnement étant estimée en utilisant les pseudo-distances corrigées calculées pour chaque satellite. 2. Procédé selon la revendication précédente dans lequel ledit procédé comprend une étape (Etp5) de correction des retards troposphériques. 3. Procédé selon une des revendications précédentes dans lequel ledit procédé comprend une étape (Etp6) de correction des erreurs système. 3049354 24 4. Procédé selon une des revendications précédentes dans lequel ledit procédé comprend une étape de compensation (Etp3) d'un biais inter-fréquences entre les signaux bi-fréquences, en amont de l'étape (Etp4) de lissage code porteuse. 5 5. Procédé selon une des revendications précédentes dans lequel le filtre de Kalman est configuré pour propager un vecteur d'état comprenant la valeur de la pseudo-distance (PD) sans bruit de mesure, la valeur de l'erreur ionosphérique (Biano) sur une fréquence après correction par le modèle de propagation de l'erreur 10 ionosphérique et les valeurs des l'ambiguïtés flottantes (sa, 6b) des mesures de phase de porteuse sur chacune des deux fréquences. 6. Procédé selon une des revendications précédentes dans lequel le filtre de Kalman est configuré pour propager un vecteur d'état comprenant en outre la valeur du biais inter-fréquences (BRF) entre les 15 signaux bi-fréquences. 7. Procédé selon une des revendications précédentes dans lequel le modèle de propagation de l'erreur ionosphérique est un modèle de Kobuchar. 8. Procédé selon une des revendications 1 à 6 dans lequel le modèle de 20 propagation de l'erreur ionosphérique est un modèle de Nequick. 9. Procédé selon une des revendications précédentes dans lequel ledit procédé comprend une étape de mémorisation de l'erreur ionosphérique. 10.Terminal de géopositionnement caractérisé en ce qu'il comprend au 25 moins un module de réception configuré pour recevoir des signaux électromagnétiques de géopositionnement, dit bi-fréquences, émis par au moins quatre satellites sur deux fréquences différentes et au moins un module de calcul configuré pour traiter lesdits signaux de géopositionnement et mettre en oeuvre le procédé de 30 géopositionnement selon l'une des revendications précédentes.