FR 2772133 A1 19990611 FR 8717042 A 19871208 La présente invention a tout d'abord pour objet un procédé de compensation des effets des perturbations magnétiques dues au porteur d'un magnétomètre, affectant la mesure de l'intensité d'un champ magnétique à mesurer par ledit magnétomètre, dans lequel on mesure à tout instant les trois composantes, dans un repère lié audit porteur, dudit champ magnétique dont on veut connaitre 1'intensité avec précision, on estime l'erreur due auxdites perturbations, et on compense lesdits effets. L'invention s'applique en particulier à la compensation des perturbations magnétiques dues à l'avion porteur d'un magnétomètre à champ total, par exemple un magnétomètre à résonance magnétique, qui affectent la mesure de l'inten- sité du champ magnétique terrestre, ou des variations de cette intensité, à l'aide de ce magnétomètre. De telles mesures sont utilisées à des fins de prospection, ou encore de détection de masses magnétiques en déplacement relatif par rapport au magnétomètre. Comme cela est connu, le fait qu'un avion soit réalisé en grande partie à l'aide de matériaux magnétiques et conducteurs entraîne que cet avion se comporte comme un ensemble de sources magnétiques qui perturbent la mesure de l'intensité du champ magnétique terrestre à l'aide d'un magnétomètre porté par cet avion. I1 est donc important d'évaluer précisément les effets des perturbations magnétiques dues à l'avion, c'est-à-dire de son aimantation, pour pouvoir compenser la mesure. L'article de P. LELIAK "Identification and evaluation of magnetic field sources of magnetic airborne detector equipped aircraft" dans IRE Transactions on aerospace and navigationnal electronics, September 1961 pages 95-105 montre que 1 'aiman- tation d'un avion peut être décomposée en trois aimantations distinctes, qui sont l'aimantation permanente, l'aimantation induite, et l'aimantation par courants de Foucault. Cet article montre aussi que l'erreur due aux effets de ces trois aimantations sur l'intensité du champ magnétique à mesurer est décomposable en une somme de termes élémentaires, chaque terme élémentaire apparaissant comme le produit d'un terme propre à l'avion par un terme lié à la position et à la vitesse de rotation de l'avion relativement à la direction du champ magnétique dont le magnétomètre à champ total mesure l'intensité. Chacun des termes propres à l'avion dépend de la répartition et des valeurs des sources magnétiques dues à l'aimantation du porteur. Chacun des termes lié à la position et à la vitesse de rotation de l'avion relativement à la direction du champ magnétique est accessible en particulier à partir de la connaissance des trois composantes, et de leurs dérivées temporelles, du champ magnétique dans un repère lié à l'avion. L'évaluation théorique de chacun des termes propres à l'avion étant très difficile, sinon impossible, on connalt déjà des procédés du type défini ci-dessus, dans lesquels, pour évaluer l'effet de l'aimantation de l'avion, on fait effectuer par celui-ci, dans une zone où l'on sait que le champ magnétique à mesurer est uniforme, une série de mouvements spécifiques, au cours d'un parcours d'identification, encore appelé passe d'identification". Par exemple, l'avion se dirigeant vers le nord, il effectue d'abord des oscillations de roulis, puis des oscillations de tangage, et enfin des oscillations de lacet. I1 change alors de cap pour se diriger vers l'est, et recommence la série d'oscillations précédentes, et ainsi de suite pour chacun des quatre caps. Les oscillations effectuées ont une fréquence déterminée, de l'ordre de quelques dixièmes de Hertz, et typiquement de 0,15 Hz, et une amplitude déterminée, par exemple + 100 en roulis et + 50 en tangage et lacet. A partir des variations d'intensité du champ magnétique, mesurées par le magnétomètre à champ total au cours de la passe d'identification, on estime les lois de variation de l'erreur dans le cas général. On prévoit alors de compenser cette erreur en appliquant au magnétomètre à champ total, par exemple, des champs magnétiques de compensation, créés par des sources commandées. Comme terreur dépend de la position et de la vitesse relatives du champ magnétique et de l'avion, on mesure à tout instant les trois composantes du champ dans un repère lié à l'avion, pour commander en conséquence les sources de compensation. Toutefois, il s'avère qu'un tel procédé ne donne pas de bons résultats, en particulier car il ne prend pas en compte le fait que les sources magnétiques dues à l'aimantation ne sont pas constantes et peuvent se trouver modifiées entre le moment où l'on fait le parcours d'identification, et le moment où l'on fait les parcours de mesure, en particulier à cause des phénomènes d'hystérésis. De plus, les oscillations à effectuer au cours du parcours d'identification sont à l'origine de contraintes mécaniques pour l'avion, et d'inconfort pour son équipage. Par ailleurs, la nécessité d'effectuer des parcours d'identification relativement souvent est une autre contrainte, puisque, pendant ces parcours, il n'est pas possible de faire des mesures. Enfin, à cause du bruit accompagnant la mesure, l'estjmation de l'erreur ne conduit pas toujours à des résultats satisfaisants. La présente invention vise à pallier cet inconvénient. A cet effet, elle a pour objet un procédé du type défini ci-dessus, caractérisé par le fait que - on échantillonne dans le temps les trois composantes mesurées et l'intensité mesurée, pour obtenir à chaque instant de rang K, K étant un entier naturel, trois composantes échantillonnées et une intensité échantillonnée, - on calcule, à partir des trois composantes échantillonnées une pluralité de variables formant une base dans l'espace des erreurs dues à l'aimantation du porteur, - on effectue un filtrage numérique desdites variables et de ladite intensité échantillonnée, au moins pour en éliminer les composantes basses frequences, et pour obtenir une pluralité de variables filtrées et une intensité filtrée, respec tivement, - on calcule une somme pondérée desdites variables filtrées, dans laquelle chaque dite variable filtrée est multipliée par un coefficient de pondération lié à l'aimantation du porteur, - on soustrait de ladite intensité filtrée ladite somme pondérée, pour obtenir une différence, et, - on calcule lesdits coefficients de pondération en fonction de ladite différence et desdites variables filtrées selon un algorithme récursif d'estimation, ladite somme pondérée représentant lorsque K est grand, la valeur de ladite erreur, et ladite différence représentant alors une mesure compensée de l'intensité dudit champ magnétique. Avec le procédé de l'invention, il n'est pas nécessaire d'effectuer des mouvements oscillatoires spécifiques avant d'effectuer un parcours de mesure, puisque, de toutes façons, lorsque K est suffisamment grand, c'est-à-dire au bout d'un certain temps de vol, les coefficients de pondération représentent la contribution des sources magnétiques dues à l'aimantation effective du porteur. I1 suffit alors, avant de procéder à un parcours de mesure, que l'avion effectue un certain nombre de mouvements proches de ceux qu'il effectuera au cours de ce parcours de mesure pour que chaque coefficient de pondération atteigne sa valeur correcte. De plus, lorsque les sources magnétiques dues à l'aimantation changent, du fait que l'algorithme de calcul des coefficients de pondération est adaptatif, ceux-ci évoluent durant les parcours, de façon à prendre en compte le changement des sources magnétiques. La mesure compensée est donc toujours correcte. Avantageusement, ledit porteur étant un aéronef, - on effectue ledit filtrage numérique pour éliminer les composantes desdites variables et de ladite intensité échantillonnée ayant une fréquence inférieure à sensible ment 0,04 Hz et une fréquence supérieure à sensiblement 0,6 Hz, et, - on fait effectuer audit aéronef, avant chaque parcours de mesure dans une zone où les variations dudit champ magné tique sont à mesurer, au moins un parcours d'identificatinn dans une zone où ledit champ magnétique est uniforme, ledit parcours d'identification comprenant une pluralité de vira ges séparés entre eux par des lignes droites, chaque virage étant relatif à un changement de cap d'une des valeurs 90 , 180 ou 270 , effectué avec un angle d'as siette latérale compris entre sensiblement 400 et sensible ment 600, avec une vitesse en cap de sensiblement 80/s, et chaque ligne droite étant de durée au moins égale à sensiblement 20 s. Dans ce cas, on obtient une compensation particulièrement rigoureuse des perturbations magnétiques dues au porteur, car, d'une part, les mouvements de l'avion au cours du parcours dtidentification sont proches de ceux qu'il effectuera au cours de la mesure, et d'autre part, la bande passante du filtrage est déterminée pour éliminer le plus possible les bruits perturbateurs. La présente invention a également pour objet un dispositif pour la mise en oeuvre du procédé précédent, et comprenant des moyens pour mesurer, à tout instant, les trois composantes, dans un repère lié audit porteur, dudit champ magnétique dont on veut connaître l'intensité avec précision, et des moyens, reliés auxdits moyens de mesure, pour estimer l'erreur due aux perturbations magnétiques dues au porteur et compenser leurs effets, dispositif caractérisé par le fait qu'il comprend - des moyens pour échantillonner dans le temps les trois composantes mesurées et l'intensité mesurée, et pour délivrer, à chaque instant de rang K, K étant un entier naturel, trois composantes échantillonnées et une intensité échantillonnée, - des premiers moyens de calcul, recevant les trois compo santes échantillonnées et agencés pour calculer une pluralité de variables formant une base dans l'espace des erreurs dues à l'aimantation du porteur, - des moyens de filtrage numérique desdites variables et de ladite intensité échantillonnée, au moins pour en éli miner les composantes basses fréquences, et délivrer une plura lité de variables filtrées et une intensité filtrée, respec tivement, - des deuxièmes moyens de calcul d'une somme pondérée desdites variables filtrées, dans laquelle chaque dite variable filtrée est multipliée par un coefficient de pondération lié à l'aimantation du porteur, - des moyens pour soustraire de ladite intensité filtrée ladite somme pondérée, et délivrer une différence, et, - des troisièmes moyens de calcul desdits coefficients de pondération en fonction de ladite différence et desdites variables filtrées selon un algorithme récursif d'estimation, ladite somme pondérée représentant, lorsque K est grand, la valeur de ladite erreur, et ladite différence représentant une mesure compensée de l'intensité dudit champ magnétique. La présente invention sera mieux comprise grâce à la description suivante de la mise en oeuvre préférée du procédé de l'invention et de la forme de réalisation préférée du dispositif de l'invention, faite en se référant aux dessins annexés, sur lesquels - la figure 1 représente une vue en perspective de l'avion porteur, du repère qui lui est lié, et du champ magnétique terrestre dans ce repère, - la figure 2 représente un schéma par blocs du dispositif de l'invention, - la figure 3 représente, de façon plus détaillée, le circuit de calcul des variables, ainsi que le circuit de filtrage et de retard du dispositif de la figure 2, - la figure 4 représente, de façon plus détaillée, le circuit de calcul des coefficients de pondération du dispositif de la figure 2, et, - la figure 5 représente, de façon schématique, plusieurs types de parcours d'identification pouvant autre effectués par l'avion porteur de la figure 1. Un procédé et un dispositif de compensation des effets de l'aimantation du porteur d'un magnétomètre, qui affectent la mesure de l'intensité H d'un champ magnétique H par ce magnétomètre sont donc maintenant décrits. Dans le but de simplifier l'écriture, et comme il est légitime de le faire dans les applications considérées, on utilisera, pour désigner un champ magnétique, la notation r réservée aux vecteurs ordinaires, alors qu'en toute rigueur, ce champ magnétique doit être noté H , puisqu'il s'agit en fait d'un vecteur axial, représentant un tenseur antisymétrique d'ordre 2. Dans l'exemple décrit ici, le champ magnétique à mesurer est le champ magnétique terrestre , et le magnétomètre est un magnétomètre à résonance magnétique nucléaire, porté par un avion 10 représenté sur la figure 1. Sur cette figure, on a représenté également le repère trirectangle O x y z lié à l'avion 10, l'axe Ox étant parallèle à l'axe longitudinal de l'avion et positif vers l'avant, l'axe Oy étant parallèle à l'axe transversal et positif vers la droite du pilote, et l'axe Oz étant orthogonal au plan Ox Oy et positif vers le bas. Les trois composantes du champ magnétique terrestre ff dans le repère O x y z sont Hxw H et H y z De façon connue, on s'intéresse aux variations de l'intensité H du champ terrestre If , et on utilise les mesures soit à des fins de prospection, pour une meilleure connaissance du champ terrestre, soit à des fins de détection, pour localiser une ou plusieurs masses magnétiques en déplacement relatif par rapport au magnétomètre, à partir de la mesure des perturbations du champ terrestre H liées à la présence de ces masses magnétiques. Compte-tenu de ordre de grandeur des variations précédentes, il importe que les mesures de l'intensité du champ soient faites avec une grande précision relative, de l'ordre de Le magnétomètre utilisé à cet effet est donc un magnétomètre de type connu, à résonance magnétique nucléaire, agencé pour délivrer un signal S qui représente l'intensité du champ magnétique dans lequel il se trouve. Or, comme cela est connu, en particulier par l'article déjà cité de P. LELIAE, lorsque l'on arrive à ce degré de précision, l'aimantation de l'avion 10 est à llorigine d'une erreur dans la mesure effectuée par le magnétomètre, du fait que le signal S représente en fait l'intensité d'un champ qui est la somme du champ magnétique terrestre et d'un champ perturbateur g produit par des sources magnétiques dues à l'aimantation de l'avion 10. On peut donc écrire S = |H + #H| (1) Comme l'intensité du champ perturbateur H est faible, relativement à l'intensité du champ , on admet que S = # H # + h. #H (2) le vecteur h représentant le vecteur unitaire de m8me direction que le champ magnétique terrestre H. Les trois composantes du vecteur h sont donc L'erreur e due à l'aimantation de l'avion 10 et qu'il faut donc évaluer vaut donc e = h . # H L'aimantation de l'avion peut être décomposée en trois aimantations distinctes qui sont : l'aimantation permanente, l'aimantation induite, et l'aimantation par courants de Foucault. L'aimantation permanente est due au fait que la structure de l'avion comprend des parties qui se comportent comme des aimants permanents ("fers durs"). Le champ pertur bateur A EI du à l'aimantation permanente est donc un p champ permanent fixe dans le repère O x y z, et dont les trois composantes dans ce repère sont Px, Py et Pz, respectivement. L'erreur dur à l'aimantation permanente vaut donc #H (6) p p ep = hx.Px + hy.Py + hz.Pz (7) L'aimantation induite est due à la perméabilité des matériaux ferromagnétiques ("fers doux") et elle dépend linéairement du champ inducteur, qui est ici le champ terrestre H . Le champ perturbateur #Hi dû à l'aimanta- tion induite est donc un champ qui est mobile dans le repère O x y z et il vaut #Hi = [S] . H (8) la matrice [S] représentant le tenseur de susceptibilité de l'avion, de la forme L'erreur due à l'aimantation induite vaut donc ei = r. t (10) soit # # # # # ei = H . [S] . H - h. [S] . # H # . h (11) 2 Soit, en remplaçant hz par : hz2 = 1 - hx2 - hy2 (12) L'aimantation par courants de Foucault apparait dans les parties conductrices de la structure et les effets Foucault sont proportionnels à la dérivée temporelle du champ magnétique terrestre H dans le repère O x y z. Le champ perturbateur #Hf du à l'effet Foucault est donc un champ mobile dans le repère O x y z et il vaut # # #Hf = [F] . H (14) H étant la dérivée temporelle du champ H dans le repère O x y z, [F] étant le tenseur des effets Foucault, de la forme L'erreur due à l'effet Foucault vaut donc # # ef = h . #Hf (16) # # ef = h. [F] . H (17) Soit en utilisant la relation hz hz# = - hx hx# - hy hy# (18) où hx#, hy# et hz# représentent les dérivées temporelles de hx, hy et hz, respectivement, Finalement, l'erreur totale e vaut e = ep + ei + ef (20) En utilisant les relations (7), (13), (19) et (20), il apparait que l'erreur e est une somme pondérée, chaque terme de cette somme étant le produit d'une variable, liée à l'attitude et à la vitesse de rotation de l'avion relativement à la direction du champ magnétique, par un coefficient de pondération propre à l'avion, et qui dépend des sources magnétiques dues à l'aimantation de cet avion. Le tableau I de la page 14 résume cette situation. Sur ce tableau, on n'a pas fait figurer le terme "constant", c'est-à-dire qui ne dépend pas de h , dans l'équation (13). Ce terme vaut # # H # Szz Ceci est légitime car, en pratique, on ne s'intéresse qu'aux variations de l'intensité du champ magnétique terrestre, et on filtre les composantes basses fréquences et donc la composante continue du signal S, ce qui fait que l'influence du terme constant ci-dessus ne se fait pas sentir sur la mesure. On peut dire en fait que les variables du tableau I forment une base dans l'espace des erreurs due à l'aimantation du porteur. En effet, on peut démontrer que Origine Variables Coefficients de pondération hx Px Aimantation hy Py permanente h~ P~ z z hx2 Sxx- Szz hy2 Syy- Szz Aimantation hx hy Sxy+ Syx induite hx hz Sxz+ Szx hy hz Syz+ Szy hx hx# Fxx- Fzz hh F hx hy Fyx Aimantation hx hz F~~ due aux hy hx# Fxy courants de Foucault hy hy# Fyy- Fzz hy hz Fzy hz hx# Fxz hz hy# Fyz TABLEAU I ces variables, au nombre de seize, sont indépendantes, et que, quelles que soient les sources magnétiques dues à l'aimantation de l'avion porteur, on peut exprimer l'erreur e comme une somme pondérée, c'est-à-dire comme une combinaison linéaire, de ces seize variables. Pour un jeu de valeurs donné de ces variables, l'influence des sources magnétiques propres à l'avion se fait sentir par l'intermédiaire des coefficients de pondération. Ainsi il apparat que l'erreur due à l'aimantation de l'avion peut être estimée facilement si l'on connaît les variables et les coefficients de pondération du tableau I. Pour déterminer les variables, on mesure les composantes Hx, H et Hz à l'aide d'un magnétomètre triaxial embar y qué sur l'avion 10, puis on procède i un échantillonnage temporel de ces composantes et enfin on calcule les variables. Pour déterminer les coefficients de pondération, on considère ces coefficients comme les coefficients a. inconnus d'un système linéaire que l'on résout par un algorithme recursif d'estimation, en l'occurrence l'algorithme des moindres carrés, grâce à un calculateur embarqué sur l'avion. Sur les figures 2,3 et 4ton a représenté un dispositif embarqué, pour la compensation des effets de l'aimantation de l'avion porteur 10. Sur ces figures figurent des blocs qui représentent des dispositifs de mesure, ou des circuits de traitement, ou encore des circuits de calcul. Les flèches qui sortent de ces blocs ou qui y pénètrent symbolisent des signaux échangés par ces différents blocs, selon la convention suivante : une flèche mince symbolise un signal, analogique ou numérique, représentatif d'une grandeur scalaire, et une flèche épaisse symbolise un signal, analogique ou numérique, représentatif d'une grandeur ayant une pluralité de composantes scalaires, c'est-àdire d'une grandeur vectorielle ou matricielle. En référance à la figure 2, un magnétomètre triaxial 4, de type connu, par exemple "Fluxgate", mesure les trois composantes Hx, Hy et Hz du champ magnétique terrestre dans le repère 0 x y z. Le magnétomètre 4 délivre trois signaux analogiques proportionnels aux trois composantes Hx, Hy et Hz du champ magnétique H dans le repère o x y z. Naturellement, il n'est pas nécessaire que le magnétomètre 4 soit très précis, une précision relative de 10-3 étant suffisante, ce qui est naturellement beaucoup moins que la précision relative de 10-7 que l'on recherche avec le magnétomètre à résonance magnétique nucléaire, tout au moins en ce qui conceren l'amplitude # H # du champ magnétique H . Compte tenu des remarques précédentes, on peut, bien sur, considérer que les effets de l'aimantation du porteur sur la mesure faite par le magnétomètre 4 sont négligeables. Trois convertisseurs analogique-numérique 5, identiques, sont commandés pour effectuer un échantillonnage temporel, et une numérisation, des signaux H, H et H, à la fréquence y d'échantillonnage F e d'un signal engendré par un circuit d'horloge non représenté. Les trois signaux en sortie des trois convertisseurs 5 sont appliqués à un filtre 12 d'harmonisation des voies, numérique, dont le rôle sera mieux compris dans la suite. Le filtre 12 est pourvu de trois sorties délivrant des signaux HxK, HyK et HzK, représentant les composantes échantillonnées à l'instant tK de rang K, K étant un entier naturel. Les signaux H , H et H sont appliqués à trois entrées xK YK z K d'un circuit de calcul numérique 6. Le circuit 6 est agencé pour calculer ici seize quantités h1 , ..., h. , ... et h16 @K @K @@K i étant un entier naturel compris entre 1 et 16, selon les formules du tableau II page 1s. h1K = HxK / MK h9K = h1K H1K# / MK h2K = HyK / MK h10K = h1K H2K# / MK h3K = HzK / MK h11K = h1K H3K# / MK h4K = h1K2 h12K = h2K H1K# / MK h5K = h2K2 h13K = h2K H2K# / MK h6K = h1K h2K h14K = h2K H3K# / MK h7K = h1K h3K h15K = h3K H1K# / MK h8K = h2K h3K h15K = h3K H2K# / MK avec : MK = module du champ de composantes HxK, HyK, H3K HiK# = dérivée temporelle de HiK TABLEAU II Comme il est facile de le remarquer en comparant les tableaux I et II, les quantités h1 , ..., hj , ... et @K @K h 6 du tableau II sont les valeurs, à l'instant t K de @@K rang K, des variables du tableau I. Un circuit de filtrage numérique 3', élimine une éventuelle composante continue de chaque variable hiK et, ici, leurs ocmposantes de fréquence inférieure à sensiblement 0,04 Hz et supérieure à sensiblement 0,6 Hz. De plus, le circuit de filtrage 3' comprend un filtre numérique du type connu dit "de blanchiment, agencé pour rendre uniforme le spectre fréquentiel du bruit qui accompagne chaque variable. Le circuit de filtrage numérique 3' délivre donc seize variables filtrées x ..., x ... et x16K' à des circuits de calcul numérique 7 et 9. Le circuit de calcul numérique 7 reçoit également seize coefficients de pondération a1K-1, ..., aiK-1, ... et a16K-1 en provenance du circuit de calcul numérique 9, et il est agencé pour calculer une somme pondérée eK des variables filtrées x1K, ..., xiK, ... et x16K dans laquelle chaque variable filtrée x. est multipliée par le coefficient de 1K pondération aiK-1# On peut donc écrire : L'intensité mesurée S, en sortie du magnétomètre 1 à résonance magnétique nucléaire, est échantillonnée à la fréquence d'échantillonnage F dans un circuit d'échantil e lonnage 2. Sur la figure 2, dans un souci de simplicité, on a représenté le circuit d'échantillonnage 2 extérieur au magnétomètre 1. En pratique, et de façon connue, le signal de sortie S du magnétomètre est un signal qui résulte d'un comptage, puisqu'il mesure la fréquence d'un signal analogique. I1 s'agit donc d'un signal numérique déjà échantillonné, et il suffit de commander le magnétomètre 1 à résonance magnétique nucléaire pour que la fréquence des échantillons soit Fe. Dans ce cas, le circuit d'échantillonage 2 est donc inclus dans le magnétomètre 1 à résonance magnétique nucléaire. Le circuit d'échantillonnage 2 est relié à un filtre 45 d'harmonisation de voies, numérique, dont le rôle sera mieux compris dans la suite. Le filtre 45 délivre en sortie une intensité échantillonnée SK, qui est appliquéeà l'entrée d'un circuit 6' de correction et de retard, dont le rôle sera mieux compris dans la suite. Le circuit 6' est relié à un circuit de filtrage numérique 3 analogue au circuit de filtrage numérique 3', et qui délivre une intensité filtrée y. Ainsi, la composante continue de cette intensité filtrée yK est nulle, ainsi que ses composantes de fréquence inférieure à sensiblement 0,04 HZ et supérieure à sensiblement 0,6 HZ, et le bruit qui l'accompagne est blanc, c'est-à-dire de spectre fréquentiel uniforme. L'intensité filtrée yE et la somme pondérée e K sont appliquées à un circuit de soustraction 8 qui délivre une différence iK valant = = YK K eK (22) La différence #K est appliquée au circuit de calcul numérique 9 qui, comme cela a été indiqué, calcule les coefficients de pondération a1K-1, ..., aiK-1, ... et a16K-1# Le rôle des filtres 12 et 45 d'harmonisation des voies est de rendre aussi semblables que possible la fonction de transfert de la voie de mesure, et la fonction de transfert de la voie de compensation. La fonction de transfert de la voie mesure est le rapport des fonctions complexes de Fourier de la valeur SK et de la valeur vraie de l'intensité, tandis que la fonction de transfert de la voie de compensation est le rapport des fonctions complexes de Fourier de la valeur de chaque signal H , H et H et de la valeur vraie de la composante correspondante K YK K du champ magnétique. Pour parvenir à une bonne harmonisation des deux voies, on détermine chacun des filtres 12 ou 45 pour que sa fonction de transfert soit identique à celle des circuits disposés en amont de l'autre filtre 45 ou 12. Ainsi le filtre du type connu passe-bas, du premier ordre. La valeur sN de son signal de sortie, à l'instant de rang N vaut sN= (eN + eN-1) /2 (23) eNet eN-1 étant les valeurs, à l'instant de rang N et à l'instant de rang (N-1), respectivement, de son signal d'entrée. Le circuit 66 de retard introduit ici un retard égal à une période d'échantillonnage, c'est-à-dire que la valeur de son signal de sortie, à l'instant de rang N, est égale à la valeur de son signal d'entrée, à l'instant de rang N-i. Le rôle des circuits 61 et 66 sera mieux compris dans la suite. Dans le circuit 6, il est prévu trois circuits 61 identiques au circuit 61 du circuit 6', chacun de ces circuits recevant un des signaux H g H et H . Chaque circuit 61 xK yK zK du circuit 6 est suivi d'un circuit 66 identique au circuit 66 du circuit 6'. Les trois signaux H' ,H' et H' en sortie des trois circuits 66 sont XK K K appliqués à trois entrées d'un circuit de multiplication 60, ainsi que, ici, aux trois entrées d'un circuit 67 de détermination du module MK du champ magnétique. Par ailleurs, les trois signaux H , H et H sont xK,HyK et HzK appliqués aux trois entrées d'un circuit 62 de calcul des dérivées. Le circuit 67 de détermination du module M K est ici agencé pour calculer le module M K selon la relation Le signal M K est appliqué à une entrée du circuit de multiplication 60 Le circuit 62 de calcul des dérivées comprend trois sorties, et trois voies indépendantes entre chacune de ses trois entrées et chacune de ses trois sorties. Sur chacune des voies, la valeur sN du signal de sortie, à l'instant de rang N, vaut sN = (eN - eN-1) / Te (25) eN et eN-1 étant les valeurs, à l'isntant de rang n et à l'instant de rang (N-1), respectivement, du signal d'entrée de la voie, et T étant la période d'échantillonnage, inver e se de la fréquence F e Les trois signaux de sortie du circuit 62 sont appliqués aux trois entrées d'un circuit 66 identique aux circuits 66 déjà rencontrés, mais comprenant trois voies indépendantes. Les trois signaux de sortie de ce circuit 66 sont appliqués à trois premières entrées d'un circuit 65 de pondération. Par ailleurs, les trois signaux de sortie des trois circuits 61 sont appliqués aux trois entrées d'un circuit 63 de calcul des dérivées, ici identique au circuit 62. Les trois signaux de sortie du circuit 63 sont appliqués aux trois entrées d'un circuit 64 de retard, pourvu de trois sorties et comprenant trois voies indépendantes entre chacune de ses trois entrées et chacune de ses trois sorties. Sur chacune des voies la valeur sN du signal de sortie, à l'instant de rang N, vaut eN et e-N-1 étant les valeurs, à l'instant de rang N et à l'instant de rang (N-1), respectivement, du signal d'entrée de la voie, et s-N-1 étant la valeur du signal de sortie à l'instant de rang (N-l). Comme cela est connu, le circuit 64 introduit en fait, sur chaque voie, un retard d'une demipériode d'échantillonnage T e Les trois signaux de sortie du circuit 64 sont appliqués à trois deuxièmes entrées du circuit 65 de pondération. Ce circuit 65 est pourvu de trois sorties, et comprend entre d'une part chacune de ses trois premières entrées et chacune de ses trois deuxièmes entrées et, d'autre part, chacune de ses trois sorties, trois voies indépendantes. Sur chacune des voies, la valeur sN du signal de sortie, à l'instant de rang N, vaut ici 1 2 eN- et eN- étant les valeurs, à l'isntant de rang N, du signal sur la première entrée de la voie et du signal sur la deuxième entrée de la voie, respectivement. Les trois signaux de sortie du circuit 65 qui, comme cela sera mieux compris dans la suite, représentent les dérivées temporelles H , H et Hz, sont appliquées au circuit de @K @K @K multiplication 60. Le circuit de multiplication 60, qui reçoit donc les signaux Hx , H' , H' , après leur passage dans les circuits 61 K YK z K et 66 et les signaux MK, Hx , H et H z est agencé pour, @K yK zk après inversion de MK, effectuer les multiplications nécessaires pour délivrer les signaux h1 , ..., hi , ... et h16, K K K conformément au tableau II de la page 17. Avant d'aborder la description détaillée du circuit de calcul numérique 9, on notera que les circuits 61 à 66 ont principalement pour fonction de calculer chaque dérivée H, H ou H comme moyenne pondérée entre une dérivée XK YAK z K aux différences finies et une dérivée homographique. Ceci sera bien compris en imaginant une variante des circuits 6 et 6' dans lesquels les circuits 61 seraient supprimés, la fonction de transfert du circuit 62 étant divisée par la fonction de transfert d'un circuit 61. Dans ce cas les signaux en sortie du circuit 62 représenteraient les dérivées homographiques des signaux appliqués à ses entrées, dérivées valables à l'instant de rang N, comme cela est connu, et les signaux en sortie du circuit 63 les dérivées aux différences finies, dérivées valables à l'instant milieu de l'intervalle entre l'instant de rang N et l'instant de rang N-l, comme cela est également connu. Le circuit de retard 64 ajoutant un retard d'une demi-pério de T /2 de façon à ce que ces dérivées soient valables à e l'instant de rang N-l, et donc en même temps que les dérivées homographiques, retardées d'une période To par les circuits 66, les signaux en sortie du circuit 65 résulteraient bien d'une moyenne pondérée des deux types de dérivées. I1 en va de même dans les circuits 6 et 6', dont la structure permet d'éviter les risques dtinstabilité liés à l'emploi d'un circuit 62 qui calculerait véritablement la dérivée homographique. La méthode de calcul des dérivées des composantes du champ magnétique, qui se trouve ainsi mise en oeuvre, sans être obligatoire, permet d'améliorer la qualité de la compensation par rapport aux méthodes qui n'utiliseraient que les dérivées aux différences finies, ou que les dérivées homographiques. En référence à la figure 4, le circuit de calcul numérique 9 comprend lui-même quatre circuits auxiliaires de calcul numérique 91, 92, 93 et 94. Le premier de ces circuits auxiliaires, le circuit 91 reçoit les seize variables filtrées x1 e.., xi , ... et @K @K x16K et une quantité matricielle PK-1, et il est agencé pour calculer un vecteur KK selon la formule formule dans laquelle - XKT désigne le vecteur ligne transposé d'un vecteur colonne dont les composantes sont les variables filtrées x1K, ..., xiK, ... et x16K, et - est un coefficient inférieur ou égal à 0,5, et dont le rôle sera précisé dans la suite. Le deuxième des circuits auxiliaires du circuit 9, le circuit 92 reçoit le vecteur KK et les seize variables filtrées x1K, ..., xiK, et x16K, et il est agencé pour calculer la matrice PK selon la formule Cette matrice PK étant mémorisée pour le cycle suivant où elle prendra le sens de Ainsi le circuit 92 est, de façon connue, pourvu d'une mémoire dans laquelle est mémorisée la valeur de la matrice pK calculée à l'instant tK. Un tel circuit, ,qui effectue donc un calcul itératif, nécessite donc la connaissance de la valeur initiale, à l'instant to, de la quantité calculée. On adopte ici Po = Po I (30) I étant la matrice unité, et po étant compris entre sensiblement 104 et sensiblement 107. Le troisième des circuits auxiliaires du circuit 9, le circuit 93, reçoit la différence #K et le vecteur KK de composantes k1K, ..., kiK, ... et k16K et calcule seize produits k1. #K,..., ki@. #K,... et k16 #K. K - K K Le quatrième des circuits auxiliaires du circuit 9, le circuit 94, reçoit ces seize produits et calcule chacun des coefficients de pondération a. selon la formule 1K a. = ai@ + ki EK (31) 1K 1K-1 1K Le circuit 94 est donc pourvu d'une mémoire dans laquelle sont mémorisées les valeurs, lors du précédent calcul, des coefficients calculés, ici les seize coefficients a1K-1, ..., aiK-1, ... et a16K-1, calculés à l'instant tK-1. Les valeurs initiales adoptées sont ici ai = O (32) O Le dispositif qui vient d'être décrit fonctionne de la façon suivante. Si l'on appelle AK-1 le vecteur dont les composantes sont les seize coefficients a1@, ..., a. K-l ... et a16K-1, on peut écrire les relations (21) et (31) sous les formes respectives : eK = XKT. AK-1 (33) AK = AK-1 + KK. #K (34) En remplaçant, dans la relation (22), eK par son expres K sion dans la relation (33) on obtient : #K = yK - XKT AK-1 (35) Et en remplaçant, dans la relation (34), #K par son expression dans la relation (35), on obtient AK = AK-1 + KIK (yK - XK@ AK-1) (36) La relation (36), rapprochée des relations (28) et (29) montre que les coefficients de pondération a. , calculés par le circuit de calcul numérique 9 en fonction de la différence #K et des variables filtrées x1K, ..., xiK, ... et x16K, sont calculés selon l'algorithme des moindres carrés récursifs, tel que cet algorithme est décrit par exemple dans l'ouvrage de M. LABARRERE, J.P. KRIEF et B. GIMONET "Filtrage et optimisation stochastique", Ecole Nationale Supérieure de l'Aéronautique et de l'espace, Paris, 1974, pages 180-181. En conséquence, lorsque K est assez grand, chaque coefficient de pondération a est une bonne estimation du 1K coefficient a. du système linéaire suivant 1K où nK représente un bruit. En conséquence, lorsque l'avion 10 se déplace dans un champ magnétique terrestre H rigoureusement constant, il est évident que les variations observées du signal S, à savoir la suite des échantillons YE, ne dépendent que de l'aimantation de l'avion 10. Lorsque K est assez grand, les coefficients de pondération a1K, ..., aiK, ... et a16K sont donc une bonne estimation des coefficients du tableau I. Ainsi, lorsque K est grand, ces coefficients de pondération sont bien liés à l'aimantation de l'avion 10. Alors la somme pondérée e K représente une bonne estimation de l'erreur due aux effets de l'aimantation de l'avion 10. Lorsque, pour les passes de mesure, l'avion se déplace dans un champ magnétique terrestre dont les variations sont à mesurer, les coefficients de pondération a1K, ..., a i K K et a restent les mêmes, car les variations à mesurer, K sans corrélation avec les variables x1~~, ... x. , ... et K x16K, apparaissent comme du bruit nK pour le circuit de calcul 9. De plus, les variations des variables x1K, ..., xiK, ... et x16K étant faibles dans la bande de fréquences utile pendant les passes de mesure, KK est très faible, et les coefficients a1, ..., aj, ... et a16 évoluent très K K K peu. La différence #K représente alors une mesure compensée de l'intensité du champ magnétique, puisque l'erreur eK due à l'aimantation de l'avion porteur 10 a été soustraite des échantillons Le coefficient y correspond, de façon connue, à une pondération exponentielle des observations passées qui permet de tenir compte du fait que, comme cela a déjà été signalé, les coefficients a1--, ..., ai~~, ... et a16 ne sont pas K 1K 10K rigoureusement constants, mais variables à cause, notamment, des phénomènes d'hystérésis, ou de l'évolution magnétique du porteur en cours de vol, par exemple lorsqu'il largue un objet magnétique. Le coefficient &gamma; est à déterminer de façon expérimentale dans un cas pratique donné. I1 dépend de la vitesse de variation des coefficients a1~~, ..., a. , .. et a16 , et K K K reste généralement faible. Sur la figure 5, on a représenté, à titre d'exemple, trois parcours d'identification P P, Pb et P, qui peuvent avantageusement être effectués par l'avion porteur 10, avant un parcours de mesure, dans une zone où ledit champ magnétique est sensiblement uniforme, afin de permettre une bonne identification des coefficients a1~~, ..., ai , et K 1K a16 K Le parcours P comprend une ligne droite de départ D, a de cap indifférent, suivie d'un virage V2, ici à droite, relatif à un changement de cap de 1800, et effectué avec un angle # d'assiette latérale de sensiblement 55. Le parcours comprend ensuite une ligne droite D, un virage V1, ici à gauche, relatif à un changement de cap de 90 , et d'angle # égal à 550, une ligne droite D, un virage analogue au virage V2, une ligne droite D, un virage analogue au virage V1 et enfin une ligne droite D. Le parcours Pb comprend une ligne droite de départ D, de cap indifférent, suivi d'un virage V3, ici à droite, relatif à un changement de cap de 2700, et d'angle # égal à 40 , une ligne droite D, un virage analogue au virage V3, mais d'angle 9 égal à 500, une ligne droite D, un virage analogue au virage V3, mais d'angle q = 600 et enfin une ligne droite D. Le parcours P comprend une ligne droite de départ D, de c cap indifférent, suivi d'un virage analogue au virage V1, mais à droite, relatif à un changement de cap de 900, d'angle q égal à 55 , une ligne droite D, un virage analogue au virage V1, une ligne droite D, un virage analogue au virage V1, une ligne droite D, un virage analogue au virage V1, mais à droite, et enfin une ligne droite D. Dans les parcours précédents, les virages correspondent à une vitesse en cap de sensiblement 8 /s et chaque ligne droite dure au moins sensiblement 20 s. Naturellement, la liste des parcours qui viennent d'entre décrits n'est pas limitative et il est suffisant qu'un parcours d'identification présente une pluralité de virages, chaque virage présentant les caractéristiques de changement de cap, assiette latérale et vitesse en cap, qui ont été indiquées, et étant séparé du suivant par une ligne droite durant au moins 20 s. I1 est à noter que les caractéristiques ci-dessus des parcours d'identification correspondent à des variations des signaux utiles situées dans la bande passante utile des filtres 3 et 3'. Naturellement, la portée de la présente demande n'est pas limitée à la description qui vient dire faite. C'est ainsi que, par exemple, les trois convertisseurs 5 de la figure 2 peuvent être remplacés par trois échantillonneurs-bloqueurs, suivis d'un multiplexeur et d'un unique convertisseur analogique-numérique. De même, on a constaté qu'il était avantageux, du point de vue de la précision du résultat, de remplacer le circuit 67 de détermination du module MK par un circuit 67'de filtrage passe-bas recevant le signal S'K représentatif de l'intensite mesurée SV après passage dans le circuit 6' de correction et de retard, comme cela est symboliquement représenté sur la figure 3. De même, les circuits 61 et 62 de la figure 6, dont les fonctions de transfert décrites sont ici particulièrement simples, peuvent être modifiés, de façon connue, pour que la valeur sN de leur signal de sortie à l'instant de rang N dépende aussi de sa valeur SN-1 à l'instant de rang N-l. De même, on peut prévoir la possibilité de réinitialiser, avec ou sans remise à zéro des coefficients a1K, ..., aiK, ... et a16K, l'algorithme d'estimation de ces coefficients. De même, le signal S en sortie du magnétomètre 1 à résonance magnétique nucléaire n'est pas nécessairement le signal délivré directement par ce dernier, et il peut avoir été déjà compensé pour tenir compte des effets gyromagnétiques par tout dispositif connu, et en particulier celui décrit dans la demande parallèle française nO 87 14 094 au nom de la demanderesse. De même, il est évident pour l'homme de métier que la décomposition des circuits de calcul du dispositif de l'invention selon les blocs des figures 2 à 4 n'est pas obligatoire, et qu'il est évidemment possible, sans sortir du cadre: de l'invention, qu'un unique calculateur embarqué effectue les différentes tâches dévolues, dans la description qui vient d'être faite, aux différents blocs. De mQme le magnétomètre utilisé n1est pas obligatoirement un magnétomètre à résonance magnétique nucléaire. Ainsi, il peut s'agir d'un magnétomètre à résonance magnétique électronique, et de façon plus générale, de tout magnétomètre du type "à champ total". Enfin, le porteur du magnétomètre n'est pas nécessairement un aéronef. C'est ainsi qu'il peut s'agir, en particulier, d'un véhicule terrestre, ou encore d'un dispositif flottant. Revendications 1 - Procédé de compensation des effets des perturbations magnétiques dues au porteur d'un magnétomètre, affectant la mesure de l'intensité d'un champ magnétique à mesurer par ledit magnétomètre, dans lequel on mesure à tout instant les trois composantes, dans un repère lié audit porteur, dudit champ magnétique dont on veut connaître l'intensité avec précision, on estime l'erreur due auxdites perturbations, et on compense lesdits effets, procédé caractérisé par le fait que - on échantillonne dans le temps les trois composantes (Hx, Hy, Hz) mesurées et l'intensité (S) mesurée, pour obtenir à chaque instant (tK) de rang K, K étant un entier K naturel, trois compsantes échantillonnée (HxK, HyK, HzK) et une intensité échantillonnée (SK), - on calcule, à partir des trois composantes échantillonnes (HxK, HyK, HzK) une pluralité de variable (h1K, ..., hj, ..., h16) formant une base dans l'espace des erreurs 1K K dues à l'aimantation du porteur (10), - on effectue un filtrage numérique desdites variables (h1K, ..., hiK, ... h16K) et de ladite intensité échantil lonnée (SE), au moins pour en éliminer les composantes bas ses fréquences, et pour obtenir une pluralité de variables filtrées (x1' XiK .... xl6) et une intensité fil K trée (YK), respectivement, - on calcule une somme pondérée (eK) desdites variables filtrées, dans laquelle chaque dite variable filtrée (si) est multipliée par un coefficient de pondération (aiK-11K) lié à l'aimantation du porteur (10), - on soustrait de ladite intensité filtrée (yK) ladite somme pondérée (eK), pour obtenir une différence (#K), et, - on calcule lesdits coefficients de pondération (a1, ..., K aiK, ... a16K) en fonction de ladite différence (#K) et desdites variables filtrées (x1K, ..., xiK, ..., x16K) selon un algorithme récursif d'estimation, ladite somme pondérée (eK) représentant lorsque K est grand, la valeur de ladite erreur, et ladite différence (# ) (i) représentant alors une mesure compensée de l'intensité dudit champ magnétique. 2 - Procédé selon la revendication 1, dans lequel, ledit porteur étant un aéronef, - on effectue ledit filtrage numérique pour éliminer les composantes desdites variables et de ladite intensité échantillonnée ayant une fréquence inférieure à sensible ment 0, 04 Hz et une fréquence supérieure à sensiblement 0,6 Hz, et, - on fait effectuer audit aéronef, avant chaque parcours de mesure dans une zone où les variations dudit champ magné tique sont à mesurer, au moins un parcours d'identification dans une zone où ledit champ magnétique est uniforme, ledit parcours d'identification (Pa ; Pb ; Pc) comprenant une pluralité de virages (V1, V2, V3) séparés entre eux par des lignes droites (D), chaque virage (V1, V2, V3) étant relatif à un changement de cap d'une des valeurs 90 , 1800 ou 2700, effectué avec un angle (q) d'assiette latérale compris entre sensiblement 400 et sensiblement 600, avec une vitesse en cap de sensiblement 8 /s, et chaque ligne droite (D) étant de durée au moins égale à sensiblement 20 s. 3 - Procédé selon l'une des revendications 1 ou 2, dans lequel, les trois composantes échantillonnées étant HXK, H et H , on calcule seize variables h1@@, ..., h. ... et yK 2K 1K 1K h16 , i étant un entier naturel variant entre 1 et 16, selon @@K les formules : h1K = HxK /MK h9K = h1K H1K# /MK h2K = HyK /MK h10K = h1K H2K# /MK h3K = HzK /MK h11K = h1K H3K# /MK h4 K = 12E K 2E = h2K H1K /MK h5K = h22K h13K = h2K H2K# /MK h6K = H1K h2K h14K = h2K H3K# /MK h7K = H1K h3K h15K = h3K H1K# /MK h8K = H2K h3K h16K = h3K H2K# /MK avec MK = module du champ de composantes HxK , HyK, HzK HiK# = dérivée temporelle de HiK 4 - Procédé selon la revendication 3, dans lequel on effectue un filtrage passe-bas de ladite intensité échan tillonnée (SK) pour obtenir ledit module MK nécessaire au K K calcul desdites variable h1K, ..., hiK' ... et h16K . 5 - Procédé selon l'une des revendications 3 ou 4, dans lequel on calcule chaque dite dérivée H. par sommation 1K pondérée d'une dérivée aux différences finies et d'une dérivée homographique. 6 - Procédé selon l'une des revendications 1 à 5, dans lequel lesdits coefficients de pondération étant a1K, ..., aiK, ... et a16K, ladite différence étnat #K et lesdites variables filtrées étant Z1 ... Z. ... et x16K' on 1K 1K 18K calcule chaque coefficient de pondération a. selon les 1K formules : aiK = aiK-1 + kiK #K aj = o o formules dans lesquelles chaque quantité k. est la compo iK sante d'un vecteur colonne KK calculé selon la formule formule dans laquelle : - XTK désigne le vecteur ligne tansposé d'un vecteur colonne XK dont les composantes sont les variables filtrées x1, ..., xj , ... et x 16 K K K - est un coefficient de pondération des effets passés de l'aimantation, et, P PK-1 est calculé à l'instant précédent selon les formules Po = po I I étant la matrice unité, et Po étant compris entre sensiblement 104 et sensiblement 107. 7 - Procédé selon l'une des revendications 1 à 6, dans lequel chaque dite variable (hiK) et ladite intensité echanK tillonnée (SK) étant accompagnés de bruit, on les filtre numériquement pour, en outre, obtenir un bruit filtré de spectre fréquentiel uniforme. 8 - Dispositif pour la mise en oeuvre du procédé selon la revendication 1, de compensation des effets des perturbations magnétiques dues au porteur d'un magnétomètre, affectant la mesure de l'intensité d'un champ magnétique à mesurer par ledit magnétomètre, comprenant des moyens pour mesurer, à tout instant, les trois composantes, dans un repère lié audit porteur, dudit champ magnétique dont on veut connaître l'intensité avec précision, et des moyens, reliés auxdits moyens de mesure, pour estimer l'erreur due auxdites perturbations et compenser lesdits effets, dispositif caractérisé par le fait qu'il comprend - des moyens (2, 45, 5, 12) pour échantillonner dans le temps les trois composantes (Hx, Hy, Hz) mesurées et l'intensité (S) mesurée, et pour délivrer, à chaque instant (tK) de rang K, K étant un entier naturel, trois composantes échan tillonnées (H , Hy, H et une intensité échantillonnée (SK), @K @K @K - des premiers moyens de calcul (6), recevant les trois composantes échantillonnées (H , HyK, HzK), et agencés @K @K @K pour calculer une pluralité de variables (h1@@, ..., hj ..., h16K) formant une base dans l'expace des erreurs dues à l'aimantation du porteur (10), - des moyens (3', 3) de filtrage numérique desdites varia bles (h1K, ..., hiK, ..., h16K) et de ladite intensité échantillonnée (SK), , au moins pour et éliminer les composantes basses fréquences , et délivrer une pluralité de variables filtrées (x1K, ..., xiK, ..., x16K) et une intensité filtrée (yK), respectivement, - des deuxièmes moyens de calcul (7) dune somme pondérée (eK) desdites variables filtrées, dans laquelle chaque dite variable filtrée (si) est multipliée par nn coef K ficient de pondération (aiK-1) lié à l'aimantation du porteur (10), - des moyens (8) pour soustraire de ladite intensité filtrée (yK) ladite somme pondérée (eK), et délivrer une diffé rence (#K), et, - des troisièmes moyens de calcul (9) desdits coefficients de pondération (a1K, ..., aiK, ..., a16K) en fonction de ladite différence (EK) et desdites variables filtrées (x1K, ..., xiK, ..., x16KK) selon un algorithme récursif d'estimation, ladite somme pondérée (eK) représentant, lorsque K est grand, la valeur de ladite erreur, et ladite différence (K) représentant alors une mesure compensée de l'intensité dudit champ magnétique. 9 - Dispositif selon la revendication 8, dans lequel lesdits moyens (3', 3) de filtrage numérique éliminent les composantes desdites variables et de ladite intensité échantillonnée ayant une fréquence inférieure à sensiblement 0,04 Hz et une fréquence supérieure à sensiblement 0,6 Hz. 10 - Dispositif selon l'une des revendications 8 ou 9, dans lequel, les trois composantes échantillonnées étant Hx, H et H , les premiers moyens de calcul (6) sont xK yK zK agencés pour calculer seize variables h1@@, ..., h , ... et K 1K h16 , i étant un entier naturel variant entre 1 et 16, selon K les formules h1K = HxK /MK h9K = h1K H1K# /MK h2K = HyK /MK h10K = h1K H2K# /MK h3K = HzK /MK h11K = h1K H3K# /MK h4K = h1K2 h12K = h2K H1K# /MK h5K = h22K h13K = h2K H2K# /MK h6K = h1K h2K h14K = h2K H3K# /MK h7K = h1K h3K h15K = h3K H1K# /MK h8K = H2K h3K h16K = h3K H2K# /MK avec MK = module du champ de composantes HxK, HyK, HzK H. = dérivée temporelle de H. iK 1K 11 - Dispositif selon la revendication 10, dans lequel les premiers moyens de calcul (6) comprennent des moyens (67') pour effectuer un filtrage passe-bas de ladite intensité échantillonnée (SK) afin d'obtenir ledit module Mx nécessaire au calcul desdites variables h1K, ..., h. ;K... et K 1X h16 K 12 - Dispositif selon l'une des revendications 10 ou Il, dans lequel les premiers moyens de calcul (6) comprennent des moyens (61-66) pour calculer chaque dite dérivée H. 1K par sommation pondérée d'une dérivée aux différences finies et d'une dérivée homographique. 13 - Dispositif selon l'une des revendications 8 à 12, dans lequel lesdits coefficients de pondération étant a1K, ..., aiK, ... et a16K, ladite différence étant et lesdites variables filtrées étant xl 9 ... xj ... et K @K x16, lesdites troisièmes moyens de calcul (9) comprennent - des quatrièmes moyens de calcul (91, 92), agencés pour calculer des quantités kl , k. , ... et k16K K @K 10K composantes d'un vecteur colonne KK calculé selon la formule: formule dans laquelle - désigne le vecteur ligne transposé d'un vecteur colonne XK sont les composantes sont les variables filtrées K x1K, ..., xiK, ... et x16K, - y est un coefficient de pondération des effets passés de l'aimantation, et, - PK-1 est calculé selon les formules Po=poI I étant la matrice unité, et po étant compris entre sensi blement 10@ et sensiblement - des cinquièmes moyens de calcul (93) agencés pour calculer les produits k. lx et, - des sixièmes moyens de calcul (94), agencés pour calculer chaque coefficient de pondération aiK selon les formules a = a. + k. i K aj = 0 o 14 - Dispositif selon l'une des revendications 8 à 13, dans lequel chaque dite variable (hi ) et ladite intensité K échantillonnée (SK) étant accompagnés de bruit, lesdits moyens (3', 3) de filtrage numérique sont agencés pour que le bruit filtré soit de spectre fréquentiel uniforme.