Autonomous and non-autonomous dynamic model based navigation system for unmanned vehicles

Abstract

A navigation system including a vehicle dynamic model (VDM) that serves as the main process model within a navigation filter is described. When used in an unmanned aerial vehicle (UAV), the navigation system may work in communication with inertial measurement units (IMUs) and environment dependent sensors such as GNSS receivers. Particularly, the navigation system is beneficial in the case of GNSS signal reception outages, where conventional IMU coasting drifts quickly. Yet, the navigation system may also be employed in other scenarios, for example during GNSS presence for improved positioning, velocity and attitude determination, or in combination with GNSS when no IMU is available by design or due to a failure. In the navigation system, a solution to VDM equations provides an estimate of position, velocity, and attitude, which can be updated within a navigation filter based on available observations, such as IMU data or GNSS measurements.

Claims

1. A navigation platform for autonomous and non-autonomous navigation of an unmanned aerial vehicle, comprising: an inertial measurement unit providing inertial measurements data; a global navigation satellite system device providing at least one of position measurement data, velocity measurement data, or receiver to satellite distance related measurements; at least one of an autopilot or a manual control device providing data of control input to actuators of the unmanned aerial vehicle; and a navigation system including: a self-calibrating vehicle dynamic model component providing linear and rotational accelerations of the unmanned aerial vehicle utilizing control input and previously estimated velocity and attitude, and a navigation filter using the vehicle dynamic model component as a main navigation process model for navigation of the unmanned aerial vehicle, the navigation filter receiving the data of control input to actuators of the unmanned aerial vehicle from the autopilot or the manual control device, the navigation filter receiving the inertial navigation data from the inertial measurement unit and updating an output of the vehicle dynamic model component with at least a portion of the inertial navigation data, the navigation filter receiving the at least one of position or velocity measurement data from the global navigation satellite system device when available and updating the output of the vehicle dynamic model component with the at least one of position or velocity measurement data when available, the navigation filter estimating position, velocity, attitude, and attitude rates of the unmanned aerial vehicle, the navigation filter estimating vehicle dynamic model parameters to calibrate the vehicle dynamic model component, the navigation filter estimating a wind velocity vector to provide the wind velocity to the vehicle dynamic model component, and the navigation filter estimating and compensating for errors in the inertial measurement unit.

2. The navigation platform of claim 1, further comprising at least one autonomous aiding device that provides navigation related data and wherein the navigation filter receives the navigation related data from the at least one autonomous aiding device and updates the output of the vehicle dynamic model component with at least a portion of the navigation related data.

3. The navigation platform of claim 1, further comprising at least one environment dependent aiding device that provides navigation related data, and wherein the navigation filter receives the navigation related data from the at least one environment dependent aiding device and updates the output of the vehicle dynamic model component with at least a portion of the navigation related data.

4. The navigation platform of claim 2, wherein the at least one autonomous aiding device includes at least one of a barometric altimeter, a magnetometer, or a speed sensor.

5. The navigation platform of claim 1, wherein the navigation system further includes at least one of a wind model and an inertial measurement unit error model, and wherein the navigation filter estimates the wind velocity using data provided by at least one of the wind model and the inertial measurement unit error model.

6. The navigation platform of claim 3, wherein the at least one environment dependent aiding device includes at least one of a radio frequency ranging device, a vision based sensor, an ultrasound sensor, or an optical flow sensor.

7. The navigation platform of claim 1, wherein the unmanned aerial vehicle is selected from the group consisting of a fixed-wing unmanned aerial vehicle and a rotary-wing unmanned aerial vehicle.

8. A navigation platform for autonomous and non-autonomous navigation of an unmanned aerial vehicle, comprising: an inertial measurement unit providing inertial navigation data; at least one of an autopilot or a manual control device providing data of control input to actuators of the unmanned aerial vehicle; and a navigation system including: a self-calibrating vehicle dynamic model component providing linear and rotational accelerations of the unmanned aerial vehicle utilizing control input and previously estimated velocity and attitude, and a navigation filter using the vehicle dynamic model component as a main navigation process model for navigation of the unmanned aerial vehicle, the navigation filter receiving the data of control input to actuators of the unmanned aerial vehicle from at least one of the autopilot or the manual control device, the navigation filter receiving the inertial navigation data from the inertial measurement unit and updating an output of the vehicle dynamic model component with at least a portion of the inertial navigation data, the navigation filter estimating position, velocity, attitude, and attitude rates of the unmanned aerial vehicle, the navigation filter estimating vehicle dynamic model parameters to calibrate the vehicle dynamic model component, the navigation filter estimating a wind velocity vector to provide the wind velocity to the vehicle dynamic model component, and the navigation filter estimating and compensating for errors in the inertial measurement unit.

9. The navigation platform of claim 8, further comprising at least one autonomous aiding device that provides navigation related data, and wherein the navigation filter receives the navigation related data from the at least one autonomous aiding device and updates the output of the vehicle dynamic model component with at least a portion of the navigation related data.

10. The navigation platform of claim 8, further comprising at least one environment dependent aiding device that provides navigation related data, and wherein the navigation filter receives the navigation related data from the at least one environment dependent aiding device and updates the output of the vehicle dynamic model component with at least a portion of the navigation related data.

11. The navigation platform of claim 9, wherein the at least one autonomous aiding device includes at least one of a barometric altimeter, a magnetometer, or a speed sensor.

12. The navigation platform of claim 8, wherein the navigation system further includes at least one of a wind model and an inertial measurement unit error model, and wherein the navigation filter estimates the wind velocity using data provided by at least one of the wind model and the inertial measurement unit error model.

13. The navigation platform of claim 10, wherein the at least one environment dependent aiding device including at least one of a radio frequency ranging device, a vision based sensor, an ultrasound sensor, or an optical flow sensor.

14. The navigation platform of claim 8, wherein the unmanned aerial vehicle is selected from the group consisting of a fixed-wing unmanned aerial vehicle and a rotary-wing unmanned aerial vehicle.

15. A navigation platform for autonomous and non-autonomous navigation of an unmanned aerial vehicle, comprising: a global navigation satellite system device providing at least one of position measurement data, velocity measurement data, or receiver to satellite distance related measurements; at least one of an autopilot or a manual control device providing data of control input to actuators of the unmanned aerial vehicle, the data of control input to actuators of the unmanned aerial vehicle; and a navigation system including: a self-calibrating vehicle dynamic model component providing linear and rotational accelerations of the unmanned aerial vehicle utilizing control input and previously estimated velocity and attitude, and a navigation filter using the vehicle dynamic model component as a main navigation process model for navigation of the unmanned aerial vehicle, the navigation filter receiving the data of control input to actuators of the unmanned aerial vehicle from at least one of the autopilot or the manual control device, the navigation filter receiving the at least one of position or velocity measurement data from the global navigation satellite system device when available and updating the output of the vehicle dynamic model component with the at least one of position or velocity measurement data when available, the navigation filter estimating position, velocity, attitude, and attitude rates of the unmanned aerial vehicle, the navigation filter estimating vehicle dynamic model parameters to calibrate the vehicle dynamic model component, and the navigation filter estimating a wind velocity vector to provide the wind velocity to the vehicle dynamic model component.

16. The navigation platform of claim 15, further comprising at least one autonomous aiding device that provides navigation related data, and wherein the navigation filter receives the navigation related data from the at least one autonomous aiding device and updates the output of the vehicle dynamic model component with at least a portion of the navigation related data.

17. The navigation platform of claim 15, further comprising at least one environment dependent aiding device that provides navigation related data and wherein the navigation filter receives the navigation related data from the at least one environment dependent aiding device and updates the output of the vehicle dynamic model component with at least a portion of the navigation related data.

18. The navigation platform of claim 17, wherein the at least one autonomous aiding device includes at least one of a barometric altimeter, a magnetometer, or a speed sensor the radio frequency ranging device is a global navigation satellite system.

19. The navigation platform of claim 15, wherein the navigation system further includes at least one of a wind model and an inertial measurement unit error model, and wherein the navigation filter estimates the wind velocity using data provided by at least one of the wind model and the inertial measurement unit error model.

20. The navigation platform of claim 17, wherein the at least one environment dependent aiding device including at least one of a radio frequency ranging device, a vision based sensor, an ultrasound sensor, or an optical flow sensor.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) Embodiments of devices, systems, and methods are illustrated in the figures of the accompanying drawings, which are meant to be exemplary and non-limiting, in which like references are intended to refer to like or corresponding parts, and in which:

(2) FIG. 1 illustrates a navigation platform of a UAV;

(3) FIG. 2 illustrates the integration of information of a vehicle dynamic model within a navigation platform of a UAV;

(4) FIG. 3 illustrates a local, body, and wind navigation system frame with airspeed (V), wind velocity (w), and UAV velocity (v);

(5) FIG. 4 illustrates a fixed-wing, single propeller UAV having aileron, elevator, and rudder control surfaces;

(6) FIG. 5 illustrates a navigation system architecture according to the present disclosure;

(7) FIG. 6 illustrates a reference trajectory and a navigation solution from a sample run with GNSS signals available during the first 100 seconds of flight in a Monte Carlo simulation;

(8) FIG. 7 illustrates a comparison between INS/GNSS and VDM/INS/GNSS with Root Mean Squared (RMS) position and yaw errors from 100 Monte Carlo runs;

(9) FIG. 8 illustrates position errors for 100 Monte Carlo runs;

(10) FIG. 9 illustrates the RMS of attitude errors from 100 Monte Carlo runs;

(11) FIGS. 10A and 10B illustrate the RMS of estimation error for IMU error terms from 100 Monte Carlo runs, together with the amplitude of simulated errors (realization values);

(12) FIG. 11 illustrates the RMS of estimation error for VDM parameters from 100 Monte Carlo runs;

(13) FIG. 12 illustrates the RMS of estimation error for wind speed from 100 Monte Carlo runs, together with real wind speed (realization value);

(14) FIGS. 13A and 13B illustrate ratios of uncertainties on states during GNSS presence and GNSS outage periods (navigation starts with GNSS available at t=0, GNSS outage begins at t=100 s, and navigation ends at t=400 s);

(15) FIG. 14 illustrates a correlation matrix at the end of navigation;

(16) FIG. 15 illustrates a reference trajectory and navigation solution from a navigation system of the present disclosure during a test flight in a 3D view;

(17) FIG. 16 illustrates a reference trajectory and navigation solution from a navigation system of the present disclosure during a test flight in a 2D view;

(18) FIG. 17 illustrates position error for a navigation system of the present disclosure during a test flight;

(19) FIG. 18 illustrates position error for an INS/GNSS navigation system during a test flight; and

(20) FIGS. 19A through 19C contain nomenclature for components of the present disclosure.

DETAILED DESCRIPTION

(21) The present disclosure provides systems, methods, and devices for increasing the accuracy and reliability of autonomous navigation of UAVs by employing a VDM as a main process model for navigation. The VDM may use input from autonomous sensors and environment dependent sensors. Autonomous sensors provide data regardless of environmental features and conditions such as texture, visibility, or signal reception condition. Illustrative autonomous sensors include IMUs, pressure sensors, and the like. On the other hand, environment dependent sensors provide data dependent upon environment features and conditions such as texture, visibility, or signal reception condition. Illustrative environment dependent sensors include radio frequency (RF) ranging sensors, GNSS receivers, and vision based sensors, and the like.

(22) FIGS. 19A through 19C contain nomenclature for components of the present disclosure.

(23) FIG. 1 illustrates a UAV navigation platform 100 wherein data from at least one autonomous sensor 102 and at least one environment dependent sensor 104 are interconnected in a navigation system 108 to provide a navigation solution 110. Time stamping of data gathered by the autonomous sensor(s) 102 and environment dependent sensor(s) 104 may be handled inside the sensors or externally (i.e., by a time stamper 106 interspersed between the sensors and the navigation system 108, as illustrated in FIG. 1). While one time stamper 106 is illustrated, it should be appreciated that more than one time stamper 106 may be used. For example, each sensor may have its own dedicated, external time stamper, autonomous sensor(s) 102 may be dedicated to a first external time stamper and environment dependent sensor(s) 104 may be dedicated to a second external time stamper, etc.

(24) FIG. 2 illustrates a UAV navigation platform 200 including a vehicle dynamic model 212 (VDM) integrated within the navigation system 108. The VDM 212 may be a process model in a Kalman filter. The navigation system 108 is fed with a control input to the UAV, which is the same data that is fed to one or more actuators 208. As illustrated, an autopilot/control unit 202, a manual control 210, and the actuator(s) 208 are in communication via switches. The actuator(s) 208 may be commanded either by the autopilot 202 or the manual control 210 via a double throw switch. The manual control 210 may communicate with the controller 206 via a single throw switch. Using this illustrative configuration of switches results in the potential for various sources of data to be sent to the navigation system 108. In an example, only data from the autopilot 202 may be fed to the actuator(s) 208 and navigation system 108, regardless of the manual control 210 data. In another example, only the manual control 210 data may be fed to the actuator(s) 208 and navigation system 108, regardless of the autopilot 202 data. In a further example, the manual control 210 data may be fed to the actuator(s) 208 and navigation system 108 through the controller 206.

(25) The autopilot 202 may include a guidance system 204 and a controller 206. In an example, data from the guidance system 204 passes through the controller 206 and is fed to the actuator(s) 208 and navigation system 108. In another example, the manual control 210 data is directly fed to the actuator(s) 208 and navigation system 108. In a further example, the manual control 210 data passes through the controller 206 (e.g., for stabilization of the UAV) and is fed to the actuator(s) 208 and navigation system 108, regardless of the guidance system 204 data.

(26) Time stamping of the autopilot 202 or manual control 210 data may be handled inside the respective devices or externally (e.g., by the time stamper 106). While it is illustrated that a single external time stamper 106 stamps all data transmitted to the navigation system 108, it should be appreciated that more than one external time stamper 106 may be used. For example, each data source may have its own dedicated, external time stamper, the sensors 102/104 may have one dedicated time stamper while all other data sources may have a separate dedicated time stamper, etc.

(27) Depending on the specific sensors used by the navigation platform 200 (i.e., implemented on the UAV), there may be certain requirements on the hardware usable within the navigation platform 200. For example, the autopilot/control unit 202 data output to the actuator(s) 208 (i.e., the control input to the navigation platform 200) may need to be accessible by the navigation system 108 (since this information is used by the VDM 212) and time-stamped with respect to other sensory data. For autonomous navigation, an IMU autonomous sensor 102 is usually a minimum sensor setup concern. However, it is possible for the navigation platform 200 to work without an IMU autonomous sensor 102, but a GNSS receiver environment dependent sensor 104 or other sensor setup allowing 3D position fixes may be required for desired results. Since an IMU autonomous sensor 102 is almost always used within a navigation platform for UAVs, omission of the IMU autonomous sensor 102 from the platform may be most beneficial as a backup solution in case of IMU autonomous sensor 102 failure, thereby increasing the reliability of the whole navigation platform. It may be desirable to properly reference the control input and the data from all the sensors to a common time scale in order for the control input and sensor data to be used within the navigation system 108. This may require a hardware/software time stamper 106 to perform data tagging with respect to a common time reference. Depending on the configuration of the sensors 102/104, a partial synchronization may be performed in other manners. For example, in the situation where both an IMU autonomous sensor 102 and a GNSS receiver environment dependent sensor 104 are available, there may be a connection between them in which the GNSS receiver environment dependent sensor 104 time frame is used to tag IMU autonomous sensor 102 data. A synchronization unit may also be available within the control unit 202.

Navigation System Frames

(28) Five frames may be used within the navigation system 108.

(29) An inertial frame (i-frame) is a non-accelerating and non-rotating reference frame that is at rest or subject to a uniform translational motion. In such a frame, the laws of Newtonian mechanics are valid. An i-frame is a celestial frame with an origin at the center of mass of Earth, such that an x.sub.i-axis of the i-frame points towards the vernal equinox and, thus, is the intersection line between the equatorial and the ecliptic plane. A z.sub.i-axis of the i-frame points towards the mean celestial pole, and a y.sub.i-axis of the i-frame completes the system to a 3D right-handed Cartesian system.

(30) An earth frame (e-frame) is an equatorial frame with an origin at the center of mass of Earth, such that an x.sub.e-axis of the e-frame points towards the Greenwich meridian, a z.sub.e-axis of the e-frame points towards the mean direction of the rotation axis of Earth, and a y.sub.e-axis of the e-frame completes the system to a 3D right-handed Cartesian system. The e-frame is therefore an Earth Centered Earth Fixed (ECEF) frame. Examples of e-frames include the International Terrestrial Reference Frame (ITRF) and the World Geodetic System WGS-84, which is the reference frame of the American Global Positioning System (GPS). When an ellipsoid of reference is associated with the frame, as in the WGS-84 datum, any point in the e-frame can be expressed either in Cartesian coordinates (x.sub.e,y.sub.e,z.sub.e) or ellipsoidal coordinates (,,h) with being the latitude, being the longitude, and h being the ellipsoidal height.

(31) A local-level frame (l-frame) is a local geodetic frame with an arbitrary origin, for example a point on Earth's surface (i.e., topocenter). An illustrative l-frame is a North-East-Down (NED) frame, which has a x.sub.l-axis that points North, a y.sub.l-axis that points East, and a z.sub.l-axis that points to the local nadir (i.e., down).

(32) A NED frame is right-handed. The nadir direction of a point is defined by its astronomical latitude and astronomical longitude, respectively noted and . If the astronomical coordinates are replaced by the ellipsoidal latitude and longitude , the z.sub.l-axis corresponds to the ellipsoidal normal, while the north direction is the tangent of the ellipsoid meridian. For the purpose of navigation, the difference between astronomical and ellipsoidal local-level frames can be neglected. The resulting frame may be referred to as an ellipsoidal tangential frame rather than a local-level frame. For purposes of the present disclosure, local-level frame, local frame, or l-frame refer to the ellipsoidal NED frame, and latitude and longitude refer to ellipsoidal quantities rather than astronomical quantities.

(33) A body frame (b-frame) is a frame attached to the body of the UAV. The b-frame may have Cartesian coordinates (x.sub.b,y.sub.b,z.sub.b). The forward or longitudinal direction of the UAV may be denoted the x.sub.b-axis, while the y.sub.b-axis points towards right and the z.sub.b-axis is downward. A rotation matrix to transform vectors from a l-frame to a b-frame is represented by Equation 1 below when Euler angles are used to express the attitude.

(34) $\begin{array}{cc}\begin{array}{c}{C}_l^b=C_1\left(r\right)C_2\left(p\right)C_3\left(y\right)\\ =\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos r& \sin r\\ 0& -\sin r& \cos r\end{array}\right]\left[\begin{array}{ccc}\cos p& 0& -\sin p\\ 0& 1& 0\\ \sin p& 0& \cos p\end{array}\right]\left[\begin{array}{ccc}\cos y& \sin y& 0\\ -\sin y& \cos y& 0\\ 0& 0& 1\end{array}\right]\end{array}& \left(\mathrm{Equation}1\right)\end{array}$

(35) The Euler angles within Equation 1 are roll (r), pitch (p), and yaw (y). Euler angles exhibit a singularity at p=/2, therefore alternative ways may be employed to express the attitude, such as quaternions.

(36) A l-frame, b-frame, and wind frame (w-frame) are described with reference to FIG. 3. An airflow, which is due to the UAV's inertial velocity v, and wind velocity w are described by an airspeed vector V. The w-frame has a first axis (i.e., x.sub.w-axis) in the direction of V, and is defined by two angles with respect to the b-frame. These angles are called the angle of attack () and the sideslip angle ().

(37) Equation 2 represents a UAV's inertial velocity and wind velocity via airspeed.
v=V+w(Equation 2)

(38) Equation 3 results when Equation 2 is applied to a body frame.

(39) $\begin{array}{cc}{v}^b=V^b+C_l^b{w}^l{V}^b=v^b-C_l^b{w}^l\left[\begin{array}{c}{V}_x^b\\ V_y^b\\ V_z^b\end{array}\right]=\left[\begin{array}{c}{v}_x^b\\ v_y^b\\ v_z^b\end{array}\right]-C_l^b\left[\begin{array}{c}{w}_x^l\\ w_y^l\\ w_z^l\end{array}\right]& \left(\mathrm{Equation}3\right)\end{array}$

(40) Equation 4 represents a rotation matrix to from a body frame to a wind frame.

(41) $\begin{array}{cc}\begin{array}{c}{C}_b^w=C_3\left(\right)C_2^T\left(\right)\\ =\left[\begin{array}{ccc}\cos & \sin & 0\\ -\sin & \cos & 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{ccc}\cos & 0& \sin \\ 0& 1& 0\\ -\sin & 0& \cos \end{array}\right].\end{array}& \left(\mathrm{Equation}4\right)\end{array}$

(42) Equation 5 represents the airspeed magnitude, the attack angle attack, the sideslip angle, and the dynamic pressure. The latter may be used to express aerodynamic forces and moments.

(43) $\begin{array}{cc}V=\sqrt{V_x^{b2}+V_y^{b2}+V_z^{b2}},=\arctan \left(\frac{V_z^b}{V_x^b}\right),=\arcsin \left(\frac{V_y^b}{V}\right),\stackrel{\_}{q}=\frac{V^2}{2}& \left(\mathrm{Equation}5\right)\end{array}$

Vehicle Dynamic Model

(44) The Vehicle Dynamic Model (VDM) of the present disclosure is a mathematical model that describes the dynamics of motion for the UAV. In order to establish the VDM, the frames in which the dynamic states are measured and the frames in which the dynamic states are resolved should be chosen. The dynamic states include position, velocity, attitude, and rotation rate (angular velocity) vectors. There are many possibilities for these choices. However, all of them carry the same amount of information and are just different representations of the same quantity. Dynamic states are represented in Equation 6.

(45) $\begin{array}{cc}{X}_n={\left[r_e^{l^T},v_e^{l^T},q_b^{l^T},{}_{\mathrm{ib}}^{b^T}\right]}^T& \left(\mathrm{Equation}6\right)\end{array}$

(46) The position vector (r.sub.e.sup.l=[,,h].sup.T) represents a position of the UAV (i.e., position of body frame origin) in the earth frame and in ellipsoidal coordinates. The velocity vector (v.sub.e.sup.l=[v.sub.N,v.sub.E,v.sub.D].sup.T) represents a velocity of the UAV as observed in the earth frame and expressed in the local (NED) frame. The q.sub.b.sup.l quaternion (q.sub.b.sup.l=[q.sub.0,q.sub.1,q.sub.2,q.sub.3].sup.T) represents an orientation of the UAV body frame with respect to the local frame. The angular velocity vector (.sub.ib.sup.b=[.sub.x,.sub.y,.sub.z].sup.T) represents an angular velocity (i.e., rotation rate) of the body frame with respect to the inertial frame, expressed in the body frame. The attitude can be represented in many different ways, such as Euler angles, rotation matrix elements, and quaternions.

(47) The VDM can be represented by Equations 7-10.

(48) $\begin{array}{cc}{\stackrel{.}{r}}_e^l=D^{-1}{v}_e^l& \left(\mathrm{Equation}7\right)\\ {\stackrel{.}{v}}_e^l=C_b^l{f}^b+g^l-\left({}_{\mathrm{el}}^l+2{}_{\mathrm{ie}}^l\right)v_e^l=C_b^l{f}^b+g^l-\left({}_{\mathrm{il}}^l+{}_{\mathrm{ie}}^l\right)v_e^l& \left(\mathrm{Equation}8\right)\\ {\stackrel{.}{q}}_b^l=\frac{1}{2}{q}_b^l.Math.{\left[{}_{\mathrm{lb}}^l\right]}_q& \left(\mathrm{Equation}9\right)\\ {\stackrel{.}{}}_{\mathrm{ib}}^b={\left(I^b\right)}^{-1}\left[M^b-{}_{\mathrm{ib}}^b\left(I^b{}_{\mathrm{ib}}^b\right)\right]& \left(\mathrm{Equation}10\right)\end{array}$

(49) In Equation 7, the D.sup.1 matrix may be represented by Equation 11.

(50) $\begin{array}{cc}{D}^{-1}=\left[\begin{array}{ccc}\frac{1}{R_M+h}& 0& 0\\ 0& \frac{1}{\left(R_p+h\right)\cos }& 0\\ 0& 0& -1\end{array}\right],& \left(\mathrm{Equation}11\right)\end{array}$

(51) R.sub.M and R.sub.P represent a meridian radius of curvature and a prime vertical radius of curvature, respectively. The meridian radius of curvature R.sub.M is represented by Equation 12 and the prime vertical radius of curvature R.sub.P is represented by Equation 13.

(52) $\begin{array}{cc}{R}_M=\frac{a\left(1-e^2\right)}{{\left(1-e^2{\sin }^2\left(\right)\right)}^{3/2}}& \left(\mathrm{Equation}12\right)\\ R_P=\frac{a}{{\left(1-e^2{\sin }^2\left(\right)\right)}^{1/2}}& \mathrm{Equation}13)\end{array}$

(53) In Equation 8, the rotation matrix C.sub.b.sup.l can be expressed as a function of quaternion q.sub.b.sup.l, as represented by Equation 14.

(54) $\begin{array}{cc}{C}_b^l=f\left(q_b^l\right)=[\begin{array}{ccc}{q}_0^2+q_1^2-q_2^2-q_3^2& 2\left(q_1{q}_2-q_0{q}_3\right)& 2\left(q_1{q}_3+q_0{q}_2\right)\\ 2\left(q_1{q}_2+q_0{q}_3\right)& q_0^2-q_1^2+q_2^2-q_3^2& 2\left(q_2{q}_3-q_0{q}_1\right)\\ 2\left(q_1{q}_3-q_0{q}_2\right)& 2\left(q_2{q}_3+q_0{q}_1\right)& q_0^2-q_1^2-q_2^2+q_3^2\end{array}]& \left(\mathrm{Equation}14\right)\end{array}$

(55) Angular velocity vectors .sub.ie.sup.l and .sub.el.sup.l are represented by Equations 15 and 16 respectively, with .sub.il.sup.l simply being the sum of Equations 15 and 16 (i.e., .sub.il.sup.l=.sub.ie.sup.l+.sub.el.sup.l)).
.sub.ie.sup.l=[.sub.ie cos 0 .sub.ie sin ].sup.T(Equation 15)
.sub.el.sup.l=[{dot over ()} cos {dot over ()} {dot over ()} sin ].sup.T(Equation 16)

(56) The gravity vector g.sup.l of Equation 8 can be determined using a gravity model for the ellipsoidal local frame. The gravity vector g.sup.l can even be considered as a constant, depending on the application and expected accuracy.

(57) In Equation 9, the symbol {circle around (x)} represents a quaternion product, which is represented by Equation 17.

(58) 0$\left(\mathrm{Equation}17\right)$$q=\left[\begin{array}{c}{q}_0\\ q_1\\ q_2\\ q_3\end{array}\right],p=\left[\begin{array}{c}{p}_0\\ p_1\\ p_2\\ p_3\end{array}\right].Math.q.Math.p=\left[\begin{array}{c}{q}_0{p}_0-q_1{p}_1-q_2{p}_2-q_3{p}_3\\ q_0{p}_1+q_1{p}_0+q_2{p}_3-q_3{p}_2\\ q_0{p}_2+q_2{p}_0-q_1{p}_3+q_3{p}_1\\ q_0{p}_3+q_3{p}_0+q_1{p}_2-q_2{p}_1\end{array}\right]$

(59) The Vehicle Dynamic Model (VDM) (as represented by Equations 7-10) may be used for any rigid body in arbitrary motion in 3D space affected by Earth gravity, external specific force f, and external moment M. Completing the model for a specific type of vehicle (e.g., a fixed wing UAV) may require proper aerodynamic models for f and M, which can be a function of dynamic states, control inputs, and physical parameters of the UAV. Physical parameters can be different for each UAV.

Navigation Filter

(60) The navigation system 108 can employ at least one of various types of estimators, such as Kalman filters, particle filters, etc. As described, the navigation system 108 includes the VDM 212 as its main model to predict navigation states of the platform. This is done, for example, by using the VDM 212 as the process model in an extended Kalman filter (EKF). As also described herein above, the VDM 212 received data from a control input (i.e., from the autopilot 202 and/or the manual control 210). It may also be beneficial to feed wind velocity data to the VDM, which can be provided by either sensory data (e.g., autonomous pressure sensors), or estimated within the navigation system 108 itself. It may further be beneficial to provide values of some physical parameters of the platform to the VDM, such as mass, moments of inertia, aerodynamic coefficients or the like. The navigation system 108 may be capable of calibrating these values online (i.e., during navigation). An adequate prior knowledge on these values can make the calibration process faster and more accurate. The navigation system 108 may fuse available sensory data with VDM prediction data to produce the final navigation solution 110.

(61) Employing the VDM as described herein brings an extra source of information into the navigation system 108 as compared to conventional navigation methods that rely only on the autonomous sensor(s) 102 and/or environment dependent sensor(s) 104. Thus, a better navigation performance results from the teachings of the present disclosure as compared to traditional methodology. The VDM reflects physical constraints on dynamic behavior of the platform, thereby removing, to some extent, physically impossible motions suggested by sensors due to their errors.

(62) Some auxiliary states may be estimated within the navigation system 108. While the auxiliary states may not be of interest from a navigation point of view, they may be needed for the navigation system 108 to function properly or with higher performance. Illustrative auxiliary states include wind velocity and some physical parameters of the platform, for example. Some part of the sensor errors (e.g., IMU bias terms) may be estimated within the navigation filter. Actuator(s) error is another possible auxiliary state to be estimated within the navigation filter. The estimated auxiliary states can be employed within the navigation filter to remove a portion of the sensor errors, thereby improving navigation performance.

Example: Fixed-Wing, Single Propeller UAV

(63) The following description involves the aforementioned teachings of the present disclosure being applied to a fixed-wing UAV embodiment. For purposes of clarity only, it may be assumed that the fixed-wing UAV includes an IMU autonomous sensor 102 and a stand-alone GNSS receiver environment dependent sensor 104. Since such limitation of sensors is only for clarity, it should be appreciated that other sensors may be present, even though not explicitly described. An Extended Kalman Filter (EKF) may be in charge of data fusion for the navigation system 108 to estimate navigation states, as well as dynamic model parameters, IMU error, and wind velocity.

Aerodynamic Forces

(64) FIG. 4 illustrates a fixed-wing, single propeller UAV 400 having aileron 402a/402b, elevator 404a/404b, and rudder 406 control surfaces. Several aerodynamic forces are experienced by the UAV 400. A thrust force is expressed in the body frame along the x.sub.b-axis while lift, lateral, and drag forces are expressed in the wind frame.

(65) The thrust force may be represented by Equation 18, wherein is the air density, n is the propeller speed, D is the propeller diameter, and C.sub.F.sub.T.sub.i is the thrust force derivatives for the particular UAV being used. (J) within Equation 18 is represented by Equation 19.

(66) $\begin{array}{cc}{F}_T=n^2{D}^4{C}_{F_T}\left(J\right),& \left(\mathrm{Equation}18\right)\\ J=\frac{V}{Dn},C_{F_T}\left(J\right)=C_{F_{T}1}+C_{F_{T}2}J+C_{F_{T}3}{J}^2.& \left(\mathrm{Equation}19\right)\end{array}$

(67) The lift force may be represented by Equation 20, wherein q is the dynamic pressure introduced in Equation 5 and S is the wing surface. C.sub.F.sub.z .sub.. . . within Equation 20 is represented by Equation 21, wherein C.sub.F.sub.z .sub.. . . the lift force derivatives for the particular UAV being used.
F.sub.z.sup.w=qSC.sub.F.sub.z(),(Equation 20)
C.sub.F.sub.z()=C.sub.F.sub.z.sub.1+C.sub.F.sub.z.sub..(Equation 21)

(68) The lateral force may be represented by Equation 22, with C.sub.F.sub.y.sub.1 of Equation 22 being represented by Equation 23, wherein C.sub.F.sub.y.sub.1 is the lateral force derivatives for the particular UAV being used.
F.sub.y.sup.w=qSC.sub.F.sub.y(),(Equation 22)
C.sub.F.sub.y()=C.sub.F.sub.y.sub.1.(Equation 23)

(69) The drag force may be represented by Equation 24, with C.sub.F.sub.x .sub.. . . of Equation 24 being represented by Equation 25, wherein C.sub.F.sub.x .sub.. . . is the drag force derivatives for the particular UAV being used.
F.sub.x.sup.w=qSC.sub.F.sub.x(,),(Equation 24)
C.sub.F.sub.x(,)=C.sub.F.sub.x.sub.1+C.sub.F.sub.x.sub.+C.sub.F.sub.x.sub.2.sup.2+C.sub.F.sub.x.sub.2.sup.2.(Equation 25)

(70) The specific force vector f, represented by Equation 26 below, is composed of these four components (i.e., the thrust force, lift force, lateral force, and drag force) divided by the mass of the UAV being used.

(71) $\begin{array}{cc}{f}^b=\frac{1}{m}\left(\left[\begin{array}{c}{F}_T\\ 0\\ 0\end{array}\right]+C_w^b\left[\begin{array}{c}{F}_x^w\\ F_y^w\\ F_z^w\end{array}\right]\right)& \left(\mathrm{Equation}26\right)\end{array}$

Aerodynamic Moments

(72) In addition to aerodynamic forces, three components of aerodynamic moments are also recognizable in the VDM of the present disclosure. These components of aerodynamic moments include a roll moment, pitch moment, and yaw moment.

(73) The roll moment is represented by Equation 27, wherein b is the wing span of the UAV being used. {tilde over ()}.sub.x within Equation 27 is represented by Equation 28. C.sub.M.sub.x .sub.. . . within Equation 28 is represented by Equation 29, wherein C.sub.M.sub.x .sub.. . . the roll moment derivatives for the particular UAV being used.

(74) $\begin{array}{cc}{M}_x^b=\stackrel{\_}{q}{\mathrm{SbC}}_{M_x}\left({}_a,,{\stackrel{~}{}}_x,{\stackrel{~}{}}_z\right),& \left(\mathrm{Equation}27\right)\\ {\stackrel{~}{}}_x=\frac{b{}_x}{2V},{\stackrel{~}{}}_y=\frac{\stackrel{\_}{c}{}_y}{2V},{\stackrel{~}{}}_z=\frac{b{}_z}{2V},& \left(\mathrm{Equation}28\right)\\ C_{M_x}\left({}_a,,{\stackrel{~}{}}_x,{\stackrel{~}{}}_z\right)=C_{M_{x}a}{}_a+C_{M_x}+C_{M_x{\stackrel{~}{}}_x}{\stackrel{~}{}}_x+C_{M_x{\stackrel{~}{}}_z}{\stackrel{~}{}}_z.& \left(\mathrm{Equation}29\right)\end{array}$

(75) The Pitch moment is represented by Equation 30, wherein c is the mean aerodynamic chord. C.sub.M.sub.y .sub.. . . within Equation 30 is represented by Equation 31, wherein C.sub.M.sub.y .sub.. . . is the pitch moment derivatives for the particular UAV being used.
M.sub.y.sup.b=qScC.sub.M.sub.y(.sub.e,,{tilde over ()}.sub.y),(Equation 30)
C.sub.M.sub.y(.sub.e,,{tilde over ()}.sub.y)=C.sub.M.sub.y.sub.1+C.sub.M.sub.y.sub.e.sub.e+C.sub.M.sub.y.sub.{tilde over ()}.sub.y{tilde over ()}.sub.y+C.sub.M.sub.y.sub.. (Equation 31)

(76) The Yaw moment is represented by Equation 32. C.sub.M.sub.z .sub.. . . within Equation 32 is represented by Equation 33, wherein C.sub.M.sub.z .sub.. . . is the yaw moment derivatives for the particular UAV being used.
M.sub.z.sup.b=qSbC.sub.M.sub.z(.sub.r,{tilde over ()}.sub.z,),(Equation 32)
C.sub.M.sub.z(.sub.r,{tilde over ()}.sub.z,)=C.sub.M.sub.z.sub..sub.r.sub.r+C.sub.M.sub.z.sub.{tilde over ()}.sub.z{tilde over ()}.sub.z+C.sub.M.sub.z.sub..(Equation 33)

(77) A moments vector M is composed of these three components (i.e., the roll moment, the pitch moment, and the yaw moment), and is represented by Equation 34.

(78) $\begin{array}{cc}{M}^b=\left[\begin{array}{c}{M}_x^b\\ M_y^b\\ M_z^b\end{array}\right]& \left(\mathrm{Equation}34\right)\end{array}$

(79) Equations 7-10 together with Equations 26-34 form the VDM for the sample fixed-wing, single propeller UAV 400 of FIG. 4. If the dynamic of the actuators is also available and intended to be considered within the navigation platform, it can be considered as a part of the VDM. In this case, for example, the rotation speed of the propeller is modeled by a first order linear differential equation, such as that represented by Equation 35. This equation may be considered as the last part of the VDM.

(80) $\begin{array}{cc}\stackrel{.}{n}=-\frac{1}{{}_n}n+\frac{1}{{}_n}{n}_c& \left(\mathrm{Equation}35\right)\end{array}$

Air Density Model

(81) Depending on accuracy requirements and fidelity of modeling, the air density may be considered as a constant or as a function of local pressure and temperature according to a suitable model for the particular UAV being used. Local pressure and temperature can be expressed as functions of the altitude. According to the International Standard Atmosphere model for low altitude, for example, the air density can be represented by Equation 36.

(82) $\begin{array}{cc}=\frac{{p_0\left[1+\mathrm{ah}/T_0\right]}^{4.2561}}{R_a{T}_0}& \left(\mathrm{Equation}36\right)\end{array}$

Navigation System

(83) An extended Kalman filter may be chosen to serve as the navigation filter implemented within the navigation system 108. The navigation system 108 of the present disclosure utilizes the VDM 212 as the main process model within a differential navigation filter. As illustrated in FIG. 5, the VDM 212 provides the navigation solution 110, which is updated by the navigation filter 502 based on available measurements. Hence, autonomous IMU sensor 102 data is treated as a measurement within the navigation filter 502, just as environment dependent GNSS receiver 104 observations are, when available. Autonomous IMU sensor 102 observations are related to system states via the VDM 212. Any other available sensory data, such as altimeter or magnetometer output, can also be integrated within the navigation filter 502 as an additional measurement. The VDM 212 is fed with the control input of the UAV, which is commanded by the autopilot 202. Another beneficial input is the wind velocity, which can be estimated by either airspeed sensors or within the navigation system 108 with no additional sensors needed. The latter approach is described herein through use of a wind model 504. The parameters of the VDM 212, which may be similar to the physical properties of the UAV, may be required within the navigation filter 502. These parameters may be pre-calibrated as fixed values. Alternatively, to increase the flexibility as well as the accuracy of the teachings of the present disclosure while minimizing the design effort, an online (i.e., during navigation) parameter estimation/refinement may be used. Calibration of such parameters may involve the parameters model 506. Autonomous IMU sensor 102 errors may also be modeled and estimated within the navigation system 108 as additional filter states.

State Space Augmentation

(84) Referring to FIG. 5, the navigation filter 502 (e.g., an EKF) may be used to estimate parameters or variables in addition to the navigation states. An augmented state vector may include navigation states (X.sub.n) (generated by the VDM 212), UAV dynamic model parameters (X.sub.p) (generated by a parameters model 506), IMU error terms (X.sub.e) generated by an IMU error model 508), and wind velocity components (X.sub.w) (generated by a wind model 504). For each of these augmentation parts, a process model may be needed. The process model may be dynamic ({dot over (X)}0) in general, but it can be static ({dot over (X)}=0) if a constant variable is being estimated within the navigation filter 502.

(85) The parameters used in the VDM 212 reflect the physical properties of the UAV. The parameters are typically constant, or change very slowly. However, the parameters' constant values may not be known initially (e.g., for a new UAV) with sufficient certainty. A static model may be used for the parameters within the navigation filter 502. The parameters vector related to the VDM contains most of the parameters used in the overall UAV dynamic modeling. The parameters used in the VDM, but not included in the augmented states, are related to the environment rather than the platform, which are known with a sufficient level of certainty. Gravity acceleration is an example of one of those parameters. Also, the mass and the moments of inertia of the UAV are not included in this vector, since they appear as scaling factors in the equations of motion and therefore they are completely correlated with the already included coefficients of aerodynamic forces and moments. The dynamic model parameters are ordered in a 261 vector and associated with a static process model as represented by Equations 37 and 38.

(86) $\left(\mathrm{Equation}37\right)$$X_p={\left[\begin{array}{cccccccccccc}& S,& \stackrel{\_}{c},& b,& D,& C_{F_{T}1},& C_{F_{T}2},& C_{F_{T}3},& C_{F_{z}1},& C_{F_z},& C_{F_{x}1},& \mathrm{.Math.}\\ \mathrm{.Math.}& C_{F_x},& C_{F_{x}2},& C_{F_{x}2},& C_{F_{y}1},& C_{M_{x}a},& C_{M_x},& C_{M_x{\stackrel{~}{}}_x},& C_{M_x{\stackrel{~}{}}_z},& C_{M_{y}1},& C_{M_{y}e},& \mathrm{.Math.}\\ \mathrm{.Math.}& C_{M_y},& C_{M_y{\stackrel{~}{}}_y},& C_{M_z{}_r},& C_{M_z},& C_{M_z{\stackrel{~}{}}_z}& {}_n& & & & & \end{array}\right]}^T$
{dot over (X)}.sub.p=[0].sub.261.(Equation 38)

(87) White noise is not considered for purposes of this disclosure, due to the fact that the parameters are expected to be constant over time, hence the static model holds with high certainty. However, since the initial values of the parameters generally have imperfections, an initial uncertainty is considered for them.

(88) The error in each accelerometer or gyroscope inside the autonomous IMU sensor 102 may be modeled as a random walk (b.sub.rw .sup.. . . ) process. Therefore, an IMU error terms vector may be represented by Equation 39, wherein the ai and gi superscripts denote the i-th accelerometer and gyroscope, respectively. The process model for X.sub.e is a linear dynamic model by definition and contains only white noise, as represented by Equation 40. G.sub.ee within Equation 40 may be represented by Equation 41. Moreover, W.sub.e within Equation 40 may be represented by Equation 42. This model is sufficient for a low-cost IMU, but can be extended as needed.
X.sub.e=[b.sub.rw.sup.a1,b.sub.rw.sup.a2,b.sub.rw.sup.a3,b.sub.rw.sup.g1,b.sub.rw.sup.g2,b.sub.rw.sup.g3].sup.T,(Equation 39)
{dot over (X)}.sub.e=G.sub.eeW.sub.e,(Equation 40)
G.sub.ee=[I].sub.66,(Equation 41)
W.sub.e=[.sub.b.sub.rw.sub.a1,.sub.b.sub.rw.sub.a2,.sub.b.sub.rw.sub.a3,.sub.b.sub.rw.sub.g1,.sub.b.sub.rw.sub.g2,.sub.b.sub.rw.sub.g3].sup.T.(Equation 42)

(89) The wind velocity, represented by Equation 43, may be stated as a vector in the local (navigation) frame consisting of the three components of wind velocity in the north, east, and down directions. While in usual cases wind velocity is expected to change slowly and smoothly from one moment to another, in general there is no information on such transition. Therefore, no deterministic part is considered in the process model of the present disclosure with respect to wind velocity, and only a white noise is assumed to rule the transition in time. Thus, the process model may be represented by Equation 44. G.sub.ww within Equation 44 may be represented by Equation 45. Moreover, W.sub.w within Equation 44 may be represented by Equation 46.

(90) $\begin{array}{cc}{X}_w={\left[w_x^l,w_y^l,w_z^l\right]}^T& \left(\mathrm{Equation}43\right)\\ X_w=G_{\mathrm{ww}}{W}_w,& \left(\mathrm{Equation}44\right)\\ G_{\mathrm{ww}}=\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right],& \left(\mathrm{Equation}45\right)\\ W_w={\left[{}_{w_x^l},{}_{w_y^l},{}_{w_y^l}\right]}^T.& \left(\mathrm{Equation}46\right)\end{array}$

Process Models Linearization

(91) Process models should be linearized to provide the dynamic matrices used in the EKF. In the most general case, the augmented state vector may be represented by Equation 47. Process models for all sub-state vectors of the augmented state vector have been presented herein above. The linearized model as represented by Equation 48 is equivalent to Equation 49.

(92) $\begin{array}{cc}X={\left[X_n^T,X_p^T,X_e^T,X_w^T\right]}^T,& \left(\mathrm{Equation}47\right)\\ X=\mathrm{FX}+\mathrm{GW},& \left(\mathrm{Equation}48\right)\\ \left[\begin{array}{c}{\stackrel{.}{X}}_n\\ {\stackrel{.}{X}}_p\\ {\stackrel{.}{X}}_e\\ {\stackrel{.}{X}}_w\end{array}\right]==\left[\begin{array}{cccc}{F}_{\mathrm{nn}}& F_{\mathrm{np}}& F_{\mathrm{ne}}& F_{\mathrm{nw}}\\ F_{\mathrm{pn}}& F_{\mathrm{pp}}& F_{\mathrm{pe}}& F_{\mathrm{pw}}\\ F_{\mathrm{en}}& F_{\mathrm{ep}}& F_{\mathrm{ee}}& F_{\mathrm{ew}}\\ F_{\mathrm{wn}}& F_{\mathrm{wp}}& F_{\mathrm{we}}& F_{\mathrm{ww}}\end{array}\right][\begin{array}{c}{X}_n\\ X_p\\ X_e\\ X_w\end{array}]+[\begin{array}{cccc}{G}_{\mathrm{nn}}& 0& 0& 0\\ 0& G_{\mathrm{pp}}& 0& 0\\ 0& 0& G_{\mathrm{ee}}& 0\\ 0& 0& 0& G_{\mathrm{ww}}\end{array}]\left[\begin{array}{c}{W}_n\\ W_p\\ W_e\\ W_w\end{array}\right].& \left(\mathrm{Equation}49\right)\end{array}$

(93) Linearizing each of the process models with respect to all state vectors, results in the determination of the corresponding row of F matrix being represented by Equation 50. Equation 50 may be evaluated either analytically or numerically. While analytical calculation is tedious (i.e., leads to very long equations for many elements), it may be desirable due to its accuracy, numerical stability, and computational efficiency.

(94) 0$\begin{array}{cc}{F}_{\mathrm{ij}}=\frac{\stackrel{.}{X_i}}{X_j},& \left(\mathrm{Equation}50\right)\end{array}$

(95) To determine the G matrix and W vector, proper uncertainty quantization should to be defined. This is beneficial, for example, when the uncertainty quantization it is not already included explicitly in the model. Assuming the sub-state vectors X.sub.n, X.sub.p, X.sub.e, and X.sub.w are uncorrelated, the off-diagonal blocks in the G matrix are all zero. Examining Equations 7-10 reveals that Equations 7 and 9 are purely kinematic equations (i.e., they are error-free). Therefore, the uncertainty terms apply only on linear and rotational accelerations of the UAV. The propeller speed is also added to the end of the X.sub.n vector, so the uncertainty on the dynamic model also appears in W.sub.n. Hence, W.sub.n may be represented by Equation 51. Additionally, G.sub.nn may thus be represented by Equation 52.

(96) $\begin{array}{cc}{W}_n={\left[{}_{{\stackrel{.}{v}}_N},{}_{{\stackrel{.}{v}}_E},{}_{{\stackrel{.}{v}}_D},{}_{{\stackrel{.}{}}_x},{}_{{\stackrel{.}{}}_y},{}_{{\stackrel{.}{}}_z},{}_{\stackrel{.}{n}}\right]}^T,& \left(\mathrm{Equation}51\right)\\ G_{\mathrm{nn}}=\left[\begin{array}{ccccccc}0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0\\ 1& 0& 0& 0& 0& 0& 0\\ 0& 1& 0& 0& 0& 0& 0\\ 0& 0& 1& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 0& 1\end{array}\right].& \left(\mathrm{Equation}52\right)\end{array}$

(97) By modeling the dynamic model parameters as constants in Equation 38, the G.sub.pp matrix and W.sub.pp vector are both empty (as represented by Equations 53 and 54).
G.sub.pp=[ ](Equation 53)
W.sub.p=[ ](Equation 54)

(98) As seen in Equation 40, the process model for IMU error terms is a linear model. Therefore, the G.sub.ee matrix and W.sub.e vector are the same ones as in Equations 41 and 42, respectively.

(99) The process model for wind velocity components is also a linear model, as represented in Equation 44. Therefore, the G.sub.ww matrix and W.sub.w vector are also the same ones seen in Equations 45 and 46, respectively.

Observation Models Linearization

(100) The linearized observation model may be derived from IMU 102 and GNSS 104 position measurements. When both IMU 102 and GNSS 104 measurements are available, the measurement vector Z may be represented by Equation 55. The linearized measurement model of vector Z may be represented by Equation 56, wherein r is the measurement noise vector with covariance matrix R. H within Equation 56 may be represented by the matrix of Equation 57.

(101) $\begin{array}{cc}Z={\left[Z_{\mathrm{IMU}}{Z}_{\mathrm{GNSS}}\right]}^{T}n& \left(\mathrm{Equation}55\right)\\ Z=\mathrm{HX}+r,& \left(\mathrm{Equation}56\right)\\ H=\frac{Z}{X}=\frac{{\left[Z_{\mathrm{IMU}}{Z}_{\mathrm{GNSS}}\right]}^T}{{\left[X_n{X}_p{X}_e{X}_w\right]}^T}=[\begin{array}{cccc}{H}_{\mathrm{IMU},n}& H_{\mathrm{IMU},p}& H_{\mathrm{IMU},e}& H_{\mathrm{IMU},w}\\ H_{\mathrm{GNSS},n}& H_{\mathrm{GNSS},p}& H_{\mathrm{GNSS},e}& H_{\mathrm{GNSS},w}\end{array}]& \left(\mathrm{Equation}57\right)\end{array}$

(102) The first row of the H matrix is derived using the IMU measurement equation (i.e., Equation 58). The IMU measures specific force components and rotation rates and includes IMU error terms.

(103) $\begin{array}{cc}{Z}_{\mathrm{IMU}}=\left[\begin{array}{c}{f}^b\\ {}_{\mathrm{ib}}^b\end{array}\right]+X_e.& \left(\mathrm{Equation}58\right)\end{array}$

(104) The second row of the H matrix is derived using the GNSS position observations in 3-D space, and is represented by Equation 59. If the GNSS position observations are completed with the velocity observations, the Z.sub.GNSS vector of Equation 59 may be extended to include the velocity observations.

(105) $\begin{array}{cc}{Z}_{\mathrm{GNSS}}=\left[\begin{array}{c}\\ \\ h\end{array}\right].& \left(\mathrm{Equation}59\right)\end{array}$

Hardware

(106) The UAV tested using the teachings herein was a custom made UAV equipped with the open-source autopilot PIXHAWK equipped with a time-stamping module. The maximum payload capacity was around 800 g, and the operational weight varied around 2500 g. Flying endurance was around 45 minutes. Thanks to its lightweight construction, the launching could be done from hand and it required only a small place for landing. The navigation sensors used included a micro-electro-mechanical system (MEMS) IMU called NavChip from INTERSENSE, and a standalone GNSS receiver inside the PIXHAWK autopilot. A geodetic grade GNSS receiver called Alpha from JAVAD was also on the UAV, working in real time kinematic (RTK) mode and providing the reference for evaluating the navigation solution 110. The sampling frequencies were 100 Hz for the IMU 102 and 1 Hz for the standalone GNSS receiver 104.

Sample Results

(107) The below description provides some simulation and experimental results as an example to show the performance of the navigation system of the present disclosure. The position results are presented in a Cartesian East, North, Up (ENU) frame and the attitude results are presented in Euler angles form (i.e., roll, pitch, and yaw).

(108) A Monte Carlo simulation was performed with 100 runs, using real 3D wind velocity data. While the trajectory and the wind was kept the same in each realization, the error in observations, initialization, and VDM parameters changed randomly for each individual run. FIG. 6 illustrates the reference trajectory, as well as the navigation solution from a sample run. The trajectory had a footprint bigger than 2 km2 km on the ground and a 1 km change in altitude.

(109) For purposes of the following description, sensor errors were modeled and added to nominal measurements. The exact values of model parameters used in simulating real sensors (i.e., the stochastic parameters) were unknown to the navigation system (i.e., the EKF), similar to real flight situations.

(110) IMU error for each accelerometer and gyroscope may be represented as the sum of a bias, a white noise, and a first order Gauss-Markov process. The numerical values of the model parameters corresponding to low-cost MEMS inertial sensors are provided in Table 1 below. The IMU sampling frequency was 100 Hz.

(111) TABLE-US-00001 TABLE 1 MEMS IMU error model parameters Error type Parameter Value Unit Accelerometers Bias 8 mg White Noise 50 g/{square root over (Hz)} 1.sup.st order Gauss-Markov 0.05 mg T 200 s Gyroscopes Bias 720 /hr White Noise 0.003 /s/{square root over (Hz)} 1.sup.st order Gauss-Markov 10 /hr T 200 s

(112) GNSS error is modeled as independent white noise signals for each channel (north, east, down). The standard deviation of the white noise is considered 1 m and the sampling frequency was 1 Hz.

(113) To start the filtering process, the initial values for all the states were required. These values normally contain some errors. The following paragraphs describe the initial errors and the initial uncertainties that were considered in simulations or within the designed EKF.

(114) For the navigation states (X.sub.n), the considered errors reflect the initialization error of the navigation system. An error of 1 m was considered for each component of the position, and the errors considered for the north, east, and down components of the velocity were 1 m/s, 0.2 m/s, and 0.5 m/s, respectively. The initial alignment error was considered to be 3 for roll and pitch and 5 for yaw, while an error of 1/s was assumed for roll, pitch, and yaw rates. An error of 15 rad/s was also considered for the propeller speed. The initial uncertainties in the state covariance matrix P were in the range from 1 to 2 times the initial errors described herein. Simulations revealed little sensitivity to the selection of these values.

(115) The errors considered for the VDM parameters (X.sub.p) were randomly distributed percentages of the real values with a standard deviation of 10%. The initial uncertainties were set in the range of the mentioned initial errors, and again, simulations showed little sensitivity to reasonable changes in these values.

(116) The initial values for the IMU error states (X.sub.e) were set to zero. The initial uncertainties were in general agreement with the IMU specifications used in the simulations, with a minor scaling with respect to the values used in simulating sensor errors, in order to preserve an approach close to reality, where even the exact error statistics of the sensors may not be known. Simulations verified little sensitivity to reasonable changes in these values.

(117) The initial values for wind velocity components (X.sub.w) were also set to zero, with an initial uncertainty of 1.5 m/s for each component. The airspeed sensor that may be on board of many (if not all) UAVs was not considered in the investigations described herein in order to further explore the capabilities of the navigation filter and avoid adding additional sensors to the minimal setup of an INS/GNSS system.

(118) Generally, a Kalman filter is provided with uncertainty values that reflect imperfections in modeling and observations. These uncertainties are discussed in the following paragraphs.

(119) The process model uncertainty is reflected in the GW term in Equations 48 and 49. Vector W contains white noise as the source of uncertainty in the process model, while matrix G is the noise shaping matrix. The four parts of W, namely and W.sub.n, W.sub.p, W.sub.e, and W.sub.w are represented in Equations 51, 54, 42, and 46, respectively. In W.sub.n, the value of 10.sup.4 was considered for the standard deviation of all terms. The small value of W.sub.n standard deviations does not necessarily mean a high confidence in the VDM, because the major portion of the uncertainty of the model is reflected in the initial uncertainty of the VDM parameters. However, when implementing the algorithm in a real scenario, a higher uncertainty may be needed. This is because there are always some unmodeled dynamics and disturbances that contribute to the deviation of model from reality, and these need to be accounted for in the error propagation. The process model for the VDM parameters as described herein is a static model with no uncertainty on the parameters being constant, so W.sub.p contains no value. The values of W.sub.e standard deviations were in general agreement with the IMU specifications used in the simulations. A standard deviation of 0.1 m/s.sup.2 was considered for each component of the wind velocity noise vector W.sub.w, to allow capturing moderate changes in wind velocity.

(120) The noise vector r in Equation 56 reflects the measurement uncertainty. The variance values for IMU measurements in the covariance matrix R were considered within the range of the IMU error statistics, neglecting the bias terms. This is because the estimated bias terms were directly subtracted from the IMU measurements before being fed to the EKF. For each channel of the GNSS measurements, the variance term in the covariance matrix was considered in the range of the error statistics of the device.

(121) FIG. 7 illustrates a comparison of the root mean square (RMS) of position and yaw errors for all the 100 Monte Carlo runs between the VDM/INS/GNSS approach described herein and a traditional INS/GNSS navigation algorithm over the whole interval. While the RMS of position error was 11.7 km for classical INS coasting after 5 minutes of autonomous navigation during GNSS outage, this was reduced to less than 110 m with the proposed VDM/INS/GNSS navigation platform under the exact same situations. That is, an improvement of two orders of magnitude in position accuracy was achieved through use of the herein disclosed navigation platform. The attitude determination also showed an improvement of 1 to 2 orders of magnitude. The improved performance of the navigation filter disclosed herein is also noticeable during the availability of GNSS in estimating yaw.

(122) The position error for all the 100 Monte Carlo runs is illustrated in FIG. 8. FIG. 8 illustrates how the error grew as time passed after GNSS outage started, and how the overall behavior was similar for individual runs. An empirical RMS was calculated from these individual errors and plotted against the predicted confidence level (1). Correctness of the navigation filter setup is demonstrated by how closely (and slightly conservatively) the error was predicted within the navigation filter.

(123) FIG. 9 illustrates an empirical RMS of attitude errors for all the 100 Monte Carlo runs, with associated predicted values of confidence (1). The results were satisfactory in terms of preserved navigation accuracy, with the RMS of error being only 0.0072 for roll, 0.020 for pitch, and 0.38 for yaw after 5 minutes of autonomous navigation during GNSS outage. In comparison, traditional INS coasting resulted in errors of 2.6 for roll, 1.5 for pitch, and 16.6 for yaw under the exact same conditions.

(124) FIGS. 10A and 10B illustrate the empirical RMS of errors in estimation of IMU error terms for all the 100 Monte Carlo runs. The predicted confidence values (1) were again very close to empirical RMS, and the accuracy of the estimation was satisfactory in comparison to the amplitude of simulated errors plotted on the same graph.

(125) The empirical RMS of mean error in estimation of VDM parameters for all the 100 Monte Carlo runs was plotted, as illustrated in FIG. 11. As FIG. 11 shows, the navigation filter was slightly optimistic in predicting the accuracy of the estimation, but the difference did not exceed 24%. The evolution of this error reveals that there was a sharp decrease in the mean error of parameter estimation from the initial 10% error to 6% during the first 40 seconds with GNSS available, which was followed by a slowly decreasing trend. The reason behind the second regime is the correlation between some parameters within the set. In such situations, groups of parameters are estimated better than individual parameters, and individual errors contribute to increasing the mean error for the whole set. The mean error was calculated in an RMS sense for all the parameters at any given time to be stricter and to cover both the mean value and the standard deviation.

(126) The RMS of error in estimation of wind speed for all the 100 Monte Carlo runs is illustrated in FIG. 12. Again, the navigation filter was able to provide a close and slightly conservative prediction on the estimation error, and the accuracy of the estimation is satisfactory for navigation purposes, especially considering that no airspeed sensor was considered in the navigation platform, so the estimation was provided simply as a courtesy of the VDM. There were no measurements on the air speed considered within the navigation filter, although they are available on most (if not all) UAVs, and will of course improve the estimation. As revealed by the slow growth of predicted error, the wind speed estimation accuracy was expected to decrease as time passed during GNSS outage. However, the rate of this growth was very well controlled, reaching only 0.7 m/s after 5 minutes of GNSS outage. All three components of the wind velocity were estimated within the navigation filter, which were converted to the wind speed when presenting results in favor of the plot being more comprehensive.

(127) Below, some aspects of the observability of navigation filter states are discussed. The discussion is based on an analysis of a covariance matrix throughout time. The close agreement between the empirical errors and the predicted confidence levels through the covariance matrix provides validity of stochastic parameters. Although the FIGS. discussed below are related to a single sample Monte Carlo run, the results were similar for all the 100 Monte Carlo runs.

(128) FIG. 13A illustrates ratios of uncertainties on states during GNSS presence periods and FIG. 13B illustrates ratios of uncertainties on states during GNSS outage periods. During the GNSS presence period (as illustrated in FIG. 13A), the uncertainties on all the states except a few VDM parameters were reduced significantly, indicating their observability, even during this relatively short period (i.e., 100 seconds) of stand-alone satellite positioning availability. The uncertainties on the VDM parameters were not increased either, which indicates no divergence is expected on their estimation. During GNSS outage, the uncertainty on position and velocity grew (as illustrated in FIG. 13B). The uncertainty on attitude was reduced by a factor of 0.3 for roll, and increased by factors of only 1.2 and 2.3 for pitch and yaw, respectively. These results are considered beneficial with respect to an attitude determination after 5 minutes of autonomous navigation. The uncertainty on rotation rates was reduced, thanks to the direct measurements from the IMU. In addition, the uncertainty on all VDM parameters was slightly reduced during GNSS outage, which means the filter was able to maintain the estimation on them even in the absence of GNSS measurements, and even slightly improve it. Conceptually, this was achieved by the redundancy of information on linear accelerations and rotation rates provided by the VDM and IMU. The uncertainty on IMU error terms was almost unchanged, which indicates the ability of the presently disclosed navigation filter to keep track of these errors even during GNSS data absence, which is not the case for inertial coasting. The uncertainty of wind velocity experienced an increase of almost 4 times, which is a useful result for navigation purposes. Further improvement can also be expected with the use of air speed sensor observations.

(129) FIG. 14 illustrates correlation coefficients among all the states within the herein disclosed navigation filter at the end of the navigation. As can be seen, the states were well decorrelated in general, which indicates good observability. The VDM parameters block showed the highest level of correlations, especially among different aerodynamic coefficients for roll, pitch, or yaw moments. This high level of correlation is distinguishable as almost all white small square blocks on the diagonal illustrated in FIG. 14. This means that rather than individual VDM parameters, groups of VDM parameters were getting estimated within the navigation filter. In this scenario, decorrelating individual VDM parameters demands longer time, more dynamic maneuvers, and possibly longer periods of GNSS availability. Considering these are practically constants for a particular navigation platform, their better estimation can be carried out only once during a calibration scenario. This results in lower uncertainties on their initial values for subsequent scenarios. For the purpose of navigation, and not necessarily parameter estimation, Monte Carlo simulations revealed that the VDM parameters were estimated adequately enough to enable the herein disclosed navigation platform to achieve a much better performance as compared to inertial coasting. Another point to notice is how well the VDM parameters decorrelated from other states of the navigation filter. There is also a considerable level of correlation between wind velocity components and position and velocity of the UAV. As depicted in FIG. 12, the uncertainty on wind velocity estimation grows steadily, yet slowly as time passes during GNSS outage. This is related to the level of process noise on the wind model. When gustier conditions are expected, higher values of process noise can be used to allow the capture of more sudden changes in wind velocity. This may come at the price of increasing the uncertainty on the navigation solution during GNSS outage. However, including air speed observations appears to be a good answer to such compromise.

Experimental Results

(130) The following discusses results of a test flight using a navigation system as described herein. The UAV used to conduct the test flight is illustrated in and described with reference to FIG. 4. The navigation sensors used to conduct the test flight are detailed herein above. FIGS. 15 and 16 illustrate the reference trajectory obtained from a RTK GNSS solution not used within the navigation system, as well as a navigation solution provided by the herein described navigation system. Available standalone GNSS measurements that were used during navigation are also depicted in FIGS. 15 and 16, as well as a last available measurement after which an artificial GNSS outage had begun. The duration of the GNSS outage was 3 minutes. As can be seen in FIGS. 15 and 16, the navigation solution according to the navigation platform of the present disclosure, after 3 minutes of autonomous navigation during GNSS outage, was still following the reference trajectory with enough accuracy for the purpose of safe guidance (e.g., returning home, landing, etc.).

(131) To get quantitative measurements on navigation quality, the 3D position error over the whole trajectory was mapped (as illustrated in FIG. 17) for the navigation system of the present disclosure. FIG. 18 illustrates quantitative measurements on navigation quality for a conventional INS/GNSS navigation system. As shown in FIG. 17, the maximum position error after 3 minutes of autonomous navigation was 69.8 m for the navigation system of the present disclosure. For the conventional INS/GNSS navigation system, this error was almost 270 times greater, exceeding 17 km. This demonstrates more than two orders of magnitude improvement in positioning accuracy. For both conventional and the herein described navigation systems, the predicted accuracy within the navigation filter is in agreement with actual error.

(132) Although illustrative embodiments of the present disclosure have been described herein with reference to the accompanying drawings, it is to be understood that the present disclosure is not limited to those precise embodiments, and that various other changes and modifications may be made by one skilled in the art without departing from the scope or spirit of the disclosure.