Pose tracking by measuring a perturbed oscillating magnetic field

Abstract

The invention relates to a method for tracking the pose of a receiver (24) relative to a transmitter (22), the method comprising obtaining a raw magnetic measurement, estimating a contribution of a perturbing element (90) to the raw magnetic measurement, calculating a refined magnetic measurement by subtracting the contribution of the perturbing element (90) from the raw magnetic measurement, and inferring a refined pose of the receiver (24) relative to the transmitter (22) from the refined magnetic measurement. The contribution of the perturbing element (90) is estimated from secondary current densities of eddy currents induced by the transmitter (22) in the perturbing element (90), the secondary current densities resulting from an interpolation of primary-current-density components of the eddy currents along different directions and at different points on a surface (94) of the perturbing element (90).

Description

[0001] DESCRIPTION

[0002] TITLE: INSTALLATION TRACKING BY MEASURING A DISTURBED OSCILLATING MAGNETIC FIELD

[0003] FIELD OF INVENTION

[0004] The present invention relates to tracking the pose, that is, the position and orientation, of a first object relative to a second object, by measuring, at the level of the first object, an oscillating magnetic field produced by a transmitter attached to the second object. More particularly, the present invention relates to implementing this tracking when this oscillating magnetic field is disturbed by the presence, in the measurement environment, of a conductive element excited by the transmitter, such that it emits a disturbing magnetic field that distorts the lines of the oscillating magnetic field.

[0005] TECHNOLOGICAL BACKGROUND

[0006] To track the pose of an object—that is, its position and orientation—within a frame of reference, magnetic pose tracking systems, commonly called EMTs (Electromagnetic Trackers), are known to be used. These systems generally comprise a transmitter fixed in the frame of reference, consisting of three coils arranged to form a fixed, substantially orthogonal trihedron, and a receiver attached to the object whose pose is to be determined. This receiver consists of three magnetic sensors arranged to form a movable, substantially orthogonal trihedron. Time-dependent electrical currents flow through the transmitter coils, generating three magnetic fields that are detected by the receiver sensors. Each sensor on the receiver measures, for each emitted magnetic field, the projection of that field onto the direction in which the sensor is located.This results in a total of nine components that allow the transition from a moving to a fixed trihedron. These nine components depend on the position and orientation of the receiver relative to the transmitter. These nine components form what will be referred to hereafter as a "magnetic matrix," in which the components of the same field are grouped in the same column of the matrix, and the components measured by the same sensor are grouped in the same row of the matrix. To allow the identification, on the receiver side, of the contribution of each coil to the detected magnetic field, it is known to excite the transmitter coils using alternating currents at their natural frequencies. The generated magnetic fields are thus alternating magnetic fields that can be distinguished from one another by demodulation using the knowledge of the coil excitation frequency.Sensitivity to constant magnetic fields can therefore be neglected and the sensors can be formed from simple coils across whose terminals the induced voltage is measured.

[0007] One drawback of magnetic technology pose tracking systems is their sensitivity to the presence of conductive objects in the environment. Consider a conductive object subjected to a single-frequency incident magnetic field (generated by a single coil) that varies over time and is not uniform across the object's volume. Eddy currents at the frequency of the incident field flow through this object, generating a reflected field that is superimposed on the incident field. This distorts the lines of the total field at that frequency, disrupting the pose measurement by the magnetic technology. The magnetic matrix M t The actual pose (r, R) can therefore be decomposed into M t = Me + M cwhere Me is the matrix measured in the absence of a disturbance (which can easily be predicted at a given exposure using models of the source and sensor defects) and M c is the contribution of the conducting object to the magnetic matrix M t after demodulation (that is, after projecting the measurement of the total magnetic field onto a "demodulation direction").

[0008] To circumvent this difficulty, various solutions have been proposed. One type of solution relies on exploiting the variation in the response of disturbances to incident fields of different frequencies. For example, methods are known from US6172499 and US6762600 in which the contribution of a disturbance is identified and removed by taking advantage of the phase variation of the disturbance's response at different frequencies. A method is also known from US11187823 for obtaining a reference signal or pose by occasionally using very low frequencies (e.g., around 100 Hz) at which conductive objects have a negligible impact. This reference signal / pose is then used to identify a correction for the fields at the main frequencies (e.g., around 30 kHz), and this correction can be extrapolated to correct the fields at the main frequencies during subsequent poses.

[0009] Solutions of this first type have the advantage of being generic (they do not require dedicated work to study each individual disturbance whose impact must be compensated for). However, they present several drawbacks. First, frequency changes are restrictive. Second, they require the sensor to be almost static when low frequencies are used.

[0010] A second type of solution involves methods that perform corrections without a detailed physical model and employ a large number of correction parameters calibrated through the acquisition of a (long and complex) calibration trajectory that meshes the target space. These parameters are typically either polynomial coefficients (or spherical harmonic development moments) or lookup tables (well known by the acronym LUT, for "Look-up Table") that associate a correction with an approximate exposure or field measurement, or a combination of both (a table that provides coefficients for a local correction calculation). Such an approach is notably described in FR2664044, which outlines the correction of disturbances caused by an extended conducting object (an aircraft cockpit) fixed relative to the source. A polynomial correction of the field measurement is proposed therein.The correction parameters are constructed from a calibration trajectory acquired in the final operating environment. Several improvements have been proposed since, notably by US6400139, which suggests using "reference sensors" fixed in the source reference frame to make the correction model more adaptable by eliminating the need for impact mapping steps and robust to potential movements. A similar approach is also described by H. Himberg in "Latency and Distortion Compensation in Augmented Environments Using Electromagnetic Trackers," PhD thesis, Virginia Commonwealth University, 2010.

[0011] However, solutions of this second type are only suitable for correcting perturbations in a very small and relatively stable space. Furthermore, the calibration phases are cumbersome: they require significant acquisition times and most often necessitate a reference system capable of providing the exact sensor position at which the measurements are acquired by the EMT. Finally, such methods are not viable when the perturbation is significant and moves with the sensor in the source frame (indeed, whereas when the perturbation is fixed in the source frame, the perturbation depends only on the sensor's position in a subspace V of ]R 3 defining the desired range of the EMT, when the perturbant moves with the sensor, the disturbance changes at each pose, which increases the space to be covered: we must then map all SO(3)xV, SO(3) being the set of rotation matrices in dimension 3).

[0012] A third type of solution involves methods that attempt a detailed modeling of conductive objects to reconstruct their impact on the fields. When the perturbation is simple (with geometric properties that trivialize the problem), exact theoretical calculations of the perturbation's impact can be used once its position is known or estimated. This is what is proposed for infinite planes in US2023083524. In other cases, a detailed modeling of the perturbation using methods approaching finite elements can be performed, the main limitation of this type of approach being the computation time, which is incompatible with real-time correction. (See H. Gietler, Object localization for autonomous systems based on electromagnetic fields.)PhD thesis, University of Klagenfurt, 2021 presents a method in which computation times are presented as acceptable thanks to a discretization of the perturbator volume and a resolution by integral calculus with simplifying assumptions on the elementary volumes, this method is not functional in practice (indeed the calculation of the necessary field integrals is not straightforward and is not even explained in the document) and is only effective in the case of a fixed perturbator in source frame (its adaptation to a fixed perturbator in sensor frame risks encountering an execution speed incompatible with a real-time application).

[0013] Thus, solutions of this third type are either limited to correcting disturbances produced by perturbants with a very simple geometric shape, or incompatible with real-time application.

[0014] It should be noted that solutions of this third type can be combined with solutions of the second type, as described for example in FR2722299. Indeed, this document describes a process using a fine modeling of small perturbators (metallic parts of a pilot's helmet) fixed with respect to the sensor, combined with a calibration of multipolar moments of a spherical harmonic development centered on other perturbators.

[0015] Finally, a fourth type of solution combines EMT measurements with inertial measurements. Such solutions are known, in particular, from US20160356601, US10746819, and K. O'Donoghue, A. Jaeger, Herman, and C. Murphy, Padraig. Sensor fusion hardware platform for robust electromagnetic navigation. IEEE International Instrumentation and Measurement Technology Conference (I2MTC), 2022. However, while these solutions offer high robustness to the appearance of disturbances or noise, they do not improve the static accuracy of the system, which is defined by EMT alone.

[0016] Thus, there is a need for a solution to correct the disturbance produced by a conductive object of relatively complex shape present in the environment of an EMT and to improve the static accuracy of this EMT, when this conductive object is moving with the sensor or when the desired range of the EMT is extended, said solution being suitable for a real-time application and not requiring quasi-static phases of the sensor.

[0017] DESCRIPTION OF THE INVENTION

[0018] One objective of the invention is to correct the perturbation produced by a relatively complex-shaped disturbing element present in the environment of an EMT, when this disturbing element is moving with the sensor and / or when the desired high-precision area of ​​the EMT is extended. Other objectives are to enable this correction in real time, to improve the static accuracy of the EMT, and to eliminate the need for quasi-static phases of the EMT sensor.

[0019] To this end, the invention relates, according to a first aspect, to a method for tracking the positioning of a receiver relative to a transmitter, the transmitter comprising at least two magnetic generators, each oriented along its own direction and capable of generating an alternating magnetic field having a distinctive characteristic enabling the alternating magnetic field generated to be distinguished from the alternating magnetic field generated by each other magnetic generator, the magnetic generators being arranged relative to each other so that their directions are non-coplanar, the method being implemented by a data processing unit and comprising the following step:

[0020] determination of a pose of the receiver relative to the transmitter at a determination instant, said determination comprising the following sub-steps:

[0021] • a) obtaining a raw magnetic measurement composed of measurements, carried out at the time of determination, of local components, along at least two non-coplanar measurement directions attached to the receiver, of at least two resulting alternating magnetic fields constituting an ambient magnetic field, each of said resulting alternating magnetic fields having a distinctive characteristic identical to that of one of the generated alternating magnetic fields,

[0022] • b) determination of an imprecise position of the receiver relative to the emitter at the time of determination,

[0023] • c) estimation, from an approximate pose based on the imprecise pose, of the contribution of a disturbing element to the raw magnetic measurement, • d) calculation of a refined magnetic measurement by subtracting the contribution of the disturbing element from the raw magnetic measurement, and

[0024] • e) deduction of a refined pose of the receiver relative to the emitter at the time of determination from the refined magnetic measurement,

[0025] in which the disturbing element is modeled as a surface mesh composed of a set of mesh points regularly distributed over a surface of the disturbing element and connected to each other by straight segments delimiting between them two-dimensional and parallelogrammatic mesh cells, and substep c) of estimating the contribution of the disturbing element includes, for each of the generated alternating magnetic fields, the estimation of local components, at the receiver level, for a receiver pose equal to the approximate pose, of a disturbing magnetic field produced by eddy currents induced in the disturbing element by said generated alternating magnetic field, by applying the Biot-Savart law to secondary current densities of these eddy currents,said secondary current densities resulting from an interpolation of primary current density components of said eddy currents, said primary current density components being a function of a pose of the disturbing element with respect to the emitter and being constituted,:,

[0026] of first components, along a first direction of the surface mesh, of first primary current densities in first internal segments of the surface mesh oriented along a second direction of the surface mesh, in particular in the middle of said first internal segments, and

[0027] of second components, along the second direction of the surface mesh, of second primary current densities in second internal segments of the surface mesh oriented along the first direction of the surface mesh, in particular in the middle of said second internal segments.

[0028] According to particular embodiments of the invention, the tracking method also has one or more of the following characteristics, taken individually or in any technically feasible combination(s):

[0029] the measurement directions are greater than or equal to three;

[0030] The primary current density components are obtained using the following formula or an equivalent formula:

[0031]

[0032] where / is a vertical matrix of the set of first components of the first primary current densities, / is a vertical matrix of the set of second components of the second secondary current densities, b is a vertical matrix of normal components of the generated alternating magnetic field, orthogonal to the surface of the disturbing element, at each of the mesh points not belonging to a contour of the surface mesh, and Ky 1 is the inverse matrix of a matrix constructed from the local Maxwell equations discretized between the primary current density components and the normal components of the generated alternating magnetic field;

[0033] the disturbing element is fixed relative to the receiver and substep c) of estimating the contribution of the disturbing element includes determining a pose of the disturbing element relative to the emitter, estimating the components of the primary current densities as a function of said pose of the disturbing element relative to the emitter, and determining the secondary current densities by interpolating the components of the primary current densities;

[0034] the Kj matrix 1 is determined prior to the implementation of the installation determination step;

[0035] the disturbing element is fixed relative to the emitter and the secondary current densities are determined prior to the implementation of the pose determination step;

[0036] the pose determination step includes several successive iterations of substeps c), d) and e), the approximate pose being constituted, for the first iteration, by the imprecise pose and, for each subsequent iteration, by the refined pose deduced at the end of the previous iteration;

[0037] the tracking process includes several iterations of the pose determination step to determine receiver poses at different determination times, the first iteration constituting a step for determining an initial receiver pose and subsequent iterations constituting steps for determining an updated pose, the inaccurate pose determined during the initial pose determination step being deduced from the raw magnetic measurement; the inaccurate pose determined during at least one of the steps for determining an updated pose is extrapolated from at least one earlier pose determined during a previous iteration of the pose determination step;

[0038] secondary current densities are determined at evaluation points regularly distributed within each mesh cell;

[0039] The local components of the disturbing magnetic field are estimated using the following formula or an equivalent formula:

[0040] b (r) = ^5's e KS,) X (r ~ S '' )

[0041]

[0042] ' ' taÂ. ' ' |rS||s

[0043] where bj(r) is a vector of said local components of the disturbing magnetic field, p0 is the magnetic permeability of free space, / (Sj) is a 3-dimensional vector of the components, in a frame attached to the disturbing element, of one of the secondary current densities,

[0044]

[0045] is the surface area of ​​an elementary cell associated with said estimated surface current density, e t is an average thickness of a portion of the surface of the disturbing element, within which the eddy currents are concentrated, at the level of the mesh cell comprising the elementary cell,

[0046]

[0047] is a vector of the coordinates of a center of said elementary cell in a frame attached to the perturbing element, r is a vector of the coordinates of the approximate pose in said frame attached to the perturbing element, x is the cross product operator and |»| is the norm 2,

[0048] the elementary cell being constituted by one of the mesh cells or by a subdivision of a subdivided mesh cell, each subdivided mesh cell being made up of elementary cells, each associated with one of the secondary current densities, regularly distributed within said subdivided mesh cell and covering the whole thereof; and the mesh cells are rectangular.

[0049] The invention also relates, according to a second aspect, to a method for calibrating an electromagnetic tracking system comprising a step of tracking a pose of a receiver of the electromagnetic tracking system relative to a transmitter of the electromagnetic tracking system by means of a tracking method according to the first aspect, on the basis of hypothetical values ​​of calibration parameters of the electromagnetic tracking system, and updating the hypothetical values ​​of the calibration parameters so as to minimize a sum of error criteria each representative of an error on a pose of the receiver relative to the transmitter determined during the tracking step.

[0050] The invention also relates, according to a third aspect, to a computer program product comprising code instructions for the implementation, by a processor, of a tracking method according to the first aspect or of a calibration method according to the second aspect.

[0051] The invention also relates, according to a fourth aspect, to a computer-readable recording medium on which a computer program product is stored according to the third aspect.

[0052] Finally, according to a fifth aspect, the invention relates to an assembly comprising a reference object, a moving object, and an electromagnetic tracking system for tracking the position of the moving object relative to the reference object, said electromagnetic tracking system comprising:

[0053] an emitter attached to the reference object, said emitter comprising at least two magnetic generators, each oriented along its own direction and capable of generating an alternating magnetic field having a distinctive characteristic enabling the alternating magnetic field to be distinguished from the alternating magnetic field generated by each other magnetic generator, the magnetic generators being arranged relative to each other so that their directions are non-coplanar, and

[0054] a receiving device comprising a receiver attached to the moving object, said receiving device being capable of measuring local components, along at least two non-coplanar directions, of at least two resultant alternating magnetic fields constituting an ambient magnetic field, each of said resultant alternating magnetic fields having a distinctive characteristic identical to that of one of the generated alternating magnetic fields, so as to obtain a raw magnetic measurement,

[0055] the assembly also including a disturbing element capable of impacting the measurement of the local components of the resulting alternating magnetic fields by the receiving device when the disturbing element is subjected to the alternating magnetic fields generated by the transmitter,

[0056] in which the electromagnetic tracking system includes a data processing unit configured for implementing a tracking method according to the first aspect or a calibration method according to the second aspect. BRIEF DESCRIPTION OF THE FIGURES

[0057] Other features and advantages of the invention will become apparent from the following description, given solely by way of example and with reference to the accompanying drawings, in which:

[0058] Figure 1 is a schematic view of an assembly consisting of a first object and a second object comprising a tracking system according to an example embodiment of the invention,

[0059] Figure 2 is a schematic view of an example of a surface mesh of a portion of the surface of a disturbing element from the whole of Figure 1, Figure 3 is a diagram of an example of a determination process that can be implemented by the monitoring systems of Figure 1, Figure 4 is a diagram illustrating a first variant of a step for providing calculation parameters for the process of Figure 3,

[0060] Figure 5 is a diagram illustrating a second variant of the calculation parameter provision step of Figure 4,

[0061] Figure 6 is a diagram illustrating a step in determining an initial pose of the process shown in Figure 3.

[0062] Figure 7 is a diagram illustrating a first variant of a substep of the initial pose determination step of Figure 6 for estimating the contribution of a disturbing element to a raw magnetic measurement produced by a receiver of the tracking system,

[0063] Figure 8 is a diagram illustrating a second variant of the sub-step in Figure 7.

[0064] Figure 8 is a diagram illustrating a step in determining an updated pose of the process shown in Figure 3, and

[0065] Figure 10 is a diagram illustrating a calibration process that can be implemented by the tracking systems in Figure 1.

[0066] DETAILED DESCRIPTION OF A PROJECT EXAMPLE

[0067] The sets 10, 11 shown in Figures 1 and 2 each include a first object 12 and a second object 14, the second object 14 being mobile relative to the first object 12.

[0068] The first object 12 consists of a reference object whose pose in a reference frame, for example an inertial frame, is known. Here and subsequently, the pose is defined as consisting of all the position and orientation data of an object in space. In three-dimensional space, the pose therefore comprises six components, consisting of:

[0069] three position components, and

[0070] three orientation components.

[0071] The first object 12, for example, consists of a fixed structure in the reference frame. Alternatively (not shown), the first object is itself mobile in the reference frame but is equipped with sensors to track its movement within the reference frame.

[0072] The second object 14 consists of a mobile object in the reference frame and whose movement we wish to follow relative to this reference frame.

[0073] The first object 12 and the second object 14 are typically arranged relative to each other so that there is no permanent line of sight between the two objects 12, 14. In other words, the second object 14 is likely to be hidden from the view of the first object 12 during its movement.

[0074] The speed of movement of the second object 14 relative to the first object 12 is typically less than 30 m / s and 10 rad / s.

[0075] For example, the second object 14 is a pen and the first object 12 is a case that is stationary relative to the writing surface: this allows the writing to be recorded. Alternatively, the second object 14 is a helmet and the first object 12 is a structure fixed in the helmet's movement environment.

[0076] A primary frame (R^) is attached to the first object 12. This primary frame (R^) is a direct orthonormal frame with origin

[0077]

[0078] formed of a triplet of axes, represented in Figure 1, comprising:

[0079] a first primary axis e*,

[0080] a second primary axis e orthogonal to the first axis e*, and

[0081] a third primary axis e orthogonal to the first and second axes e and e yA secondary frame (R2) is attached to the second object 12. This secondary frame (R2) is a direct orthonormal frame with origin O2, formed by a triplet of axes, represented in Figure 1, comprising:

[0082] a first secondary axis e,

[0083] a second secondary axis e y orthogonal to the first axis e, and

[0084] a third secondary axis e z 2 orthogonal to the first and second axes

[0085]

[0086] summer y 2. The position of the second object 14 relative to the first object 12 can be characterized by the set formed by the rotation matrix R 1 2from the primary frame (R1) to the secondary frame (R2) and the vector going from the origin of the primary frame (R1) to the origin O2 of the secondary frame (R2). Here and subsequently, we will use the notation convention for rotation matrices from a frame A to a frame B as matrices whose columns are the unit vectors of frame B expressed in frame A. Thus, in this specific case, the rotation matrix R from the primary frame (R1) to the secondary frame (R2) is the matrix whose columns are the coordinates of the axes

[0087]

[0088] e 2, e 2 of the secondary frame (R2) expressed in the coordinate system of the primary frame (R^.

[0089] The assembly 10 also includes a tracking system 20 for monitoring the placement of the second object 14 relative to the first object 12. This tracking system 20 consists of an electromagnetic tracking system, or EMT. It comprises a transmitting device 21 which includes a transmitter 22 attached to the first object 12 and a receiving device 23 which includes a receiver 24 attached to the second object 14.

[0090] The transmitter 22 and receiver 24 are typically configured to be functional when separated by a distance of less than 10 m. This helps to limit the size of the transmitter 22 and the receiver 24.

[0091] The emitter 22 comprises at least two, here three, magnetic generators 30, 32, 34. Each magnetic generator 30, 32, 34 behaves essentially as a magnetic dipole having an alternating magnetic dipole moment of amplitude m and frequency f. It is oriented along a proper direction, respectively g1, g2, g3, defined here and thereafter as the direction of said dipole moment.

[0092] The amplitude m of the dipole moment is preferably substantially equal for all generators 30, 32, 34. The frequency f of the dipole moment, however, is specific to each generator 30, 32, 34, so that each generator 30, 32, 34 is thus capable of generating an alternating magnetic field with a distinctive characteristic that allows this alternating magnetic field to be distinguished from the alternating magnetic field generated by each other magnetic generator 30, 32, 34, this distinctive characteristic being the frequency of the magnetic field. Alternatively, the frequency f of the dipole moment is substantially equal for all generators 30, 32, 34, the generators 30, 32, 34 then exhibiting another distinctive characteristic that allows the production of magnetic fields that can be distinguished from one another. The frequency f of the dipole moment is typically between 500 Hz and 500 kHz, preferably between 5 kHz and 20 kHz.

[0093] For this purpose, each magnetic generator 30, 32, 34 comprises a coil, respectively 36, 38, 40, connected to a coil excitation current generator, respectively 42, 44, 46. Each coil 36, 38, 40 is oriented along a proper direction, defined by the axis around which the coil extends, this direction constituting the orientation direction of the magnetic generator 30, 32, 34 to which it belongs.

[0094] Each current generator 42, 44, 46 belongs to the transmitting device 21. Each current generator 42, 44, 46 is, for example, as shown, integrated into the first object 12. Alternatively (not shown), at least part of the current generators 42, 44, 46 is offset relative to the first object 12.

[0095] Each current generator 42, 44, 46 is configured to generate the excitation current of the corresponding coil 36, 38, 40 at a natural frequency, this frequency constituting the frequency of the alternating magnetic field generated by the magnetic generator 30, 32, 34 to which the coil 36, 38, 40 connected to this current generator 42, 44, 46 belongs.

[0096] Each current generator 42, 44, 46 is specifically configured to generate the excitation current of the corresponding coil 36, 38, 40 at a frequency between 500 Hz and 500 kHz, preferably between 5 kHz and 20 kHz.

[0097] The magnetic generators 30, 32, 34 are arranged relative to each other such that each is oriented along a proper direction g1, g2, g3, said directions g1, g2, g3 being non-coplanar. In particular, the magnetic generators 30, 32, 34 are arranged such that their directions g1, g2, g3 are such that:

[0098] a first direction of generator g! is coincident with the first primary axis e*,

[0099] a second direction of generator g2 either included in a plane defined by the first and second primary axes e*, e and forms a first angle 0 between 0 and^ radian (inclusive) with the second primary axis e, and

[0100] a third direction of generator g3 either included in a plane passing through the third primary axis e and forming a second angle <p compris entre 0 et —Q radian (exclus) avec le plan défini par les premier et troisième axes primaires e*, e, ladite troisième direction de générateur g3formant un troisième angle i compris entre 0 et^ radian (inclus) avec le troisième axe primaire e zFor example, the magnetic generators 30, 32, 34 are arranged so that their directions g1, g2, g3 are substantially orthogonal to each other, each direction g1, g2, g3 being typically, as shown, substantially collinear with one of the axes e1, e2, e3 of the primary frame (R^, that is to say that the first and second angles 0, <p sont inférieurs ou égaux à 0,05 radian.

[0101] Preferably, the magnetic generators 30, 32, 34 have substantially coincident centers, meaning they can be modeled as magnetic dipoles located at approximately the same point. To this end, the coils 36, 38, 40 have substantially coincident centers, meaning the centers of the coils are separated by a distance less than the average radius of the coils 36, 38, 40. The center of the magnetic generators 30, 32, 34 constitutes the origin

[0102]

[0103] of the primary frame (R^.

[0104] The transmitting device 21 also includes a clock 48 providing a reference signal representative of time. This reference signal is used by the magnetic generators 30, 32, 34 to vary over time the alternating magnetic fields they generate.

[0105] The receiving device 23 is capable of measuring the local components, at the level of the receiver 24, along the three axes e2, e y 2, e z 2 of the secondary frame (R2), of three resulting alternating magnetic fields constituting an ambient magnetic field, each of said resulting alternating magnetic fields having a distinctive characteristic, allowing said resulting alternating magnetic field to be distinguished from the other resulting alternating magnetic fields, identical to that of one of the generated alternating magnetic fields

[0106] For this purpose, the receiver 24 comprises three magnetic sensors 50, 52, 54 each capable of measuring a component of the ambient magnetic field along a measurement direction c1, c2, c3 specific to the sensor 50, 52, 54, these measurement directions c1, c2, c3 being non-coplanar. These sensors 50, 52, 54 are preferably arranged so that their measurement directions c1, c2, c3 are substantially orthogonal to each other, each measurement direction c1, c2, c3 being typically, as shown, substantially collinear with one of the axes e3, e3, e of the primary frame (R^.

[0107] Each magnetic sensor 50, 52, 54 is here formed of a coil, respectively 56, 58, 60. Each coil 56, 58, 60 is oriented along a proper direction, defined by the axis around which the coil extends, this direction constituting the measurement direction of the magnetic sensor 50, 52, 54 to which it belongs.

[0108] Preferably, the magnetic sensors 50, 52, 54 have substantially coincident centers. To this end, the coils 56, 58, 60 have substantially coincident centers, that is to say, the centers of the coils are separated in pairs by a distance less than the average radius of the coils 56, 58, 60. The center of the magnetic sensors 50, 52, 54 constitutes the origin O2 of the secondary frame (R2).

[0109] The receiving device 23 also includes tensiometers, 62, 64, 66 respectively, each connected to one of the coils, 56, 58, 60 respectively, to measure the voltage across said coil. Each tensiometer 62, 64, 66 is, for example, as shown, integrated into the first object 12. Alternatively (not shown), at least some of the tensiometers 62, 64, 66 are located remotely from the first object 12. Each tensiometer 62, 64, 66 is designed to produce a voltage signal representative of the voltage across the associated coil 56, 58, 60 and therefore representative of the component of the ambient magnetic field along one of the measurement directions c1, c2, c3.

[0110] For this purpose, each tensiometer 62, 64, 66 is typically made in the form of an amplifier coupled to an analog-to-digital converter (ADC).

[0111] The voltage signals produced by tensiometers 62, 64, 66 together form a measurement signal. Since the measurement directions c1, c2, c3 are not coplanar, this measurement signal is representative of the ambient magnetic field at the receiver 24.

[0112] The receiving device 23 further includes a processing unit 70 to process the measurement signal in order to deduce the local components, at the level of the receiver 24, along the three axes e, e,

[0113]

[0114] of the secondary frame (R2), resulting alternating magnetic fields. Such processing is known in itself and will not be detailed further. It typically involves demodulating the measurement signal with demodulation signals at frequencies equal to those of the resulting alternating magnetic fields.

[0115] The determination system 20 further includes a data processing unit 80. This data processing unit 80 is configured to deduce from the measurements made by the receiving device 23 the position of the receiver 24 relative to the transmitter 22. It is also configured to deduce from this position the position of the first object 12 relative to the second object 14.

[0116] For this purpose, the data processing unit 80 is typically a computer. It comprises a processor or CPU (Central Processing Unit) 82 and a memory 84 of the RAM (Random Access Memory) and / or ROM (Read Only Memory) type. The processor 82 is configured to execute instructions loaded into the memory 84. When the processing unit 80 is powered on, the processor 82 is able to read instructions from the memory 84 and execute them. These instructions form a computer program causing the processor 82 to implement processes 1000 (Figure 3) and 2000 (Figure 10), which will be detailed below.

[0117] The data processing unit 80 also includes a buffer memory 86 for the temporary storage of information necessary for the implementation of processes 1000 and 2000.

[0118] Alternatively, the data processing unit 80 is made up of a microcontroller.

[0119] In the example shown, the data processing unit 80 is integrated into the second object 14. Alternatively (not shown), at least part of the data processing unit 80 is remote, for example, in a mobile terminal (not shown) and / or a remote server (not shown). In other words, at least part of the steps of process 1000 and / or 2000 is performed by a mobile terminal and / or a remote server. The second object 14 then includes a communication system, typically a wireless communication system, configured to send the data from the receiving device 23 to the mobile terminal and / or the remote server.

[0120] The assembly 10 also includes a disturbing element 90 designed to impact the measurement of the local components of the resulting alternating magnetic fields by the receiving device 23 when the disturbing element 90 is subjected to the alternating magnetic fields generated by the transmitter 22.

[0121] This disturbing element 90 consists of an object made of current-conducting material. It is configured such that the eddy currents induced by the emitter 22 in the disturbing element 90 are concentrated in a portion of the surface 92 (Figure 2) of the disturbing element 90, delimiting a surface 94 of the disturbing element 90. In other words, any eddy currents induced by the emitter 22 in the disturbing element 90 that do not propagate in said portion of the surface 92 are negligible. This portion of the surface 92 preferably has a thickness at every point less than 5 mm, for example less than 2.5 mm, while being less than one-tenth, for example less than one-hundredth, of each other dimension of said portion of the surface 92.For this purpose, the disturbing element 90 is typically made up of a so-called thin plate, that is to say a plate whose thickness at any point is less than 5 mm, for example less than 2.5 mm, while being less than one-tenth, for example less than one-hundredth, of each other dimension of said plate; the surface portion 92 is then made up of the whole plate. Alternatively, the disturbing element 90 has a conductivity such that the skin thickness of the disturbing element 90 at any point at the frequencies of the alternating magnetic fields generated by the emitter 22 is less than 1.67 mm, for example less than 0.85 mm, while being less than one-third of one-tenth, for example one-third of one-hundredth, of each dimension of a surface of the disturbing element 90 oriented towards the emitter 22; the surface portion 92 is then constituted by a portion of the disturbing element 90 extending over three times the skin thickness.

[0122] The surface portion 92 can be of any shape. In particular, it can be solid or perforated, flat or curved, and have edges of any shape. To keep the explanation simple and understandable, we will henceforth consider a surface portion 92 that is flat, solid (i.e., without a hole), parallelepiped-shaped, and of uniform thickness and conductivity (i.e., identical at every point), the disturbing element 90 being typically constituted by said surface portion 92 (i.e., consisting of a thin, flat, solid, parallelepiped-shaped plate). This surface portion 92 has a length l x a width l y a thickness l z and a conductivity of o.

[0123] A tertiary frame (R3) is attached to the perturbing element 90. This tertiary frame (R3) is a direct orthonormal frame with origin O3 formed by a triplet of axes, represented in Figure 1, comprising:

[0124] a first tertiary axis e3,

[0125] a second tertiary axis e y 3orthogonal to the first axis e3, and

[0126] a third tertiary axis e3 orthogonal to the first and second axes e3 and e3. Here, the third tertiary axis e3 is chosen so as to be orthogonal to the surface 94 of the perturbing element 90 delimited by the portion of surface 92. The first tertiary axis e3 is preferably chosen so as to be parallel to a first edge 96 of said surface 94 and the second tertiary axis e3 is advantageously chosen so as to be parallel to a second edge 98 of said surface 94.

[0127] The position of the disturbing element 90 relative to the first object 12 can thus be characterized by the set formed by the rotation matrix R 1 3 from the primary frame (R^ to the tertiary frame (R3) and the vector

[0128]

[0129] (expressed in the primary frame (R^) going from the origin of the primary frame (R^) to the origin O3 of the tertiary frame (R3) and, conversely, the pose of the first object 12 relative to the perturbing element 90 can be characterized by the set formed by the rotation matrix R 3 1 from the tertiary frame (R3) to the primary frame (R^ (equal to the transpose of the matrix R) 1J3 ) and the vector

[0130]

[0131] (expressed in the tertiary frame (R3)) going from the origin O3 of the tertiary frame (R3) to the origin O1 of the primary frame (R1). Similarly, the pose of the perturbing element 90 relative to the second object 14 can be characterized by the set formed by the rotation matrix R 23 from the secondary frame (R2) to the tertiary frame (R3) and the vector

[0132]

[0133] (expressed in the secondary frame (R2)) going from the origin O2 of the secondary frame (R2) to the origin O3 of the tertiary frame (R3) and, conversely, the pose of the second object 14 relative to the perturbing element 90 can be characterized by the set formed by the rotation matrix R 32 from the tertiary frame (R3) to the secondary frame (R2) (equal to the transpose of the matrix R 23 ) and the vector

[0134]

[0135] (expressed in the tertiary frame (R3)) going from the origin O3 of the tertiary frame (R3) to the origin O2 of the secondary frame (R^.

[0136] In a first embodiment of the assembly 10, the disturbing element 90 is fixed to the first object 12. The position of the disturbing element 90 relative to the first object 12 (and conversely, the position of the first object 12 relative to the disturbing element 90) is then constant. In other words, the rotation matrices R 1 3 and R3 1 and the vectors

[0137]

[0138] and are constant. In a second embodiment of the assembly 10, the disturbing element 90 is fixed to the second object 14. The pose of the disturbing element 90 relative to the second object 14 (and conversely, the pose of the second object 14 relative to the disturbing element 90) is then constant. In other words, the rotation matrices R 2,3 and R 3,2 and the vectors

[0139]

[0140] and are constant.

[0141] With reference to Figure 2, the perturbing element 90 is modeled as a surface mesh 100. This surface mesh 100 is composed of a set of mesh points 102 regularly distributed over the surface 94 of the perturbing element 90 delimited by the portion of surface 92. These mesh points 102 are connected to each other by straight segments 104 delimiting between them two-dimensional mesh cells 106.

[0142] The straight segments 104 include external segments 107, each bordering one of the mesh cells 106 and thus defining a boundary of the mesh 100. The straight segments 104 also include internal segments 108, each bordering two of the mesh cells 106 and thus forming the boundary between these cells 106. The external segments 107 and internal segments 108 together constitute the entirety of the straight segments 104; in other words, there is no straight segment 104 that is not an external segment 107 or an internal segment 108.

[0143] For the surface mesh 100, we identify a first direction d1 and a second direction d2, different from the first direction d. The internal segments 108 comprise first internal segments 109 oriented along this second direction d2 and second internal segments 110 oriented along this first direction d1. The first internal segments 109 are thus all parallel to each other, and the second internal segments 110 are all parallel to each other. The first and second internal segments 109 and 110 together constitute the entirety of the internal segments 108; in other words, there is no internal segment 108 that is not a first internal segment 109 or a second internal segment 110.

[0144] Advantageously, the first direction dï is included in the plane defined by the first and third tertiary axes e3, ef. Here, the first direction dï is in particular substantially parallel to the first tertiary axis e3.

[0145] Preferably, the second direction d2 is contained within the plane defined by the second and third tertiary axes e y 3, e3. Here, the second direction d2 is in particular substantially parallel to the second tertiary axis e3.

[0146] The mesh cells 106 are parallelogrammatic, specifically rectangular, here square. Each is bounded along the first direction d1 by a pair of second segments 110 of the surface mesh 100, and along the second direction d2 by a pair of first segments 109 of the surface mesh 100. Each mesh cell 106 has a first dimension along the first direction d1 and a second dimension along the second direction d2. Here, since the mesh cells 106 are square, these first and second dimensions are equal. Furthermore, the mesh cells 106 are identical in size, so the first and second dimensions of each cell 106 are equal to the first and second dimensions of every other cell 106. We will denote h the value of these first and second dimensions in what follows.

[0147] The 106 mesh cells together cover substantially the entirety of a relevant region of the surface 94, that is to say that at least 90% of the relevant region of the surface 94 is covered by the 106 mesh cells. By "relevant region", we mean here and thereafter a region which affects the shape of the global eddy currents in the portion of surface 92.

[0148] To keep the explanation simple and understandable, the 106 mesh cells in the detailed example here are nine in number. However, in practice, the number of 106 mesh cells will be much higher and will typically range from 100 to 10,000.

[0149] This model of the disruptive element 90 is typically stored in memory 84 of the data processing unit 80.

[0150] A method 1000 for tracking the placement of the receiver 24 relative to the transmitter 22, implemented by the data processing unit 80, will now be described, with reference to Figure 3.

[0151] This process 1000 includes a first step 1100 of providing calculation parameters, followed by a step 1200 of determining an initial pose at an initial determination time, itself followed by a repeated step 1300 of determining an updated pose at a later determination time.

[0152] Step 1100, which involves providing calculation parameters, has two variants, depending on whether one considers the first or second embodiment described above (disruptive element 90 attached to the first or second object 12, 14). The first of these variants, corresponding to the first embodiment (disruptive element 90 attached to the first object 12), will now be described, with reference to Figure 4.

[0153] In this first variant, step 1100 begins with a substep 1110 of providing parameters for the perturbing element 90. These parameters of the perturbing element 90 typically include the surface mesh 100, the dimensions l x , L y and z of the surface portion 92, and the conductivity 0. They also include the position of the disturbing element 90 relative to the object, here the first object 12, to which it is attached. These parameters are typically stored in the memory 84 of the data processing unit 80.

[0154] Substep 1110 is followed by substep 1120 of determining normal components of each of the alternating magnetic fields generated by the emitter 22, orthogonal to the surface 94 of the disturbing element 90 delimited by the portion of surface 92, at several of the mesh points 102, in particular at each of the mesh points 102 not belonging to the contour of the mesh 100. These normal components are typically determined from the pose of each mesh point 102 relative to the emitter 22 (known because the pose of the disturbing element 90 relative to the first object 12 is predetermined), from a theoretical model of the alternating magnetic fields generated by the emitter 22, and from an algorithm for correcting this theoretical model aimed at correcting, for example, the defects of the emitter 22, by techniques known to the person skilled in the art.

[0155] This substep 1120 is itself followed by a substep 1130 of estimation of primary current density components in the surface portion 92 of the disturbing element 90 at the level of the internal segments 108 of the surface mesh 100, the said primary current densities being eddy current densities induced in the disturbing element 90 by the magnetic fields generated by the emitter 22.

[0156] The primary current density components thus estimated consist of first components, along the first direction dï of the 100 mesh, of first primary current densities at first primary evaluation points o x ,..., o x b(Figure 2), and second components, along the second direction d2 of the mesh 100, of second primary current densities n equal to the number of second internal segments 110) at second primary evaluation points oo (Figure 2). They are estimated for each of the alternating magnetic fields generated by the emitter 22.

[0157] The first primary assessment points o x ,..., o x b belong in particular to the first internal segments 109 and are typically constituted by the midpoints of these first internal segments 109. The second primary evaluation points o,..., o belong in particular to the second internal segments 110 and are typically constituted by the midpoints of these second internal segments 110.

[0158] Primary current densities are typically volume current densities.

[0159] These primary current density components are typically obtained using the following formula or an equivalent formula:

[0160]

[0161] where / is a vertical matrix of the set of first components of the first primary current densities j*,..., j x m For an alternating magnetic field j generated by emitter 22, / is a vertical matrix of the set of second components of the second primary current densities j,..., j for the alternating magnetic field j, b is a vertical matrix of the normal components of the alternating magnetic field j determined during substep 1120, and Kj 1 is the inverse matrix of a matrix K y constructed from local Maxwell equations discretized between the estimated primary current density components and the normal components of the alternating magnetic field j.

[0162] In the detailed (and simplified) example given here, the matrix K y is given by the following formula:

[0163] / I ji o a)j 1 \

[0164] i Kz _ I — C x + iz N X M X — C v + iz N V M V ]

[0165]

[0166] J >- — — r I < J

[0167] h X 4n

[0168] n x x oy fai yyi

[0169] iM Jj x \ ^xn / /

[0170] Or:

[0171] i is imaginary such that i 2 = — 1;

[0172] (JÔJ is the angular frequency (i.e., angular frequency) of the alternating magnetic field j

[0173] D x and D v are sparse matrices, each containing the factors +1 and -1 such that the following equality is satisfied: D x y

[0174]

[0175] + = 0 so as to reflect the result of applying the divergence operator to the two terms of the Maxwell-Ampère equation in quasi-stationary regime on each of the mesh cells 106;

[0176] M x is a square matrix whose:

[0177] o the diagonal terms (M x ) kk are equal to f vk ^d 3 r', where y x is a prismatic volume of the surface portion 92 which extends over the entire thickness of the surface portion 92 and whose base is a parallelogram, centered on one of the first primary evaluation points x ..., o, whose edges are parallel to the first and second directions of 1; d2et which extends, following the first direction d 1; from the center of one of the mesh cells 106 framing said evaluation point o kup to the center of the other mesh cell 106 framing said evaluation point and, following the

[0178]

[0179] second direction d2, from one of the segments 104 framing said evaluation point o k up to the other segment 104 framing said evaluation point o k , and o the other terms (M x ) ife (with k ≠ l) are equal to ||7 x ||i 1 -i where

[0180]

[0181] PfcOfl 1 1 is the norm 2 of a vector going from one of the first primary evaluation points of,..., of to another of the first primary evaluation points of,..., of and ||7 X || is the volume of the prismatic volume Vf centered on the evaluation point of;

[0182] M y is a square matrix whose:

[0183] o the diagonal terms (M y ) fefe are equal to f vk ^d 3r', where Vf is a prismatic volume of the surface portion 92 that extends over the entire thickness of the surface portion 92 and whose base is a parallelogram, centered on one of the second primary evaluation points of,..., of, whose edges are parallel to the first and second directions d 1; d2et which extends, following the first direction d 1; from one of the segments 104 framing said evaluation point of to the other segment 104 framing said evaluation point of, and, along the second direction d2, from the center of one of the mesh cells 106 framing said evaluation point of to the center of the other mesh cell 106 framing said evaluation point of, and

[0184] o the other terms (M y ) ife (with k ≠ l) are equal to ||yy||i where \ofof I

[0185]

[0186] is the norm 2 of a vector going from one of the second primary evaluation points of,..., of to another of the second primary evaluation points of,..., of and ||7y|| is the volume of the prismatic volume Vf centered on the evaluation point of;

[0187] C x , Cy, N x and Ny are sparse matrices, each containing the factors +1 and -1, such that the following equality is satisfied:

[0188]

[0189] -^(c x jf + C y jf ) = — + ^

[0190]

[0191] (N X M X / J + NyMy / J}) so as to reflect the result of applying the Maxwell-Faraday law in quasi-stationary regime at each of the mesh points 102, the term

[0192]

[0193] i(c x / J + C yjf) representing the (discreted) curl of the current density at each of the mesh points 102 and the term

[0194]

[0195]

[0196] (

[0197]

[0198] N X M X / J + NyMy / j) representing the self-induced magnetic field by the disturbing element at each of the mesh points 102, resulting from the application of the rotational operator to the vector potentials of this self-induced magnetic field obtained by application of the Biot and Savart law at each of the evaluation points of a primary current density.

[0199] In particular, by taking the ordering of the first primary evaluation points o*,..., o x b and the second primary assessment points y ,..., o y b shown in Figure 2, the D matrices x , D y , Cx , C y , N x and N y can be written:

[0200] Z 01

[0201] 0

[0202] -1

[0203] — o

[0204] 0

[0205] \ ° o

[0206] 0 0 0 0\

[0207] 0 0 oo\

[0208] 0 0 0 0

[0209] 1 0 0 0

[0210] -1 1 0 0

[0211] 0 -1 0 0

[0212] 0 0 1 0 /

[0213] 0 0 -1 1 /

[0214] / I -1 0 0 0 0 \

[0215] r > 0 1 -1 0 0 0

[0216] - I 0 0 0 1 -1 0 /

[0217] \0 0 0 0 1 -1 /

[0218] / -I 0 1 0 0 0\

[0219] r _ 0 -1 0 1 0 0 1

[0220] 4 - 0 0 -1 0 1 0 /

[0221] \ 0 0 0 -1 0 1 /

[0222] / I -1 0 0 0 0 \

[0223] - T > / 0 1 -1 0 0 0

[0224] - I 0 0 0 1 -1 0 /

[0225] \0 0 0 0 1 -1 /

[0226] / —I 0 1 0 0 0\

[0227] - T > I 0 -1 0 1 0 0 1

[0228] 4 “ I 0 0 -1 0 1 0 /

[0229]

[0230] \ 0 0 0 -1 0 1 /

[0231] Note that the number of rows of the matrices D x and D y is not equal to the number of mesh cells 106 (8 rows instead of 9) because applying the divergence operator to the two terms of the Maxwell-Ampère equation in quasi-stationary regime on each of the mesh cells 106 leads to redundant equations.

[0232] For further information on how to determine the K matrix and the underlying theory, the reader may refer to the articles J. Nagel, "Fast finite-difference calculation of eddy currents in thin metal sheets," *Applied Computational Electromagnetics Society Journal*, 33:575–584, June 2018, and J. Nagel, "Finite-difference simulation of eddy currents in nonmagnetic sheets via electric vector potential," *IEEE Transactions on Magnetics*, 55:1–8, December 2019, the contents of which are incorporated here by reference.

[0233] Note that an example of a formula equivalent to the formula

[0234] = K j -1

[0235]

[0236] The formula given above is as follows:

[0237]

[0238] where Tj is the vector of a potential such that:

[0239] the first components of the first primary current densities constitute the variation along the second direction d2 of the mesh 100 at the mesh points 102 where the normal components of the alternating magnetic fields generated by the emitter 22 are determined, and

[0240] the opposites of the second components of the second primary current densities constitute the variation along the first direction dï at the said mesh points 102.

[0241] Substep 1130 is followed by substep 1140 for determining secondary current densities in the surface portion 92 at secondary evaluation points s1, ..., s 12 (Figure 2), the said secondary current densities remaining densities of the eddy currents induced in the disturbing element 90 by the magnetic fields generated by the emitter 22.

[0242] These secondary assessment points s1, ..., s 12are regularly distributed within each 106 mesh cell. In other words, each 106 mesh cell includes at least one secondary evaluation point s1, ..., s 12 and, for each mesh cell 106, the secondary evaluation points s1, ..., s 12 included by said mesh cell 106 are regularly distributed within it. The distribution of said secondary evaluation points s1, ..., s 12 can nevertheless vary from one mesh cell 106 to another, that is to say that the number of secondary evaluation points s1, ..., s 12 within a mesh cell 106 may be different from the number of secondary evaluation points s1, ..., s 12 inside another 106 mesh cell. Thus, there can be 106 mesh cells with a single secondary evaluation point s1, ..., s 12, consisting of the center of mesh cell 106, and other mesh cells 106 containing several secondary evaluation points s1, ..., s 12 A particular example of this is given in Figure 2, where all the 106 mesh cells include a single secondary evaluation point s 1; s 12 , with the exception of the central mesh cell 106, which includes four secondary evaluation points s5, s6, s7, s8.

[0243] Preferably, a density of secondary evaluation points s1, ..., s will be chosen. 12 in the mesh cell(s) 106 closest to receiver 24 and / or in the mesh cell(s) 106 where the current densities are highest.

[0244] At each secondary evaluation point s1, ..., s 12 is associated an elementary cell 120 centered on said secondary evaluation point s1, ..., s 12Each mesh cell 106 consists of the elementary cell(s) associated with the secondary evaluation point(s) s1, ..., s 12 included by said mesh cell 106. Thus, for each mesh cell 106 comprising a single secondary evaluation point s1, ..., s 12 , the elementary cell 120 associated with this secondary evaluation point s1, ..., s 12 is constituted by the mesh cell 106 and, for each mesh cell 106 comprising several secondary evaluation points s1, ..., s 12 , the 120 elementary cells associated with these secondary evaluation points s1, ..., s 12 are each made up of a subdivision of the mesh cell 106 and cover the entirety of it.

[0245] These secondary current densities are volume current densities represented by a 3-dimensional vector. They are obtained by interpolating the components of the primary current densities. The interpolation is, for example, a linear interpolation. We then have:

[0246] j s = pr3jf + ^4jf +i

[0247]

[0248] \ o

[0249] Or:

[0250] j s is a secondary current density, measured at a secondary evaluation point s1, ..., s 12 of a 106 mesh cell, expressed in the tertiary frame (R3),

[0251] Jk et Jfc+i are ^ es first components of first primary current densities in segments 104 of the mesh 100 bordering the mesh cell 106 along the first direction of 1;these first components being constituted by the first components determined during substep 1130 for the first primary current densities evaluated in internal segments 108, and being zero for the first primary current densities evaluated in external segments 107,

[0252] jf and jf +1 are the second components of second primary current densities in segments 104 of the mesh 100 bordering the mesh cell 106 along the second direction d2, these second components being constituted by the second components determined during substep 1130 for the second primary current densities evaluated in internal segments 108, and being zero for the second primary current densities evaluated in external segments 107, and

[0253] K 1; K2, K3 and K4 are interpolation coefficients.

[0254] For example:

[0255] fKjx k + Çl -K)j x k+1

[0256] is = \ «'iï + (i -or +i

[0257]

[0258] \ o

[0259] with K and K' being interpolation coefficients, these interpolation coefficients K, K' being preferably barycentric weights such as:

[0260] K = | where |nifcSi| is the norm of the vector going from point m k

[0261]

[0262] evaluation of one j k of the first components of first primary current densities in segments 104 of the mesh 100 bordering the mesh cell 106 along the first direction dï at the secondary evaluation point Si and |

[0263]

[0264] m£ +1 Si| is the norm of the vector going from point m k+1 evaluation of the other j k+1of the first components of the first primary current densities in segments 104 of the mesh 100 bordering the mesh cell 106 along the first direction di at the secondary evaluation point s b and K = \ H, where m sj is the norm of the vector going from the point

[0265]

[0266]

[0267] evaluation of one of the second components of the second primary current densities in segments 104 of the mesh 100 bordering the mesh cell 106 along the second direction d2 at the secondary evaluation point Sj and |

[0268]

[0269] mf +1 Si| is the norm of the vector going from point m +1 evaluation of the other j +1of the second components of second primary current densities in segments 104 of the mesh 100 bordering the mesh cell 106 along the second direction d2 at the secondary evaluation point s,.

[0270] The secondary current densities are thus determined prior to the implementation of steps 1200, 1300 of determining an initial pose and of determining an updated pose.

[0271] A second variant of step 1100, corresponding to the second embodiment (disruptive element 90 attached to the second object 14), will now be described, with reference to Figure 5.

[0272] In this second variant, step 1100 always begins with step 1110, which provides parameters for the 'perturbing element 90'. As in the first variant, these parameters of the perturbing element 90 typically include the surface mesh 100, the dimensions l x , L y and zof the surface portion 92, and the conductivity 0, and the position of the disturbing element 90 relative to the object to which it is attached. However, unlike the first variant, this position is not relative to the first object 12, since the disturbing element 90 is not attached to the latter, but relative to the second object 14.

[0273] This second variant differs further from the first variant in that step 1100 does not include substeps 1120, 1130 and 1140. Indeed, since the position of the disturbing element 90 relative to the first object 12 is not known, it is not possible to determine the normal components of each of the alternating magnetic fields generated by the emitter 22 at several of the mesh points 102. It is a fortiori not possible to determine the current densities of the eddy currents induced by these magnetic fields.

[0274] In this second variant, step 1100 nevertheless includes a substep 1150 for determining each matrix Kj 1 (or of each matrix K'y) 1 ) for calculating the primary current density components from the normal components of each of the alternating magnetic fields generated by the emitter 22 at each of the mesh points 102 not belonging to the contour of the mesh 100, and a substep 1160 for determining an interpolation matrix grouping the interpolation coefficients for calculating the secondary current densities from the primary current density components. It is indeed clear that the interpolation equations presented in substep 1140 can be written in matrix form so as to give all the secondary current densities from all the primary current density components in a single operation. It is also clear that:

[0275] the K matrices y depend only on the frequencies of the magnetic fields generated by the emitter 22 and on the parameters of the disturbing element 90 (in particular its surface mesh 100) and can therefore be calculated (as well as their inverses Kj 1 ) without knowing the position of the interfering element 90 relative to the emitter 22,

[0276] The interpolation coefficients depend exclusively on parameters of the surface mesh 100 (in particular the distances between the primary and secondary evaluation points) and can therefore be calculated (as well as the matrix grouping them) without knowing the pose of the disturbing element 90 relative to the emitter 22.

[0277] It also includes a substep 1170 for determining a calculation matrix for disturbing magnetic fields at the receiver 24. This matrix is ​​configured to allow the calculation, by a simple matrix operation, of local components, at the receiver 24, of disturbing magnetic fields produced by eddy currents induced in the disturbing element 90 by the magnetic fields generated by the emitter 22, when the secondary densities of these eddy currents are known. Indeed, it is known that eddy currents propagating in a current-conducting material such as that of the disturbing element 90 produce magnetic fields at the same frequency as said currents.In application of the Biot and Savart law, the local components of such magnetic fields at a measurement point are a function of these secondary current densities, of the surfaces of the elementary cells 120 associated with the secondary evaluation points s1, ..., s. 12 of these secondary current densities, the average thicknesses of the surface portion 92 of the disturbing element 90 at the level of these elementary cells 120 and the vectors linking said measurement point to the secondary evaluation points s1, ..., s 12. The position of the disturbing element 90 relative to the second object 14 (and therefore relative to the receiver 24) being known, as well as the dimensions of the disturbing element 90 and the parameters of the surface mesh 100, it is therefore possible to establish a matrix function linking the local components of such a magnetic field at the receiver 24 to the secondary densities of eddy currents induced in the disturbing element 90.

[0278] Thus, in substep 1170, a calculation matrix Q is determined such that:

[0279] B(r) = £IJ S

[0280] Or:

[0281] J s is the following matrix:

[0282] /

[0283]

[0284] y 1 ^) y 2 (*i) y 3 (*i)\

[0285] I y 1 ^) y 2 (s f ) y 3 (s f ) I

[0286] V^SJV) y2 (sjv) j 3 s

[0287] where each j (Si) is a 2-dimensional vector of components, along the first and second directions d 1; d2, of one of the secondary current densities, evaluated at a secondary evaluation point s b eddy currents induced in the disturbing element 90 by an alternating magnetic field y, and A / is the total number of secondary evaluation points s b and B(r) is a matrix of the local components, at the receiver 24, of the disturbing magnetic fields produced by these eddy currents. The matrix Q is typically given by the following formula:

[0288] Po / $i e i $N e N

[0289]

[0290] / 0 (r — Sj)3

[0291] with Sit = I — (r — S;)30

[0292]

[0293] V (r-Si)2- (r-sOi

[0294] Or:

[0295] |i0 is the magnetic permeability of free space,

[0296] If is a vector of the coordinates of a secondary evaluation point s, in the tertiary frame (R3),

[0297] If is the surface area of ​​the elementary cell 120 associated with the secondary evaluation point s b

[0298] is an average thickness of the surface portion 92 of the perturbing element 90 at the level of the mesh cell 106 comprising the elementary cell 120,

[0299] A / is the total number of secondary assessment points b

[0300] r is a vector of the coordinates of the pose of receptor 24 in the tertiary frame (R3),

[0301] |»| is standard 2, and

[0302] (•) k is the k component of the vector shown in parentheses.

[0303] These substeps 1150, 1160 and 1170 can be implemented in parallel or sequentially.

[0304] Step 1200 of determining an initial pose of the second object 14 relative to the first object 12 will now be described, with reference to Figure 6.

[0305] This step 1200 consists of a first iteration of a step 1201 of determining a precise pose of the second object 12.

[0306] This step of determining a precise pose 1201 begins with a substep 1210 of obtaining a raw magnetic measurement M b produced by the receiving device 23 at the initial determination time. This raw magnetic measurement M bis composed of measurements, performed by the receiving device 23 at the initial determination time, of the local components, at the receiver 24, along the three axes e2, e2, e2 of the secondary frame (R2), of the three resulting alternating magnetic fields constituting the ambient magnetic field. During this sub-step of acquisition 1210, said raw magnetic measurement M b is typically transmitted by the receiving device 23 and received by the data processing unit 80.

[0307] Substep 1210 is followed by substep 1220 for determining an imprecise pose p0 of the receiver 24 relative to the transmitter 22 at the initial determination time. During this substep 1220, said imprecise pose p0 is advantageously deduced, for this first iteration, from the raw magnetic measurement M b, typically by means of a theoretical model of the alternating magnetic fields generated by the emitter 22 and an algorithm for correcting this theoretical model aimed at correcting, for example, the defects of the emitter 22 and those of the receiver 24, by techniques known to the person skilled in the art.

[0308] Substep 1220 is followed by several substeps 1230, 1240, 1250, 1260, 1270 which will be repeated.

[0309] The first of these substeps 1230, 1240, 1250, 1260, 1270 is a substep 1230 of determining an approximate pose p j k , where k=1,...,m is an index incremented at each iteration of substep 1230. This approximate pose p j;k is a function of the imprecise pose p0. During the first iteration of substep 1230, this approximate pose, denoted pi, is constituted by the imprecise pose p0.

[0310] This substep 1230 is followed by a substep 1240 estimating the contribution of the disturbing element 90 to the raw magnetic measurement M b Indeed, each disturbing magnetic field produced by the eddy currents induced in the disturbing element 90 by an incident magnetic field such as one of those generated by the emitter 22 exhibits the same intrinsic characteristics, and in particular the same differentiating characteristics, as that incident magnetic field. Consequently, each of the resulting alternating magnetic fields measured at the receiver 24 is formed not only from one of the alternating magnetic fields generated by the emitter 22, but also from the corresponding disturbing magnetic field re-emitted by the disturbing element 90. Therefore, the raw magnetic measurement M b can be written M b = M e + M p, where Me is the contribution of emitter 22, i.e., the alternating magnetic fields generated by emitter 22, to the raw magnetic measurement M b , and M p is the contribution of the disturbing element 90, i.e., the alternating disturbing magnetic fields re-emitted by the disturbing element 90, to the raw magnetic measurement M b .

[0311] Substep 1240 has two variants, depending on whether one considers the first or second embodiment described above (disruptive element 90 attached to the first or second object 12, 14). The first of these variants, corresponding to the first embodiment (disruptive element 90 attached to the first object 12), will now be described, with reference to Figure 7.

[0312] In this first embodiment, substep 1240 begins with substep 1241 for estimating local components, at the receiver 24, of the disturbing magnetic fields produced by the disturbing element 90. These local components are estimated from the approximate pose p j k and secondary current densities induced in the disturbing element 90 by the magnetic fields generated by the emitter 22. The approximate pose p j k having been determined during the previous substep 1230 and the secondary current densities induced in the disturbing element 90 by the magnetic fields generated by the emitter 22 having been determined prior to the implementation of the initial pose determination step 1200, during the calculation parameter provision step 1100, substep 1241 can indeed be implemented without any further prior step.

[0313] During this substep 1241, the data processing unit 80 estimates, for each of the alternating magnetic fields generated by the transmitter 22, the local components, at the receiver 24, of the disturbing magnetic field produced by the eddy currents induced in the disturbing element 90 by said generated alternating magnetic field, taking the approximate pose p j;k as the position of receiver 24 in the primary frame (R^). For this purpose, the data processing unit 80 applies the Biot-Savart law to the secondary current densities of these eddy currents determined during sub-step 1140 of step 1100, which provides the calculation parameters. These local components of the disturbing magnetic field are thus typically estimated using the following formula or an equivalent formula:

[0314]

[0315] , I J taZ. ' ' |r— Sj| 3

[0316] Or:

[0317] bj(r) is a vector of said local components of the disturbing magnetic field, here expressed in the tertiary frame (R3),

[0318] Po is the magnetic permeability of free space.

[0319] each / (Sj) is a 3-dimensional vector of the components, in the tertiary frame (R3), of one of the secondary current densities, evaluated at a secondary evaluation point s b

[0320] If is the surface area of ​​the elementary cell 120 associated with the secondary evaluation point s b

[0321] is an average thickness of the surface portion 92 of the perturbing element 90 at the level of the mesh cell 106 comprising the elementary cell 120,

[0322] If is a vector of the coordinates of the secondary evaluation point s, in the tertiary frame (R3),

[0323] r is a vector of the coordinates of the approximate pose in the tertiary frame (R3) (such a vector can be easily determined by means of the knowledge of the approximate pose and the pose of the perturbing element 90 in the primary frame (Ri)),

[0324] x is the vector product operator, and |»| is the 2 norm.

[0325] Substep 1241 is followed by substep 1242 for deducing the contribution of the disturbing element 90 to the raw magnetic measurement M b based on the local components of the disturbing magnetic fields thus estimated.

[0326] During this substep 1242, the data processing unit 80 projects these local components onto the measurement directions c1, c2, c3 of the receiver 24, taking into account the phase shifts between these local components and the demodulation signals used by the receiving device 23 to demodulate the measurement signal produced by the receiver 24.

[0327] The projection of the local components onto the measurement directions c1, c2, c3 of receiver 24 is typically obtained by applying to said local components a rotation matrix defining the orientation of the disturbing element 90 relative to the approximate pose, constituted by the matrix R 23 (assuming the secondary frame (R2) is positioned within the primary frame (R^ equal to the approximate position). This matrix R 2J3 can be easily obtained by the product

[0328]

[0329] where is the transpose of the rotation matrix defining the approximate pose in the primary frame (R^.

[0330] The phase shifts between these local components and the demodulation signals are typically accounted for by temporal projection, onto the demodulation direction, of the spatial projection of said local components onto the measurement directions c1, c2, c3. In other words, for each spatial projection of the local components onto one of the measurement directions c1, c2, c3 of receiver 24, the term is calculated

[0331]

[0332] cos ( <p{), où:

[0333] B is the amplitude of the projection of the local components of a disturbing magnetic field B J on a measurement direction c b And

[0334] <p est le déphasage, relativement aux signaux de démodulation utilisés par le dispositif de réception 23 pour démoduler le signal de mesure produit par le récepteur 24, de la projection des composantes locales du champ magnétique perturbateur B J on the measurement direction c,.

[0335] It should be noted that the amplitude B can be intrinsically deduced from the projection of the local components of a disturbing magnetic field B J on a measurement direction Cj. As for the phase shift

[0336]

[0337] it is equal to the sum A <p- + 3<pj, où:

[0338] <pj est un premier déphasage entre, d’une part, la projection des composantes locales du champ magnétique perturbateur B J on the measurement direction c, and, on the other hand, the alternating magnetic field generated by the emitter 22 at the origin of the eddy currents inducing this field B J, which first phase shift can be intrinsically deduced from the projection of the local components of the disturbing magnetic field B J on the measurement direction c b and 8 <pj est un deuxième déphasage entre, d’une part, les signaux de démodulation et, d’autre part, le champ magnétique alternatif généré par l’émetteur 22 à l’origine des courants de Foucault induisant le champ magnétique perturbateur B J , which second phase shift is known thanks to a model of the transmitter 22 and, if necessary, thanks to a synchronization system of the receiving device 23 with the transmitter 22.

[0339] It should also be noted that, alternatively, one can first perform a time projection of the local components of the disturbing magnetic fields onto the demodulation direction, and then transfer the results thus obtained into the secondary frame (R2) by applying the matrix R 2 3 to the matrix of these results.

[0340] This yields a matrix estimating the theoretical contribution of the disturbing element 90 to the raw magnetic measurement, whose coefficients consist of the terms

[0341]

[0342] B i j cos(φ i j ). A correction is then applied to this theoretical matrix to obtain a matrix M p estimation of the actual contribution of the disturbing element 90 to the raw magnetic measurement. This correction is typically adapted to simulate measurement errors due in particular to defects in the receiver 24, using techniques known to those skilled in the art.

[0343] A second variant of substep 1240, corresponding to the second embodiment (disruptive element 90 attached to the second object 14), will now be described, with reference to Figure 8.

[0344] In this second variant, substep 1240 comprises the same substeps 1241 and 1242 as described previously. These substeps 1241 and 1242 differ from those described above only in that the local components at receiver 24 of the disturbing magnetic fields produced by the disturbing element 90 are estimated using the calculation matrix Q determined during substep 1170 of step 1100, which provides the calculation parameters. The pre-existence of this calculation matrix Q allows for a faster calculation of these local components once the secondary current densities of the eddy currents are known.

[0345] However, unlike the first variant described above, these secondary current densities were not determined prior to the implementation of substep 1240. Substep 1240 therefore also includes substeps 1243, 1244, 1245, and 1246, which precede the implementation of substep 1241 for determining these secondary current densities. These preliminary substeps 1243, 1244, 1245, and 1246 include a first substep 1243 for determining the pose of the disturbing element 90 relative to the emitter 22. This pose is determined from the approximate pose of the second object 14 relative to the first object 12 and the (known) pose of the disturbing element 90 relative to the second object 14.

[0346] This substep 1243 is followed by a substep 1244 for determining the normal components of each of the alternating magnetic fields generated by the emitter 22, orthogonal to the surface 94 of the disturbing element 90 delimited by the portion of surface 92, at several of the mesh points 102, in particular at each of the mesh points 102 not belonging to the contour of the mesh 100. These normal components are typically determined from the pose of each mesh point 102 relative to the emitter 22 (known from the fact that a pose of the disturbing element 90 relative to the first object 12 was determined in substep 1243), a theoretical model of the alternating magnetic fields generated by the emitter 22, and an algorithm for correcting this theoretical model aimed at correcting, for example, the defects of the emitter 22, using known techniques. The man of the trade.

[0347] This substep 1244 is itself followed by a substep 1245 of estimation of components of the primary current densities of the eddy current densities induced in the disturbing element 90 by the magnetic fields generated by the emitter 22.

[0348] These primary current density components are typically obtained using the following formula or an equivalent formula:

[0349]

[0350] where / is a vertical matrix of the set of first components of the first primary current densities j,...,

[0351]

[0352] for an alternating magnetic field j generated by the emitter 22,

[0353]

[0354] is a vertical matrix of the set of second components of the second primary current densities y'i,..., j for the alternating magnetic field j, b is a vertical matrix of the normal components of the alternating magnetic field j determined in substep 1244, and Kj 1 is the inverse matrix of matrix K y which was determined prior to the implementation of the initial installation determination step 1200, during sub-step 1150 of step 1100 of providing calculation parameters.

[0355] The fact that the Ky matrix 1Pre-calculating the values ​​allows for a substantial time saving in the execution of substep 1245. Substep 1245 is followed by substep 1246, which determines the secondary current densities of the eddy currents induced in the disturbing element 90 by the magnetic fields generated by the emitter 22. As described for substep 1140, these secondary current densities are obtained by interpolating the components of the primary current densities. However, unlike substep 1140, this interpolation is performed using the interpolation matrix determined prior to the implementation of the initial pose determination step 1200, during substep 1160 of step 1100, which provides the calculation parameters.

[0356] Again, the fact that this interpolation matrix is ​​pre-calculated allows for a substantial time saving in the execution of this substep 1246.

[0357] Substep 1246 is followed by substeps 1241, then 1242.

[0358] Returning to Figure 6, substep 1240 is followed by substep 1250, which calculates a refined magnetic measurement. Put simply, in this substep 1250, the processing unit 80 subtracts the contribution of the disturbing element 90 from the raw magnetic measurement. In other words, it performs the following calculation: M a = M b — M p , where M a is the refined magnetic measurement, M b is the raw magnetic measurement, and M p is the actual contribution of the disruptive element 90 estimated during substep 1240. Note that this estimated contribution M p is different from the actual contribution M p because the position of the second object 14 relative to the first object 12 is not precisely known. Thus, although constituting a better estimate of the contribution M eof transmitter 22 than is the raw measurement M b , the refined measure M a is not equal to this contribution Me and a residual M e — M a remains.

[0359] Substep 1250 is followed by substep 1260 of deduction of a refined pose p f;k , where k=1,...,m is an index incremented at each iteration of substep 1260. This refined pose p f k is deduced from the refined magnetic measurement M a , typically by means of a theoretical model of the alternating magnetic fields generated by the emitter 22 and an algorithm for correcting this theoretical model aimed at correcting, for example, the defects of the emitter 22 and those of the receiver 24, by techniques known to the person skilled in the art.

[0360] Substep 1260 is followed by substep 1270, which calculates the difference between the approximate pose p j;k and the refined pose p f;k, here by means of a predefined function f, and comparison of this deviation to a threshold, typically a predetermined threshold.

[0361] If this difference is greater than the threshold, then substeps 1230, 1240, 1250, 1260 and 1270 are repeated, the approximate pose p j)k+1 being constituted, for the new iteration of these substeps, by the refined pose p f;k deduced at the end of the previous iteration. If the difference is less than the threshold, then substep 1270 is followed by substep 1280 in which the precise pose (which, for this first iteration, constitutes the initial pose) is determined to be the refined pose p f;k deduced at the end of the last iteration of substep 1260.

[0362] Step 1300 of determining an updated pose will now be described, with reference to Figure 9.

[0363] This step 1300 consists of a new iteration of the step for determining a precise pose 1201. Step 1300 therefore includes the same sub-steps 1210, 1220, 1230, 1240, 1250, 1260, 1270, 1280. Reference is made here to the description of these sub-steps given above.

[0364] This new iteration of the precise pose determination step 1201 differs from the first iteration only in the following characteristics:

[0365] The imprecise pose p0 determined during substep 1220 is preferably extrapolated from at least one previous pose determined during a previous iteration of the precise pose determination step 1201, rather than deduced from the raw magnetic measurement M b (although such an implementation of substep 1220 remains possible, in a non-preferred variant of the invention), and

[0366] Logically, the precise position determined during sub-step 1280 constitutes the updated position and not the initial position.

[0367] A calibration procedure 2000 for the tracking system 20, implemented by the data processing unit 80, will now be described, with reference to Figure 10.

[0368] This process 2000 includes a first step 2100 of providing initial hypothetical values ​​for calibration parameters of the tracking system 20. These calibration parameters typically include one or more of the following parameters:

[0369] the positioning parameters of the disturbing element 90 relative to that of objects 12, 14 with respect to which it is fixed,

[0370] geometric parameters of the perturbing element 90, such as dimensional parameters of its surface portion 92, parameters of concavity or convexity, parameters of holes formed in the surface portion 92, etc.

[0371] physical parameters of the disturbing element 90, such as the conductivity (possibly non-uniform and anisotropic) of the surface portion 92. These hypothetical values ​​of the calibration parameters are for example stored in the memory 84 of the data processing unit 80.

[0372] Other parameters of the tracking system 20 are, however, predetermined and are not subject to calibration.

[0373] Step 2100 is followed by a pose tracking step 2200 based on hypothetical values. During this step 2200, the pose of the second object 14 relative to the first object 12 is tracked using the pose tracking method 1000, for a limited number of preselected poses, using the hypothetical values ​​of the calibration parameters as the parameters for the perturbing element 90. For this implementation of the pose tracking method 1000, it is preferable that the imprecise pose p0 determined during each iteration of substep 1220 be systematically deduced from the raw magnetic measurement M b or, alternatively, deduced from a measurement provided by an external sensor (not shown), for example an imaging motion capture system.

[0374] Step 2200 is followed by a step 2300 of calculation, for each determined installation, of a criterion representative of an error on the determined installation.

[0375] This error criterion is preferably chosen so as to be an increasing function of the likelihood of an error on the determined pose.

[0376] According to a first variant of this step 2300, this criterion is a function of the diagonal terms of a diagonal matrix S of positive real values ​​S 11; S 22 and S 33 ordered so that > ​​S 22 > S 33 This diagonal matrix results from a singular value decomposition of a corrected matrix M. c deduced from the matrix of the refined measure M a obtained during the last iteration of substep 1250, said corrected matrix M c being obtained by means of a correction algorithm aimed at correcting the errors of the refined measurement M alinked, for example, to the defects of the transmitter 22 and those of the receiver 24 by techniques known to those skilled in the art. This singular value decomposition is such that the diagonal matrix S satisfies the following relation:

[0377] M c = PSQ T

[0378] where P is a first matrix belonging to the SO(3) group of rotation matrices in dimension 3 and Q T is the transpose of a second matrix Q also belonging to the group SO(3).

[0379] The error criterion is thus determined, for example, by the following formula:

[0380] D = \left(\frac{(S_{11} - 2S_{22})^2 + (S_{11} - 2S_{33})^2}{S_{11}^2 + 4S_{22}^2 + 4S_{33}^2}\right)^{1 / 2}

[0381]

[0382] S 2 11 + 4S 2 22 + 4S 2 33 )

[0383] According to a second variant of this step 2300, the representative criterion of an error on the determined pose is the variation of the gap between the pose determined for the receiver 24 and the pose determined, using the same tracking procedure 1000 and with the same hypothetical values ​​of the calibration parameters, for another receiver (not shown), fixed relative to the receiver 24.

[0384] Step 2300 is followed by step 2400 of calculation of a sum of said error criteria on all poses determined during substep 2200.

[0385] If this is the first iteration of step 2400, it is followed by a step 2700 of updating the hypothetical values, then steps 2200, 2300 and 2400 are repeated, the placement follow-up of step 2200 being carried out on the basis of the hypothetical values ​​thus updated.

[0386] If this is a subsequent iteration of step 2400, it is followed by substep 2500, which calculates a gradient of the sum relative to the variation of the hypothetical values. This gradient is the difference between the sum and the sum calculated in the previous iteration of substep 2400, relative to the variation of the hypothetical values ​​between the two iterations. This gradient is then compared to a threshold in step 2600. If it is greater than the threshold, steps 2700 and 2200 to 2600 are repeated. If it is less than the threshold, step 2600 is followed by step 2800, which validates the hypothetical values ​​provided in the last iteration of step 2700.

[0387] The hypothetical discounted values ​​provided in step 2700 are preferably established on the basis of a gradient descent algorithm for which the function to be minimized is the sum determined in step 2400.

[0388] Thanks to the invention described above, it is possible to correct in real time, for example at a frequency of 40 Hz, the measurements of an electromagnetic tracking system 20 disturbed by a complex-shaped interfering element 90, even for poses of the receiver 24 close to the interfering element 90. The calculations required for this correction are indeed simple, allowing for rapid execution, yet still enable a relatively accurate modeling of the interfering magnetic fields re-emitted by the interfering element 90, thus giving good accuracy to the corrected measurement. Consequently, the accuracy of the electromagnetic tracking system 20 is improved.

[0389] This correction also improves the accuracy of the static measurements of the electromagnetic tracking system 20, since it directly corrects the measurements of the receiver 24. It also eliminates the need for quasi-static phases of the receiver 24, since it does not require variation in the frequency of the alternating magnetic fields generated by the transmitter 22.

[0390] Finally, the invention enables autonomous calibration of the electromagnetic tracking system 20 without requiring the use of a reference system. This calibration is also relatively fast, especially compared to that of a mapping model.

Claims

DEMANDS 1. A method (1000) for tracking the placement of a receiver (24) relative to a transmitter (22), the transmitter (22) comprising at least two magnetic generators (30, 32, 34), each oriented along a proper direction (g1, g2, g3) and capable of generating an alternating magnetic field having a distinctive characteristic enabling the alternating magnetic field to be distinguished from the alternating magnetic field generated by each other magnetic generator (30, 32, 34), the magnetic generators (30, 32, 34) being arranged relative to each other such that their directions (g1, g2, g3) are non-coplanar, the method (1000) being implemented by a data processing unit (80) and comprising the following step: determination (1201) of a pose of the receiver (24) relative to the transmitter (22) at a determination instant, said determination (1201) comprising the following substeps: • a) obtaining (1210) a raw magnetic measurement composed of measurements, carried out at the time of determination, of local components, along at least two non-coplanar measurement directions (c1, c2, c3) attached to the receiver (24), of at least two resulting alternating magnetic fields constituting an ambient magnetic field, each of said resulting alternating magnetic fields having a distinctive characteristic identical to that of one of the generated alternating magnetic fields, • b) determination (1220) of an imprecise pose of the receiver (24) relative to the emitter (22) at the time of determination, • c) estimation (1240), based on an approximate pose that is a function of the imprecise pose, of a contribution from a disturbing element (90) to the raw magnetic measurement, • d) calculation (1250) of a refined magnetic measurement by subtracting the contribution of the disturbing element (90) from the raw magnetic measurement, and e) deduction (1260) of a refined pose of the receiver (24) relative to the emitter (22) at the instant of determination from the refined magnetic measurement, in which the disturbing element (90) is modeled as a surface mesh (100) composed of a set of mesh points (102) regularly distributed over a surface (94) of the disturbing element (90) and connected to each other by straight segments (104) delimiting between them two-dimensional and parallelogrammatic mesh cells (106), and substep c) of estimation (1240) of the contribution of the disturbing element (90) comprises, for each of the generated alternating magnetic fields, the estimation (1241) of local components, at the level of the receiver (24), for a pose of the receiver (24) equal to the approximate pose, of a disturbing magnetic field produced by eddy currents induced in the disturbing element (90) by said generated alternating magnetic field, by applying the Biot-Savart law to densities secondary currents of these eddy currents,said secondary current densities resulting from an interpolation of primary current density components of said eddy currents, said primary current density components being a function of a pose of the disturbing element (90) with respect to the emitter (22) and being constituted,:, of first components, along a first direction (d1) of the surface mesh (100), of first primary current densities in first internal segments (109) of the surface mesh (100) oriented along a second direction (d2) of the surface mesh (100), in particular in the middle of said first internal segments (109), and of second components, along the second direction (d2) of the surface mesh (100), of second primary current densities in second internal segments (110) of the surface mesh (100) oriented along the first direction (d1) of the surface mesh (100), in particular in the middle of said second internal segments (110).

2. A tracking method (1000) according to claim 1, wherein the primary current density components are obtained by means of the following formula or an equivalent formula: where / is a vertical matrix of the set of first components of the first primary current densities, / is a vertical matrix of the set of second components of the second secondary current densities, b is a vertical matrix of normal components of the generated alternating magnetic field, orthogonal to the surface (94) of the disturbing element (90), at each of the mesh points (102) not belonging to a contour of the surface mesh (100), and Ky 1 is the inverse matrix of a matrix constructed from the local Maxwell equations discretized between the primary current density components and the normal components of the generated alternating magnetic field.

3. A tracking method (1000) according to claim 1 or 2, wherein the disturbing element (90) is fixed relative to the receiver (24) and substep c) of estimating (1240) the contribution of the disturbing element (90) comprises determining (1243) a pose of the disturbing element (90) relative to the emitter (22), estimating (1245) the components of the primary current densities as a function of said pose of the disturbing element (90) relative to the emitter (22), and determining (1246) the secondary current densities by interpolating the components of the primary current densities.

4. A tracking method (1000) according to claims 2 and 3 taken together, wherein the matrix Kj 1 is determined prior to the implementation of the installation determination step (1201).

5. A tracking method (1000) according to claim 1 or 2, wherein the disturbing element (90) is fixed relative to the emitter (22) and the secondary current densities are determined prior to the implementation of the pose determination step (1201).

6. A tracking method (1000) according to any one of the preceding claims, wherein the pose determination step (1201) comprises several successive iterations of substeps c), d) and e), the approximate pose being constituted, for the first iteration, by the imprecise pose and, for each subsequent iteration, by the refined pose deduced at the end of the previous iteration.

7. A tracking method (1000) according to any one of the preceding claims, comprising several iterations of the pose determination step (1201) to determine poses of the receiver (24) at different determination times, the first iteration constituting a step (1200) for determining an initial pose of the receiver (24) and subsequent iterations constituting steps (1300) for determining an updated pose, the imprecise pose determined during the initial pose determination step (1200) being deduced from the raw magnetic measurement.

8. A tracking method (1000) according to claim 7, wherein the imprecise pose determined during at least one of the steps for determining an updated pose (1300) is extrapolated from at least one prior pose determined during a prior iteration of the pose determination step (1201).

9. A monitoring method (1000) according to any one of the preceding claims, wherein the secondary current densities are determined at evaluation points (s1, ..., s 12 ) regularly distributed within each mesh cell (106).

10. A tracking method (1000) according to any one of the preceding claims, wherein the local components of the disturbing magnetic field are estimated using the following formula or an equivalent formula: ' ' taZ. ' ' |rS||3 where bj(r) is a vector of said local components of the disturbing magnetic field, p0 is the magnetic permeability of free space, / (Sj) is a 3-dimensional vector of the components, in a frame attached to the disturbing element (90), of one of the secondary current densities, is the surface area of ​​an elementary cell (120) associated with said estimated surface current density, e t is an average thickness of a portion of the surface (92) of the disturbing element (90), within which the eddy currents are concentrated, at the level of the mesh cell (106) comprising the elementary cell (120), s t is a vector of the coordinates of a center of said elementary cell (120) in a frame attached to the perturbing element (90), r is a vector of the coordinates of the approximate pose in said frame attached to the perturbing element (90), x is the cross product operator and | • | is the norm 2, the elementary cell (120) being constituted by one of the mesh cells (106) or by a subdivision of a subdivided mesh cell (106), each subdivided mesh cell (106) being made up of elementary cells (120), each associated with one of the secondary current densities, regularly distributed within said subdivided mesh cell (120) and covering the whole thereof.

11. Calibration method (2000) of an electromagnetic tracking system (20) comprising a step (2200) of tracking a pose of a receiver (24) of the electromagnetic tracking system (20) relative to a transmitter (22) of the electromagnetic tracking system (20) by means of a tracking method (1000) according to any one of the preceding claims, on the basis of hypothetical values ​​of calibration parameters of the electromagnetic tracking system (20), and updating (2700) the hypothetical values ​​of the calibration parameters so as to minimize a sum of error criteria each representative of an error on a pose of the receiver (24) relative to the transmitter (22) determined during the tracking step (2200).

12. Product computer program comprising code instructions for the implementation, by a processor, of a tracking method (1000) or a calibration method (2000) according to any one of the preceding claims.

13. Computer-readable recording medium on which a computer program product according to claim 12 is stored.

14. Assembly (10) comprising a reference object (12), a moving object (14) and an electromagnetic tracking system (20) for tracking the position of the moving object (14) relative to the reference object (12), said electromagnetic tracking system (20) comprising: an emitter (22) integral with the reference object (12), said emitter (22) comprising at least two magnetic generators (30, 32, 34), each oriented along a proper direction (g1, g2, g3) and capable of generating an alternating magnetic field having a distinctive characteristic enabling the alternating magnetic field to be distinguished from the alternating magnetic field generated by each other magnetic generator (30, 32, 34), the magnetic generators (30, 32, 34) being arranged relatively to each other such that their directions (g1, g2, g3) are non-coplanar, and a receiving device (23) comprising a receiver (24) integral with the moving object (14), said receiving device (23) being capable of measuring local components, along at least two non-coplanar directions (c1, c2, c3), of at least two resultant alternating magnetic fields constituting an ambient magnetic field, each of said resultant alternating magnetic fields having a distinctive characteristic identical to that of one of the generated alternating magnetic fields, so as to obtain a raw magnetic measurement, the assembly (10) also comprising a disturbing element (90) capable of impacting the measurement of the local components of the resulting alternating magnetic fields by the receiving device (23) when the disturbing element (90) is subjected to the alternating magnetic fields generated by the transmitter (22), in which the electromagnetic tracking system (20) includes a data processing unit (80) configured for the implementation of a tracking method (1000) or a calibration method (2000) according to any one of claims 1 to 11.

Images