Hybrid Model-Based and Data-Driven Wind Velocity Estimator for an Autonomous Robotic Airship

TL;DR

This paper introduces a hybrid wind velocity estimator for autonomous airships, combining model-based and data-driven methods to improve estimation accuracy using limited sensors.

Contribution

It proposes a novel hybrid estimator that fuses Extended Kalman Filter and Neural Network approaches for wind velocity estimation in airships.

Findings

Hybrid estimator outperforms individual methods in simulations

Fusion improves robustness and accuracy of wind estimation

Approaches use minimal sensor data for effective estimation

Abstract

In the context of autonomous airships, several works in control and guidance use wind velocity to design a control law. However, in general, this information is not directly measured in robotic airships. This paper presents three alternative versions for estimation of wind velocity. Firstly, an Extended Kalman Filter is designed as a model-based approach. Then a Neural Network is designed as a data-driven approach. Finally, a hybrid estimator is proposed by performing a fusion between the previous designed estimators: model-based and data-driven. All approaches consider only Global Positioning System (GPS), Inertial Measurement Unit (IMU) and a one dimensional Pitot tube as available sensors. Simulations in a realistic nonlinear model of the airship suggest that the cooperation between these two techniques increases the estimation performance.

Tables

Table 1. Table 1: Correlation coefficient R-value and MSE.

Dataset	R-value	MSE
Training	0.99731	0.0206
Validation	0.99732	0.0205
Test	0.99720	0.0205
Total	0.99729	0.0205

Table 2. Table 2: Sensor noise standard deviation

Data	$\sigma \left(w_y\right)$
$\varphi ,\theta$ (rad)	5.2$\cdot {10}^{-3}$
$\psi$ (rad)	0.1
${\text{V}}_g$ (m/s)	0.4
$V_{pitot}$ (m/s)	6.04$\cdot {10}^{-4}$

Table 3. Table 3: Sample frequency specification

Sensor	Frequency	Sampled data
IMU	100 Hz	$𝚽$
GPS	4 Hz	${\text{V}}_g$
Pitot tube	18 Hz	$V_{pitot}$

Table 4. Table 4: RMS value of the estimation error.

	${\tilde{V}}_{N_w}$ (m/s)	${\tilde{V}}_{E_w}$ (m/s)
Cho2011	1.01	1.74
EKF	0.58 (-42.6%)	1.42 (-18.3%)
NN	1.19 (+17.8%)	1.25 (-27.0%)
Hybrid	0.74 (-26.7%)	0.71 (-59.2%)

Table 5. Table 5: RMS value of the estimation error in second simulation.

	${\tilde{V}}_{N_w}$ (m/s)	${\tilde{V}}_{E_w}$ (m/s)
Cho2011	0.71	0.52
EKF	0.52 (-26.5%)	0.38 (-26.9%)
NN	1.01 (+42.2%)	1.21 (+132.7%)
Hybrid	0.46 (-35.2%)	0.52 (0%)

Equations

$$\dot{\mathbf{\xi }}=g(\mathbf{\xi },\text{x}),$$

$$\dot{\text{x}}=f(\text{x},\text{u},\text{d}),$$

$$f(\textbf{x},\textbf{u},\textbf{d})=\textbf{M}^{-1}\Big{(}\textbf{F}_{k}+\textbf{F}_{a}+\textbf{F}_{p}+\textbf{F}_{g}+\textbf{F}_{w}\Big{)}\mathrm{,}$$

$$g(\mathbf{\xi },\text{x})=\left[\begin{array}{cc}{\text{S}}_{\Phi }^T& 0_{3\times 3}\\ 0_{3\times 3}& {\text{R}}_{\Phi }^T\end{array}\right]\left[\begin{array}{c}{\text{V}}_g\\ \mathbf{\Omega }\end{array}\right],$$

$${\text{S}}_{\Phi }=\left[\begin{array}{ccc}{\text{c}}_{\psi }{\text{c}}_{\theta }& {\text{s}}_{\psi }{\text{c}}_{\theta }& -s_{\theta }\\ {\text{c}}_{\psi }{\text{s}}_{\theta }{\text{s}}_{\varphi }-{\text{s}}_{\psi }{\text{c}}_{\varphi }& {\text{s}}_{\psi }{\text{s}}_{\theta }{\text{s}}_{\varphi }+{\text{c}}_{\psi }{\text{c}}_{\varphi }& {\text{c}}_{\theta }{\text{s}}_{\varphi }\\ {\text{c}}_{\psi }{\text{s}}_{\theta }{\text{c}}_{\varphi }+{\text{s}}_{\psi }{\text{s}}_{\varphi }& {\text{s}}_{\psi }{\text{s}}_{\theta }{\text{c}}_{\varphi }-{\text{c}}_{\psi }{\text{s}}_{\varphi }& {\text{c}}_{\theta }{\text{c}}_{\varphi }\end{array}\right],$$

$${\text{R}}_{\Phi }=\left[\begin{array}{ccc}1& {\text{s}}_{\varphi }{\text{t}}_{\theta }& {\text{c}}_{\varphi }{\text{t}}_{\theta }\\ 0& {\text{c}}_{\varphi }& -{\text{s}}_{\varphi }\\ 0& {\text{s}}_{\varphi }/{\text{c}}_{\theta }& {\text{c}}_{\varphi }/{\text{c}}_{\theta }\end{array}\right],$$

$${\text{V}}_a={\text{V}}_g-{\text{V}}_w,$$

$$V_t=||{\text{V}}_a|{|}_2=\sqrt{u_a^2+v_a^2+w_a^2}.$$

$$\beta ={\sin }^{-1}\frac{v_a}{V_t},$$

$$\alpha ={\tan }^{-1}\frac{w_a}{u_a}.$$

$$w_a=u_a\frac{\sin \alpha }{\cos \alpha }\mathrm{and}{v}_a=V_t\sin \beta .$$

$$V_t=\frac{u_a}{\cos \alpha \cos \beta }.$$

$$\Delta P=\eta\big{(}u_{a}\big{)}^{2}\mathrm{,}$$

$$V_{pitot}=\sqrt{\Delta P},$$

$$V_t=\frac{V_{pitot}}{\sqrt{\eta }\cos \alpha \cos \beta },$$

$$c_f=\sqrt{\eta }\cos \alpha \cos \beta ,$$

$$V_t=\frac{1}{c_f}{V}_{pitot}.$$

$${\text{V}}_{NEDa}={\text{V}}_{NED}-{\text{V}}_{NEDw}=S_{\Phi }^T{\text{V}}_a.$$

$$V_{pitot}^{2}=c_{f}^{2}\big{(}(V_{N}-V_{N_{w}})^{2}+(V_{E}-V_{E_{w}})^{2}+(V_{D})^{2}\big{)}\mathrm{.}$$

$$V_{N_a}=u_a\cos \psi \cos \theta ,$$

$$V_{E_a}=u_a\sin \psi \cos \theta ,$$

$$V_N=\frac{V_{pitot}}{c_f}\cos \psi \cos \theta +V_{N_w},$$

$$V_E=\frac{V_{pitot}}{c_f}\sin \psi \cos \theta +V_{E_w}.$$

$${\mathbf{\chi }}_{k+1}=\text{F}{\mathbf{\chi }}_k+{\mathbf{\nu }}_k,$$

$${\text{z}}_k=\text{h}({\mathbf{\chi }}_k)+{\mathbf{\upsilon }}_k$$

$$\text{F}=\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right]\text{; }$$

$${\mathbf{\chi }}_{k|k-1}=\text{F}{\mathbf{\chi }}_{k-1},$$

$${\text{P}}_{k|k-1}=\text{F}{\text{P}}_{k-1}{\text{F}}^T+\text{Q}.$$

$$\textbf{H}_{k}=\left[\frac{\partial\textbf{h}(\boldsymbol{\chi})}{\partial V_{N_{w}}},~{}\frac{\partial\textbf{h}(\boldsymbol{\chi})}{\partial V_{E_{w}}},~{}\frac{\partial\textbf{h}(\boldsymbol{\chi})}{\partial c_{f}}\right]\Bigg{|}_{\boldsymbol{\chi}_{k|k-1}},$$

$$\frac{\partial \text{h}(\mathbf{\chi })}{\partial V_{N_w}}=[-2{\hat{c}}_f^2{\hat{V}}_{N_w}(V_N-{\hat{V}}_{N_w})10{]}^T,$$

$$\frac{\partial \text{h}(\mathbf{\chi })}{\partial V_{E_w}}=[-2{\hat{c}}_f^2{\hat{V}}_{E_w}(V_E-{\hat{V}}_{E_w})01{]}^T\text{and}$$

Full text

††journal: Journal of the Brazilian Society of Mechanical Sciences and Engineering

1 Introduction

Recently, Unmanned Aerial Vehicles (UAVs) become useful in several applications due to their economic efficiency and mobility. For outdoor applications, air related information such as: angle of attack, sideslip and wind velocity, helps to improve the control performance.

Outdoor airships commonly have a guidance control to track a trajectory. The first idea is that the attitude reference shall be coincident with the reference trajectory attitude. However, there are two situations when this is not desirable: in the presence of wind disturbances (an almost certainty when flying outdoors) and if the objective is ground-hover (since the desired attitude is arbitrarily defined).

Airships, generally have tail surfaces to produce moments, however they do not have an actuator to directly produce a lateral force. Thus, when the airship is to keep a stationary position with respect to a ground target, it must align against the wind and use its forward propulsion to balance the aerodynamic drag. Therefore, by knowing the wind speed we may improve the control performance.For example, in zhou2019 a path following control approach is presented. The control law designed assume that attack and sideslip angles are measured. Similarly, in wangYueying2019 ; zheng2017 the same information is necessary for the path following approaches proposed.

This work is placed in the context of project DRONI marton2019 . The project aims to develop an autonomous airship for performing environmental monitoring tasks in remote amazon rain forest areas. Such tasks include aerodynamic flights and ground-hover (i. e. keep a stationary position with respect to a ground target). Therefore, whenever there is wind, the DRONI airship must try to align itself with the relative airspeed, thus reducing the sideslip angle. This implies that guidance control depends on information about wind velocity and attitude. However, measuring such elements is not a trivial task.

The most common solution is to estimate the wind velocity in order to extract the necessary information about the vehicle motion. The Model-based techniques are the most popular strategies. As an example, in perry2008 , it is proposed an approach for estimate the angle of attack and sideslip angle by the kinematic equations of motion of an aerobatic UAV. Meanwhile, with the same kinematic equations, in cho2011 an Extended Kalman Filter (EKF) is proposed for estimating the wind heading and velocity using an aircraft with a single GPS and Pitot tube. In johansen2015 a wind velocity observer also based in the kinematics is proposed for small UAVs with experimental results. Similarly, in shen2015 is also proposed an EKF for wind velocity estimation, however applied to a Stationary stratospheric airship in simulation environment. Then in rhudy2017 are presented four Model-based solutions considering an aircraft with four different possible configurations of sensors.

Recently, the Machine Learning approach has become popular in the field of robotics ribeiro2019 . The impressive growing of computational resources and increasing acquired data over the years have increased the potential of these Data-driven techniques. These strategies were already introduced in applications such as control of aircraft chaturvedi2002 and air data estimation for a Micro-UAV samy2010 . The wind estimation problem is addressed in allison2019 using a quadrotor and with a Machine-learning approach. However, a Data-driven online estimation of wind velocity for robotic airships is still a challenge.

This paper presents an alternative version of a Model-based wind velocity estimator using the EKF technique similar to the solution presented in cho2011 , however taking the DRONI airship as a case study. Then, a Data-driven approach of estimation using a Neural Network (NN) is proposed. Finally, a hybrid version that uses both Model-based and Data driven techniques is considered. The main tool to validate the proposed estimation approaches is a dynamical realistic non-linear model of an airship in Simulink/MATLAB. This tool is a result of the research group efforts since the project AURORA elfes2002 which was improved during the projects DIVA moutinho2016 and DRONI marton2019 .

This paper is organized as follows: the airship nonlinear modeling is summarized in Sect. 2; then the kinematic equations of motion are analyzed in Sect. 3; an EKF is designed for wind velocity estimation in Sect. 4; the NN approach of wind velocity estimation is presented in Sect. 5; the hybrid version of wind velocity estimation is presented in Sect. 6; validating simulations take place in Sect. 7 by establishing a comparison between the proposed approaches and the approach presented by cho2011 ; finally some conclusions are drawn in Sect. 8.

2 Airship modeling

The DRONI airship has a nonrigid body filled with helium gas, however without a ballonet. It is composed by a hull with 11m length and 2.48m diameter equipped with: 4 vectored propellers with independent thrusters (see Figure 1) and tail surfaces (rudder, elevator and aileron).

The DRONI airship is instrumented with a Xsens Mti-G 700 which contains an Inertial Measurement Unit (IMU) and a Global Positioning System (GPS). This set of sensors is located at the gondola which is 1.43 meters below the balloon Center of Volume (CV) as shown in Figure 2. Additionally, it is equipped with a one dimensional Pitot tube located at the airship nose. The airship has a total mass of 46.9 kg and buoyancy mass of -43.8 kg, resulting in an apparent mass of 3.1 kg.

This section presents a summary of the mathematical modeling of an airship. For a detailed description, please refer to moutinho2016 . The airship nonlinear model can be expressed as a state-space model given by:

[TABLE]

where:

• •

the kinematic states $\boldsymbol{\xi}=[\textbf{P}_{NED}^{T}\;\;\boldsymbol{\Phi}^{T}]^{T}$ include the cartesian positions $\textbf{P}_{NED}=[P_{N}\;\;P_{E}\;\;P_{D}]^{T}$ and angular position $\Phi=[\phi\;\;\theta\;\;\psi]^{T}$ in the North-East-Down frame;

• •

the dynamic states $\textbf{x}=[\textbf{V}_{g}^{T}\;\;\boldsymbol{\Omega}^{T}]^{T}$ include the linear speed $\textbf{V}_{g}=[u\;\;v\;\;w]^{T}$ and angular speed $\boldsymbol{\Omega}=[p\;\;q\;\;r]^{T}$ in the body frame;

• •

the input vector $\textbf{u}=[\delta_{e}~{}\delta_{a}~{}\delta_{r}~{}\delta_{0}~{}\mu_{0}]^{T}$ includes: $\delta_{e},\delta_{a}$ and $\delta_{r}$ which are elevator, aileron and rudder deflection (rad); $\delta_{0}$ as the normalized thrusters voltage (V/V); $\mu_{0}$ as the common vectoring angle of the thrusters (rad);

• •

and, finally, the disturbance vector d that includes wind velocities and gust parameters.

Another parameter of concern is known as added mass. As stated in carichner2013book , the added mass phenomenon is not significant enough to impact airplane performance. Thus, in general, it is not included in the modeling of fixed-wing aircrafts. However, the added mass is always present for any body which is moving through a fluid (see lamb1918 ; thomasson2000 ). Basically, a body in a fluid behaves as though it has more mass than it actually does. Such behavior varies with the motion nature (accelerating, decelerating or turning). Thus, a significant effect on the dynamics may be experienced by vehicles which have a density similar to the displaced external fluid (such as airships).

As stated by moutinho2016 ; fossen2011 , when the displaced fluid mass is not negligible, as is the case for airships, balloons and submarines, the equations of motion are usually derived from the Lagrangian approach. However, the resultant dynamics can be represented by the Newton-Euler equations including five components of forces and moments, namely: $\textbf{F}_{k}$ containing the kinetic forces and moments, including Coriolis and centrifugal force terms; $\textbf{F}_{a}$ given by aerodynamic forces and moments; $\textbf{F}_{p}$ given by propulsion forces and moments; $\textbf{F}_{g}$ given by gravity forces and moments, which are function of the difference between the weight and buoyancy forces; and $\textbf{F}_{w}$ given by the wind induced forces and moments. Therefore, the linear and angular accelerations are given by:

[TABLE]

where M includes apparent mass and inertial coefficients of the airship. These equations are referenced in the body frame centered in the Center of Buoyancy (CB) that is approximately equivalent to the Center of Volume (CV) as shown in Figure 2.

Finally, the linear and angular positions are updated through kinematic equations (3).

[TABLE]

where $0_{3\times 3}$ is a $3\times 3$ null matrix, and $\textbf{S}_{\Phi}\in\mathbb{R}^{3\times 3}$ and $\textbf{R}_{\Phi}\in\mathbb{R}^{3\times 3}$ are rotational matrices from airship body to NED frame, given by:

[TABLE]

[TABLE]

where $\text{c}_{x}=\cos(x)$, $\text{s}_{x}=\sin(x)$ and $\text{t}_{x}=\tan(x)$.

3 Kinematic Equations of motion

In this section we find three equations which are correlated with the Pitot probe measurement, wind speed ($\textbf{V}_{w}$), airship groundspeed ($\textbf{V}_{g}$) and orientation ($\boldsymbol{\Phi}$).

Let the airship motion be represented by its inertial velocity $\textbf{V}_{g}$. Similarly, the wind is described by an inertial velocity $\textbf{V}_{w}$. The airship relative air velocity is called airspeed ($\textbf{V}_{a}$) and it is given by:

[TABLE]

where $\textbf{V}_{w}=[u_{w}\,v_{w}\,w_{w}]^{T}$ and $\textbf{V}_{a}=[u_{a}\,v_{a}\,w_{a}]^{T}$.

The Euclidean norm of the airspeed is called true airspeed ($V_{t}$) and it is given by:

[TABLE]

Other important definitions are the sideslip angle $\beta$ and angle of attack $\alpha$. The sideslip angle is a relative orientation between the vertical plane of the vehicle and the vector $\textbf{V}_{a}$. Moreover, the angle of attack $\alpha$ is given by the angle between the vector $\textbf{V}_{a}$ and the horizontal plane of the vehicle, as shown in Figure 3.

Threfore, we can define $\beta$ and $\alpha$ by the following statements:

[TABLE]

[TABLE]

Thus, an equivalent formula is given by:

[TABLE]

Finally, we obtain:

[TABLE]

The Pitot tube is located in the airship nose, thus it will be measuring the longitudinal dynamic pressure $\Delta P$ in the body frame. Also, it has a correlation with the airspeed as shown below:

[TABLE]

where $\eta$ is the calibrating factor that is correlated with the air density and pitot efficiency. Now consider the following variable transformation:

[TABLE]

Thus, $V_{pitot}$ is correlated with the true airspeed through the following statement:

[TABLE]

Because there are uncertainties in $\eta$ and the angles $\alpha$ and $\beta$ are unmeasurable with the actual available sensors, those values will be estimated together as a scale factor $c_{f}$ given by:

[TABLE]

therefore (12) becomes:

[TABLE]

Now, consider that the rotation of vector $\textbf{V}_{g}$ from the body frame to the NED frame is given by $\textbf{V}_{NED}=[V_{N}~{}V_{E}~{}V_{D}]^{T}$. Also, consider that such rotation applied to vector $\textbf{V}_{w}$ is given by $\textbf{V}_{NED{w}}$. Thus, the airspeed in NED frame is given by:

[TABLE]

It is known that a rotational operation does not change the vector module, thus the following statement is valid, considering the wind strictly horizontal:

[TABLE]

Supposing that, the airship starts from an initial condition where $\alpha$ and $\beta$ are negligible ($\alpha\approx\beta\approx 0$) we have $u_{a}=V_{t}$ and $v_{a}=w_{a}=0$. Therefore in the global frame we have:

[TABLE]

[TABLE]

where $\psi$ and $\theta$ are the yaw and pitch angles, respectively. Since we have (4), then:

[TABLE]

[TABLE]

The values of $V_{pitot}$, $V_{N}$ and $V_{E}$ are measurable by the Pitot tube and GPS, therefore (16), (19) and (20) can be used as observation equations, while $c_{f}$, $V_{N_{w}}$ and $V_{E_{w}}$ are estimated states in the EKF. Note that, the Euler angles ($\phi$, $\theta$ and $\psi$) can be measured by the IMU.

4 Extended Kalman Filter

In this section we propose an EKF in order to estimate the incident wind in the airship body. Assuming that the wind is strictly horizontal the main goal is to obtain an estimation of the wind velocity in the horizontal plane (North-East) and the Pitot probe scale factor. One main benefit is expected from this method. The estimator presented in cho2011 , uses only one measurement update equation which is similar to (16). By introducing (19) and (20) in the measurement update stage we expect to have a better performance than using only (16).

There is no given model to determine the wind behavior. Thus, we assume that the wind is constant with a Gaussian input with significant covariance. In addition, as long as we do not have a sideslip sensor, thus the uncertain and time varying factor $c_{f}$ will be estimated, which leads us to the following reduced system model:

[TABLE]

where: $\boldsymbol{\chi}_{k}=[V_{N_{wk}}~{}V_{E_{wk}}~{}c_{fk}]^{T}$ is the state vector in the instant $t=kt_{s}$; $\textbf{z}_{k}=[V_{pitot_{k}}^{2}~{}V_{N_{k}}~{}V_{E_{k}}]^{T}$ is the system output in the instant $t=kt_{s}$; $t_{s}$ is the sample time in seconds; $\textbf{h}(\boldsymbol{\chi}_{k})$ is the output function, which can be computed through $\eqref{eq:z1}$, $\eqref{eq:z2}$ and $\eqref{eq:z3}$;

[TABLE]

$\mathfrak{\boldsymbol{\nu}}_{k}\sim\textbf{N}(0,\textbf{Q})$ is the process noise with Gaussian distribution and covariance Q; and, finally, $\mathfrak{\boldsymbol{\upsilon}}_{k}\sim\textbf{N}(0,\textbf{R})$ is the measurement noise also with Gaussian distribution and covariance R.

Given the model described in (21a), we can update the state and covariance matrix (P) as follows:

[TABLE]

[TABLE]

For the accomplishment of the EKF final stage we obtain the Jacobian matrix of $\textbf{h}(\boldsymbol{\chi}_{k})$ evaluated in the measured values of $V_{pitot}$, $\textbf{V}_{NED}$, $\boldsymbol{\Phi}$ and $\boldsymbol{\chi}_{k|k-1}$ which is given by:

[TABLE]

where:

[TABLE]

[TABLE]

[TABLE]

Finally, the standard algorithm of EKF can be applied as follows:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where P is the covariance matrix, C is the covariance error, $\tilde{\textbf{y}}$ is the measurement error, K is the Kalman gain and I is the identity matrix with appropriate dimensions.

5 Neural Network

This section proposes a Neural Network (NN) to estimate the wind speed in NED frame and the Pitot probe scale factor $c_{f}$. The main benefit expected from this method is the NN intrinsic property of create a mapping from input parameters to the output. Because the NN is trained by measured data, it is able to detect abrupt variation in the wind velocity. However, it is also expected to be very sensitive to measurement errors.

As the analytical model (21a)–(21b) shows, the model output is composed by (16), (19) and (20), which are nonlinear equations in the model states, vehicle velocity, orientation and Pitot pressure. In order to avoid complex nonlinearities, the measured data is remapped into 8 inputs given by:

[TABLE]

Meanwhile the output vector $\boldsymbol{\chi}_{nn}$ of the NN are the estimated wind velocities in the horizontal plane and the scale factor $c_{f}$ as shown below:

[TABLE]

The NN was designed in the MATLAB Neural Network ToolboxTM. It is a three-layer fitting NN, which has three nonlinear hidden layers containing 24 neurons each and three linear outputs. The activation function of the nonlinear neurons are sigmoidal. The resulting flow chart is shown in Figure 4.

The training dataset is composed by simulations in variations of the two trajectories shown in Figure 5. The trajectories (a) and (b) are rotated of {0, 45, 90, 135, 180, 225, 270, 315} degrees around origin generating 16 scenarios (see Figure 15) in which the airship performs curves and straight lines in different directions. For each scenario were performed simulations with wind speed at $|\vec{V}_{w}|=$ {0, 1, 2, 3, 4, 5}m/s and heading $\phi_{w}=${0, 22.5, 45, 67.5, 90, 112.5, 135, 157.5, 180, 202.5, 225, 247.5, 270, 292.5, 315, 337.5} degrees, where \phi_{w}=\tan^{-1}\big{(}\frac{V_{E_{w}}}{V_{N_{w}}}\big{)}. Hence, a total of 1281 simulations were performed. In all simulations the airship performs a typical cruise flight at 7m/s airspeed and constant altitude of 50 meters.

The dataset stored has approximately 375 MB. Training phase and further evaluations were carried out on a desktop architecture which features a four-core 4.00 GHz Intel Core i7-6700K Processor, NVIDIA GeForce GTX 950, 32GB of RAM and Ubuntu 16 LTSOS. The system took 22 minutes and 17 seconds to train the afore-mentioned dataset.

Training was achieved using Matlab Neural Network Toolbox. The Scaled Conjugate Gradient (SCG) algorithm is used with 5000 epochs of training iterations using 70% of the collected data randomly taken as the training set, 15% used for validation set and 15% as the test set. The performance evaluation is made by Mean Squared Error (MSE).

The residual error of the trained NN is shown in Figure 6. Note that, about 80% of the total data are distributed around zero, in which 60% were taken as training data, 10% as validation data and about 10% as test data.

The correlation coefficient (R-value), which is a linear regression between the NN predicted values and the targets, and best MSE are shown in Table 1. Values of R closer to 1 indicate better agreement between targets and predicted values. Note that both R-value and MSE are close for Training, Validation and Test dataset, which indicates that the NN is not overfitting the data.

In order to avoid high frequency oscillations in the estimation, a low-pass filter was introduced for filtering the input values $\textbf{V}_{g}$, $\boldsymbol{\Phi}$ and $V_{pitot}$. The low-pass filter can be expressed by the following Laplace representation:

[TABLE]

where $\tau$ is the time constant (here we use $\tau=1.5$ seconds), $Z(s)$ is the Laplace transform of the input and $Y(s)$ is the Laplace transform of the filter output. After filtering the measured sensor data, the input vector $\textbf{z}_{nn}$ is computed and passed by the trained NN as shown in the block diagram from Figure 7. Finally, the estimated wind velocity $\boldsymbol{\chi}_{nn}$ is given to the airship navigation system.

6 Hybrid estimator

Here we propose a hybrid estimator that performs a fusion between both estimators, namely: from the EKF designed in Sect. 4 and NN designed in Sect. 5. This fusion is performed by changing the measure update stage of the EKF approach. The NN output $\boldsymbol{\chi}_{nn}$ is added to the measurement vector of the EKF as a redundant measure. Thus, resulting in the new measurement vector $\textbf{z}_{h_{k}}$, updating function $\textbf{h}_{h_{k}}(\boldsymbol{\chi}_{k})$ and its respective Jacobian $\textbf{H}_{h_{k}}$ shown below:

[TABLE]

where $\textbf{I}_{3}$ is the identity matrix of third order. Then the EKF standard algorithm is used by updating the dimensions of the matrices $\textbf{C}_{k}$, $\textbf{K}_{k}$ and R. The resulting estimator has a cascaded form as illustrated in Figure 8.

It is important to highlight that the NN used here is the same NN previously designed and trained in Sect. 5.

7 Simulation results

In this section two simulations are performed in order to evaluate all the three approaches presented before and establish a comparison with the traditional model-based approach proposed by Cho et. al. cho2011 .

During simulation the airship is well controlled with ideal feedback, meanwhile the estimators evaluated are receiving noisy data from the modeled sensors as shown in the block diagram from Figure 9. The sensors are modeled in the simulation environment with sample frequency as specified in Table 3. Also, each sensor has a generic modeling including a Gaussian noise as specified by the manufacturer in Table 2. The generic sensor modeling is shown in Figure 10, where $w_{z}$ is a Gaussian noise, $\hat{z}$ is the true simulated value and $z$ is the sensor output. All estimators use the same sample frequency of 16Hz ($t_{s}=0.0625$ seconds).

An online repository111https://github.com/leve-fem/airship_estimator is available containing all approaches presented here. The algorithms are implemented in C/C++ and Python inside the Robot Operating System (ROS) ros . Also, in this same repository a link to the dataset used for the NN training task and the simulations performed are available for future researches.

7.1 First scenario

The first scenario considers wind with absolute value $|\vec{V}_{w}|=2m/s$ and heading $\psi_{w}=\frac{\pi}{2}$ rad (blowing from East to West). Then in the instant $t=160s$ the wind is intensified to $|\vec{V}_{w}|=3m/s$ and its heading is changed to $\psi_{w}=\pi$ (blowing from South to North).

The airship is controlled to follow the trajectory shown in Figure 11. Five instants are highlighted with gray background in order to establish further comparisons with the results in Figure 12. Moreover, results using the estimator proposed by cho2011 were introduced as “Cho2011” in order to establish a comparison. The covariance matrices used in the Model-based approaches can be found in Appendix.

In Figure 12 is possible to note that, the “Cho2011” estimator has to acquire information about the motion in all directions before it converges to the correct values. After the instant (I), all estimators converges to values within of an acceptable error. It is important to note that, the NN has some high frequency oscillations at the trajectory curves, which deteriorate the performance. When the airship is following a straight line the designed NN has a good estimation, however sometimes with a bias from the real value.

Also in Figure 12, we can note that when the wind velocity has a significant variation in the instant (III), the two Model-based approaches (“EKF” and “Cho2011”) do not converge immediately because both depend on information (given by the Pitot tube) about the other directions to converge to the correct wind velocity. Meanwhile the NN clearly has an instantly reaction to these variations. Even though the NN converges for a biased value, such information was sufficient to correct the estimation of the Hybrid approach before the instants (IV) and (V). In these final instants, the Model-based estimators finally converges for values within a range of acceptable error.

In Table 4 are shown the RMS values of the estimation errors $\tilde{V}_{N_{w}}=\hat{V}_{N_{w}}-V_{N_{w}}$ and $\tilde{V}_{E_{w}}=\hat{V}_{E_{w}}-V_{E_{w}}$ and the percentage in relation to the approach presented by Cho et. al.cho2011 . By the RMS values we can observe that the NN had a better estimation of $V_{E_{w}}$ in comparison to the Model-based approaches. However, for the component $V_{N_{w}}$, the Model-based approaches presented a better performance. Meanwhile the “Hybrid” which has the information of both approaches had the best performance in the estimation of $V_{E_{w}}$ and acceptable estimation in the $V_{N_{w}}$ component.

Figure 13 shows the histograms of computational time necessary for each method. Although NN and Hybrid approaches have higher computational time (100 $\sim$ 160 microseconds), they are able of running at almost 1000Hz. Meanwhile, both Model-based approaches (EKF and Cho2011) are able of running at more than 1000Hz with a computational time between 25 and 35 microseconds. Note that the proposed approach by Cho et. al. cho2011 has lowest computational cost, since the EKF proposed here has two additional measurement update equations.

7.2 Second scenario

In a second simulation in the same trajectory of the previous simulation, it is considered a different variation in the wind velocity. Initially the wind has absolute value $|\vec{V}_{w}|=2m/s$ and heading $\psi_{w}=0$ (blowing from North to South), then in the instant $t=160s$ the wind is intensified to $|\vec{V}_{w}|=3m/s$ and its heading is changed to $\psi_{w}=\frac{\pi}{2}$rad (blowing from East to West). The results are shown in Figure 14 and resultant RMS error is shown in Tab. 5.

In this second simulation, the airship starts with the nose against the wind, thus a few seconds of simulation are sufficient to the Model-based approaches estimate correctly the wind velocity components. Note that as the airship do curves in the scenario in instants (I) and (II), NN and Cho2011 approaches presents some estimation errors, while the EKF and Hybrid maintain the estimation near the true value. The improvement made in the model-based approaches by including two additional measurement update equations avoids high oscillations.

In the overall simulation, Model-based approaches had better performance. It occurs because the observability is not degraded during the wind transition at instant $t=160$ seconds. After the instant (III), the airship is traveling to the East, thus the wind is against its nose. Therefore, the Pitot tube is able to capture the relative speed perfectly. By table 5 we can also note that the EKF presented similar estimation performance for both wind components ($V_{N_{w}}$ and $V_{E_{w}}$), differently from the previous simulation. Clearly, in this scenario, the NN is not too advantageous. It presented the higher RMS errors.

8 Conclusions

In this paper we presented Model-based and Data-driven approaches for estimation of wind velocity for a robotic airship. The Model-based approach uses only kinematic equations of motion of the airship for the design of an EKF. The Data-driven proposed approach is composed by a NN trained with a dataset containing 1281 simulations in different conditions. Also, a novel Hybrid approach is proposed, by performing a fusion between the designed Model-based and Data-driven approaches with a cascaded structure.

A high-fidelity and realistic airship dynamic model for Matlab/Simulink platform was used in the simulations. Two simulation scenarios were presented. In the first, the wind variation is not detected properly by the model-based approaches, thus the Data-driven is advantageous in this scenario. In the second, the Model-based are more precise in the estimation, and the Data-driven performance was the same or worse.

The performance results obtained showed that the proposed EKF has a slightly better performance in comparison to the strategy found in the literature. It occurs because of the two additional measurement update equations that we introduced. Although, the performance was improved, EKF still have the same drawbacks when in situations of low observability. On the other hand, the NN presents higher sensibility to wind variations in such situations, however with biased estimations and high frequency oscillations due to sensor noise and non-trained situations. The NN approach is independent of any mathematical model, however, requires a sufficient and representative database.

The Hybrid approach maintains satisfactory performance in the two simulation scenarios presented. This approach combines the advantages of both estimators (model-based and data-driven) and mitigates the drawbacks, once it has the information about the model and about previous missions (training dataset). However, it requires a representative dataset as well as a mathematical modeling.

As a general result, we can conclude that the cooperation between both approaches (Model-based and Data-driven) can be highly effective for solving estimation problems with observation deficiency. For the specific problem of wind estimation we obtained satisfactory results with the Hybrid approach. However, since the proposed methods are only analyzed theoretically and validated via simulation, an actual benchmark or field test is needed in the subsequent work to verify the proposed approaches. Future efforts will be made to validate these results outside of a simulation environment.

Appendix

The covariance matrices used for each Model-based approach are shown below:

[TABLE]

Acknowledgements.

The present work was sponsored by the Brazilian agencies CNPq and FAPESP, through Projects DRONI (CNPq 402112/2013-0), INCT-SAC (CNPq 465755 /2014-3; FAPESP 2014/50851-0) and scholarships (CNPq 305 600/2017-6; FAPESP BEP 2017/11423-0). Also, this work was supported by Fundação para a Ciência e a Tecnologia (FCT), through IDMEC, under LAETA, project UID/EMS/ 50022/2013 (Portugal). Moreover, the authors are grateful to Erasmus Mundus SMART2 for the financial support through project reference 552042-EM-1-2014-1-FR-ERA MUNDUS-EMA2. Furthermore, the authors appreciate all the support given by members of the Advanced Computing, Control & Embedded Systems Laboratory (ACCES-Lab) and Laboratory of Study in Exterior Vehicles (LEVE) from University of Campinas.