Signal accumulation in problems with constraints and its application to real-life magnetometry

We consider the method of signal accumulation for problems with natural constraints. We focus on the case when the constraint has the implicit form of the type ∥x∥=C with the unknown constant C to be determined as part of the problem solution. The construction and derivation are carried out based on application to the real-life magnetometry problems related to public security systems.


Introduction
Magnetic object localization and classification are an old problem that is finding its way into a multitude of applications. Passive magnetometry has been employed in fields ranging from geological exploration, biomedical treatment and wreck removal to localization of unexploded ordinance [1]. It comprises the advantages of all-weather performance, simple equipment and convenient signal processing.
This technique appears to be especially useful in security systems where there is a need to detect and classify metallic objects of moderate size and distinguish potentially hazardous objects (such as small-arms or cold arms) from benign objects. This is due to the fact that every gun has a barrel, which is usually made of ferromagnetic material which, under the influence of Earth's magnetic field, becomes a magnet itself and, hence, produces a magnetic field of its own. The easiest model of a * Author to whom any correspondence should be addressed.
Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. magnet is a point dipole, which is employed in most of the works on the matter [2][3][4].
To locate a moving magnet and measure its magnetic properties, magnetic sensor arrays are commonly used. The United States Naval Research Laboratory employed a superconducting gradiometer scheme to locate a moving magnetic dipole source as early as 1975. Wynn, instead of only measuring a magnetic induction vector, relied on the measurements of a magnetic gradient tensor to track a magnetic dipole, and managed to track a moving magnetic dipole using continuous measurement performed by a static measuring station [2,3]. Subsequently, numerous methods have been developed to locate the target with magnetic sensor arrays. In 2003, Heath devised MATLAB algorithms for the three-dimensional (3D) inversion of potential field tensor data using Monte Carlo and Downhill Simplex approaches [5]. In 2006, Nara derived a simple formula for the localization of a magnetic dipole using measurements of a magnetic induction vector at two points in close proximity and employing the magnetic gradient tensor calculated from these measurements [4]. In all these cases the authors concentrated on the localization of the magnet first, implying that given the magnetic dipole's coordinates, its magnetic induction could be easily calculated from the first principles of magnetostatics.
In 2007, Arie formulated the problem as an overdetermined nonlinear equation set using a magnetic dipole model for the stationary target and employing a synthetic aperture via a moving magnetic sensor to collect bulk measurement data. He used a simulated annealing algorithm inspired by equilibrium thermodynamics to rapidly find a good approximation to the global optimum of this equation set [6]. In 2014, Wahlstrom indicated that the sensor models could be combined with a standard motion model and a standard nonlinear filter to track metallic objects in a magnetometer network [7].
In 2017, Gao et al [8] proposed a method of target localization with the alternating magnetic field based on coherent demodulation, but the single alternating magnetic dipole should have no roll and move at a constant speed.
More recently, Liu et al [9] returned to the magnetic gradiometer technique and focused on the detection of metallic objects based on their magnetic signatures without trying to classify detected objects at the danger scale. Although this approach seems to be sensitive enough, it is still an approximation because the calculation of the magnetic gradient tensor based on measurements of magnetic induction in two distinct points of space cannot be exact due to the high nonlinearity of magnetic induction dependence on range. This is why we abandon the gradiometry approach and base our method on 3D magnetometers only. This is similar to the approach used by Chiba and Nara [10], who return to the 3D-magnetometer technique to find the position of a mobile device using magnetic signals from two rotating magnets and using a fast Fourier transform to retrieve spectral properties of the resulting magnetic signals in what may be called a cooperative environment. Our system uses absolutely non-cooperative targets and a sensor fusion approach to not only localize the magnetic target but to accurately measure its magnetic strength as well.
It is worth mentioning different empirical methods like the one described in [11]. There, the authors also used a gradiometer measuring system (see above) together with advanced statistical methods for nonlinear empirical pattern recognition and classification of magnetic signatures. More specifically, [11] considers the time dependence of signals, which are measured at different time frames (and corresponding to different positions of a moving target), and then employs a Fourier-based digital signal processing (DSP) technique (called 'Joint Time/Frequency Analysis') to process them. The inescapable problem with this approach is that these time series contain not only information on the magnetic properties of the target in question but also information on the target's movement. The item's magnetic signature is variable and does not behave linearly. The gait and speed of passage, proximity to the center of the portal and background clutter all impact and alter the concealed weapon detector(CWD) magnetic signature. This makes it difficult to reduce the response to a mathematical-based linear algorithm.
Our approach, however, gets rid of this difficulty by employing a matched filter: using the information on the target's movements obtained independently (via ToF cameras, which measure the positions of the target at every time frame) we are able to detach that information which is pertinent to the target's movement from the target's magnetic properties (i.e. the m vector) and deal only with the latter. This is described in detail in section 2 below. Therefore, we are not forced to consider any time series and to employ any threshold analysis; instead, we deal only with the bunch of raw data and perform measurements of the physical value of interest, namely the m . The latter may be used for threat discrimination as usual using either threshold analysis, neural networks, machine learning, etc. Our goal is to provide input information for the following steps of target acquisition, which we thoroughly do, introducing no additional assumptions, and to find an optimal solution (minimum of a discrepancy).
Another misfortune of the approach in [11] is that the details of the algorithm remain unknown as proprietary software is employed. This makes direct comparison via simulations impossible. We plan to extend and deepen our mathematical model and calibrate it with the use of a real inspection zone, where we will also compare it with existing empirical algorithms. At the moment, these are beyond the scope of the present work, which is mostly about the theoretical development of the novel approach designed for the constrained problems widely present in real life but essentially missing in research.
With regard to accurate measurement of vector-valued variables under heavy noise conditions, the common practice is the approach where one smooths out the measurement results with a type of filter: the Kalman filter and its variations being among the most popular. As simple as it is, i.e. an easy-toimplement recursive filter with very moderate computational requirements, it allows one to smooth the measurement results, taking into account the a priori postulated motion model for the measured variable and using up all the information contained in the measurements supplied by the measuring device.
In practice, however, we often encounter situations when the measured quantity is restricted by various constraints, which are inescapable consequences of the underlying physics, for example, the measured vector might have a constant norm invariant in time. Such constraints inherent in the problem represent additional information, which common filtering techniques (and the Kalman filter in particular) fail to use. This problem has attracted a lot of attention recently. Significant progress has been achieved in the investigation of the state-ofthe-art constrained Kalman filter (see [12][13][14][15][16][17][18] and many others). Still, such an approach can account for known constraints of the type of equations or inequalities (and even these are best treated in case they are linear/nonlinear constraints requiring approximations of the Taylor expansion type), but it fails to use the undefined constraints. For example, it can deal with the situation of the constraint of the type x = C where the constant C is known, but it fails in such a constraint when C is unknown and subject to determination.
Another well-known technique for combating measurement noise is that of signal accumulation [19][20][21], exploited in a multitude of engineering applications (Reed-Solomon codes being one of the most famous examples). The basic idea here is that if the measured quantity represents a sum of a non-random signal of interest and a random white noise (which, according to the central limit theorem, usually implied to be Gaussian), then the additive accumulation of measurements leads to coherent signal accumulation and incoherent noise accumulation, which increases the signal-to-noise (SNR) ratio, which in turn leads to smaller measurement error. In this way, if we combine N measurements, the signal rises roughly N times, the noise rises √ N times and the resulting SNR rises √ N times. This is a very general technique which finds its way into very different applications, the most common of which is a matched filter technique, which is the same approach in disguise. But as powerful as it is, it fails to account for the inherent constraints a signal may be subject to.
The main contributions of the present manuscript are the following.
• In addition to detection of an object possessing a magnetic signature we calculate its position on a person's body and its magnetic signature itself which, in our case, is a magnetic moment vector of the representing dipole. As shown in this paper this set of data greatly enhances classification of the object and helps to avoid most false alarms. • We show how to employ the latter technique for the problem of the measurement of the vector quantity subject to various constraints. • Finally, though we take as an example the problem of magnetic measurements under noise encountered in real-life applications, the scope of the approach presented is in no way limited by magnetometry alone but can be employed for any linear system under suitable constraints.
The paper outline is as follows: section 1 describes the reallife magnetometry problem in detail, laying out the approach proposed; section 2 establishes the main equation system subject to solution; section 3 explains our approach to solving it; section 4 contains the main existence and uniqueness theorem, stating the solution exists and describing a numerical algorithm for finding it; section 5 presents the numerical data obtained via computer simulation and contains a discussion of the benefits of the approach proposed.

Real-life problem in magnetometry
When designing a system for detection of concealed smallarms worn by people entering a public place, a common approach is to employ magnetometry. A point magnet is characterized by a vector of magnetic moment m, which is an intrinsic property of the magnet. Its orientation represents the orientation of the magnet itself and its norm does not depend on anything and does not change in time, being a characteristic feature of the magnet (i.e. its 'magnetic force'). Determination of this magnetic force can help us to distinguish between benign and hazardous objects, and its orientation can help to determine its position on a human's body (because if an elongated object magnetized along its longitudinal axis swings in sync with a part of the person-e.g. his arm or his leg-then our confidence in the magnet's position grows).
To start the process of zone inspection we need to know whether any people are there. For this, 3D cameras (based either on the time-of-flight principle, structured light approach or stereo-pair design) are usually employed. They not only give the system notice of the people's presence in the inspection zone but measure people's coordinates as well, which gives us the number of people present and approximate whereabouts of the magnets' positions in the zone. But due to the strong dependence of the magnet's magnetic field on the magnet's coordinates (1/r 3 dependence on distance) this kind of accuracy is clearly insufficient and subject to refinement.
Magnetic induction B produced by a magnet with the magnetic moment m in the point r relative to the magnet is given by a well-known formula from magnetostatics [22]: Here, and below, we also use ·, · for the standard dot product between two vectors. We see that as the dependence of B on the magnet's coordinates is strongly nonlinear, its dependence on m is linear and can be described by a linear operator A: Here, the vectors x and d are the magnet's and magnetometer's coordinates, respectively. This means that we deal with a linear system here: at least as long as the measurement of a magnetic moment m is concerned.
Having placed several magnetometers at the known points at the boundary of the inspection zone we can measure magnetic induction in those points of space. This means that we have a number of equations where d is the number of detectors. If we denote we can rewrite the system as a single multi-dimensional equation for a combined vector (x, m) ⊤ . We solve this system in a leastsquares sense, i.e. we find such a vector (x, m) ⊤ such that the value of the discrepancy is at a minimum. This problem is highly nonlinear in x due to nonlinearity of the magnetic induction dependence on the magnet's coordinates and is quadratic in m, which allows us to employ linear algebra methods to solve it in m and, hence, to decrease the dimensionality of the problem left to a nonlinear minimizer employed to find x. The solution of the m-problem (given the vector x) is well-known and is obtained via the Moore-Penrose pseudo-inverse matrix: which can be found using the standard QR-decomposition or singular value decomposition (SVD) techniques. Now, the problem with magnetic measurements is that we cannot shield our measuring system from outside magnetic interference and noise, of which there are many, in contrast, for example, with electric measurements (where a simple Faraday cage will do). Thus, the problem of fighting the noise is of utmost importance in this measurement scheme as well as in many other industrial applications. To deal with the situation we can employ the signal accumulation over time technique.
While a person under inspection passes through our system, the measuring device may make measurements of the person's coordinates and magnetic fields generated by it multiple times a second. Given a typical hiking speed of 1 m s −1 , a length of inspection zone of ∼3 m and a measuring frame rate of tens of frames per second, it gives us a number of measurements of the order of a hundred. If we measure the person's coordinates at each time frame via the optical device already mentioned, then we know a person's trajectory and, hence, a magnet's trajectory (provided a magnet is stationary in a person's reference frame, which is typically the case) up to an unknown translation vector and suffice it to find a position of a magnet at a single time frame only (e.g. the last frame)others will be found by applying the known trajectory shifts. In fact, we apply a type of matched filter here. This means that in this accumulation scheme, the dimensionality of the nonlinear minimization problem does not increase and it remains 3D. Figure 1 represents a schematic diagram of the analysis of a magnetic moment. Here, we use independent data for magnetic induction and the trajectory of a person. We emphasize that the minimization process is split into the explicitly solvable part related to m and nonlinear minimization in x. This reduces the number of parameters for nonlinear minimization to three.
The properties of the person's magnet do not change over time as the person passes through our system. If we make an assumption that this means that the vector m remains constant over time, we can combine equation (4) for various time frames in the large system of 3 dT equations where, as before, d is the number of detectors and T is the number of available time frames. The solution of this larger system remains exactly the same as before. We just need to tweak the A * at each time frame to account for the person's movement (which we know exactly via the optical subsystem). This is trivial because if y 1 , . . . , y T are the known positions of a fixed point on a person's body and x is a magnet's position at the last time frame, then the magnet's position at an arbitrary time frame is clearly Now, we notice that a person may put a dangerous item in his arm or his leg and swing or rotate it intentionally while on the move. This means that the m vector will not remain constant as the person passes through; only the vector's norm (which is a magnet's strength and, hence, the intrinsic property of the item) will remain invariant over time. Invariance of the m norm condition is additional information, which manifests itself in a number of additional equations which we want to account for. This brings us to the need to solve the signal accumulation problem in case we have constraints of the type of equations. As discussed in the introduction, the incorporation of any constraints into nonlinear minimization analysis is hard and, to our knowledge, there are no results for the particular case of implicit constraints of the form m = C where the constant C is unknown. The major advantage of our approach is that the constraint can now be added to the minimization in m part of the analysis (see figure 1) and still allow for the explicit solution followed by the standard nonlinear minimization in x.
One can argue that not only the direction but also the magnitude of the magnetic field is variable. It is definitely true in a general setting. Nevertheless, in arms detection applications the magnetic noise present is usually a small addition of magnetic fields from weak sources, such as magnetic fields from AC, flowing in wires inside modern buildings; radio waves emitted by various modern wireless devices; magnetic fields from power sources-electric washing machines, elevators, etc-all of which are situated outside of our inspection zone. We deem our inspection zone to be free of such sources of magnetic noise (which is usually the case), and this means that every interfering noise source (if any) is closer to magnetic detectors (which are located on the border of the inspection zone) than to the magnet of interest (which, by definition, is always inside the inspection zone). Therefore, the influence of noise on the magnet itself (which could in principle change its magnetic force and, hence, the value of m ), is much weaker than its influence on the detectors and can be neglected.
Moreover, the typical small-arms item is a ferromagnet weighing several pounds. Even for a material with relatively low coercivity (e.g. 0.16 kA m −1 for iron) it would take an enormous magnetic field to change its magnetic moment vector norm-or a very long time to do so via a realistic magnetic field. Both situations cannot occur on a timescale of several seconds while a person passes the inspection zone in a modern public place (e.g. an airport or a shopping mall). Significant re-magnetization of a weapon would require quite a substantial source of the magnetic field inside the inspection zone, and as such this source itself would be considered an alarming object for the system. This makes the idealization of the constant magnetic norm reasonable for the problem and this is why we concentrate on combating the measurement noise, assuming the estimated quantity-the m value-to be (unknown) constant.
Recalling the formulation of the previous problem, we reformulate this system in terms of the minimum search: Find the minimum of the discrepancy where the matrix A and vectors m and b are given by given that This is a problem for a conditional extremum, which is solved via the Lagrange multipliers method. First, let us deal with the constraints (7). In these equations we have T terms, so that we have T − 1 equality signs (i.e. equations), and we can write them as: or, which is the same, If we introduce the trivial notation for the scalar function of the vector argument for the multi-index and for the projector operator depending on the multi-index and if, as well, for any integer t denote then our constraints are given by the relations for t = 1, . . . , T − 1. Therefore, our problem is reformulated as Find the minimum of the function Using our notation we transform the Lagrangian as: As always, setting partial derivatives in λ 1 , λ 2 , . . . , λ T−1 to zero gives us the constraining equation (9) (in these equations we will substitute the solution m(λ 1 , λ 2 , . . . , λ T−1 ) found with the aim of determining the values of the Lagrange multipliers). Now we will deal with finding the solution vector m.

Solving the Lagrange problem
Let us calculate the gradient 1 of the Lagrangian (10). Making use of a well-known relation we deduce that a gradient of a scalar function f(x) = x 2 is given by where as ∇x we denote a matrix made by row-gradients of each of the vector's x components (i.e. in essence, this is a Jacobian of the x vector function). Thus If we transpose the equation and set it to zero we arrive at: Using the block-diagonal structure of the A matrix we write down these 3 T equations in terms of 3-components of the m vector, i.e. in terms of m t vectors: Note that this system of equations has an asymmetry between the first T − 1 equations and the last one. To symmetrize it we introduce a new miscellaneous Lagrange multiplier λ T : Then, the last equation takes the form analogous to the rest, but we get an additional relation, which follows from (11): This is the main system of equations, which is subject to solution. Note that equation (12) looks independent for different m t in the same way as this was before, in a single-frame method, i.e. it is as if we have 'uncoupled' them. However, the coupling between them remains because Lagrange multipliers λ 1 , . . . , λ T are bound by the last equation (12). Moreover, the T − 1 constraining equation (9) should hold as well.

Problems of minimization
Note that the scheme of the solution remains the same: we find the vector m as a function of the Lagrange multiplier's vector λ = (λ 1 , . . . , λ T ) ⊤ (where λ T is a function of other multipliers) and substitute it into the constraining equation (9) to determine λ 1 , . . . , λ T−1 (and, hence, λ T ). However, this is a difficult task as these equations cannot be solved analytically. Therefore, we present a way to avoid solving it.
Constraining equation (9) means just that the vectors' m t norm value does not depend on t (which is, in fact, our ultimate goal), but the value itself is not specified. If, however, we could take an arbitrary value s for the norm squared and solve every equation (12) (i.e. select the corresponding value of the Lagrange multiplier λ t (s)) against m t in such a way that m t (λ t (s)) 2 = s, then equation (9) would hold automatically! Admittedly, in this case the validity of the last equation of system (12), i.e.
is not guaranteed. But the function Λ = Λ(s) is a scalar function of the scalar argument s, and our task would come down to numerically finding its root only, which is not hard to achieve via any standard root-finding algorithm. This is a short description of the approach. Now, we discuss it in more detail and prove the existence and uniqueness theorem for the solution.

Existence and uniqueness theorem
Let us have a look at the first T equation (12). Take any of them: ( We see that the vector m t (λ) for t = 1, . . . , T (and its norm together with it, of course) depends in fact only on a single component of the λ vector, namely, the λ t ; therefore, we can write the solution in the form: Denote by 0 < µ 1 t ⩽ µ 2 t ⩽ µ 3 t the eigenvalues of the positive semi-definite matrix A ⊤ t A t and by ψ j t corresponding orthonormal eigenvectors. Then we have It is easy to see that the smooth function m t (λ t ) 2 tends to 0 as λ t → −∞ and is monotonically growing to +∞ as λ t approaches some ν t > 0 from the left. Generically, ν t = µ 1 Theorem 5.1. The Lagrange system (9) and (12) admits the unique solution with the value of λ in the region of regularity.
Proof. Due to smoothness and monotonicity, for any s > 0 there is a unique value of λ(s) in the region of regularity such that Moreover, for the function Obviously, Λ(s) is smooth and monotonically growing. Thus, there is the unique value s * > 0 such that Λ(s * ) = 0. The value λ(s * ) and the corresponding vectors (13) solve our Lagrange system (9) and (12). The solution is unique for λ in the region of regularity.

Numerical results
To illustrate the proposed technique a numerical simulation has been carried out. A single magnet with the magnetic moment's m 0 (t) norm of 1.3 A·m 2 (and, correspondingly, its norm squared that is equal to s 0 = 1.69) is being moved along a straight line during 50 frames. The m 0 vector is being swung sideways with an amplitude of 30 • in a plane, orthogonal to the trajectory. The magnet's magnetic induction is being measured by 20 3D magnetometers located around the inspection zone, whose width is about 1 m. A random Gaussian-distributed 3D noise vector with zero mean and a std of 10 nT is added to each magnetometer's measurement (noise vectors for different magnetometers are independent).
For a typical value of a magnetic moment of 1 A·m 2 the value of the magnetic induction B at a distance of 1 m from the magnet along the magnet's polar line is 200 nT; therefore, our noise value of 10 nT gives an SNR of the order of 20 or 26 dB, which is a realistic value in applications. If we decrease the distance twofold then the magnetic induction will increase eightfold, to 1600 nT, and vice versa; if we take twice the distance the induction will fall to 25.5 nT. Thus, the B varies as the magnet being worn moves along the trajectory; therefore, the SNR rises at trajectory points which are near to a given magnetometer and falls in those points which are far away. Thus, the notion of the SNR seems to be slightly superficial here.
From our experience we can say that typically benign objects (such as notebooks, magnetic clips, earphones, smartphones, etc) possess a magnetic moment value of less than 0.5 A·m 2 , and hazardous objects possess a value of over 1 A·m 2 in general; therefore, we adopt the values of 0.3 A·m 2 (benign object) and 1.3 A·m 2 (hazardous object) for our numerical example. To complete the picture we also include the results for relatively large noise with a std of 50 nT.
We process the measurement data and calculate the resulting value of m t 2 for each time frame and show the estimated square of the magnetic moment's norm in relative units (i.e. divided by the real value of s 0 ; see figure 2). We do it in several ways: • First, we apply the traditional scheme by calculating the m t vector (and, hence, its norm) independently at each frame according to equation (5). The resulting values of m t 2 /s 0 (shown as triangles) are not constant and change chaotically; • Second, we apply our technique for finding the s * value by solving system (9) and calculate the corresponding m t vectors at each time frame. These vectors turn out to have the same value of m t 2 = s * , which is shown as circles in the figure; • Finally, to emphasize the importance of the relaxed condition that m t is constant versus the stricter assumption that the vector m t itself is constant, we included the results for the fixed vector method. We accumulate all the measurements and calculate a single magnetic moment vector m, presuming it does not change over the trajectory. The resulting constant value of the (incorrect) squared magnetic moment norm, divided by s 0 , is represented by squares.
We see that as the single-frame estimation method provides the non-constant vector norm, deviating wildly from the real value, the presented method provides far better results: the norm of the resulting vectors remains constant and is much nearer to the real value. The fixed vector method provides norm estimation, which significantly underestimates the value as it anticipates the direction of the vector to be constant while in reality it is not. Figure 3 shows the same quantities for a seemingly benign object with s 0 = 0.09. Here, the fixed vector technique again provides the underestimated value of s (see the squares in   figure 3) seem to be more intense due to the lower value of the magnetic moment. Figures 4 and 5 show the graphs of the relative squared norm for the same object, but in the presence of the much more severe noise of 50 nT instead of 10 nT. Here, by the dotted line we also included the mean value of s for the one-frame solutions (represented by triangles). We see that in the case of severe noise, the simple mean tends to overestimate the true value of s, which is detrimental from the point of view of the false alarm rate.
In general, calculation of the mean value of s for a oneframe solution can be performed as follows. We denote by m 0 the true solution of the noiseless system Am 0 = b 0 . Then, in the presence of the noise δb we get m = A −1 (b 0 + δb). We assume that the noise δb has a normal symmetric distribution with zero mean and variance of each component equal to σ 2 . Then, we have E( m 2 ) = m 0 2 + E( A −1 δb 2 ).  Let ν j > 0, j = 1, 2, 3, be the eigenvalues of the positive definite matrix A. Then, by the properties of the normal distributions we get This calculation shows that such a mean value overestimates the real value of s 0 . The effect is stronger for larger noise or a relatively smaller signal, which may lead to a false alarm. The same effect can be observed if one uses a suboptimal measurement configuration, where each trajectory point's aperture is small (see the contribution from the eigenvalues of the matrix A). In practice, false alarming on benign objects is equally undesirable behavior as is missing the signal from a hazardous object. Our method shows itself to be better suited for the task and much less overestimates the real value of interest. Figure 6 shows the graph of the square of the distance between the real vector m 0 (t), which was put into the model,  and the calculated vector m t in relative units for each time frame, i.e. the ( m 0 (t) − m t 2 )/s 0 for s 0 = 1.69 (a presumably hazardous object). As before, the triangles represent the standard one-frame method and the circles are the new method introduced in the present paper. Even though the main information to be determined in the current setting is the norm of the magnetic moment, it is still important to note for further considerations that our method provides consistently better approximation of the true magnetic moment vector as well. Figure 7 shows the graph for the same object, but for a higher noise level of 50 nT. Figures 8 and 9 show the graph for the benign object with s 0 = 0.09 for the noise levels of 10 nT and 50 nT, respectively.

Conclusion
In the present manuscript we introduced a new approach to the problem of accumulation of the signal in the presence of an implicit constraint of the form x = C, where the constant C is also to be determined. Theoretical analysis as well as the theorem of existence and uniqueness of the solution are given in the setting of the real-life magnetometry problem. Numerical simulations show that our method is superior compared to both the classical approach of a single-frame capturing when there are no constraints in the problem and the stricter approach when the magnetic moment vector itself is assumed to be fixed. In addition, our approach is shown to be more robust in the sense that it is less likely to trigger a false alarm from benign objects when relatively high noise is present. Finally, we believe that the method can be widely applied to other problems where implicit constraint is the core part of the construction.

Data availability statement
All data that support the findings of this study are included within the article (and any supplementary files).