Patent Publication Number: US-2016245655-A1

Title: Adjusted navigation

Description:
BACKGROUND OF THE INVENTION 
     The present invention relates to adjusted navigation. More in particular, the present invention relates to a method and a device for adjusting the navigation state determined by navigation systems, especially but not exclusively inertial navigation systems. The present invention particularly relates to a method and a device for providing adjusted or corrected navigation states. 
     Navigation systems are designed to determine their position, and optionally other variables such as velocity, at a certain point in time. Systems such as GPS (Global Positioning System) determine their position with the aid of external references, that is, the locations of satellites. Inertial navigation systems use inertia referenced sensors like accelerometers and gyroscopes, also known as inertial sensors, to determine movements relative to a starting position. The output of an inertial navigation system (such as position, velocity etc.) is determined by numerically integrating in time the inertial sensor output. As inertial sensors always introduce errors, however small, the position determined by the system tends to increasingly deviate from the real position as time progresses. For this reason, the navigation states produced by inertial navigation systems need to be adjusted from time to time. 
     U.S. Pat. No. 8,700,326 and its European counterpart EP 2504664 (TNO), which are herewith incorporated in this document, disclose a navigation system provided with a magnetic field sensor unit arranged for providing a magnetic field signal indicative of a magnetic field strength. This known system comprises an integration unit for calculating a path-integral of the magnetic field signal so as to produce a value which is related to its position. The known system is arranged to exchange data with a similar system when they are at approximately the same location, that is, when they encounter each other. The exchanged data allow the systems to determine their locations and to correct navigation errors. The system of 
     U.S. Pat. No. 8,700,326 is very useful but requires at least two similar systems to encounter each other in order to produce a position update. 
     United States Patent Application US 2013/0179075 (Indooratlas) discloses an apparatus for determining a position inside a building. This known device first uses a non-magnetic field based location system to determine a location estimate and then uses an Earth magnetic field measurement and an Earth magnetic field map to further determine the location estimate. This known device has the disadvantage that it is less accurate in situations where the magnetic field is highly uniform. Further, the location estimate of this known device is highly correlated to the location estimate of the non-magnetic field based location system. Further use of this location to correct Inertial navigation systems is therefore not possible. 
     SUMMARY OF THE INVENTION 
     It is an object of the present invention to overcome these and other problems of the Prior Art and to provide a navigation adjustment method and system which is capable of operating independently, that is, which does not require a second system to provide a navigation update but which still has a high accuracy. 
     Accordingly, the present invention provides a navigation device, comprising:
         a navigation state estimation unit for repeatedly providing a navigation state estimate, said navigation state comprising at least one of a position, a velocity and an attitude of the device,   a prediction unit for providing, using the navigation state estimate, a predicted magnetic property value,   a storage unit for retrieving, using the navigation state estimate, a previously stored magnetic property value,   a subtraction unit for producing the difference of the predicted magnetic property value and the retrieved magnetic property value,   a filter unit for producing, using said difference, a navigation state correction value, and   a correction unit for producing, using said navigation state correction value, a corrected navigation state estimate,       

     wherein the prediction unit comprises:
         a magnetic measuring unit for repeatedly providing a measured magnetic field strength value, and   an integrating unit for integrating over distance the measured magnetic field strength values to produce the magnetic property value corresponding with the navigation state estimate.
 
By providing a prediction unit for providing, using the navigation state estimate, a predicted magnetic property value, and by further providing a storage unit for retrieving, using the navigation state estimate, a previously stored magnetic property value, it is possible to compare the predicted and the retrieved magnetic property values corresponding with the estimated navigation state. The comparison preferably results in a difference value which is in turn used to correct or adjust the navigation state estimate.
       

     The filter unit preferably comprises a Kalman filter. Such filters, which are well known per se, are very suitable to produce navigation state correction values, based on a single difference value or a multiple of difference values, in the case a multiple of sensors are used. It is for instance possible to use GPS position measurements (if available) in combination with a barometer height measurement and magnetic field measurements, or any other navigation related sensor measurement, at the same time. This seamless integration with other sensors is one of the advantages of the method. 
     The Kalman filter is preferably updated when producing each navigation state correction value. However, in some embodiments the Kalman filter is not updated if the difference of the predicted magnetic property value and the retrieved magnetic property value substantially equals zero. That is, when the input of the Kalman filter is approximately zero, the filter is preferably no longer updated to preserve its existing state. In these cases the Kalman filter merely performs a prediction calculation (coasting). 
     It is further preferred that the Kalman filter is not updated if the difference of the predicted magnetic property value and the retrieved magnetic property value exceeds a threshold value. That is, if the difference is too large, it is determined that a navigation error has been made that cannot be corrected. Exceeding this threshold may advantageously be recorded and/or be brought to the attention of a user. 
     In a first embodiment, the difference of the predicted magnetic property value and the retrieved magnetic property value is given by 
       {tilde over (ψ)}−ψ table ( p{tilde over ( p )}   end )=(   H     avg   − H     end )·   {circumflex over (p)} , 
 
     where H avg  is the average 
     magnetic field strength during a measurement trajectory,  H   end  is the magnetic field strength at the end of the measurement trajectory, and  {circumflex over (p)}  is the estimated position error. In this embodiment, the difference value is determined using the average and the final magnetic field strength, and the navigation state error. 
     In a second embodiment, the difference of the predicted magnetic property value and the retrieved magnetic property value is given by 
       {tilde over (ψ)} ψ t ( {tilde over ( p )}   e )=   {circumflex over (V)}   ·(∫ 0   T     H ·dt )| {circumflex over (α)} ·(∫ 0   T     H ·t·dt )    {circumflex over (p)}   ·( H   end ).
 
     where {tilde over ({circumflex over (V)})} is the estimated velocity error,  H  is the magnetic field strength during a measurement trajectory,  H   end   is the magnetic field strength at the end of the measurement trajectory,  {circumflex over (α)}  is the acceleration error, and  {circumflex over (p)}  is the position error. In this embodiment, the constant velocity error  {circumflex over (V)}  taken into account. 
     Although the navigation state estimation unit is preferably configured for inertial navigation (INS), it may alternatively, or additionally, be configured for other dead-reckoning navigation techniques, such as odometry, step counting systems etc. 
     The present invention further provides a vehicle, provided with a navigation device as defined above. The vehicle may be a road vehicle, such as a car or a bicycle, but may also be an airplane, a helicopter or a ship. 
     The present invention additionally provides a set of vehicles arranged for, preferably wireless, data exchange, the vehicles together comprising all constituent units of the navigation device defined above, a first vehicle comprising the magnetic measuring unit and a second vehicle comprising the navigation state estimation unit. In such a set of vehicles, the navigation device is distributed over two or more vehicles, the relatively expensive magnetic measuring unit being present in the first vehicle only. The set may comprise more than one second vehicle. 
     In an embodiment of such a set of vehicles the magnetic measuring unit of said first vehicle is part of a first navigation device comprised in the first vehicle, and the navigation state estimation unit of the second vehicle is part of a second navigation device, according to one of the embodiment as specified above, comprised in the second vehicle, wherein said magnetic measuring unit of said first vehicle determines magnetic property values for storage in said storage unit of said second navigation device, to be used by the subtraction unit of the second navigation device for producing the difference of the predicted magnetic property value and the retrieved magnetic property value. In this embodiment the navigation system of the first vehicle prepares the magnetic property values to be used as reference data by the navigation system of the at least one second vehicle. Provided that this reference data is sufficiently accurate, the navigation state estimation unit of the second navigation devices can be less accurate, as the correction unit can sufficiently correct the navigation state estimate using the reference data. Therewith the set of vehicles only needs a first vehicle with a high quality navigation state estimation unit while a plurality of second vehicles can be used that do not need a high quality navigation state estimation unit, therewith providing for a high overall quality at relatively low costs. 
     The previously stored magnetic property values may be obtained by moving a suitable apparatus, preferably using a vehicle, and accurately determining its navigation state (including its position) a large number of times, for example by using GPS, triangulation involving radar and/or infra-red light, and/or other means. Such apparatus may comprise a magnetic measuring unit for repeatedly providing a measured magnetic field strength value, and an integrating unit for integrating the measured magnetic field strength values over distance to produce the magnetic property value corresponding with the determined navigation state. The magnetic property values thus found may be transferred to a suitable storage unit, such as a ROM (Read Only Memory) unit from which they can be retrieved when required. 
     The present invention also provides a navigation method comprising the steps of:
         repeatedly providing a navigation state estimate,   providing, using the navigation state estimate, a predicted magnetic property value,   retrieving, using the navigation state estimate, a previously stored magnetic property value,   producing the difference of the predicted magnetic property value and the retrieved magnetic property value,   producing, using said difference, a navigation state correction value, and   correcting, using said navigation state correction value, a corrected navigation state estimate,       

     wherein the step of providing a predicted magnetic property value comprises:
         repeatedly providing a measured magnetic field strength value, and   integrating the measured magnetic field strength values over distance to produce the magnetic property value corresponding with the navigation state estimate.       

     The step of producing a navigation state correction value may advantageously comprise Kalman filtering, which Kalman filtering preferably comprises updating a filtering state when producing each navigation state correction value. The Kalman filtering state may not be updated if the difference of the predicted magnetic property value and the retrieved magnetic property value substantially equals zero. 
     It is preferred that the previously stored magnetic property values are obtained by repeatedly providing a measured magnetic field strength value, and integrating the measured magnetic field strength values over distance. 
     In an advantageous embodiment, the step of repeatedly providing a measured magnetic field strength value is carried out in a first vehicle, and the other steps are carried out in a second vehicle, said values being transmitted from the first to the second vehicle, and the second vehicle preferably following the first vehicle. In an embodiment, the first vehicle carries the relatively expensive apparatus required to measure the magnetic field strength while the second (and possible third, fourth, etc.) vehicle need only carry relatively inexpensive navigation apparatus. The magnetic property values may be transmitted from the first to the second vehicle wirelessly, e.g. using a mobile telecommunication network. 
     The present invention additionally provides a computer program product for carrying out the method as defined above. A computer program product may comprise a set of computer executable instructions stored on a tangible data carrier, such as a CD, a DVD or a memory stick. The set of computer executable instructions, which allow a programmable computer to carry out the method as defined above, may also be available for downloading from a remote server, for example via the Internet. 
     Reference is made to the book by Paul D. Groves: “Principles of GNSS, Inertial, and Multisensor Integrated Navigation Systems”, ISBN-13: 978-1-58053-255-6, 2007, which explains navigation systems, including Kalman filters, in more detail. 
    
    
     
       BRIEF DESCRIPTION OF THE DRAWINGS 
       The present invention will further be explained below with reference to exemplary embodiments illustrated in the accompanying drawings, in which: 
         FIG. 1  schematically shows navigation paths in accordance with the present invention. 
         FIG. 2  schematically shows a navigation device in accordance with the present invention. 
         FIG. 3  schematically shows a prediction unit of the navigation device illustrated in  FIG. 2 . 
     
    
    
     DETAILED DESCRIPTION OF EMBODIMENTS 
     The present invention is inspired on Ampere&#39;s circuital law, which states: 
     
       
      
       
       C 
         H ·d l =l 
      
     
     in which: 
       H  is the field strength of the local static magnetic field, 
     C is a closed path over which the integration is done, 
     d l  is a small element of the closed path C, 
     I is the net electric DC current flowing through the surface surrounded by the closed path. 
     In practical indoor/outdoor situations, usually no DC current flows in the space where navigation takes place. AC currents (power lines) do flow, but these do not pose a problem as their effects in general average out in time and further, power lines usually go in pairs (such that the net AC current is zero). If we apply Ampere&#39;s law to closed paths in such a space without net DC currents it then simplifies to: 
     
       
      
       
       C 
         H ·d l =O 
      
     
     In U.S. Pat. No. 8,700,326 a method is disclosed how the accelerometer biases of multiple INS (Inertial Navigation System) can be estimated by partial evaluation of the integral in Ampere&#39;s law. This is done by detecting encounters between the navigating objects. This made the invention in general only applicable to systems with a plurality of navigating objects and depending on encounters between them. 
     The present invention solves this problem, making it applicable to a single navigating system. The concept is called “Amperian Navigation Correction” (ANC), referring to Ampere&#39;s law for the magneto-static field. 
     Notation 
     Bars are used to indicate vectors (example:  ν ) 
     A tilde is used to further indicate estimated values (example: {tilde over ( ν )}) 
     A hat is used to indicate the error of a value (example: {circumflex over ( ν )}) 
     Vector Subscripts &amp; Superscripts: 
     Subscript b: used to denote body referenced vectors (example:  ν   b ) 
     Subscript e: used to denote Earth referenced vectors (example: ν e ) 
     Unless explicitly denoted otherwise, vectors are always in the Earth reference system. 
     Superscript “table”: used to indicate a “table look-up” value 
     Superscript “end”: used to indicate a final or last measured value 
     Superscript m: used to denote a measured value 
     Magnetic Field B and H 
     Historically the “magnetic induction” or “magnetic flux density” is denoted by the letter B and the “magnetic field” is denoted by the letter H. Nowadays, the name “magnetic field” is commonly used for both the B-field and the H-field. In this document this modern naming convention will be followed, using the name “magnetic field” for both the B-field and H-field. 
     In the absence of DC currents in a defined space where navigation takes place, Ampere&#39;s law reduces to    C   H ·d l =O for any closed path in that space. For this space, an unambiguous scalar potential field is defined as a function of the position  p . The scalar potential is found by integrating over any arbitrary path C that ends in  p . The origin position where integration starts is Ō. The scalar potential defined for the origin position is arbitrary and may be set to zero. 
       ψ( p )−∫ Ō     p     H ·d l +ψ(Ō)−∫ Ō     p     H ·d l 
 
     So as long as we define an origin position Ō and a scalar potential ψ(Ō) for that position, the scalar potential function is completely defined given a fixed magnetic field. 
     For the method described in this document, it is required that the magnetic field is known at locations where navigation can take place. For places where navigation will not take place, the magnetic field need not be known. Based on the known magnetic field, this scalar potential ψ( p ) can be calculated a priori. The calculated map of the scalar potential can be stored in a (remote) database or locally in a storage unit of the navigating object. 
     Reference is made to  FIG. 1 . Suppose an object is navigating over some true path S, starting in the origin Ō and taking up a travelling time T. At the end of the path, the scalar potential ψ table ({tilde over ( p )} end ) is determined by looking it up in a table, indexed by the position {tilde over ( p )} end , which is estimated by the INS. {tilde over ( p )} end  deviates from the true position  p   end  by an error  {circumflex over (p)}   end . 
     At the true position  p   end  the magnetic field strength is measured  B   b   end,meas  using magnetometers. From this measurement we can derive  H   b   end,meas  using the relation 
         B   b   end,meas =μ H   b   end,meas ,   )1)
 
     with μ the permeability of the materials of the sensor and surrounding the sensor. It is assumed that non-magnetic materials are used and that navigation takes place in a non-magnetic medium (such as air) such that μ=μ 0 , the permeability of free space. The measurement is done in the body axis system (or also called sensor axis system) and rotated to the Earth (or navigation) axis system using the body to Earth rotation matrix of the INS. 
           {tilde over (H)}     end   = R     be   · H     b   end,meas    (2)
 
     It is assumed that in the small volume around  p   end  the magnetic field strength is uniform. This volume is indicated by the shaded area in  FIG. 1 . If the quality of the IMU (Inertial Measurement Unit) or similar navigation unit is sufficient, this condition will automatically be met. The higher the IMU quality, the lower the expected error after some time. For some applications, the homogeneity of the field can be “designed” such that this condition is met with respect to the expected error after some time. The Earth magnetic field can for instance be made non-homogeneous by placing magnets or magnetic (iron) materials in or near the space where navigation takes place, or DC current sources near the space where navigation takes place. 
     It is further assumed that the attitude errors of the INS are very small such that the following approximation is valid:  H   end ≈ {tilde over (H)}   end . Then, the scalar potential at the INS position is equal to the scalar potential at the true position plus the contribution due to the INS position error: 
       ψ table ( {tilde over ( p )}   end )=ψ end   + H     end   · {circumflex over (p)}     (3)
 
     During the travelling time T, the scalar potential is also evaluated based on the INS readings. It is assumed that the sampling frequency F s  used by the INS is sufficiently high and the corresponding sampling interval T s  is sufficiently small. The numerical integration step derived from the INS readings can then be approximated by an infinitesimal integration step: Δ {tilde over (x)} ≅d {tilde over (x)} . The evaluation of the scalar potential then becomes: 
       {tilde over (ψ)}=∫ S     {tilde over (H)} ·d {tilde over (x)}   
 
     And the estimated magnetic field is the measured magnetic field rotated to the navigation frame by the estimated body to Earth rotation matrix (equation 2): 
     
       
      
         {tilde over (H)} ={tilde over (R)} 
       be 
       · H   
       b 
       meas 
      
     
     The integration step as derived from the INS readings is the true integration step plus some small error: 
     
       
      
       d{circumflex over (x)}=dx+d{circumflex over (x)} 
      
     
     This integration step error depends on the INS velocity error. As  {tilde over (H)} ≈ H , the evaluated scalar potential becomes: 
       {tilde over (ψ)}=∫ S     {tilde over (H)} ·d {tilde over (x)} ≈∫   S     H ·d x +∫   S     H ·d {circumflex over (x)}   
 
     The first part of the right hand side is the evaluation of the integral over the true path S and the true magnetic field strength and is identical to ψ end . 
       {tilde over (ψ)}=ψ end |∫ S     H ·d {circumflex over (x)}     (4)
 
     The rightmost part of the right hand side is the deviation of ψ end  due to the INS errors. d {circumflex over (x)}  is an error that only depends on the velocity error and changes during the time interval T. 
     In the absence of aiding information, d {circumflex over (x)}  changes due to a linear combination of two primary effects. Each of these effects will be analyzed individually: 
     Effect 1: a constant velocity error. Due to constant velocity errors (the starting value) the position error grows linearly. For this effect, d {circumflex over (x)}  can be written as: 
         d {circumflex over (x)} = {circumflex over (V)} ·dt.   (5)
 
     We use a capital  {circumflex over (V)} to distinguish it from the INS velocity error  {circumflex over (V)} . The INS velocity error is not constant but when other aiding information stops, it starts at  {circumflex over (V)} and then grows due to the constant attitude errors. 
     Effect 2: a constant attitude error. Constant attitude errors result in constant acceleration errors via gravity compensation. Due to constant acceleration errors, velocity errors are growing linearly. This has an exponential effect on the position error. For this effect, d {circumflex over (x)}  can be written as: 
         d {circumflex over (x)} = {circumflex over (α)} ·t·dt   (6)
 
     The combined effect is: 
         d{circumflex over (x)} =(   {circumflex over (V)} +{circumflex over (α)}·t )  (7)
 
     Substituting in equation 4: 
       {tilde over (ψ)}=ψ end +∫ O   T     H   ·(   {circumflex over (V)} + {circumflex over (α)} ·t )· dt 
 
       {tilde over (ψ)}=ψ end +{∫ O   T     H   {circumflex over (V)} ·dt }+{∫ O   T     H · {circumflex over (α)} ·t·dt }  (8)
 
     In which the two different effects are expressed between the two sets of braces. Because of the substitution, integration is now done over time instead of distance. T is the total integration time. 
     It is noted that the number of effects can in principle be extended to a third effect, the gyro bias errors  b   g , resulting in: 
         d {circumflex over (x)} = b     g ·½ t   2   ·dt 
 
     Constant gyro bias errors are however assumed to be compensated by the technology disclosed in US2014005975, which document is herewith incorporated in this document. 
     It is further noted that a constant INS position error cannot be resolved using this method as d {circumflex over (x)} =O, so  {circumflex over (p)}   end  is the position error growth during the interval T. 
     When analyzing effect 1, the first effect of a constant velocity error on the scalar potential is considered. Substituting equation 5 in equation 4, the evaluated scalar potential becomes: 
     
       
         
           
             
               
                 
                   
                     
                       
                         
                           ψ 
                           ~ 
                         
                         = 
                           
                          
                         
                           
                             ψ 
                             end 
                           
                           + 
                           
                             
                               ∫ 
                               0 
                               T 
                             
                              
                             
                               
                                 H 
                                 _ 
                               
                               · 
                               
                                 
                                   V 
                                   _ 
                                 
                                 ^ 
                               
                                
                               
                                   
                               
                               · 
                               
                                  
                                 t 
                               
                             
                           
                         
                       
                     
                   
                   
                     
                       
                         
                           = 
                             
                            
                           
                             
                               ψ 
                               end 
                             
                             + 
                             
                               
                                 V 
                                 _ 
                               
                               ^ 
                             
                           
                         
                          
                         
                           · 
                           
                             
                               ∫ 
                               0 
                               T 
                             
                              
                             
                               H 
                               _ 
                             
                           
                         
                          
                         
                           · 
                           
                               
                           
                            
                           
                              
                             t 
                           
                         
                       
                     
                   
                   
                     
                       
                         = 
                           
                          
                         
                           
                             ψ 
                             end 
                           
                           + 
                           
                             
                               
                                 
                                   V 
                                   _ 
                                 
                                 ^ 
                               
                               · 
                               T 
                               · 
                               
                                 1 
                                 T 
                               
                             
                              
                             
                               
                                 ∫ 
                                 0 
                                 T 
                               
                                
                               
                                 
                                   H 
                                   _ 
                                 
                                 · 
                                 
                                     
                                 
                                  
                                 
                                    
                                   t 
                                 
                               
                             
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   9 
                   ) 
                 
               
             
           
         
       
     
     Again T is the duration of the time interval.  {circumflex over (V)} ·T is equal to the INS error growth  {circumflex over (p)}  if constant velocity errors are the only cause of the position error growth: 
         {circumflex over (p)} − {circumflex over (V)} ·T   (10)
 
     Further, 
     
       
         
           
             
               1 
               T 
             
              
             
               
                 ∫ 
                 0 
                 T 
               
                
               
                 
                   H 
                   _ 
                 
                 · 
                 
                     
                 
                  
                 
                    
                   t 
                 
               
             
           
         
       
     
     is the time-averaged value of the magnetic field strength during the interval T. 
     
       
         
           
             
               
                 H 
                 _ 
               
               avg 
             
             = 
             
               
                 1 
                 T 
               
                
               
                 
                   ∫ 
                   0 
                   T 
                 
                  
                 
                   
                     H 
                     _ 
                   
                   · 
                   
                       
                   
                    
                   
                      
                     t 
                   
                 
               
             
           
         
       
     
     Substituting these in equation 9, the equation for the evaluated scalar potential becomes: 
       {tilde over (ψ)}=ψ end +   {circumflex over (p)} · H     avg 
 
     Subtracting equation 3 from this equation: 
       {tilde over (ψ)}−ψ table ( {tilde over ( p )}   end )=(ψ end   + {circumflex over (p)} · H     avg )−(ψ end   + H     end   · {circumflex over (p)} · H     avg   − H     end 
 
       {tilde over (ψ)}−ψ table ( {tilde over ( p )}   end )=(   H     avg   − H     end )·   {circumflex over (p)}     (11)
 
     This results in a measurement equation for the Kalman filter, that relates the INS error  {dot over (p)} , given a constant velocity error, to a measurement  r ={tilde over (ψ)}−ψ table ({tilde over ( p )} end ), and can be implemented in a standard navigation filter structure as shown in  FIG. 2 . The term between the parenthesis can be evaluated during navigation and is therefore known. 
     As an alternative, it is also possible to directly relate the constant velocity error and the INS position error to the measured scalar potential difference. Subtracting equation 3 from equation 9, results in: 
       {tilde over (ψ)}−ψ table ( {tilde over (p)}   end )=(ψ end   +{circumflex over (V)}·∫   O   T   H·dt )−(ψ end   +H   end   · {tilde over (p)}   )
 
       {tilde over (ψ)} ψ table ( {tilde over ( p )}   end )=   {circumflex over (V)} ·∫   O   T     H ·dt  H     end   · {circumflex over (p)}     (12)
 
     Again, a measurement equation is obtained. ∫ 0   T   H ·dt and  H   end  are both known and can be calculated during navigation. 
     When analyzing effect 2, a comparable equation is obtained for constant acceleration errors (due to attitude). Substituting equation 6 in equation 4 results in: 
       {tilde over (ψ)}=ψ end +∫ O   T     H · {circumflex over (α)} ·t·dt=ψ   end   + {circumflex over (α)} ·∫   O   T     H ·t·dt   (13)
 
     For constant acceleration errors, the effect on the position error is: 
         {circumflex over (p)} = 2   1   {circumflex over (α)} T 2 , or  {circumflex over (α)} = {circumflex over (p)}   T     2     2   (14)
 
     Substituting this in equation 13 results in: 
     
       
         
           
             
               
                 
                   
                     ψ 
                     ~ 
                   
                   = 
                   
                     
                       ψ 
                       end 
                     
                     + 
                     
                       
                         
                           
                             p 
                             _ 
                           
                           ^ 
                         
                         · 
                         
                           2 
                           
                             T 
                             2 
                           
                         
                       
                        
                       
                         
                           ∫ 
                           0 
                           T 
                         
                          
                         
                           
                             H 
                             _ 
                           
                           · 
                           t 
                           · 
                           
                              
                             t 
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   15 
                   ) 
                 
               
             
           
         
       
     
     Subtracting equation 3 from this: 
     
       
         
           
             
               
                 
                   
                     
                       
                         ψ 
                         ~ 
                       
                       - 
                       
                         
                           ψ 
                           table 
                         
                          
                         
                           ( 
                           
                             
                               
                                 p 
                                 ~ 
                               
                               _ 
                             
                             end 
                           
                           ) 
                         
                       
                     
                     = 
                     
                       
                         ( 
                         
                           
                             ψ 
                             end 
                           
                           + 
                           
                             
                               
                                 
                                   p 
                                   ~ 
                                 
                                 _ 
                               
                               · 
                               
                                 2 
                                 
                                   T 
                                   2 
                                 
                               
                             
                              
                             
                               
                                 ∫ 
                                 0 
                                 T 
                               
                                
                               
                                 
                                   H 
                                   _ 
                                 
                                 · 
                                 t 
                                 · 
                                 
                                     
                                 
                                  
                                 
                                    
                                   t 
                                 
                               
                             
                           
                         
                         ) 
                       
                       - 
                       
                         ( 
                         
                           
                             ψ 
                             end 
                           
                           + 
                           
                             
                               
                                 H 
                                 _ 
                               
                               end 
                             
                             · 
                             
                               
                                 p 
                                 _ 
                               
                               ^ 
                             
                           
                         
                         ) 
                       
                     
                   
                    
                   
                     
 
                   
                    
                   
                       
                   
                    
                   
                     
                       
                         ψ 
                         ~ 
                       
                       - 
                       
                         
                           ψ 
                           table 
                         
                          
                         
                           ( 
                           
                             
                               
                                 p 
                                 ~ 
                               
                               _ 
                             
                             end 
                           
                           ) 
                         
                       
                     
                     = 
                     
                       
                         ( 
                         
                           
                             
                               2 
                               
                                 T 
                                 2 
                               
                             
                              
                             
                               
                                 ∫ 
                                 0 
                                 T 
                               
                                
                               
                                 H 
                                 · 
                                 t 
                                 · 
                                 
                                     
                                 
                                  
                                 
                                    
                                   t 
                                 
                               
                             
                           
                           - 
                           
                             
                               H 
                               _ 
                             
                             end 
                           
                         
                         ) 
                       
                       · 
                       
                         
                           p 
                           _ 
                         
                         ^ 
                       
                     
                   
                 
               
               
                 
                   ( 
                   16 
                   ) 
                 
               
             
           
         
       
     
     This results in a measurement equation for the Kalman filter that relates the INS error  {tilde over (p)} , given a constant acceleration error (or corresponding attitude error), to a measurement  r ={tilde over (ψ)}−ψ table ({tilde over ( p )} end ), and can be implemented in a fairly standard navigation filter structure as shown in  FIG. 2 . The term between the parenthesis can be evaluated during navigation and is therefore known. 
     For this effect an alternative measurement equation also exists. Subtracting equation 3 from equation 13 results in: 
       {tilde over (ψ)}−ψ table ( {tilde over ( p )}   end )=(ψ end + {circumflex over (α)} ·∫ O   T     H ·t·dt )−(ψ end   + H     end   · {circumflex over (p)}   )
 
       {tilde over (ψ)}−ψ table ( {tilde over ( p )}   end )=   {circumflex over (α)} ·∫   O   T     H ·dt  H ·t·dt− {circumflex over (p)} · H     end    (17)
 
     ∫ O   T   H ·t·dt and  H   end  are known and can be evaluated during navigation. 
     Although gyro bias errors are assumed to be compensated by the technology disclosed in US2014005975, which is herewith incorporated in this document in its entirety, its effect can also be observed. Only for completeness, it is briefly discussed here. For this third effect, equation 4 likewise becomes: 
     
       
         
           
             
               
                 
                   
                     ψ 
                     ~ 
                   
                   = 
                     
                    
                   
                     
                       ψ 
                       end 
                     
                     + 
                     
                       
                         ∫ 
                         0 
                         T 
                       
                        
                       
                         
                           
                             H 
                             _ 
                           
                           · 
                           
                             
                               b 
                               _ 
                             
                             g 
                           
                           · 
                           
                             1 
                             2 
                           
                         
                          
                         
                           
                             t 
                             2 
                           
                           · 
                           
                               
                           
                            
                           
                              
                             t 
                           
                         
                       
                     
                   
                 
               
             
             
               
                 
                   = 
                     
                    
                   
                     
                       ψ 
                       end 
                     
                     + 
                     
                       
                         1 
                         2 
                       
                        
                       
                         
                           
                             b 
                             _ 
                           
                           g 
                         
                         · 
                         
                           
                             ∫ 
                             0 
                             T 
                           
                            
                           
                             
                               H 
                               _ 
                             
                             · 
                             
                               t 
                               2 
                             
                             · 
                             
                                 
                             
                              
                             
                                
                               t 
                             
                           
                         
                       
                     
                   
                 
               
             
           
         
       
       
         
           
             
               Substituting 
                
               
                   
               
                
               
                 
                   p 
                   _ 
                 
                 ^ 
               
             
             = 
             
               
                 1 
                 c 
               
                
               
                 
                   b 
                   _ 
                 
                 g 
               
                
               
                 T 
                 3 
               
                
               
                 : 
               
             
           
         
       
       
         
           
             
               ψ 
               ~ 
             
             = 
             
               
                 ψ 
                 end 
               
               + 
               
                 
                   
                     
                       p 
                       _ 
                     
                     ^ 
                   
                   · 
                   
                     
                       o 
                       ^ 
                     
                     
                       T 
                       3 
                     
                   
                 
                  
                 
                   
                     ∫ 
                     0 
                     T 
                   
                    
                   
                     
                       H 
                       _ 
                     
                     · 
                     
                       t 
                       2 
                     
                     · 
                     
                         
                     
                      
                     
                        
                       t 
                     
                   
                 
               
             
           
         
       
     
     As mentioned before, d {circumflex over (x)} is a linear combination of these effects, resulting in a general expression for the measurement equation. It can be written as a combined measurement equation that directly relates an evaluated difference in scalar potential to the INS position error growth, a constant velocity term and a constant acceleration term. 
     Repeating equation 8: 
       {tilde over (ψ)}=ψ end +{∫ O   T     H · {circumflex over (V)} ·dt }+{∫ O   T     H · {circumflex over (α)} ·t·dt }
 
     Substract equation 3 and rearranging: 
       {tilde over (ψ)}−ψ t  ( {tilde over ( p )}   e )−   {circumflex over (V)}   ·(∫ O   T     H ·dt )+ û ·(∫ O   T     H ·t·dt )−   {circumflex over (p)}   ·(   H     end )  (18)
 
     This results in a measurement equation that relates the constant velocity error, constant INS attitude error and INS position error growth to an evaluated difference in scalar potential. The terms between the parenthesis are known and can be evaluated during navigation. This expression can likewise be extended to further include the gyro bias. 
     Further substituting equations 10 and 14 results in: 
     
       
         
           
             
               
                 
                   
                     
                       ψ 
                       ~ 
                     
                     - 
                     
                       
                         ψ 
                         t 
                       
                        
                       
                         ( 
                         
                           
                             
                               p 
                               ~ 
                             
                             _ 
                           
                           e 
                         
                         ) 
                       
                     
                   
                   = 
                   
                     
                       ( 
                       
                         
                           
                             
                               1 
                               τ 
                             
                              
                             
                               
                                 ∫ 
                                 0 
                                 T 
                               
                                
                               
                                 
                                   H 
                                   _ 
                                 
                                 · 
                                 
                                     
                                 
                                  
                                 
                                    
                                   t 
                                 
                               
                             
                           
                           + 
                           
                             
                               2 
                               
                                 τ 
                                 2 
                               
                             
                              
                             
                               
                                 ∫ 
                                 0 
                                 T 
                               
                                
                               
                                 H 
                                 _ 
                               
                             
                           
                         
                          
                         
                           
                             · 
                             t 
                             · 
                             
                                 
                             
                              
                             
                                
                               t 
                             
                           
                           - 
                           
                             
                               H 
                               _ 
                             
                             end 
                           
                         
                       
                       ) 
                     
                     · 
                     
                       
                         p 
                         _ 
                       
                       ^ 
                     
                   
                 
               
               
                 
                   ( 
                   19 
                   ) 
                 
               
             
           
         
       
     
     This results in a measurement equation that relates an evaluated difference in scalar potential to the INS position error growth, given all the combined effects. The term between the parenthesis is known and can be evaluated during navigation. This expression can also be extended to further include the gyro bias. 
     In summary, a general expression (equation 19) is found for a measurement equation that relates a difference in scalar potential to the INS position error growth. This general expression has the form: 
       Δ ψ   = {circumflex over (p)} ·C   0 ( T )
 
     in which Δψ is the evaluated difference in scalar potential. C 0 (T), is a constant, given an integration time T. 
     An alternative general expression (equation 18) for a measurement equation is found that relates a difference in scalar potential to the INS position error growth, a constant velocity error and a constant acceleration error. This general expression has the form: 
       Δ ψ   = {circumflex over (p)} ·C   1 ( T )+   {circumflex over (V)} ·C   2 ( T )+ {circumflex over (α)} · C   2 ( T )
 
     The constant velocity error  {circumflex over (V)} is not the same as the (varying) INS velocity error  {circumflex over (V)} . 
     Therefore, use of this form of measurement equation requires an additional state for the constant velocity error  {circumflex over (V)} to be estimated by the Kalman filter. 
     It is noted that observability of the INS error sources depends on the amount of variation of the magnetic field. If the magnetic field is completely homogeneous, error sources are not observable. Correction may be postponed until a variation in the magnetic field is detected. 
     During a trajectory (also referred to as path), the magnetic field can be compared with the expected magnetic field by looking it up in a table. If the difference exceeds a threshold, INS correction can be omitted (if other sensors are available), or postponed until the magnetic field is as expected. 
     Because the magnetic field is known beforehand, the evaluation period T can be taken such that maximum observability is achieved. 
     The navigation device  1  shown merely by way of non-limiting example in 
       FIG. 2  comprises a navigation state estimation unit  11 , a prediction unit  12 , a storage unit  13 , a subtraction unit  14 , a filter unit  15 , and a correction unit  16 . In the embodiment shown, the navigation state estimation unit  11  is constituted by an INS (Inertial Navigation System) unit which produces an estimated navigation state comprising an estimated position, an estimated velocity and an estimated attitude. The estimated position is fed to the storage unit  13  to retrieve a stored magnetic property value (Ψ). According to the present invention, this magnetic property value is the integral over distance of the magnetic field strength from a reference point to the current position. Using the estimated position, the prediction unit  12  produces a predicted magnetic property value, as will be explained later in more detail with reference to  FIG. 3 . If the predicted and the retrieved magnetic property values are identical, then the navigation error is apparently zero and no correction is required. Only rarely, it may occur that the magnetic field over the part of the trajectory is uniform, in which cases no correction can be applied. 
     Normally, however, the predicted and the retrieved magnetic property values are not identical. Their values are compared in the subtraction unit  14 , which produces their difference r which is output to the filter unit  15 . The subtraction unit  14  is preferably further configured for detecting two special conditions: 
     a. The difference of the predicted and the retrieved magnetic property values is smaller than a first threshold, meaning that they are practically identical; and 
     b. The difference of the predicted and the retrieved magnetic property values is greater than a second threshold, meaning that these values diverge more than expected. 
     When either condition is met, special measures may be taken, as will later be discussed in more detail. 
     In the embodiment shown in  FIG. 2 , the filter unit  15  comprises a Kalman filter, which may be known per se. Using the difference of the predicted and the retrieved magnetic property values, the filter unit  15  produces a correction value which is fed to the correction unit  16 , which in turn uses the correction value to correct the estimated navigation state produced by the navigation state estimating unit and to produce a corrected navigation state estimate. In the embodiment shown, the correction unit  16  uses subtraction to correct the navigation state, but other correction mechanisms, such as multiplicative correction, may also be used, provided a suitable correction value is used. As a standard practice the INS may regularly be reinitialized to the corrected states to avoid large navigation state errors. 
     It will be understood that the storage unit  13  preferably comprises a ROM (Read Only Memory) or similar data storage. The storage unit  13  may be physically present in the vehicle in which the device is used, or may be a remote storage with data communication facilities. 
     When condition a. above is met, meaning that the predicted and the retrieved magnetic property values are identical or nearly identical, the updating process of the Kalman filter may be suspended in order to preserve the filter&#39;s internal states. This suspension ends when the first threshold is exceeded again. 
     When condition b. above is met, meaning that the predicted and the retrieved magnetic property values are too different to allow a regular correction, the updating process of the Kalman filter may also be suspended. 
     The prediction unit  12  schematically illustrated in  FIG. 3  comprises a magnetic measuring unit  121  and an integrating unit  122 . The magnetic measuring unit  121  measures the magnetic field H during the course of a trajectory. Preferably, the magnetic field is measured at suitable sampling intervals, which may be expressed in time and/or distance. The measured magnetic field values are fed to the integration unit  122 , which integrates these values over distance so as to produce the magnetic property value of the present invention, using for example {tilde over (ψ)}=∫ S   {tilde over (H)} ·d {tilde over (x)}  as explained above. The thus produced “predicted” magnetic property values are then fed to the subtraction unit  14 , where they are compared with the “expected” magnetic property values retrieved from the storage unit  14 . 
     The present invention is based upon the insight that navigation errors can be corrected using partial integrals over distance of the magnetic field, and that comparing such integrals during use with previously made integrals allows navigation corrections to be made. Use of the integrals, instead of the magnetic field measurements, allows for the seamless integration with other sensors in integrated navigation systems as described in detail in Groves. 
     The present invention can advantageously be used together with the invention disclosed in US patent application US2014005975 and its European counterpart EP2678637 (TNO file PLT 2011 015), both of which are herewith incorporated in this document. More in particular, the present invention can advantageously be combined with a measuring device and/or method for measuring a vectorial physical quantity, comprising the steps of
         measuring a component of the vectorial physical quantity with at least a first and a second sensor therewith generating a respective sense signal respectively, said sense signals each including a first part indicative for the component of the measured quantity and a second, noise part, said sensors having an Allan variance curve with a minimum value for a particular integration time, said Allan variance curve having a first tangent line to a portion of the curve for which the integration time approaches zero, and having a second, horizontal, tangent line of constant standard deviation corresponding to said Allan minimum value, said first and said second tangent line having an intersection point for a second particular integration time,   causing a relative rotation between said at least a first and a second sensor,   providing a difference signal indicative for a weighted difference between the sense signals, said weighted difference being substantially independent of the vectorial physical quantity, comprising the steps of   weighting the sense signals dependent on the relative rotation to provide mutually weighted sense signals, the signal weighting unit having at least one multiplication unit for multiplying one of the sense signals with a weighting factor,   providing a difference signal indicative for a difference between the mutually weighted sense signals,   estimating a correlated low frequency noise component from the difference signal and from information about the relative rotation between the sensors and for generating an estimated noise signal indicative for the estimated value of said noise component, said low frequency noise estimation having an associated effective integration time that is at least two times a sample time with which the sense signals are obtained, which effective integration time is less than the smallest particular integration time of the sensors comprised in the measuring device and which effective integration time is at most two times said second integration time,   providing an output signal indicative for a sensed value of the vectorial physical quantity corrected for the noise as estimated by said low frequency noise estimation.       

     It is noted that any terms used in this document should not be construed so as to limit the scope of the present invention. In particular, the words “comprise(s)” and “comprising” are not meant to exclude any elements not specifically stated. Single (circuit) elements may be substituted with multiple (circuit) elements or with their equivalents. 
     It will be understood by those skilled in the art that the present invention is not limited to the embodiments illustrated above and that many modifications and additions may be made without departing from the scope of the invention as defined in the appending claims.