Abstract:
Embodiments of the present disclosure are directed to a method and apparatus for calibration of gyro sensors by using magnetic sensor measurements and background computation during normal product operation. In one embodiment, magnetic sensor measurements are used to adjust gyro gain by comparing measured magnetic vector orientation with its expected orientation computed from gyro integration. The background process constantly compares this discrepancy for various values of gyro gain and selects the one that minimizes such error on average. In one embodiment, device orientation obtained by gyro integration is improved by using magnetic sensor measurements.

Description:
RELATED APPLICATIONS 
     This application claims the benefit of priority of U.S. Provisional Application No. 61/456,467, filed Nov. 8, 2010, the content of which is hereby incorporated by reference in its entirety. 
    
    
     BACKGROUND OF THE DISCLOSURE 
     1. Field of the Disclosure 
     The present disclosure relates to devices incorporating magnetic sensors and more particularly to systems including inertial gyroscopic sensors. 
     2. State of the Art 
     Low cost inertial and other sensors became very popular in numerous consumer electronic devices such as smart phone, game controllers, and many others. This imposes very challenging requirements on the developers. The less expensive are products, the less technically sophisticated its customer can be expected and the product should provide robust and care-free operating experience even while using low cost low precision sensors. 
     It is a common practice in the device manufacturing, to calibrate all the sensors employed by the device at the time of device manufacturing before it is shipped to the customer. Such, so called factory calibration, can be performed very accurately using special equipment, controlled environment, and trained personnel. However, by the time the product is used by the customer the sensor calibration is off due to the temperature change, different battery voltage, and other effects that are beyond control of the design engineers. 
     This issue of proper calibration is especially important for the sensors, which, in order to be used in the navigation or attitude determination solution, need to be integrated over time. With integration calibration errors grow with time and make the resulted solution very quickly unacceptable. 
     This means that to ensure proper device operation, for every sensor used in the device, there should be an appropriate sw and algorithms that allow for autonomous real time calibration of such effects as signal drift, and effective gain change which happen during normal device operation due to the environmental effects listed above. In other words, a practical design must incorporate internal algorithms that autonomously monitor and adjust sensor calibration parameters in real time. 
     This also means that there should be methods that allow reducing the effect of integrated sensor error. 
     It is, therefore, an object of the present disclosure to resolve one or more of the foregoing difficulties by providing a method that automatically adjusts calibration parameters of the gyro sensors in real time during normal device operation. 
     It is also an object of this disclosure to provide a method that allows reduction of error in gyro integration. 
     A typical system today has a full set of sensors such as 3-axis accelerometer, 3-axis magnetic sensor, and 3-axis gyro sensor. All these sensors need to be calibrated to allow effective use. The most common way of accelerometer calibration is to use time moments when system is stationary and to use the gravity vector as the natural standard to determine accelerometer bias and gain. These stationary situations can be set up by the user during user calibration or detected during real time operation, e.g. via “zero motion detector” as in U.S. Pat. No. 5,991,692 since in a stationary position only gravity vector effects accelerometer readings. Using multiple measurements with device in different orientation, a full accelerometer calibration can be performed. As an example see U.S. Pat. No. 6,729,176. 
     Magnetic sensors are commonly calibrated using the natural Earth magnetic field. The calibration can be done by recording measurements in several predetermined device orientations, or as a background process using operational or user calibration device 3D motion, e.g., US Patent Application Publication No. 20090070056. Such procedure creates new calibration parameters every time the calibration is performed. 
     Since GPS sensor became widely used and very common in navigation systems where magnetic sensors are employed, it is often used for the calibration of magnetic sensors by comparing GPS velocity vector direction with the direction obtained from the magnetic sensor. Then direct computation of Kalman Filter procedure is used to derive the required calibration parameters for the magnetic sensors. 
     There is a significant body of work that allows a person who is of ordinary skills in the art, to derive an appropriate procedure to accurately calibrate magnetic and accelerometer sensors used in the device by using natural and ever present forces of Earth Gravity and Earth Magnetic field which are well known and tabulated with high accuracy at any location on the Earth, and by using such additional sensors as GPS and Temperature. 
     As anyone skillful in the arts is aware, when device orientation is derived from the measurements of the well calibrated magnetic and gravity sensors, each attitude determination is independent. Therefore, unavoidable errors in determination of the magnetic and gravity vector directions are not accumulated in time. Quite opposite, by smoothing the resulted trajectory the orientation error in each trajectory point can be reduced. 
     However, in the situations when device experiences an unknown acceleration, the internal accelerometers cannot be used to determine device orientation relevant to the Earth gravity vector. Magnetic vector alone is not sufficient to restore device orientation, nor a single antenna GPS can be used to restore body orientation. In such cases one has no alternative but to employ gyro sensors to determine device 3D orientation. 
     The Gyro sensors are measuring a rate of turn around its local coordinate axis. When properly integrated, this provides a total 3D body rotation from its initial orientation. However, since gyro signal must be integrated to obtain orientation any error in gyro bias or gain rapidly grows in time. 
     Gyro bias can be determined by observing gyro signals at the moment when there is no rotation is present, assuming that such moments can be detected. However, to calibrate gyro gain one needs to perform very accurate and well known rotation around each axis which is difficult to perform even in factory setting and all but impossible for the in-field real time calibration. 
     Tekawy, et al., U.S. Pat. No. 7,667,645, U.S. Pat. No. 7,393,422 is using power difference from different GPS satellites in multiple GPS antennas to calibrate gyro drift. 
     Drag, et al., U.S. Pat. No. 7,657,183, uses optical sensors as a primary method of orientation. This patent mentions that gyro can be calibrated using these optical signals but does not explain how it can be done. 
     Achkar, et al., U.S. Pat. No. 5,562,266, is also using Sun and North Pole detectors to estimate gyro drift (bias) via a Kalman Filter, which is also similar to the method used by Basuthakur, et al, U.S. Pat. No. 5,452,869 
     Therefore, at the present, the gain calibration methods need to rely on the external determination of the device trajectory orientation, usually GPS for vehicle, star finder for satellites, or sonar for underwater apparatus. However, such methods don&#39;t work when device experiences rotation which is independent on its trajectory which is a common case for the hand held devices. 
     It is an object of the present disclosure to provide a practical and efficient method that automatically adjusts calibration parameters of the gyro sensors in real time during normal device operation. It is also an object of this disclosure to provide a practical and efficient method of improving accuracy of the device orientation determination. 
     SUMMARY OF THE DISCLOSURE 
     The present disclosure instructs anyone who is proficient in the art, how to perform the appropriate real time computations which employs measurements of the magnetic sensors, prior calibrated by any of the methods that are known to those who are proficient in the arts, and integration of the gyro sensors, to determine correction to the gain of the gyro sensors. 
     This disclosure observes that if gyro sensors are properly calibrated, then the device rotation computed by a proper integration of the gyro signals, should rotate the magnetic vector components from their orientation at the moment when the integration starts to their observed orientation at the moment when the gyro integration ends. Then any measure of mismatched between the observed and rotated magnetic vectors can be used to improve calibration of the gyro sensors. 
     With such measure of closeness established, the current disclosure teaches that the larger the gyro calibration error, the larger this measure of mismatch to become and that the minimum of this measure corresponds to the best calibration of the gyro sensors. Therefore, if one builds dependence between the employed gyro parameters and the error between the measured and estimated orientation of the magnetic vector, the correct value of the gyro calibration parameters will yield a minimum to this surface. 
     This disclosure also instructs those who are of ordinary skills, that device orientation, computed by gyro integration can be improved by making sure that the computed rotation is consistent with the observed orientation of magnetic vector. 
     One embodiment of the present disclosure is a device that includes a microprocessor having an associated memory; at least two sets of independent sensors in the device capable of sensing device 3D orientation, and/or 3D motion where at least one sensor set is not well calibrated at the beginning of device motion; and a microprocessor program operable in the microprocessor. This program performs steps of a) reading, storing, and processing sensor data from the sensors at a sequence of time moments t=0, 1, . . . k; b) estimating device 3D orientation at time t=k using sensor data at t=0, 1, . . . , k; c) forecasting a next time moment t=(k+1) value of one set of the above mentioned sensors by using sensor measurements for the previous time moments: t=m, . . . , k 0&lt;=m&lt;=k; d) constructing a difference between forecasted and observed sensor values; and e) using the difference to continuously improve determination of the above estimated device 3D orientation and/or for the calibration of the set of sensors that was not well calibrated at the beginning of the device motion. 
     In one embodiment of the device, the sensor sets are a magnetic set and a gyro set. In one embodiment the gyro sensor set is not well calibrated at the beginning of the device motion. 
     In one embodiment a method of calibration of the gyro sensor set comprises: observing gyro sensor signals over a predetermined time period; correcting gyro signals by applying one of a first guess of the unknown gyro calibration parameters or their nominal values; determining a perceived change of the device orientation during this time period based on the observed and such corrected gyro signals; comparing magnetic sensor values at the beginning and the end of an integration period, appropriately corrected by the derived change of device orientation; and correcting gyro calibration parameters initial values by applying correction derived from a measure of difference between the observed and nominal rotation corrected magnetic signals. 
     In one exemplary embodiment of the device, only one calibration parameter is determined for all gyro channels. The calibration parameter may be an unknown gyro gain correction to the nominal values. In one embodiment the calibration parameter may be an unknown gyro bias correction to the nominal values. 
     In another embodiment, multiple calibration parameters may be determined. In such an embodiment the calibration parameters may be unknown gain corrections for each of the gyro sensors. Alternatively the calibration parameters may be unknown bias corrections for each of the gyro sensors. In another embodiment, the calibration parameters are a combination of unknown bias corrections and unknown gain corrections for all gyro sensors. 
     In another embodiment the measure of the difference between the magnetic signal at the end of the time period and the rotation corrected value of the magnetic vector at the beginning of the time period is a vector norm of the vector difference of these two vectors. In an embodiment, at the same time interval multiple alternative rotations may be determined with different assumed values of the gyro parameters under calibration. Alternatively, at different time intervals different values of the gyro calibration parameters may be used for the rotation determination. In such cases, correction parameter values may be selected by the method of Optimal Experiment Design. 
     In one embodiment calibration parameters may be determined by finding such parameter combination that yields a minimum to the computed measure of difference between the observed and rotation corrected magnetic signals. Alternatively all computations may be done off line with the data that was recorded for one set of calibration parameters and then recomputed with artificial change in the calibration parameters. 
     In one embodiment a method of determination of the device rotation in accordance with the present disclosure includes observing gyro sensor signals over a predetermined time period; correcting the gyro signals by applying gyro calibration parameters; determining a perceived change of the device orientation during this time period based on the observed and such corrected gyro signals. The method further includes observing magnetic sensor signals at the beginning and the end of this time period; correcting magnetic sensor signal by applying magnetic calibration parameters; determining a perceived change of the device orientation during this time which is consistent with the observed change in the calibration corrected magnetic sensor signals and is the closest to the rotation determined from the gyro signals; and determining the device rotation during this time interval as a combination of the rotation derived from the gyro sensors and a rotation that is consistent with the magnetic sensor observations. 
     In such embodiments, the rotation closeness measure is a norm of quaternion difference between the quaternion that is computed from the gyro signal integration and a quaternion that provides device rotation that is consistent with the magnetic sensor observations. The computation may be done off line using recording of the sensor values during an actual motion. 
    
    
     
       BRIEF DESCRIPTION OF THE DRAWINGS 
         FIG. 1  depicts a device that consists of a microcomputer ( 100 ) which includes Processor ( 101 ), Memory Bus ( 102 ), Memory ( 103 ), ( 104 ), ( 105 ), and I/O interface ( 106 ). Microcomputer has Display ( 111 ) and Input Device ( 112 ). I/O Interface is connected to a bank of sensors with their ND converters ( 121 ): 3 axis Gyro sensors ( 122 ), 3 axis Magnetic Sensors ( 123 ) and 3 axis Accelerometers ( 124 ). CPU is sampling all the sensors via I/O ( 106 ) and stores the corresponding values in the Storage Device ( 105 ) or in the Main Memory ( 103 ). The appropriate computer program, that processes such obtained sensory data is stored in ROM ( 104 ) and is executed by CPU ( 101 ). 
         FIG. 2  shows a trajectory restoration via gyro integration when nominal gyro gain provided by the manufacturer was used. Dash line shows the true trajectory. 
         FIG. 3 . provides the same trajectory restoration when an optimal gain found with one implementation of the presented method was applied. 
         FIG. 4  depicts a graph of the experimental data on how Error in Magnetic vector restoration depends on the Gyro Gain variation from its manufacturer nominal value. 
         FIG. 5  provides algorithm flow chart on how Error of magnetic vector restoration is computed in real time. 
         FIG. 6  depicts how optimal value of the Gyro Gain is computed from multiple Error values obtained for different values of Gyro Gain. 
         FIG. 7  shows that multiple Error values can be computed in parallel or in serious depending on the CPU capabilities. The actual trajectory computation are first performed using a nominal manufacturer gain value G0 and then switched to the optimal gain Gopt as soon as it becomes available. 
         FIG. 8  illustrates how rotation obtained from the Gyro Integration—quaternion Le, can be improved by computing rotation Lf which is consistent with the observed resulted orientation of the magnetic vector M and is the closest to the Gyro rotation Le. Quaternion Lm provides a rotation that rotate magnetic vector M from its observed position M(k) to the shown orientation M(k+1) which is its observed orientation at t=k+1. Quaternion Lf is then constructed of two rotations—Lm and an additional rotation around M. 
     
    
    
     DETAILED DESCRIPTION 
     In the following description, numerous specific details are set forth in order to provide a more thorough disclosure. It will be apparent, however, to one skilled in the art, that the art disclosed may be practiced without these specific details. In some instances, well-known features may have not been described in detail so as not to obscure the art disclosed. 
     In a preferred implementation we use quaternion arithmetic for rotation determination. However, this selection does not constrain the disclosure and is used for the illustration purposes only. Those who are proficient in the art know that there are other rotation representations such as Rotational Matrix, Euler angles, or other methods. 
     Quaternion is an ordered set of four real numbers q=(a,b,c,d) that has some particular mathematical properties. A detailed description of quaternion can be found in numerous textbooks on the subject. For our purpose we use unit quaternion.
 
norm( q )=sqrt( a*a+b*b+c*c+d*d )=1  (1)
 
     Union quaternion represents rotation around a unit vector e=(e1,e2,e3) by an angle α
 
 q =[cos(α/2),sin(α/2)* e]   (2)
 
     So if there a vector m=[mx,my,mz] represented in the coordinate system (x,y,z) is rotated around vector e by angle α, then the new coordinates of vector m are
 
 m′=q*m*q′   (3)
 
     In quaternion representation a rigid body rotation can be represented by the following differential equation
 
 dq ( t )/ dt= 0.5*ω( t )* q ( t )  (4)
 
     Where ω=[ωx, ωy, ωz] is a vector of angular rate of rotation in the body coordinate system which is what is measured by the properly calibrated gyro sensors. 
     The solution of this differential equation can be obtained numerically using such methods as Runge-Kutta or any other method known to the knowledgeable in the arts. At the same time one needs to take into consideration that when a rigid body is rotated first by a quaternion q1 and then by a quaternion q2 then the total rotation can be represented by their ordered product
 
 q 12 =q 2 *q 1  (5)
 
     Suppose that at the time moment T=t0 the magnetic vector measurements in the unit coordinate frame are
 
 m 0 =[mx 0 ,my 0 ,mz 0]  (6)
 
     and at the time T=tk the measured magnetic vector has the following internal components
 
 mk=[mxk,myk,mzk]   (7)
 
     The rotation from T=t0 to T=tk can be computed by integrating equation (4) under constraints (5). Suppose that this integration produced quaternion qk. Then, using eq. (3) the expected composition of the magnetic vector in the rotated coordinate system should be
 
 mk _est= qk*m 0* qk′   (8)
 
     At this point one can compare the observed and estimated value of magnetic vector and to create a measure of closeness. As a method of illustration we will use a vector norm of the difference of these two vectors
 
 em ( k )=norm( mk−mk _est)  (9)
 
     This operation can be performed at the next time interval t(k)−t(k+1), etc. As a result, one can compute accumulated error measure
 
 Em=em ( k )+ em ( k+ 1), . . . , em ( k+n )  (10)
 
     The raw gyro sensor measurements are digital values of the voltage or current generated by the sensor. These values are denoted as g=[gx, gy, gz]. These values are related to the angular rotation rate w that is used in eq (4). It is commonly known by the experienced in the art, that the vectors w and g are related by a linear relationship for each of their components as the following:
 
 gk=αk*ωk+bk;k=x,y,z   (11)
 
     Therefore, the angular rate ω, which is used in eq (4) can be derived from the gyro raw signals as
 
ω k=ck*gk−bk;k=x,y,z   (12)
 
     In the equation (12) ck represents gyro gain, and bk is gyro bias. The values of ck and bk should be accurately known at the time of device operation in order to allow correct determination of the device orientation. While these parameters are often determined during factory calibration their actual values during operation are slightly different. Suppose the actual value of ω differs from the value computed as eq (12) by a small gain correction μ which is unknown.
 
ω′=(1+μ)*ω  (13)
 
     To determine the unknown parameter μ we propose to run several gyro integrations (4) in parallel or in sequence, with several different values of p. In the preferred implementation, as an illustration but not as a limitation, three integrations are run in parallel with μ0=0, μ1=+d, μ2=−d, where d is a small positive number 0&lt;d&lt;&lt;1. The results of the integration with μ=0 is used for the device operation while the other integrations are run on the background and their results are used for the determination of the true value of μ as described below. 
     The closeness measure Em—eq (10), is computed for each of these integration. As a result, the set {Em(k), m(k)} constitutes a dependence curve Em vs μ as illustrated in  FIG. 4 . Then an optimal value of μ that yields the minimum to this function can be found using quadratic fit or any other function minimum determination method that is known to those proficient in the arts.  FIG. 4  shows a quadratic fit and its minimum as found on the data collected during operation of some actual device. 
     After the optimal value of μ=μ* is found, it is used for the nominal integration: μ0=μ*, μ1=μ*+d, μ2=μ*−d, and the whole process continues so any possible real time change in gyro gain get detected immediately. 
     The presented illustration determines only one parameter as a min of a single variable function. However, a similar procedure can be performed if more parameters need to be estimated. Indeed, suppose that one wants to find correction parameters μ for each of three gyro axis individually—μx, μy, μz. Then, instead of repeating the same patent of 3 parallel integrations with corrections [0,0,0], [d,d,d], and [−d,−d,−d] one can change one parameter for each run, e.g. [0,0,0], [d,d,d], [−d,−d,−d] for the first run, then after a an optimal common factor μ* is determined as described above, the next run can be done with parameters [0,0,0], [d,0,0], [−d,0,0], so additional correction μx* is determined, then μy* and so on. 
     Alternatively, one can design an experiment where all three parameters are determined in each run by employing extra parallel integrations. Since the computational load to resolve eq (4) is not significant, a parallel computation of the multiple threads at the same time is quite achievable with today microprocessors typically used in such systems. 
     In another implementation only one integration thread for one value of μ is performed each time and then results of several integrations are compared even if these integration were performed at different time. 
     Even when all calibration parameters are optimized, the gyro signals have unavoidable noise which integration lead to orientation errors that grow with time proportionally to the square root of time of integration. These errors can be reduced but not eliminated by using Dead Zone technique, where any gyro signal which is smaller than some predetermined threshold is forced to zeros. The dead Zone technique improves the overall SNR (Signal to Noise Ratio) but does not eliminate the full effect of the noise integration. 
     Therefore, to alleviate the above mentioned difficulties, this disclosure teaches how to improve orientation determination using magnetic vector in addition to the gyro integration. 
     If the magnetic sensors are properly calibrated, then one expects that the rotation computed from the gyro integration according to the eq (4)-(5) should be consistent with a rotation that brings vector m(t0) into m(t). Lg is a quaternion computed from the eq (4)-(5). 
     A rotation that rotates vector m from m0 to m1 can be performed with a quaternion Lm:
 
 m 0 e=m 0/norm( m 0);m1 e=m 1/norm( m 1);
 
 e=m 0 e×m 1 e;  
 
α= a  cos( m 0 e  dot  m 1 e )
 
 Lm =[cos(α/2),sin(α/2)* e]   (14)
 
     An additional rotation around new vector m1 is possible without changing vector m1. This additional rotation can be represented by a quaternion Lr
 
 Lr =[cos(β/2),sin(β/2)* m 1 e]   (15)
 
     Where β is an angle of additional rotation, not defined at this moment. This disclosure teaches, that the angle β should be selected such, that the total rotation by quaternion Lm and Lr, (denote it as Le(β)) should be the closest to the rotation Lg
 
 Le (β)= Lr (β)* Lm   (16)
 
     In the preferred implementation of this disclosure for the measure of quaternion closeness the norm of a quaternion difference is selected. Therefore, an additional rotation β is determined from minimization of norm(Lg−Le(β))
 
β-&gt;min(norm( Lg−Le (β))  (17)
 
     Eq. 17 can be solved by many methods familiar to those proficient in the arts. Let the minimal value of (17) is achieved at β=β*. 
     Therefore, at the end of each gyro integration step, there are two quaternion that potentially describe rotation: Lg and Le(β*). There are multiple ways how these two quaternion can be used—e.g., as a weighted sum with the weights inverse proportional to the estimated accuracy of each quaternion. The preferred implementation is using an average of these two quaternion.  FIG. 8  shows improvement that the preferred implementation of this disclosure achieved for a test case rotation in controlled conditions. 
     It will be readily seen by one of ordinary skill in the art that the disclosed embodiments fulfill one or more of the advantages set forth above. After reading the foregoing specification, one of ordinary skill will be able to affect various changes, substitutions of equivalents and various other embodiments as broadly disclosed herein. It is therefore intended that the protection granted hereon be limited only by the definition contained in the appended claims and equivalents thereof.