Automotive navigation system and method to utilize internal geometry of sensor position with respect to rear wheel axis

[0001] Embodiments disclosed here relate to a method involving a vehicle navigation system, and more particularly, to a navigation method and system utilizing the internal geometry of the sensor position with respect to the vehicle's rear-wheel axis for maintaining high positioning accuracy even when GPS signals are lost for a long period of time.

[0002] The inertial navigation system (INS) is a widely used technology for guidance and navigation. The INS is composed of an inertial measurement unit (IMU) and a processor wherein an IMU contains accelerometers and gyroscopes which are inertial sensors detecting platform motion with respect to an inertial coordinate system. An important advantage of the INS is independence from external support, such as positional signals from artificial satellites, however, it cannot maintain high accuracy for long distance by itself because of accumulating sensor errors over time.

[0003] More recent development in global positioning system (GPS) has enabled low-cost navigation without growing error. The GPS, however, involves occasional large multipath errors in urban canyons (i.e., urban areas surrounded by high rise buildings) and signal dropouts inside buildings or tunnels. Therefore, efforts have been made to develop integrated INS/GPS navigation systems by combining the GPS and INS using a Kalman filter algorithm to remedy the performance problems in both systems.

[0004] Inertial sensors (accelerometers and gyroscopes) for an IMU used to be expensive and large, thus only used in high precision applications, for example, aerospace and military navigation. To establish an IMU with compact packaging and an inexpensive manner, efforts have been made to develop micro-electro mechanical system (MEMS) sensors, resulting in commercialization of low-cost and small inertial sensors. However, MEMS sensors involve large bias and noise. Low cost MEMS sensors have been largely adopted by cost sensitive navigation products such as automotive and portable navigation systems. In integrated MEMS IMU/GPS navigation systems, however, errors quickly accumulate into a large amount as soon as GPS signals drop out due to buildings, tunnels, etc.

[0005] Bye et al. suggested, in U.S. Pat. No. 6,859,727 entitled “ATTITUDE CHANGE KALMAN FILTER MEASUREMENT APPARATUS AND METHOD”, in Col. 4, lines 61-62, a Kalman filter based calibration method as “for a non-rotating IMU, the externally observed attitude or heading change (at the aiding source) is taken to be zero”. When this concept is applied as a condition of zero side velocity of a ground vehicle, divergent navigation solutions due to erroneous MEMS sensors are largely improved to be non-divergent solutions with reasonable positioning accuracy. However, as discovered by the inventor of this application, vehicle's side velocity has an analytical non-zero term. Suppressing the non-zero velocity into zero will cause unfavorable side effects such as a shortage of the total velocity and erroneous increment of pitch angle estimate, which results in a large positioning error when GPS signals are lost.

[0006] Therefore, there is a need for a new navigation system and method using low-cost MEMS IMU with capability of maintaining high accuracy even when GPS is lost for a long period of time by evaluating not only constant values but also non-constant and non-zero analytical conditions.

[0007] It is, therefore, an object of disclosure to provide an embodiment which is an integrated INS/GPS navigation system and method incorporating low-cost MEMS in an IMU (inertial measurement unit) to utilize the internal geometry of the sensor position with respect to the vehicle's rear-wheel axis for maintaining high positioning accuracy even when GPS signals are lost for a long period of time.

[0008] One aspect of the embodiment is that the proposed navigation method uses an analytical condition derived from a vehicle's mechanical condition so-called Ackermann Steering Geometry (see Genta, G., “MOTOR VEHICLE DYNAMICS Modeling and Simulation”, World Scientific Publishing Co. Pte. Ltd., 1997, 5 Toh Tuck Link, Singapore, pp. 206-207) for enhancement in navigation accuracy. The analytical condition is a relationship between the vehicle's lateral directional velocity, the distance of the sensor position with respect to the rear wheel axis, and the angular rate with respect to the vehicle's z-axis.

[0009] Another aspect of the embodiment is a Kalman filter based navigation method to utilize the aforementioned analytical condition by incorporating the distance of the sensor position with respect to the rear wheel axis into the Kalman filter's states as an auxiliary parameter so that the system can automatically estimate the distance between the sensor position and the rear wheel axis without need of manually measuring the distance upon installation of the navigation system.

[0010] Another aspect of the embodiment is the Kalman filter based navigation method to continuously utilize the aforementioned analytical condition as an auxiliary measurement at a high frequency executed independently of GPS measurement.

[0011] Another aspect of the embodiment is that the system with zero distance of the sensor position with respect to the rear-wheel axis will achieve the highest positioning accuracy, thus such sensor position is suggested as “the best sensor position” for automotive navigation. This is because the analytical condition reduces to “zero lateral velocity” without need of evaluating gyro outputs including bias estimates. The system with the best sensor position is practically achieved by placing a sensor IMU at the bottom-center of the rear trunk which is approximately above the center of the rear wheel axis for majority of vehicles.

[0012] Another aspect of the embodiment is to show a vehicle contour image and the navigation system's position on a display with proper geometry between the vehicle contour and navigation system's position in which the distance between the navigation system and the rear wheel axis is automatically estimated.

[0013] Another aspect of the embodiments is to align the navigation system's absolute position (latitude and longitude) with respect to objects surrounding the vehicle by placing the GPS antenna above the navigation system in which the distance between the navigation system and the rear wheel axis is automatically estimated.

[0014] According to the embodiments: (1) regardless of the sensor position, the distance of the sensor position with respect to the rear wheel axis will be automatically estimated without need of measuring the distance by hand, which will be utilized to enhance navigation accuracy; (2) high positioning accuracy is maintained even when GPS signals are lost for a long period of time using a low-cost MEMS IMU; (3) the best sensor position to achieve the highest navigation accuracy is the center of the rear-wheel axis which is practically available by placing a sensor IMU at the bottom-center of the trunk; (4) a driver's safety consciousness is enhanced by the visual aid from the display showing the vehicle contour with proper geometry with respect to the navigation system as well as to the surrounding objects.

[0035] Various embodiments will be described in detail with reference to the accompanying drawings. The embodiments described here are related to an integrated INS/GPS navigation system and method incorporating low-cost micro-electro mechanical system (MEMS) sensors in an IMU (inertial measurement unit) to utilize the internal geometry of the sensor position with respect to the vehicle's rear-wheel axis for maintaining high positioning accuracy even when GPS signals are lost for a long period of time.

[0036] FIGS. 1A-1E show examples of structure related to the embodiments of the integrated INS/GPS navigation system and method. FIG. 1A is a block diagram showing an example of comprehensive system architecture of the embodiment of the integrated INS/GPS navigation system with the input-output relationship, FIG. 1B is a schematic diagram showing an example of structure of an IMU (inertial measurement unit) with low-cost MEMS sensors, FIG. 1C is a schematic diagram showing an example of structure of an auxiliary measurement unit, FIG. 1D is a schematic diagram showing an example of structure of a GPS measurement unit, and FIG. 1E is a block diagram showing an example of structure of a navigation operation unit.

[0037] With reference to FIG. 1A, the integrated INS/GPS navigation system includes an IMU (inertial measurement unit) 10 containing MEMS inertial sensors of a three-axis accelerometer and three-axis gyro, an inertial navigation system (INS) computational unit 20, an iterated extended Kalman filter (IEKF) unit 30, a calibration unit (Cal) 40 having an auxiliary measurement unit (Aux) 50 and a global positioning system (GPS) unit 60, a navigation operation unit 70, and a display 80. In the example of FIG. 1A, the IMU 10 and the Cal 40 are configured to input various parameters, the INS 20, the IEKF 30 and the navigation operation unit 70 are configured to process the input parameters to produce state estimates of the vehicle including position estimates, and the display 80 is configured to output the resultant position estimates, etc.

[0038] The IMU 10 measures vehicle's accelerations and angular rates. In the example of FIG. 1B, the IMU 10 includes a processor 12, and the inertial sensors consisting of three (three-axis) accelerometers Acc x-z and three (three-axis) gyroscopes Gyro x-z. The accelerometers Acc x-z detect accelerations in the three (X,Y,Z) coordinates with respect to the sensor fixed coordinate system described in FIG. 3, and the gyroscopes Gyro x-z detect angular rates about the three (X,Y,Z) coordinates with respect to the sensor fixed coordinate system described in FIG. 3. As noted above, the inertial sensors are established by low-cost MEMS (micro-electro mechanical system) sensors. The processor 12 calculates the accelerations and angular rates of the vehicle based on the measured signals from the inertial sensors Acc x-z and Gyro x-z. The IMU 10 produces the measured data at a rate of, for example, 25 times per second (25 Hz), which is supplied to the INS 20.

[0039] The INS 20 is, for example, so called a Six Degrees of Freedom (6DOF) INS which executes the INS computation to update navigation state estimates (position, velocity, orientation, and sensor bias) and their covariances, i.e., uncertainties of estimates, upon the sensor measurement. The INS 20 updates the navigation state estimates at a rate of, for example, 25 times per second (25 Hz). The navigation state estimates are periodically calibrated by a Kalman filter, for example, Iterated Extended Kalman Filter (IEKF) 30 shown in FIG. 1A. Basic structure and operation of an INS and a Kalman filter will be described in detail with mathematical equations in the later sections of “INS Technology” and “Kalman Filter Technology”.

[0040] The Aux (auxiliary measurement unit) 50 is unique to the embodiments of the integrated INS/GPS navigation system and method. Typically, the Aux 50 is configured by a processor as shown in FIG. 1C to conduct an operation prescribed by a program. More specifically, the Aux 50 produces the auxiliary measurement data (reference data) involving a distance (d) between the navigation system (ex. inertial sensors of IMU) and a vehicle's rear wheel axis. In the example of FIG. 1C, based on vehicle side velocity v_by and vehicle's z-axis angular rate ω_bz from the INS 20 as shown in FIG. 1A, the Aux 50 produces the auxiliary measurement data which is expressed by 0=v_by−dω_bz. The auxiliary measurement data (reference data) will be described in more detail later with respect to the vehicle side velocity in the description of analytical conditions. The Aux 50 sends the auxiliary measurement data to the Kalman filter 30 at a rate of, for example, 5 times per second (5 Hz).

[0041] In the example of FIG. 1D, the GPS 60 is configured by a GPS antenna, a GPS receiver 60, and a processor 62. Through the GPS antenna, the GPS receiver 61 receives GPS signals from a plurality of artificial satellites, and the processor 62 calculates the estimated location of the vehicle by comparing clock signals and position data included in the GPS signals. More specifically, the GPS 60 measures the GPS antenna position and velocity based on range and range-rate information between the GPS antenna and multiple satellites in the field of view to send the measurements to the Kalman filter 30. Typically, the GPS 60 produces the position and velocity data every one second (1 Hz).

[0042] The state estimates from the INS 20, the measurement data from the GPS 60, and the auxiliary measurement data from the Aux 50 are combined by the Kalman filter 30 which optimally estimates, in real time, the states of the navigation system based on such noisy measurement data. Namely, the navigation state estimates from the INS 20 are periodically calibrated by the Kalman filter 30 by taking the differences between the INS state estimates and the calibration measurements (auxiliary measurement and GPS measurement) obtained from the Aux 50 and the GPS 60, respectively. As noted above, the auxiliary measurement from the Aux 50 involves the analytical condition which is the relationship between the vehicle's lateral directional (side) velocity, the distance of the sensor position with respect to the rear wheel axis, and the angular rate with respect to the vehicle's z-axis. By incorporating this analytical condition in the calibration process by the Kalman filter 30, a positioning error caused by the vehicle side velocity can be canceled or minimized, a theory of which is described later.

[0043] The integrated INS/GPS navigation system of FIG. 1A estimates the following parameters: three position parameters (latitude, longitude, and altitude), three velocity parameters (velocities along the sensor fixed x, y, and z-axes), three orientation parameters (roll, pitch, and yaw angles of the sensor fixed coordinate system with respect to the North-East-Down coordinate system), six sensor biases (accelerometers and gyro biases along the sensor fixed x, y, and z-axes, respectively), two sensor attachment angles (pitch and yaw angles with respect to the vehicle body fixed coordinate system), and the position of the sensor IMU with respect to the vehicle's rear wheel axis.

[0044] In one embodiment, the position information including latitude, longitude, altitude, and the position of the sensor IMU with respect to the rear wheel axis will be used to display the vehicle contour and the navigation system with proper geometry between the navigation system's position and the vehicle contour on the display 80 in which the distance between the navigation system and the rear wheel axis is automatically estimated by the navigation system. The navigation system's absolute position (latitude, longitude, (and altitude if necessary)) and position information from a map database are aligned by placing the GPS antenna above the navigation system.

[0045] The navigation operation unit 70 is provided to conduct an overall operation of the navigation system for specifying a destination, searching and calculating an optimum route to the destination, conducting the route guidance operation to the destination, displaying a vehicle contour image with respect to images of surrounding objects, etc. In the example of FIG. 1E, the navigation operation unit 70 includes an input device 71 for selecting a menu, specifying a destination, executing a command, etc., a processor 72 for controlling an overall operation of the navigation system, and a display controller 73 for controlling the operation of the display 80. The input device 71 can be various hard keys, a touch screen formed on the display 80, a remote controller, a voice interface, etc.

[0046] In the block diagram of FIG. 1E, the navigation operation unit 70 further includes a map database (data storage device) 74 such as a hard disc, CD-ROM, DVD, flash memory, etc., for storing the map data (position data of links, nodes, polygons, etc.), a ROM 75 for storing various programs for navigation operations, and a RAM 76 for storing operational data or a processing result such as a guidance route. The map data is used to calculate a route to the destination, to produce a map image on the navigation screen such as on the display 80, and to provide various information on points of interest (POI), etc. An example of programs stored in ROM 75 includes a route search program to search and calculate possible routes to the destination, and a map matching program for matching the position estimates from the INS 20 with link, node and polygon data from the map database.

[0047] In FIG. 1E, the navigation operation unit 70 further includes a sensor unit 77 for detecting distances form other vehicles, pedestrians, structures, etc., a wireless communication device 78 for wireless communication with a remote server such as a traffic information server, a local event server, an internet server, a social network server, etc., and a contour database 79 for storing information on the vehicle contour. The sensor unit 77 may be configured by a plurality of radar sensors, cameras, etc. to measure shapes and distances from the outer objects. The vehicle contour data for the contour database 79 may be available from vehicle manufacturers or from data books in the industry. The vehicle contour data is used for showing images on the display 80 so that the user can easily comprehend the relationship among the vehicle contour, a location of the navigation system in the vehicle, and the objects surrounding the vehicle.

[0048] FIG. 2 is a flowchart showing an example of process in the integrated INS/GPS navigation system and method involving the preferred embodiment of FIGS. 1A-1E. At step 111, the IMU sensor measures vehicle three-axis accelerations and three-axis angular rates. This step is conducted by the IMU 10 having the accelerometers x-z and gyroscopes x-z shown in FIGS. 1A and 1B. As noted above, this step is conducted at a relatively high frequency, e.g., 25 Hz. The measured accelerations and angular rates are sent to the INS (inertial navigation system) in step 112 to update the navigation state estimates (position, velocity, orientation, and sensor bias) and their covariances. An example of INS is a “Six Degrees of Freedom (6DOF) INS’ as indicated by the INS 20 in FIG. 1A. As also noted above, this step is conducted at relatively high frequency, for example, 25 Hz.

[0049] At step 113, the Aux (auxiliary measurement unit) 50 shown in FIGS. 1A and 1C produces the auxiliary measurement data (reference data) involving a distance (d) between the navigation system (ex. inertial sensors of IMU) and a vehicle's rear wheel axis and a vehicle side velocity. More specifically, based on the vehicle side velocity v_by and the vehicle's z-axis angular rate ω_bz, the Aux 50 produces the auxiliary measurement data expressed by, for example, 0=v_by−dω_bz, which is sent to the Kalman filter 30 at a rate of, for example, 5 Hz.

[0050] The parameters from the INS are periodically calibrated at lower frequencies in step 114 according to the IEKF (iterated extended Kalman filter) technology by incorporating the auxiliary measurement data (reference data) from the Auxiliary measurement unit. In the example of FIG. 1A, this process is conducted by the Kalman filter 30 which receives the state estimates from the INS 20 and the auxiliary measurement data from the Aux 50.

[0051] At step 115, the process incorporates the GPS measurements. In the example of Figure D, the GPS (GPS Measurement Unit) 60 receives GPS signals from a plurality of artificial satellites and calculates the estimated position and velocity of the vehicle by comparing clock signals and position data included in the GPS signals. Typically, the GPS 60 produces the position and velocity data every one second (1 Hz), which is sent to the Kalman filter 30 as shown in FIG. 1A.

[0052] In step 116, the navigation state estimates (position, velocity, orientation, and sensor bias) from the INS 20 are calibrated by the Kalman filter 30 by taking the differences between the INS state estimates and the GPS measurements obtained from the GPS 60.

[0053] At step 117, the process repeats the above steps 111-116 as an integrated INS/GPS navigation system to continuously optimize the navigation position estimates. The GPS signals may not be available for a long period of time if a vehicle is in a tunnel, building, or valley of high-rise buildings, thus in such a case, the calibration data based on the GPS measurement cannot be used by the Kalman filter 30. Even in such a situation, the integrated INS/GPS navigation system and method disclosed here is able to maintain the high positioning accuracy since it utilizes the auxiliary measurement data from the Aux 50 to minimize the position error associated with the vehicle side velocity.

[0054] In step 118, the process combines the optimized position estimates obtained in the foregoing steps with the position information retrieved from the map database and the vehicle contour data. This process is conducted to match the optimized position estimates of the navigation system in the vehicle with the objects surrounding the vehicle and the vehicle contour. In the example of FIGS. 1A and 1E, this step is conducted by the navigation operation unit 70 which includes the map database 74 and the contour database 79.

[0055] In step 119, the navigation system displays the vehicle contour and the position of the navigation system in the vehicle with proper geometry between the navigation system's position and the vehicle contour. As noted above, based on the auxiliary measurement data (reference data) produced by the Aux 50, the distance between the navigation system and the rear wheel axis is automatically estimated by the navigation system. The navigation system's absolute position (latitude, longitude, (and altitude if necessary)) and the surrounding objects produced from the map database are aligned by placing the GPS antenna above the navigation system.

[0056] Based upon the system architecture and flowchart discussed so far, descriptions will be made regarding the detailed and theoretical navigation method to utilize the internal geometry of the sensor position with respect to the rear wheel axis for navigation accuracy enhancement. As preparation, brief explanation is present regarding the three major coordinate systems that the embodiments of the integrated INS/GPS navigation system will use in the following theoretical description.

[0057] FIG. 3 defines the sensor fixed coordinate system in which each of x, y, and z-axis represents a sensor's sensitive direction. A parameter with respect to the sensor coordinate system is indicated by a subscript of “s”.

[0058] FIG. 4 defines the vehicle body fixed coordinate system with a top view and a bird view in which the x-axis is defined toward the vehicle's forward direction, the y-axis is defined toward the vehicle's right-hand direction, and the z-axis is defined in the vehicle's downward direction from the ceiling to floor. A parameter with respect to the vehicle body fixed coordinate system is indicated by a subscript of “b”. The direction of the z_b-axis in the top view is represented by “the forward rotation of a right screw”, which is a common notation in physics. In other words, the clockwise rotation represents that we see the rear view of an arrow; the counter-clockwise rotation represents that we see the front view of an arrow.

[0059] FIG. 5 defines the North-East-Down (NED) coordinate system in which the x-axis is defined in the northerly direction, the y-axis is defined in the easterly direction, and the z-axis is defined in the vertically downward direction. A parameter with respect to the NED coordinate system is indicated by a subscript of “n”.

Analytical Condition

[0060] Here, the analytical condition incorporated in the embodiments of the integrated INS/GPS navigation system and method is described in detail. The analytical condition is derived from the internal geometry of the sensor position with respect to the rear wheel axis. This section corresponds to the function of the Aux (Auxiliary measurement unit) 50 in the block diagram of FIG. 1A and the step 113 in the flowchart of FIG. 2.

[0061] FIG. 6 shows so-called Ackermann Steering Geometry (see Genta, G., “MOTOR VEHICLE DYNAMICS Modeling and Simulation”, World Scientific Publishing Co. Pte. Ltd., 1997, 5 Toh Tuck Link, Singapore, pp. 206-207) in which a tangential line from the center of each front wheel passes a single center of cornering. FIG. 6 also shows that the rear wheels do not tilt in cornering for majority of vehicles.

[0062] Because of this basic mechanism, a point in the vehicle in front of the rear wheel axis has non-zero side velocity in the y_b direction. FIG. 7 illustrates that the navigation system attached to the vehicle has the velocity vector “V” in the direction tangential to the line toward the center of cornering which is not exactly aligned in the x_b direction. The existence of the side velocity, or the y_b component of V, or v_by, is largely noticeable for a driver of a van or pickup which has a long body from the front to end.

[0063] In FIG. 8, the analytical value of v_by is derived. Define the following parameters:

[0064] ω_bz: vehicle's directional angular rate with respect to the z_b-axis

[0065] @: sideslip angle between the x_b direction and V direction

[0066] d: distance of the sensor position with respect to the rear wheel axis

[0067] R: the radius of cornering

[0068] V: magnitude of V

Using these parameters defined, the following theoretical derivations hold true

[0069] $\begin{array}{cc}\begin{array}{c}{v}_{\mathrm{by}}=V\sin @\\ =\left(R{\omega }_{\mathrm{bz}}\right)\sin @\end{array}\mathrm{Meanwhile}\sin @=\frac{d}{R}\mathrm{Substituting}\mathrm{this}\mathrm{into}{v}_{\mathrm{by}}=\left(R{\omega }_{\mathrm{bz}}\right)\sin @,\mathrm{we}\mathrm{have}{v}_{\mathrm{by}}=\left(R{\omega }_{\mathrm{bz}}\right)\frac{d}{R}\mathrm{or}{v}_{\mathrm{by}}=d{\omega }_{\mathrm{bz}}& \left(A\right)\end{array}$

[0070] To incorporate this analytical condition into the Kalman filter algorithm for avoiding need of measuring the distance “d” by hand, it is possible to

[0071]

1. incorporate the constant parameter “d” into the INS navigation states and the Kalman
filter's states		METHOD 1
2. use the following auxiliary measurement in the Kalman filter's calibration process in
addition to GPS measurement
	0 = vby − dωbz	METHOD 2

These methods will be further described in the following sections.

[0072] FIG. 9 shows that only when the sensor IMU is placed on the rear wheel axis, zero distance of the sensor position is achieved, which makes the analytical condition to reduce to simplified “zero lateral velocity”:
0=v_by  (A)′
This special condition is within the scope of the analytical condition of Equation (A) which can be achieved by d=0. FIG. 10 shows that, in many cases, placing a sensor IMU at the bottom of the trunk of the vehicle will yield zero distance (d=0) of the sensor position with respect to the rear wheel axis.
INS Technology

[0073] Next, description will be made regarding the INS technology to update navigation state estimates including platform position, velocity, and orientation. This section corresponds to the function of the INS (inertial navigation system) 20 in the block diagram of FIG. 1A and the step 112 in the flowchart of FIG. 2. In the following equations, a dot above a variable, e.g., if {dot over (v)}_s, represents the time-rate of the variable v_s.
{dot over (v)}_s=ω_s×v_s+a_s+g₅  (1) Velocity Rate Equation
in which

[0074] $v_s=\left[\begin{array}{c}{v}_{\mathrm{sx}}\\ v_{\mathrm{sy}}\\ v_{\mathrm{sz}}\end{array}\right]\text{:}\mathrm{velocity}\mathrm{vector}\mathrm{with}\mathrm{respect}\mathrm{to}\mathrm{the}\mathrm{sensor}\text{-}\mathrm{fixed}$ $\mathrm{coorinated}\mathrm{system}$

[0075] ${\omega }_s=\left[\begin{array}{c}{\omega }_{\mathrm{sx}}\\ {\omega }_{\mathrm{sy}}\\ {\omega }_{\mathrm{sz}}\end{array}\right]\text{:}\mathrm{three}\text{-}\mathrm{axis}\mathrm{gyro}\mathrm{output}\mathrm{vector}\mathrm{with}\mathrm{respect}\mathrm{to}\mathrm{the}$ $\mathrm{sensor}\text{-}\mathrm{fixed}\mathrm{coordinate}\mathrm{system}$

[0076] $a_s=\left[\begin{array}{c}{a}_{\mathrm{sx}}\\ a_{\mathrm{sy}}\\ a_{\mathrm{sz}}\end{array}\right]\text{:}\mathrm{three}\text{-}\mathrm{axis}\mathrm{accelerometer}\mathrm{output}\mathrm{vector}\mathrm{with}\mathrm{respect}$ $\mathrm{to}\mathrm{the}\mathrm{sensor}\text{-}\mathrm{fixed}\mathrm{coordinate}\mathrm{system}$
g_s: gravity vector transformed into the sensor-fixed coordinate system
{dot over (c)}₀₀=c₀₁ω_sz−c₀₂ω_sy
{dot over (c)}₁₀=c₁₁ω_sz−c₁₂ω_sy
{dot over (c)}₂₀=c₂₁ω_sz−c₂₂ω_sy
{dot over (c)}₂₁=−c₂₀ω_sz+c₂₂ω_sx  (2) Orientation Rate Equation
in which

[0077] $\left[\begin{array}{ccc}{c}_{00}& c_{01}& c_{02}\\ c_{10}& c_{11}& c_{12}\\ c_{20}& c_{21}& c_{22}\end{array}\right]=T_{\mathrm{ns}}\text{:}\mathrm{transformation}\mathrm{matrix}\mathrm{from}\mathrm{the}\mathrm{sensor}$ $\mathrm{coordinate}\mathrm{system}\mathrm{to}\mathrm{the}\mathrm{NED}\mathrm{coordinate}\mathrm{system}.$
Equation (2) is derived in the conventional INS (Inertial Navigation System) technology as shown below:

[0078] $\left[\begin{array}{ccc}{\stackrel{.}{c}}_{00}& {\stackrel{.}{c}}_{01}& {\stackrel{.}{c}}_{02}\\ {\stackrel{.}{c}}_{10}& {\stackrel{.}{c}}_{11}& {\stackrel{.}{c}}_{12}\\ {\stackrel{.}{c}}_{20}& {\stackrel{.}{c}}_{21}& {\stackrel{.}{c}}_{22}\end{array}\right]=\left[\begin{array}{ccc}{c}_{00}& c_{01}& c_{02}\\ c_{10}& c_{11}& c_{12}\\ c_{20}& c_{21}& c_{22}\end{array}\right]\left[\begin{array}{ccc}0& -{\omega }_{\mathrm{sz}}& {\omega }_{\mathrm{sy}}\\ {\omega }_{\mathrm{sz}}& 0& -{\omega }_{\mathrm{sx}}\\ -{\omega }_{\mathrm{sy}}& {\omega }_{\mathrm{sx}}& 0\end{array}\right]\text{:}$ $\mathrm{orientation}\mathrm{rate}\mathrm{equation}$

[0079] The transformation matrix, T_ns, carries the information of angles between the sensor fixed coordinate system and the NED coordinate system. When it is necessary to convert the transformation-matrix representation into Eularian angles between the sensor and NED coordinate systems, it is possible to execute the following steps (see also U.S. Pat. No. 7,957,898 “Computational Scheme for MEMS Inertial Navigation Systems” issued to Hoshizaki, T.):

[0080] $\left[\begin{array}{ccc}{C}_{E2}{C}_{E3}& -C_{E1}{S}_{E3}+S_{E1}{S}_{E2}{C}_{E3}& S_{E1}{S}_{E3}+C_{E1}{S}_{E2}{C}_{E3}\\ C_{E2}{S}_{E3}& C_{E1}{C}_{E3}+S_{E1}{S}_{E2}{S}_{E3}& -S_{E1}{C}_{E3}+C_{E1}{S}_{E2}{S}_{E3}\\ -S_{E2}& S_{E1}{C}_{E2}& C_{E1}{C}_{E2}\end{array}\right]=\hspace{1em}\left[\begin{array}{ccc}{c}_{00}& c_{01}& c_{02}\\ c_{10}& c_{11}& c_{12}\\ c_{20}& c_{21}& c_{22}\end{array}\right]\begin{array}{ccc}\mathrm{Step}1.& \sin \left(E_2\right)=-c_{20}& \cos \left(E_2\right)=\sqrt[]{1-{\sin }^2\left(E_2\right)}>0\\ \mathrm{Step}2.& \sin \left(E_1\right)=\frac{c_{21}}{\cos \left(E_2\right)}& \cos \left(E_1\right)=\sqrt[]{1-{\sin }^2\left(E_1\right)}>0\\ \mathrm{Step}3.& \sin \left(E_3\right)=\frac{c_{10}}{\cos \left(E_2\right)}& \cos \left(E_3\right)=\frac{c_{00}}{\cos \left(E_2\right)}\end{array}$
in which
E₃: yaw angle of the sensor coordinate system with respect to the NED coordinate system
E₂: pitch angle of the sensor coordinate system with respect to the NED coordinate system
E₁: roll angle of the sensor coordinate system with respect to the NED coordinate system
C_E1, S_E1 . . . : cos(E₁), sin(E₁), and so on
In the above steps, the following automotive platform conditions are assumed:

• • [0081] −90 deg<E₂<+90 deg −90 deg<E₁<+90 deg
    {dot over (N)}=v_nx
    {dot over (E)}=v_ny
    {dot over (D)}=v_nz  (3) Position Rate Equation
    in which
    N; northerly displacement
    E: easterly displacement
    D: downward displacement

[0082] $\left[\begin{array}{c}{v}_{\mathrm{nx}}\\ v_{\mathrm{ny}}\\ v_{\mathrm{nz}}\end{array}\right]=v_n=T_{\mathrm{ns}}{v}_s=\hspace{1em}\left[\begin{array}{ccc}{c}_{00}& c_{01}& c_{02}\\ c_{10}& c_{11}& c_{12}\\ c_{20}& c_{21}& c_{22}\end{array}\right]\left[\begin{array}{c}{v}_{\mathrm{sx}}\\ v_{\mathrm{sy}}\\ v_{\mathrm{sz}}\end{array}\right]\text{:}\mathrm{velocity}\mathrm{vector}\mathrm{transformed}\mathrm{into}\mathrm{NED}\mathrm{coordinate}\mathrm{system}.$
The following equations represent constant dynamics. Sensor biases are assumed constant as follows, although they drift slowly according to the temperature change in reality.
{dot over (b)}_ωx=0
{dot over (b)}_ωy=0
{dot over (b)}_ωz=0  (4-1)
{dot over (b)}_ax=0
{dot over (b)}_ay=0
{dot over (b)}_az=0  (4-2)
{dot over (p)}₀₀=0
{dot over (p)}₁₀=0
{dot over (p)}₂₀=0  (4-3)
{dot over (d)}=0  (4-4)
in which

[0083] $\left[\begin{array}{c}{b}_{\omega x}\\ b_{\omega y}\\ b_{\omega z}\end{array}\right]=B_{\omega }\text{:}\mathrm{gyro}\mathrm{bias}\mathrm{vector}\mathrm{with}\mathrm{respect}\mathrm{to}\mathrm{sensor}\mathrm{fixed}$ $\mathrm{coordinate}\mathrm{system}$

[0084] $\left[\begin{array}{c}{b}_{\mathrm{ax}}\\ b_{\mathrm{ay}}\\ b_{\mathrm{az}}\end{array}\right]=B_a\text{:}\mathrm{accelerometer}\mathrm{bias}\mathrm{vector}\mathrm{with}\mathrm{respect}\mathrm{to}\mathrm{sensor}$ $\mathrm{fixed}\mathrm{coordinate}\mathrm{system}$

[0085] $\hspace{1em}\left[\begin{array}{ccc}{p}_{00}& p_{01}& p_{02}\\ p_{10}& p_{11}& p_{12}\\ p_{20}& p_{21}& p_{22}\end{array}\right]=T_{\mathrm{bs}}\text{:}\mathrm{transformation}\mathrm{matrix}\mathrm{from}\mathrm{the}\mathrm{sensor}\mathrm{coordinated}\text{-}\mathrm{fixed}\mathrm{system}\mathrm{to}\mathrm{vehicle}\mathrm{body}\mathrm{fixed}\text{}\mathrm{coordinate}\mathrm{system}$
d: distance between the sensor IMU position and the rear wheel axis

[0086] The transformation matrix, T_hs, carries the information of angles between the sensor fixed coordinate system and the vehicle body fixed coordinate system. When it is necessary to convert the transformation-matrix representation into Eulerian angles between the sensor and vehicle body fixed coordinate systems, it is possible to execute the following steps (see also U.S. Pat. No. 7,957,898 “Computational Scheme for MEMS Inertial Navigation Systems” issued to Hoshizaki, T.):

[0087] $\left[\begin{array}{ccc}{C}_{A2}{C}_{A3}& -S_{A3}& S_{A2}{C}_{A3}\\ C_{A2}{S}_{A3}& C_{A3}& S_{A2}{S}_{A3}\\ -S_{A2}& 0& C_{A2}\end{array}\right]=\hspace{1em}\left[\begin{array}{ccc}{p}_{00}& p_{01}& p_{02}\\ p_{10}& p_{11}& p_{12}\\ p_{20}& 0& p_{22}\end{array}\right]\begin{array}{ccc}\mathrm{Step}1.& \sin \left(A_2\right)=-p_{20}& \cos \left(E_2\right)=\sqrt[]{1-{\sin }^2\left(A_2\right)}>0\\ \mathrm{Step}2.& \sin \left(A_3\right)=\frac{p_{10}}{\cos \left(A_2\right)}& \cos \left(E_1\right)=\sqrt[]{1-{\sin }^2\left(A_3\right)>0}\end{array}$
in which
A₃: yaw angle of the sensor coordinate system with respect to the vehicle body fixed coordinate system
A₂: pitch angle of the sensor coordinate system with respect to the vehicle body fixed coordinate system
A₁: roll angle of the sensor coordinate system with respect to the vehicle body fixed coordinate system
C_A2, S_A2, . . . : cos(A₂), sin(A₂), and so on
In the above steps, the following practical conditions are assumed:
A₁=0 −90 deg<A₂<+90 deg −90 deg<A₃<+90 deg

[0088] Summarizing Equations (1) through (4) reduces to a vector representation:
{dot over (x)}=f(x,ω_s,a_s)  (5)
in which the non-linear state vector is defined by
x=[v_sx,v_sy,v_sz,N,E,D,c₀₀,c₁₀,c₂₀,c₂₁,b_ωx,b_ωy,b_ωz,b_ax,b_ay,b_ax,p₀₀,p₁₀,p₂₀,d]
where incorporation of “d” into the navigation states is one of the unique methods of the navigation system as mentioned earlier as METHOD 1.

[0089] To time-integrate Equation (5) on a processor, the following well-known Runge-Kutta 4th order equation is used (see Kreyszig, E., “Advanced Engineering Mathematics”, John Wiley & Sons, 1999, New York, N.Y., pp. 947-948):

[0090] $\begin{array}{cc}A=T\times f\left(x_k,{\omega }_{s,k},a_{s,k}\right)B=T\times f\left(x_k+\frac{1}{2}A,{\omega }_{s,k},a_{s,k}\right)C=T\times f\left(x_k+\frac{1}{2}B,{\omega }_{s,k},a_{s,k}\right)D=T\times f\left(x_k+C,{\omega }_{s,k},a_{s,k}\right)x_{k+1}=x_k+\frac{1}{6}\left(A+2B+2C+D\right)& \left(6\right)\end{array}$
in which
T: sampling time, e.g., 0.04 sec for 25 Hz
x_k value of x at the k-th time epoch of t=t_k=T×k
The INS technology described above is MEMS based simplified INS method based on conventional INS technology with simplification of small terms such as Earth rotation and Earth curvature (see U.S. Pat. No. 7,957,898 “Computational Scheme for MEMS Inertial Navigation Systems” issued to Hoshizaki, T.).
Kalman Filter Technology

[0091] The Kalman filter technology is described which is used to calibrate the INS estimates based upon reference measurement. This section corresponds to the functions of the Kalman filter 30 in the block diagram of FIG. 1A and the step 114 in the flowchart of FIG. 2.

Kalman Filter State Equation:

[0092] The first step of Kalman filter implementation is to linearize Equation (5) around a set of particular estimates of {circumflex over (x)}_k to approximate the dynamics of the small error δx with respect to the currently known estimates, {circumflex over (x)}_k. Here, a parameter with a hat represents that it is an estimate of the parameter, e.g., {circumflex over (x)}_k is the estimated amount of the parameter, x_k. Since the linearized equation remains accurate only for small value of δx around {circumflex over (x)}_k, it is called “small perturbation equation”.
δx_k+1=F({circumflex over (x)}_k)δx_kΓ_k({circumflex over (x)}_kw_k  (7)
in which the Kalman filter's state vector (small perturbation vector) is given by
δx=[δv_sx,δv_sy,δv_sz,δN,δE,δD,δα,δβ,δγ,b_ωx,b_ωy,b_ωz,b_axb_ayb_az,δb,δc,δd]
where incorporation of “d” into the Kalman filter states is one of the unique methods of this navigation system as mentioned earlier as METHOD 1.
{circumflex over (x)}_k: estimated value of X_k at t=t_k
x_k={circumflex over (x)}_k+δx: relationship between the estimated value, {circumflex over (x)}_k, and the exact value, x_k

[0093] $w_k=\left[\begin{array}{c}{w}_{\omega x}\\ w_{\omega y}\\ w_{\omega z}\\ w_{\mathrm{ax}}\\ w_{\mathrm{ay}}\\ w_{\mathrm{az}}\end{array}\right]\text{:}\mathrm{input}\mathrm{noise}\mathrm{vector}\mathrm{regarding}\mathrm{three}\text{-}\mathrm{axis}\mathrm{gyro}$ $\mathrm{output},{\omega }_s,\mathrm{and}\mathrm{three}\text{-}\mathrm{axis}\mathrm{accelerometer}\mathrm{output},\text{}{a}_s,\mathrm{with}\mathrm{white}\mathrm{model}\mathrm{in}\mathrm{the}\mathrm{discrete}\mathrm{time}\mathrm{space}$
Standard deviation (σ) of each white noise is defined as follows:

[0094] σ_ωx=N_ωx

[0095] σ_ωy=N_ωy

[0096] σ_ωz=N_ωz

[0097] σ_ax=N_ax

[0098] σ_ay=N_ay

[0099] σ_ax=N_ax

[0100] Note that these statistical noise specifications are defined in the discrete time space to be used for Equation (7). Statistical noise specification can be obtained by taking sensor measurements at the designated frequency, at 25 Hz for the above example, for enough time in the static condition.

[0101] Here, “δα, δβ, δγ” are small perturbations of “c₀₀, c₁₀, c₂₀, c₂₁” where “δα, δβ, δγ” and “δc₀₀, δc₁₀, δc₂₀, δC₂₁” have the following relationship:

[0102] $\begin{array}{cc}\hspace{1em}\left[\begin{array}{ccc}\delta c_{00}& \delta c_{01}& \delta c_{02}\\ \delta c_{10}& \delta c_{11}& \delta c_{12}\\ \delta c_{20}& \delta c_{21}& \delta c_{22}\end{array}\right]=-\left[\begin{array}{ccc}0& -\mathrm{\delta \gamma }& \mathrm{\delta \beta }\\ \mathrm{\delta \gamma }& 0& -\mathrm{\delta \alpha }\\ -\mathrm{\delta \beta }& \mathrm{\delta \alpha }& 0\end{array}\right]\hspace{1em}\left[\begin{array}{ccc}{c}_{00}& c_{01}& c_{02}\\ c_{10}& c_{11}& c_{12}\\ c_{20}& c_{21}& c_{22}\end{array}\right]& \left(8\right)\end{array}$
Similarly, “δb, δc” are small perturbations of “p₀₀, p₁₀, p₂₀,” where “δb, δc” and “δp₀₀, δp₁₀, δp₂₀” have the following relationship:

[0103] $\begin{array}{cc}\left[\begin{array}{ccc}\delta p_{00}& \delta p_{01}& \delta p_{02}\\ \delta p_{10}& \delta p_{11}& \delta p_{12}\\ \delta p_{20}& \delta p_{21}& \delta p_{22}\end{array}\right]=-\left[\begin{array}{ccc}0& -\delta c& \delta b\\ \delta c& 0& -\delta a\\ -\delta b& \delta a& 0\end{array}\right]\hspace{1em}\left[\begin{array}{ccc}{p}_{00}& p_{01}& p_{02}\\ p_{10}& p_{11}& p_{12}\\ p_{20}& p_{21}& p_{22}\end{array}\right]& \left(9\right)\end{array}$
Here, it is assumed that there is no roll angle of the sensor-fixed coordinate system with respect to the vehicle-fixed coordinate system (A₁=0), so as its small perturbation δa=0).

[0104] Using the first order approximation, the matrices of F({circumflex over (x)}_k) and Γ_k({circumflex over (x)}_k) are given by the following Equations (10) and (11). The matrix element with no indication means 0. The following notations are used in the matrix representations:

[0105] $I=\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right]\text{:}3\mathrm{by}3\mathrm{identity}\mathrm{matrix}$

[0106] $g_n=\left[\begin{array}{c}0\\ 0\\ 9.8m\text{/}{s}^2\end{array}\right]\text{:}\mathrm{gravity}\mathrm{vector}\mathrm{with}\mathrm{respect}\mathrm{to}\mathrm{the}\mathrm{NED}$ $\mathrm{coordinate}\mathrm{system}$

[0107] $\mathrm{Rot}\left({\omega }_s\right)=\left[\begin{array}{ccc}0& -{\omega }_{\mathrm{sz}}& {\omega }_{\mathrm{sy}}\\ {\omega }_{\mathrm{sz}}& 0& -{\omega }_{\mathrm{sx}}\\ -{\omega }_{\mathrm{sy}}& {\omega }_{\mathrm{sx}}& 0\end{array}\right],\mathrm{where}{\omega }_s=\left[\begin{array}{c}{\omega }_{\mathrm{sx}}\\ {\omega }_{\mathrm{sy}}\\ {\omega }_{\mathrm{sz}}\end{array}\right]\text{:}3\mathrm{by}3\mathrm{matrix}\mathrm{formation}\mathrm{using}a$ $3\mathrm{by}1\mathrm{vector}\mathrm{of}{\omega }_s$
T: sampling time, e.g., T=0.04 sec for 25 Hz

[0108]

(10)
F({circumflex over (x)}k)
	δvsx	δvsy	δvsz	δN	δE	δD	δα	δβ	δγ	bωx	bωy	bωz	bax	bay	baz	δb	δc	δd
δvsx	I-Rot(ωs)T				−{circumflex over (T)}snRot(gn)T	Rot({circumflex over (v)}s)T	IT
δvsy
δvsz
δN	{circumflex over (T)}nsT	I	Rot({circumflex over (v)}n)T
δE
δD
δα							I	-{circumflex over (T)}nsT
δβ
δγ
bωx										I
bωy
bωz
bax													I
bay
baz
δb																I
δc
δd

[0109]

Γ({circumflex over (x)}k)
(11)
	wωx, wωy, wωz	wax, way, waz
	δvxb	Rot({circumflex over (v)}s)	I
	δvyb
	δvzb
	δN
	δE
	δD
	δα	−{circumflex over (T)}ns
	δβ
	δγ
	bωxb
	bωyb
	bωzb
	baxb
	bayb
	bazb
	δb
	δc
	δd

[0110] While Equation (6) updates navigation states at a high frequency, uncertainties of the navigation estimates accumulate over time. The uncertainties of the navigation states can be mathematically represented by a covariance matrix as follows:

[0111] $P=E\left[{\begin{array}{cc}\delta x& \delta x\end{array}}^T\right]=\left[\begin{array}{cccc}E\left[\delta v_{\mathrm{sx}}\delta v_{\mathrm{sx}}\right]& E\left[\delta v_{\mathrm{sx}}\delta v_{\mathrm{sy}}\right]& E\left[\delta v_{\mathrm{sx}}\delta v_{\mathrm{sz}}\right]& \\ E\left[\delta v_{\mathrm{sy}}\delta v_{\mathrm{sx}}\right]& E\left[\delta v_{\mathrm{sy}}\delta v_{\mathrm{sy}}\right]& E\left[\delta v_{\mathrm{sy}}\delta v_{\mathrm{sz}}\right]& \dots \\ E\left[\delta v_{\mathrm{sz}}\delta v_{\mathrm{sx}}\right]& E\left[\delta v_{\mathrm{sz}}\delta v_{\mathrm{sy}}\right]& E\left[\delta v_{\mathrm{sz}}\delta v_{\mathrm{sz}}\right]& \\ & ⋮& & \ddots \end{array}\right]\text{:}$ $\mathrm{covariance}\mathrm{matrix}\mathrm{defined}\mathrm{for}\mathrm{small}\mathrm{perturbation}\mathrm{states}$
The covariance matrix must be also updated along with Equation (6) at a high frequency according to the following equation so that the Kalman filter method can calibrate the navigation state estimates:
P_k+1⁻F({circumflex over (x)}_k)P_k⁻F({circumflex over (x)}_k)^T+Γ({circumflex over (x)}_k)Q_kΓ({circumflex over (x)}_k)^T  (12)
in which
superscript “T”: transpose of the matrix
superscript “−”: before the Kalmanf filter calibration at the time-epoch of t_k
superscript “+”: after the Kalmanf filter calibration at the time-epoch of t_k

[0112] $Q_k=E\left[\begin{array}{cc}{w}_k& w_k^T\end{array}\right]=\left[\begin{array}{cccccc}{N}_{\omega x}^2& 0& 0& 0& 0& 0\\ 0& N_{\omega y}^2& 0& 0& 0& 0\\ 0& 0& N_{\omega z}^2& 0& 0& 0\\ 0& 0& 0& N_{ax}^2& 0& 0\\ 0& 0& 0& 0& N_{ay}^2& 0\\ 0& 0& 0& 0& 0& N_{\mathrm{az}}^2\end{array}\right]$
covariance matrix of w_k
Kalman Filter Measurement Equation:

[0113] The second step of Kalman filter implementation is to find the relationship between available reference measurements and the Kalman filter states. Such a relationship is called measurement equation. A measurement equation is derived for each measurement in the following.

[0114] Here, further analysis is present with respect to the analytical condition that has been derived from the vehicle mechanism noted above. This subsection will complete the details of function of the Aux (auxiliary measurement unit) 50 in the block diagram of FIG. 1A and the step 113 in the flowchart of FIG. 2. In the foregoing discussion, the following analytical condition has been achieved:
V_bydω_bz
To incorporate this condition into the Kalman filter algorithm, an auxiliary measurement of z₁ is incorporated as
z₁=v_by−dω_bz whose reference value is always z₁=0  (13-1)
where incorporation of auxiliary measurement equation (13-1) into the Kalman filter measurement is one of the unique methods of this navigation system as mentioned earlier as METHOD 2.

[0115] This is a non-linear measurement equation. To incorporate this into the Kalman filter algorithm, the equation is linearized in terms of a set of known estimates of {circumflex over (x)}_k. After a rigorous mathematical derivation, the linear perturbation equation is found as:

[0116] $\delta v_s=\left[\begin{array}{c}\delta v_{\mathrm{sx}}\\ \delta v_{\mathrm{sy}}\\ \delta v_{\mathrm{sz}}\end{array}\right]$ $B_{\omega }=\left[\begin{array}{c}{b}_{\omega x}\\ b_{\omega y}\\ b_{\omega z}\end{array}\right]$ $\Delta =[\begin{array}{c}\delta a\\ \delta b\\ \delta c\end{array}]$ $w_{\omega }=\left[\begin{array}{c}{w}_{\omega x}\\ w_{\omega y}\\ w_{\omega z}\end{array}\right]$
Similarly, since a vehicle is predominantly attached to the road surface, it is also true that:
v_bz=0
To incorporate this condition into the Kalman filter algorithm, an auxiliary measurement of z₂ is incorporated as
z₂=v_bz whose reference value is always z₂=0  (14-1)
This is the 3rd row of the following vector equation in terms of the navigation states:
z=T_bsv_x
This is a non-linear measurement equation. To incorporate this into the Kalman filter algorithm, the non-linear equation is linearized in terms of a set of particular estimates of {circumflex over (x)}_k. After a short derivation, the linear perturbation equation is found as the 3rd row of:
δz=Rot({circumflex over (v)}_b)Δ+{circumflex over (T)}_bsδv_s
or,
δz₂=[{circumflex over (p)}₂₀0{circumflex over (p)}₂₂]δv_s+[−{circumflex over (v)}_by{circumflex over (v)}_bx0]Δ  (14-2)

[0117] Next is to analyze the correlation between GPS measurements and the Kalman filter states. This subsection corresponds to the functions of the GPS 60 in the block diagram of FIG. 1A and the step 115 in the flowchart of FIG. 2.

[0118] GPS position measurement gives us absolute measurement (i.e., estimates with bounded error) of N, E, and D which can be converted from latitude, longitude, and altitude (see U.S. Pat. No. 7,957,898 “Computational Scheme for MEMS Inertial Navigation Systems” issued to Hoshizaki, T.). Therefore, the measurement equation will be
z_p=p_n  (15-1)
δz_p=δp_n  (15-2)
in which

[0119] $p_n=\left[\begin{array}{c}N\\ E\\ D\end{array}\right]$ $\delta p_n=\left[\begin{array}{c}\delta N\\ \delta E\\ \delta D\end{array}\right]$

[0120] GPS velocity measurement gives us absolute measurement of

[0121] $v_n=\left[\begin{array}{c}{v}_{\mathrm{nx}}\\ v_{\mathrm{ny}}\\ v_{\mathrm{nz}}\end{array}\right].$
Therefore, the measurement equation will be
z_v=v_n=T_nsv_s  (16-1)
This is a non-linear measurement equation. To incorporate this into the Kalman filter algorithm, the non-linear equation is linearlized in terms of a set of particular estimates of {circumflex over (x)}_k to obtain the following equation.
δz_v=Rot({circumflex over (v)}_n)ε+{circumflex over (T)}_nsδv_x  (16-2)

[0122] Summarizing Equations (13-1), (14-1), (15-1), and (16-1) reduces to a vector representation of:
z_k=h_k(x_k)  (17)
in which

[0123] $z_k=\left[\begin{array}{c}{z}_1\\ z_2\\ z_p\\ z_v\end{array}\right],h_k\left(x_k\right)=\left[\begin{array}{c}{v}_{\mathrm{by}}-d{\omega }_{\mathrm{bz}}\\ v_{\mathrm{bz}}\\ p_n\\ v_n\end{array}\right]$
Summarizing Equations of (13-2), (14-2), (15-2), and (16-2) reduces to a vector representation of:
δz_k=H({circumflex over (x)}_k)δx_k+n_k  (17-2)
in which

[0124] $\delta z_k=\left[\begin{array}{c}\delta z_1\\ \delta z_2\\ \delta z_p\\ \delta z_v\end{array}\right]$ $n_k=\left[\begin{array}{c}{n}_{\mathrm{by}}\\ n_{\mathrm{bz}}\\ n_p\\ n_v\end{array}\right]$
n_k is a measurement error vector which is assumed to be white noise. The size of each measurement error is described in the following:
σ_nby√{square root over ({circumflex over (d)}²{circumflex over (p)}₂₀²N_ωx²+{circumflex over (d)}²{circumflex over (p)}₂₂²N_ωz²)}: standard deviation for n_by derived from Equation (13-2);

[0125] This is a design parameter and can be adjusted by investigation of measurement residuals (i.e., differences between reference measurements and state estimates).

[0126] A constant parameter of σ_nby=0.05 (m/s) is also a good candidate.

σ_nbz: standard deviation for n_bz; This is a design parameter and can be adjusted by investigation of measurement residuals. σ_nbz=0.05 (m/s) is a good candidate.

[0127] $\left[\begin{array}{c}{\sigma }_N\\ {\sigma }_E\\ {\sigma }_D\end{array}\right]\text{:}\mathrm{standard}\mathrm{deviations}\mathrm{for}\mathrm{GPS}\mathrm{velocity}\mathrm{measurement}$ $\mathrm{errors}\left(n_p\right)\mathrm{to}\mathrm{be}\mathrm{given}\mathrm{by}\mathrm{the}\mathrm{GPS}\mathrm{receiver}\mathrm{in}\mathrm{real}\text{-}\mathrm{time}$

[0128] $\left[\begin{array}{c}{\sigma }_{\mathrm{vnx}}\\ {\sigma }_{\mathrm{vny}}\\ {\sigma }_{\mathrm{vnz}}\end{array}\right]\text{:}\mathrm{standard}\mathrm{deviations}\mathrm{for}\mathrm{GPS}\mathrm{velocity}\mathrm{measurement}$ $\mathrm{errors}\left(n_v\right)\mathrm{to}\mathrm{be}\mathrm{given}\mathrm{by}\mathrm{the}\mathrm{GPS}\mathrm{receiver}\mathrm{in}\mathrm{real}\text{-}\mathrm{time}$

[0129] The measurement matrix H({circumflex over (x)}_k) is given in the following Equation (18). The matrix element with no indication means 0. Notice that:

• [0130] (i) H₁ is always available no matter if GPS signals are available or not. IEKF calibration is always executed based on H₁ at an intermediate frequency, e.g., 5 Hz.
• [0131] (ii) IEKF calibration based on GPS measurement is executed at 1 Hz using H₂ only when GPS signals are available.

[0132] Execution of IEKF calibration based on auxiliary measurement equation for navigation accuracy enhancement is another unique feature of the embodiments.

[0133]

	(18)
	H({circumflex over (x)}k)
		δvxb	δvyb	δvzb	δN	δE	δD	δα	δβ	δγ	bωx	bωy	bωz	bax	bay	baz	δb	δc	δd
H1	δzvby	{circumflex over (p)}10	{circumflex over (p)}11	{circumflex over (p)}11							{circumflex over (d)}{circumflex over (p)}20	0	{circumflex over (d)}{circumflex over (p)}22				-{circumflex over (d)}{circumflex over (ω)}bx	-{circumflex over (V)}bx	_{circumflex over (ω)}bz
	δvvbz	{circumflex over (p)}20	0	{circumflex over (p)}22													{circumflex over (V)}bx	0
H2	δzp				1
						1
							1
	δzv	ĉ00	ĉ01	ĉ02				0	-{circumflex over (V)}nz	-{circumflex over (V)}ny
		ĉ10	ĉ11	ĉ12				-{circumflex over (V)}nz	0	-{circumflex over (V)}nx
		ĉ20	ĉ21	ĉ22				-{circumflex over (V)}ny	-{circumflex over (V)}nx	0

[0134] When measurements are available, the following iterative calibration steps are executed for i=0 to “n” (see Gelb, A., Applied Optimal Estimation, THE M.I.T. PRESS, 1974, Cambridge, Mass., pp. 190-191). In many cases, n=1 or 2 gives great improvement compared to no iteration (n=0). There are few cases in which more that 10-time iterations are required.
k_k,i=P_k⁻H_k^T({circumflex over (x)}_k,i⁺)(H_k({circumflex over (x)}_k,i⁺)P_k⁻H_k^T({circumflex over (x)}_k,i⁺)+R_k)⁻¹  (19) Computation of Kalman Gain, K
{circumflex over (x)}_k,i+1⁺={circumflex over (x)}_k⁻+K_k,i[z_k−h_k({circumflex over (x)}_k,i⁺)−H_k({circumflex over (x)}_k,i⁺)({circumflex over (x)}_k⁻−{circumflex over (x)}_k,i⁺)]  (20) Calibration of State Estimates, x
P_k,i+1⁺=(I−K_k,iH_k({circumflex over (x)}_k,i⁺))P_k⁻  (21) Calibration of Covariance, P
in which

• [0135] {circumflex over (x)}_k,0⁺={circumflex over (x)}_k⁻ sign in the superscript represents that the parameter is calibrated, the “−” sign in the superscript represents that the parameter is not calibrated yet
• [0136] R_k: covariance matrix regarding measurements
  Since the size of H({circumflex over (x)}_k) dynamically changes, the associated R_k also changes according to the following:

[0137] $R_k=\left[\begin{array}{cc}{\sigma }_{\mathrm{nby}}^2& 0\\ 0& {\sigma }_{\mathrm{nbz}}^2\end{array}\right]\text{:}\mathrm{for}\mathrm{Auxiliary}\mathrm{Measurement}$

[0138] $R_k=\left[\begin{array}{cccccc}{\sigma }_N^2& 0& 0& 0& 0& 0\\ 0& {\sigma }_E^2& 0& 0& 0& 0\\ 0& 0& {\sigma }_D^2& 0& 0& 0\\ 0& 0& 0& {\sigma }_{\mathrm{vnx}}^2& 0& 0\\ 0& 0& 0& 0& {\sigma }_{\mathrm{vny}}^2& 0\\ 0& 0& 0& 0& 0& {\sigma }_{\mathrm{vnz}}^2\end{array}\right]\text{:}\mathrm{for}\mathrm{GPS}\mathrm{Measurement}$

[0139] Calibration is made not only for the navigation states but also continuously made for inertial sensor outputs in the INS computation with the latest bias estimates according to the following manner:
ω_sx⁺=ω_sx⁻+b_ωx
ω_sy⁺=ω_sy⁻+b_ωy
ω_sz⁺=ω_sz⁻+b_ωz
a_sx⁺=a_sx⁻+b_ax
a_sy⁺=a_sy⁻+b_ay
a_sz⁺=a_sz⁻+b_az
in which

[0140] ${\omega }_s=\left[\begin{array}{c}{\omega }_{\mathrm{sx}}\\ {\omega }_{\mathrm{sy}}\\ {\omega }_{\mathrm{sz}}\end{array}\right]\text{:}\mathrm{three}\text{-}\mathrm{axis}\mathrm{gyro}\mathrm{output}\mathrm{vector}\mathrm{with}\mathrm{respect}\mathrm{to}\mathrm{the}$ $\mathrm{sensor}\text{-}\mathrm{fixed}\mathrm{coordinate}\mathrm{system}$

[0141] $a_s=\left[\begin{array}{c}{a}_{\mathrm{sx}}\\ a_{\mathrm{sy}}\\ a_{\mathrm{sz}}\end{array}\right]\text{:}\mathrm{three}\text{-}\mathrm{axis}\mathrm{accelerometer}\mathrm{output}\mathrm{vector}\mathrm{with}\mathrm{respect}$ $\mathrm{to}\mathrm{the}\mathrm{sensor}\text{-}\mathrm{fixed}\mathrm{coordinate}$
Simulation Results

[0142] FIG. 11A through FIG. 18B show simulation results of the embodiments of the INS/GPS navigation system and method described above using the sensor data collected along an actual on-road drive. A mid-sized van was used for this test drive, and the distance of the sensor position with respect to the rear wheel axis is known to be about 1.7 m. Therefore, the actual v_by component (side velocity) is definitely non-zero.

[0143] FIG. 11A shows the top view of the navigation solutions with the embodiments applied to the sensor data in which a vehicle drives through a spiral passage of a three-dimensional parking garage building in Yokohama, Japan. The drive inside the parking garage lasts about 4 minutes in which there is no GPS signal available. The spiral trajectory makes circles in FIG. 11A which are desirably located approximately in the same place. FIG. 11B shows the same navigation solutions in the bird view to visualize the vehicle's three-dimensional motion. Although FIG. 11B shows the almost successful result, there is a discrepancy between the entrance height and the exit height for about 10 m which represents navigation inaccuracy.

[0144] FIGS. 12A and 12B show the top view and bird view, respectively, of the navigation solutions with the conventional technology applying the measurement of “v_by=0” without consideration of the internal geometry of the sensor position. Degradation in positioning accuracy is noticeable in FIGS. 12A and 12B compared to the positioning accuracy of the embodiments of the new navigation system and method shown in FIGS. 11A and 11B.

[0145] The primary reason of the degradation in navigation accuracy of the conventional technology is that suppressing the analytical non-zero value of “V_by=dω_bz” into 0 has caused the following unfavorable side effects: (I) shortage of the total amount in the estimate of the forward velocity; (II) erroneous increment of pitch-angle estimate; (III) large error in velocity estimation during a straight drive in either of the forward or backward direction which directly results in large positioning error especially when backing; (IV) large error in velocity estimation during cornering which results in erroneous shift of a circular path. These facts of (I) through (IV) are noticeable in the following figures.

[0146] (I) FIG. 13A shows the estimated velocity components in the vehicle body fixed axes made by the embodiments of the new navigation system and method, and FIG. 13B shows the same parameters made by the conventional technology of “v_by=0”. Notice that v_by component is much smaller in FIG. 13B than FIG. 13A because of erroneous suppression of “v_by=0” in FIG. 13B. Since v_by is a component of the total velocity, if v_by is suppressed to be smaller than the actual value, so does the total velocity, resulting in that v_bx is also smaller in FIG. 13B than FIG. 13A up to 2 m/s for the entire data period. FIG. 13C shows the difference of v_bx (vehicle's forward velocity) in FIG. 13B from the data in FIG. 13A which ensures the smaller velocity estimates in FIG. 13B than FIG. 13A.

[0147] (II) FIG. 14A shows the vehicle's pitch-angle estimate, i.e., E₂ (sensor pitch angle to the NED surface)−A₂ (sensor attachment pitch angle to vehicle) made by the embodiments of the new navigation system and method. FIG. 14B shows the same parameter made by the conventional technology of “v_by=0”. Comparison of FIGS. 14A and 14B tells that pitch angle in FIG. 14B is larger than FIG. 14A up to 2 degrees as the cost of suppressing the total velocity as described in (I). FIG. 14C shows the difference of the data in FIG. 14B from the data in FIG. 14A. It is clear that the conventional technology makes larger estimates of pitch angle up to two degrees over the new navigation system and method.

[0148] (III) When pitch-angle estimate is larger than the true value, unnecessary deceleration happens in the navigation computation not only when a vehicle is cornering but also when a vehicle is straightly driving. This often results in excessive speed estimate in backing or erroneously estimated backward motion when a vehicle is actually going forward. FIG. 15A shows a magnified view of the navigation trajectory made by the embodiments (the same as FIGS. 11A and 11B) around the vehicle's proper backing motion in the top floor of the spiral parking garage. FIG. 15B shows a magnified view of the navigation trajectory made by the conventional technology of “v_by=0” (the same as FIGS. 12A and 12) around the vehicle's excessive backing motion in the top floor of the spiral parking garage. The excessive backing motion in FIG. 15B made by the conventional technology is the primary source of large positioning error. Theoretically speaking, 2-degree error in pitch-angle estimation results in 17 m of transitional positioning error after 10 seconds in a forward or backward straight drive as calculated in the following:

[0149] $\begin{array}{c}\mathrm{Straight}\mathrm{Path}\mathrm{Positioning}\mathrm{Error}=9.8\left(m\text{/}{s}^2\right)*\sin \left(2*\mathrm{pi}\text{/}180\right)*0.5*10*10\\ =17.1\left(m\right)\\ =0.171*t^2;t=\mathrm{time}\mathrm{in}\left(\sec \right)\end{array}$
FIG. 16 shows the theoretical amount of the positioning error in a straight path due to pitch angle-estimation error of 2 degrees over time based on the equation above.

[0150] (IV) A shortage in forward speed estimation happens often cyclically in cornering as already shown in FIG. 13C which results in erroneous shift of a circular path. FIG. 17A shows a plotting of the theoretical positioning error in cornering due to cyclic velocity-estimation error. The plotting with dots represents the exact path emulating the circular motion in FIG. 11A with the cornering radius (R) of 10 m and cornering angular rate (ω) of 45 deg/sec. The plotting with “x” represents the erroneous path emulating the conventional technology of “v_by=0” in which the direction is always inclined to the left to the true angle with @=9.8 deg (equivalent to the amount with d=1.7 m and R=10 m in FIG. 8). The angle inclination distorts the circular path into an oval shape, and the cyclic error in velocity estimation causes a transitional shift of the circular path. The theoretical amount of positioning error particularly to this example is 1.54 m/cycle representing continuous shift of the circular path to go down. FIG. 17B shows the velocity profiles used to create FIG. 17A in which the plotting with dots (top line) represents the true velocity used to make the circular path and the plotting with “x” (two broken lines) represents the velocity with cyclic error used to create the moving oval trajectory. Aggregation of the errors of FIGS. 16 and 17A appears in FIGS. 12A and 12B.

[0151] FIG. 18A shows the time history of the estimated distance between the position of the sensor IMU and the vehicle's rear-wheel axis made by the Kalman filter process. FIG. 18B shows the associated uncertainty (sigma value) obtained from the covariance matrix in the Kalman filtering process according to the following equation.
σ of d=√{square root over (P[18,18])}(m)
in which P is the covariance matrix. These figures show that, the more cornering, the more calibrated the distance estimation. After undergoing the intensive cornering in a spiral parking garage, the estimation of the distance between the position of the sensor IMU and the vehicle's rear-wheel axis is well converged.
Display

[0152] In addition to the normal navigation functionalities, the integrated INS/GPS navigation system and method provides a unique display method utilizing the internal geometry of the sensor position with respect to the vehicle's rear wheel axis. As described above, the embodiments of the integrated INS/GPS navigation system and method will produce the position estimates including the internal geometry of the sensor IMU position with respect to the vehicle's rear wheel axis with high accuracy. Such position estimates will be aligned with vehicle contour information and position information of links, nodes, polygons, etc., i.e., objects (map image) surrounding the vehicle, derived from the map database. This section corresponds to the functions of the navigation operation unit 70 and the display 80 in the block diagram of FIG. 1A and the steps 118 and 119 in the flowchart of FIG. 2. As noted above, such a matching process can be conducted by the navigation operation unit 70 with use of the information from the map database 74 and the contour database 79.

[0153] FIG. 19 shows a display image with the vehicle contour and navigation system's position with proper geometry between the vehicle contour and navigation system's position in which the distance between the navigation system and the rear wheel axis is automatically estimated by the method of present invention. Showing the vehicle contour with proper geometry with respect to the navigation system as well as to the objects surrounding the vehicle visually aids a driver's safety consciousness. The navigation system's absolute position (latitude and longitude) and map database are aligned by placing the GPS antenna above the navigation system. FIG. 20 shows an instruction to a user or an installer indicating to place the GPS antenna above the navigation system.

[0154] As has been described above, the embodiments of the integrated INS/GPS navigation system and method achieve the following advantageous effects: (1) regardless of the sensor position, the distance of the sensor position with respect to the rear wheel axis will be automatically estimated without need of measuring the distance by hand, which will be utilized to enhance navigation accuracy; (2) high positioning accuracy is maintained even when GPS signals are lost for a long period of time using a low-cost MEMS IMU; (3) the best sensor position to achieve the highest navigation accuracy is the center of the rear-wheel axis which is practically available by placing the sensor IMU at the bottom-center of the trunk of the vehicle; (4) the driver's safety consciousness is enhanced by the visual aid from the display showing the vehicle contour with proper geometry with respect to the navigation system as well as to the surrounding objects.

[0155] The detailed description in the foregoing is intended as a description on examples of apparatus, method, mathematical expressions, etc., in accordance with aspects of the present invention and is not intended to represent the only form in which the present invention may be prepared or utilized. Further, although the invention is described herein with reference to the preferred embodiments, one skilled in the art will readily appreciate that various modifications and variations may be made without departing from the spirit and scope of the present invention. Such modifications and variations are considered to be within the purview and scope of the appended claims and their equivalents. For example, the present invention can also be applied to the embodiments with a reduced-axis sensor configuration, such as one-axis accelerometer and one-axis tyro, and speed information obtained from vehicle's Controller-Area Network (CAN) bus.