METHOD FOR LOCATING A MAGNETIC OBJECT
20180306872 ยท 2018-10-25
Inventors
Cpc classification
G01R33/0064
PHYSICS
International classification
Abstract
A method allowing the determination of the location and/or the orientation of a magnetic field source in space, and more particularly a method for determining the relative position in space of at least one magnetic field source in relation to at least one magnetic field sensor comprises the steps of acquiring measurements of the magnetic field, computing a solution of the expression of the magnetic field generated by at least one magnetic field source by modeling each magnetic field source by an element chosen from among a superposition of solenoids and a superposition of charged planes, then by estimating the value of a complete elliptic integral linked to the model by an algorithm using a Landen transformation and computing at least one element chosen from among the position and the orientation of each magnetic field source.
Claims
1. A method for determining the relative position in space of at least one magnetic field source in relation to at least one magnetic field sensor comprising the steps of: a) acquiring measurements of the magnetic field generated by at least one said magnetic field source by means of at least one magnetic field sensor; b) computing, by means of at least one processor, a solution of the expression of the magnetic field generated by at least one said magnetic field source by modeling each said magnetic field source by an element chosen from among a superposition of solenoids and a superposition of charged planes, then by estimating the value of a complete elliptic integral linked to said model by an iterative Bulirsch algorithm using a Landen transformation; c) computing at least one element chosen from among the position and the orientation of each said magnetic field source by using at least one method chosen from among: a minimization of the cost between said measurements of said magnetic field sensors and said solutions computed in the step b); a Bayesian filtering, at successive instants, of each said measurement of each said magnetic field sensor, said filtering being a function of each said solution computed in the step b).
2. The method as claimed in claim 1, wherein at least the Bayesian filtering is chosen as method for the step c), said Bayesian filter being a Kalman filter.
3. The method as claimed in claim 1, wherein at least one said magnetic field source is in motion in relation to at least one said magnetic field sensor and wherein at least one element chosen from among the position and the orientation of at least one said magnetic field source and the position and the orientation of at least one said magnetic field sensor is computed from a Kalman filtering of said measurements of at least one said magnetic field sensor, said Kalman filtering being a function of each said computed solution, at each of said successive instants.
4. The method as claimed in claim 1, wherein at least one said magnetic field source can be placed at less than ten centimeters from at least one said magnetic field sensor.
5. A human/machine interface comprising at least one magnetic field source, at least one magnetic field sensor, and at least one processing unit programmed appropriately to implement each step of claim 1.
6. The human/machine interface as claimed in claim 5, comprising a plurality of magnetic field sensors forming overall a non-planar array.
7. The human/machine interface as claimed in claim 5, wherein at least one said magnetic field source is chosen from among a single-layer solenoid and a solid cylindrical magnet.
8. The human/machine interface as claimed in claim 5, comprising at least two magnetic field sources, wherein at least one of said magnetic field sources is a solenoid, and wherein the processing unit is programmed to command a periodic emission of the intensity of the magnetic field of said solenoid at a different frequency from at least one other of said magnetic field sources.
9. The human/machine interface as claimed in claim 5, wherein at least one said magnetic field source is a permanent magnet composed of an alloy of elements chosen from among iron, vanadium, neodyme, boron, aluminum, nickel, titanium, copper, samarium and cobalt.
Description
[0025] The invention will be better understood and other advantages, details and features thereof will become apparent from the following explanatory description, given by way of example with reference to the attached drawings in which:
[0026]
[0027]
[0028]
[0029]
[0030]
[0031]
[0032]
[0033] Hereinafter in the text, the term magnetic field sensor is defined by a triaxial magnetic field sensor: the measurement of the magnetic field is performed on three different axes of space at the point of the magnetic field sensor.
[0034]
[0035] In a particular embodiment of the invention, the magnetic field source 3 can be a permanent magnet, that is to say a magnet whose magnetic moment is non-nil in the absence of external magnetic field. The form of the magnetic field source 3 can be a solid cylinder or a holed cylinder. The magnetic field source 3 can be manufactured in an alloy of elements chosen from among iron, vanadium, neodyme, boron, aluminum, nickel, titanium, copper, samarium and cobalt.
[0036] In a particular embodiment of the invention, the magnetic field sensors 4 of one and the same element of reference 2, when the device uses a number greater than or equal to four thereof, are not coplanar. This arrangement makes it possible to obtain information on the magnetic field gradient for each of the axes defined in
[0037] The magnetic field source 3 can also be a solenoid passed through by an alternating electrical current.
[0038] The human/machine interface 1 can comprise several magnetic field sources 3. In the case where several of these magnetic field sources 3 are solenoids passed through by alternating currents, it is possible to impose a frequency of variation of the electrical current particular to each of the solenoids, in order to differentiate the magnetic field sources 3 upon their detection.
[0039] More generally, in embodiments of the invention, the human/machine interface can comprise at least one field source 3 produced by a solenoid; the processing unit 7 can be programmed to command a periodic emission of the intensity of the magnetic field of said solenoid, at a different frequency from at least one other frequency of emission of the magnetic field of another magnetic field source 3. This feature makes it possible to distinguish the origin of the magnetic field received by the magnetic field sensor or sensors (4) by analyzing the different frequencies that make up one or more signals. It allows the simultaneous location of several magnetic field sources 3.
[0040] In another particular embodiment of the invention, the magnetic field source 3 is produced by the passage of electrical current in a single-layer solenoid.
[0041] The measurement element 2 can, in a particular embodiment of the invention, comprise an array of magnetic field sensors 4. The magnetic field sensors 4 can be placed in rows and in columns to form a matrix.
[0042] Each magnetic field sensor 4 can be fixed with no degree of freedom to the other magnetic field sensors 4. To this end, the magnetic field sensors 4 are fixed onto the top face of a rigid substrate. Each magnetic field sensor 4 measures the direction and the intensity of the magnetic field generated by one or more magnetic field sources 3 at the location of the sensor. For that, each magnetic field sensor 4 measures the norm of the orthogonal projection of the magnetic field generated by the magnetic field source or sources 3, at the location of this magnetic field sensor 4, on the three axes of this magnetic field sensor 4. In an exemplary embodiment, these three axes are mutually orthogonal. For example, a magnetic field sensor 4 will measure the components of the magnetic field on axes colinear to i, j and k.
[0043] Each magnetic field sensor 4 is connected to the processing unit 7 via an information transmission bus 8.
[0044] The processing unit 7 is capable of determining the position and/or the orientation of the field source or sources in the reference frame (i,j,k) of the measurement element 2 from the measurements of the magnetic field sensor or sensors 4. The processing unit 7 comprises, among other things, one or more processors and memories. The processor or processors are programmed to execute instructions stored in the memory or memories.
[0045] In particular, the processing unit 7 associates a mathematical model with all or some of the measurements of the magnetic field sensors 4 to deduce therefrom the position and/or the orientation of the magnetic field source or sources 3.
[0046] The processing unit 7 is also capable of restoring the location and/or the orientation of one or more magnetic field sources 3 to digital platforms via the output 9.
[0047]
[0048] Trmolet de Lacheisserie, . (Magntisme: Fondements, tome l de Collection Grenoble Sciences. [Magnetism: The Basics, volume 1 of Collection Grenoble Sciences.] EDP Sciences; 2000) demonstrates that induction B generated by the magnetized substance and/or by imposed electrical currents of density j.sub.0, at any point of space, is reduced to a problem of electromagnetism of vacuum in which two types of currents are to be considered: [0049] the real (or free) currents of volume density j.sub.0 [0050] the currents associated with the magnetized substance (or linked currents), of volume density j.sub.m and surface density j.sub.ms given by:
j.sub.m=rotM(1)
j.sub.ms=Mn(2)
in which n is the unitary vector normal to the surface of the magnetized substance illustrated in
[0051] The relationship (1) indicates that, in the case of a uniform magnetization, j.sub.m=0. Thus, a uniformly magnetized cylinder is equivalent to a solenoid passed through by a surface current density |j.sub.ms|=|M|.
[0052]
[0053]
.sub.m=div M(3)
.sub.m=M.Math.n(4)
[0054] The magnetic field considered is, here, the magnetic field due to the substance (sometimes called demagnetizing field of the magnetized substance, and field of the substance outside of the magnetized substance). The magnetic induction B is obtained by taking into account the ambient magnetic field H.sub.0 and the magnetization M by the relationship B=.sub.0(H.sub.0+H.sub.m+M). In the case of a magnet of FeNdB type, the ambient field, assumed to be of the order of the Earth's magnetic field in this embodiment of the invention, will be disregarded with respect to the remnant magnetization. When the magnetization is uniform, according to (3), .sub.m=0; the magnetization can be characterized by surface charge densities .sub.m. For example, the field H.sub.m created (in the substance and outside) by a cylinder uniformly magnetized along its axis is the same as that generated by two disks of uniform surface densities .sub.m=|M| situated on each of the planar faces of the cylinder.
[0055]
[0056] Durand E. (Electrostatique et Magntostatique, Tome l Les distributions [Electrostatics and Magnetostatics, Volume I: Distributions], Edition Masson et Cie, p. 246-250) has presented the computations of the magnetic inductions generated for hypothetical distributions of magnetic currents or masses, in particular for a solenoid passed through by a surface current j.sub.ms. As described previously, a solenoid of main axis z passed through by a surface current j.sub.s is the amperian equivalent of a cylindrical permanent magnet of axis z, and of uniform magnetization M=|js|z. Thus, the knowledge of the magnetic induction generated at any point of space by a solenoid passed through by a surface current j.sub.s makes it possible to compute the magnetic induction generated at any point of space by a cylindrical magnet of uniform magnetization M=|js|z.
[0057] In cylindrical coordinates, the center of the reference frame being the geometrical center of the solenoid O, a point M being identified by its coordinates (r,,z), the induction produced at any point of space by the solenoid is independent of the variable because of the symmetry of the problem. The components B.sub.r and B.sub.z are given by:
and K, E and the complete elliptic integrals of the first, second and third species:
[0058] One of the features of the invention is combining the use of the amperian and/or coulombian models presented above with a resolution of the elliptical integrals (9), (10) and (11) by so-called Bulirsch algorithms (Bulirsch R., Numerical Calculation of Elliptic Integrals and Elliptic Functions, Numerische Mathematik 7, 78-90, 1965). These are iterative algorithms allowing the evaluation of common mathematical functions, in this case of the elliptical integrals, by using the fewest possible iterations.
[0059] An advantageous implementation of this algorithm, in terms of efficiency, is presented in Derby, N., & Olbert, S. (2010), Cylindrical magnets and ideal solenoids, American Journal of Physics, 78(3), 229-235. This implementation makes it possible to compute the function C (for Complete Elliptic Integral), defined by:
[0060] It is then possible to rewrite the equations (5) and (6) according to the notations of Derby et al.:
[0061] The function C can be computed efficiently by using the Bulirsch algorithms, numerically translating a Landen transformation.
[0062] The Bulirsch publication considers three distinct cases. In the case of the computation of a complete integral of first species, which can be written in the form:
the algorithm translates a numerical integral computation using a Landen transformation. In the case of a complete integral of second species (including the case of the integrals of first species), which can be written in the form:
the algorithm also translates a numerical integral computation using a Landen transformation. It is also possible to note that, for a complex argument of the elliptic integral, the resolution of the elliptic integrals of first and second species use a Gauss transformation because the Landen transformation is not numerically stable in this case, but most of the embodiments of the invention use Landen transformations.
[0063] Finally, a third algorithm is proposed in the case of a complete elliptic integral of third specifies, that is to say of the form:
[0064] This third algorithm also translates a numerical integral computation using a Landen transformation (whereas the more general form, in the case of a non-complete elliptic integral is computed by using a Gauss transformation). In the embodiments of the invention, the computation of the integrals K(m), E(n) and (n,m) or C(k.sub.c, p, c, s) corresponds to the computation of complete integrals: the embodiments of the invention use Landen transformations to advantageously compute these integrals.
[0065] The computation of the function C requires on average five iterations to converge to a first accuracy of 10.sup.12 which can make it possible to use this method to model, for example, the magnetic field generated by solenoids of different sizes without using the dipole approximation. The inventors have discovered that the use of this method for computing the magnetic field could make it possible to make the computation load acceptable for a standard processing unit, or a standard processor, while modeling the magnetic field at distances very close to the magnetic field source or sources. For example, in the case of a magnetic source of a characteristic dimension of a centimeter, an embodiment of the invention makes it possible to model the magnetic field, by virtue of the computation of the complete elliptic integrals by a Landen transformation, at a distance less than ten centimeters, which is impossible by using a dipole model of the magnetic field sources. When measuring the position and/or the relative orientation of a magnetic field source 3 in relation to a sensor 4, it is possible to bring the source 3 and the sensor 4 to within a distance less than 10 cm. In other embodiments of the invention, the modeling of the magnetic field 3 is done by a coulombian model. It is possible to compute equivalent complete elliptic integrals equivalent to K, E and by the Bulirsch algorithms.
[0066]
[0067] The formulae (13) to (19) are expressed in the local cylindrical reference frame of the magnetic field source 3.
[0068] To switch from the global reference frame (O,i,j,k) to a local reference frame linked to the solenoid S(S,i,j,k), it is necessary to perform a translation and two rotations. The translation is performed with the vector p.sub.s then a rotation about the axis k makes it possible to place the solenoid in the plane (ik). Finally, a rotation about j makes it possible to straighten the solenoid such that its axis coincides with the axis k. The position of the point M in the reference frame S is given by Ms:
M.sub.S=R.sub.y()R.sub.z()(OMOS)(23)
[0069] The point M is described in the Cartesian reference frame linked to the magnetic source.
[0070] The panel B of
M(,0,z)=R.sub.z()M.sub.S(24)
[0071] The effect of this rotation is presented by the panel C of
B.sub.=R.sub.z()R.sub.y()R.sub.z()(B.sub.,0,B.sub.z)(25)
[0072] In embodiments of the invention, the rotation matrices are generated by geometrical operations not involving trigonometric functions of arc cosine or arc tangent type for reasons of numerical stability and speed of material hardware execution.
[0073] It is advantageously possible to use a second method for switching from the global reference frame w to the local reference frame S linked to the magnetic field source 3. P.sub.m/s=OMOS is used to define the relative position of the measurement point M in relation to the center of the solenoid S. The vector M is defined as the moment of the magnetic field source 3 (colinear to the axis of the magnetic field source 3). The normalized moment of the magnetic field source 3 is defined by M.sub.n=M/M. The unitary vector from S to M is defined by P.sub.m/s,n=P.sub.m/s/P. In the cylindrical reference frame linked to the magnetic field source 3, the coordinates of the point M are (x,y,z)=(0,0,P.Math.M.sub.n). The magnetic induction vector B.sub. in the global reference frame is then defined after the computation of the magnetic induction vector B.sub.S in the local reference frame S by B.sub.=B.sub.S.Math.M.sub.n+B.sub.S.Math.(M.sub.nP.sub.m/s,nM.sub.n).
[0074] The invention uses estimation methods which use, on the one hand, the analytical models of induction introduced previously, and, on the other hand, measurements from sensors. The combined use of these methods then makes it possible to compute the relative positions and/or orientations of one or more magnetic field source(s) 3 in relation to one or more magnetic field sensors 4. It can also be used to estimate the position and/or the orientation of one or more magnetic sensors 4 in relation to one or more magnetic field sources 3.
[0075] In the embodiments of the invention, a step occurs after the acquisition of the measurements from the magnetic field sensor or sensors 4 and after the computation by the processor unit 7 of the expression of the magnetic field 3 generated by at least one magnetic field source 3 by modeling the sources 3 by a superposition of solenoids and/or a superposition of uniformly charged planes, as described previously. This step consists in estimating the location and/or the orientation of a magnetic field source 3 and/or of a magnetic field sensor 4 according to a method chosen between an optimization, a Bayesian filtering or a method that is a hybrid of an optimization and a Bayesian filtering.
[0076] In the case of the optimization, the aim is to find a vector of parameters x (for example the position and/or the orientation of the magnetic field source 3, the magnetic moment of the magnetic field source 3, or the position and/or the orientation of the magnetic field sensor 4) which minimize a given cost. In this case, the optimization consists in minimizing the costs between the measurements of said sensors and the solutions computed in the modeling of the magnetic field. For example, it is possible to use the least squares method:
x=argmin.sub.x((zh(x)).sup.T(zh(x)))(26)
[0077] Here z is the measurement performed by the magnetic field sensor or sensors 4. h(x)=B.sub.(x) is the model proposed in the preceding section. To resolve this problem, it is for example possible to use the Gauss-Newton method or the Levenberg-Marquardt method.
[0078] In the case of the Bayesian filter, it is possible to use the 1st order Markovienne approximation to estimate the position of one or more magnetic field sources 3 in relation to the magnetic field sensor(s) or vice-versa. This operation can be performed in two steps: [0079] a prediction using the state at the preceding instant using a dynamic model of state evolution to predict the current state, [0080] an update based on the innovation supported by the measurements given by the magnetic field sensors 4 in relation to the predicted measurements, computed on the basis of the amperian model explained previously.
[0081] Generally, to estimate the a posteriori state from an a priori state and from a measurement, the Bayes' theorem can be used:
[0082] The computation of the likelihood involves a function h which links the state to the measurement. Here, the state is the position and/or the orientation of the magnetic source and the measurement is the magnetic field measured by the magnetic field sensor or sensors 4. This principle applies to all the Bayesian estimators. In particular, a Kalman filter can be used to resolve this problem. In this computation, it is possible to estimate the magnetic moment vector instead of the orientation to eliminate certain trigonometric computations.
[0083]
[0084] In a Kalman filter, the state is described by the first two moments of its probability distribution. Since the latter is assumed Gaussian, the two moments are sufficient to describe it completely. The expectation of the state is called
[0085] A Kalman filter is said to be extended if it reuses the same equations used in a conventional Kalman filter by locally linearizing the measurement function:
h(x.sub.t+1)=H.sub.t+1x.sub.t+1(28)
[0086] The function is computed from (25). H.sub.t+1 is the jacobian matrix evaluated at
P.sub.t+1|t=FP.sub.t|tF.sup.T+Q.sub.t(29)
where F is an identity matrix and Q is the covariance matrix of the state noise. It describes the uncertainty on the evolution of the state. The update phase makes it possible to use the proposed model to correct the prediction on the basis of the information provided by the current measurement, called innovation {tilde over (y)}.sub.t+1. The following applies:
P.sub.t+1|t+1=P.sub.t+1|tK.sub.t+1H.sub.t+1P.sub.t+1|t(30)
[0087] Here, H.sub.t+1=.sub.x[h(x.sub.t+1)] is the jacobian matrix of the measurement function h evaluated at
{tilde over (y)}.sub.t+1=z.sub.t+1h(
K.sub.t+1=P.sub.t+1|tH.sub.t+1.sup.T.sub.{tilde over (y)}.sub.
.sub.{tilde over (y)}.sub.
[0088] In embodiments of the invention, the position and/or the relative orientation of a field source 3 in relation to the position and/or the orientation of a sensor 4 are computed from a Kalman filtering as described previously. The same principle applies for the Kalman filter of UKF type, as well as the particle filters and all the other Bayesian filters.
[0089] It is advantageously possible to use an algorithm called LMA (described in Aloui, S., Villien, C., & Lesecq, S. (2014), A new approach for motion capture using magnetic field: models, algorithms and first results. International Journal of Adaptive Control and Signal Processing) in place of the Kalman filter. This approach is based on the use of a hybrid Bayesian/optimization algorithm of Levenberg-Marquard type. When using the LMA algorithm, a prediction step is initially performed, similar to the prediction step described previously when using a Kalman filter. Secondly, an update step is performed, similar to the Kalman method described previously, but by replacing the cost optimization equations with:
[0090] Each iteration i{1 . . . n.sub.i} is computed by:
x.sub.t+1|t+1.sup.i=x.sub.t+1|t+1.sup.i1+.sub.t+1|t+1.sup.i1(35)
in which:
.sub.t+1|t+1.sup.i1=[(P.sub.t+1|t+1.sup.i1).sup.1+.sup.idiag((P.sub.t+1|t+1.sup.i1).sup.1)].sup.1e.sup.i1(36)
and:
e.sup.i1=[(H.sub.t+1.sup.i1).sup.TR.sub.t+1.sup.1(z.sub.th(x.sub.t+1|t+1.sup.i1))+(P.sub.t+1|t).sup.1(x.sub.t+1|tx.sub.t+1|t1.sup.i1)](37)
(P.sub.t+1|t+1.sup.i1).sup.1=H.sub.t+1.sup.i1TR.sub.t+1.sup.1H.sub.t+1.sup.i1+(P.sub.t+1|t).sup.1(38)
where is the damping factor, non-constant. The damping factor is initialized at .sup.0=1. For each iteration i, the error e.sup.i=e.sup.i.sub.2 is computed. If e.sup.i>e.sup.i1, x.sub.t+1| t.sup.i is left as such and is reduced (for example .sup.i+1=.sup.i/10). Otherwise, x.sub.t+1|t.sup.i is rejected (x.sub.t+1|t.sup.i=x.sub.t+1|t.sup.i1) and is increased (for example .sup.i+1=10.sup.i).
[0091] The covariance matrix of the a posteriori estimated distribution is P.sub.t+1|t+1=P.sub.t+1|t+1.sup.n.sup.