# A Novel Estimation Method of Wind Speed and Direction Based on GNSS/INS/ADS Integrated Navigation

Authors: Jingjuan Zhang; Hao Cong; Xueyun Wang
Journal: Journal of Sensors
 (2021)

---

## Abstract

Arriving on time is of great importance for flight management and passenger experience. One of the essential factors that impacts on-time arrival is the wind condition. Accurate information of wind speed and direction around the fuselage helps to improve the performance of on-time arrival and four-dimensional trajectory (4DT) planning. To determine accurate wind information in real-time, a novel airborne estimation method of wind speed and direction is proposed in this paper. Inertial Navigation System (INS), Global Satellite Navigation System (GNSS), and Air Data System (ADS) are fused in an Unscented Kalman Filter (UKF), which provides great accuracy and robustness in nonlinearity conditions. The dynamic models of wind are established, and implementations of the UKF are detailed. Finally, simulations are designed and the effectiveness of the proposed method is verified through the comparison with the traditional direct measurement method. Results demonstrate that the accuracy of wind speed and direction obtained by our method is nearly two times higher than the traditional direct measurement method.


---

## Body

## 1. Introduction

The estimation of instantaneous wind speed and direction at the location of the aircraft can be a significant problem due to the inevitable presence of wind in the flying. This problem is particularly important in four-dimensional trajectory (4DT) planning [1, 2]. Four-dimensional trajectory planning is aimed at conducting trajectory planning in four-dimensional space [3, 4]. The resulting flight trajectory includes all the three-dimensional positions where the aircraft is flying and the flight time when it is in that position, which means both the three-dimensional coordinates of the trajectory and the time to reach the three-dimensional coordinate point are planned. Since wind speed and direction can have a great impact on the time dimension of four-dimensional trajectory planning, accurate real-time wind information is needed in 4DT. Besides, the ambient wind speed and direction are also very useful information for various objectives such as flight control [5], micro aerial vehicles [6], atmospheric energy harvesting [7–9], air traffic control [10], aerial refueling [11], and collision avoidance. A popular approach to the wind estimation problem involves applying the aircraft dynamic model to predict the information of wind through Kalman filter [12, 13].A problem with the existing work using Extended Kalman Filter (EKF) [14–17] for wind estimation [18] is that the wind dynamic model is nonlinear so the linearization accuracy of EKF may not be enough or even diffuse. Lefas developed a simple filter for wind estimation using magnetic heading, true air speed, and radar measurements; however, the data transmitted through the air-ground data link of the radar system can result in wildly erroneous measurements [10]. Langelaan et al. conducted a thorough simulation research on a method of directly estimating wind field; however, the measurement error of wind speed increases with the increase of airspeed [1]. Cho et al. proposed a method of horizontal wind estimation based on the EKF, which combines simulation and experimental data; however, it is usually not feasible in practical application that the aircraft needs to make banking turns or circle maneuvers, especially for civilian flights [19]. The horizontal wind speed and direction are predicted by the assumption of random moving noise, and then, these states are adjusted by the wind triangle comparison of air speed and ground speed provided by GPS [20, 21]. There is also some work in wind estimation using the Unscented Kalman Filter (UKF), though this work is limited in assuming a known a priori wind direction [22]. Taylor presented a simulation study of parameter estimation techniques for air data calibration and wind estimation [23]. This method assumes that the inertial velocities and Euler angles are accurate, that the calibration models are correct, and that the steady-state component of wind is constant throughout the maneuver. Liu et al. [24] comprehensively evaluated the sideslip angle estimation methods of vehicles in a comparative study. The results show that under the dynamic model, the UKF has higher performance than EKF to a certain extent. To overcome the time-varying and complexity of environmental noise information in practical application, Wang et al. [25] derived an AEKF method from a maximum-likelihood estimation for rolling state estimation. The results show that AEKF has better performance and feasibility than EKF. Under the condition of multisensor measurement, Ding et al. [26] proposed a weighted measurement fusion Kalman estimator based on centralized fusion to estimate the longitudinal vehicle speed. Through the combination of multiple sensor signals, three virtual sensors are formed, each virtual sensor generates a corresponding speed information, and the speed information in a long time is obtained by weighting coefficient estimation. It has high precision and reliability. This shows that the combination of data fusion and Kalman filter is of great significance.The work introduced herein is to deduce the estimation formulas of three-dimensional wind speed and direction based on the measurements from airborne Global Satellite Navigation System (GNSS), the output results of Inertial Navigation System (INS), and the information of airspeed, angle of attack, and sideslip angle provided by an Air Data System (ADS). The work also extends Cho’s work including vertical wind and wind direction, as well as angle of attack and sideslip information. In addition, the three-dimensional position is also provided as measurements, which is of great significance to improve the accuracy of wind estimation. The contribution of this work is that it provides an innovative method based on UKF for real-time and accurate estimation of wind field information around aircraft, as well as navigation, position, and flow angle information.The organization of this paper proceeds as follows. Section2 describes the aircraft dynamic equations. Section 3 introduces the Unscented Kalman Filter used in wind speed and direction estimation in this paper. Section 4 presents various simulations and their results. Finally, Section 5 presents the concluding remarks.


## 2. Problem Formulations

The estimation algorithm uses information from an Air Data System (ADS), Global Satellite Navigation System (GNSS), and Inertial Navigation System (INS) to estimate the three-dimensional components of wind speed and direction.


### 2.1. Introduction of Coordinates

First, the coordinates that are used in the paper are defined.On−XnYnZn denotes the navigation frame which coincides with the geographic frame (Og−XgYgZg), north east and down (NED).Ob−XbYbZb denotes the aircraft body coordinate frame. Ob is at the center of gravity, Xb is in the plane of symmetry of the aircraft and parallel to the design axis of the aircraft and points to the nose, Yb is perpendicular to the plane of symmetry of the aircraft and points to the right of the fuselage, and Zb is in the plane of symmetry of the aircraft, perpendicular to the Xb and pointing below the fuselage. The rotation from NED to the aircraft body coordinate system is defined by Euler angles, as shown in Figure 1.Figure 1
Relationship between the aircraft body coordinate frame and the geographic system.Ow−XwYwZw denotes the wind coordinate frame. Ow is at the center of gravity and fixed connection with aircraft body; Xw is pointing at the direction of airspeed Va; Zw is in the plane of symmetry of the aircraft, perpendicular to the Xw and pointing below the fuselage; and Yb is defined by the right-hand rule. Rotation from the aircraft body coordinate frame to the wind coordinate frame is defined by the angle of attack α and sideslip β, as shown in Figure 2.Figure 2
Relationship between the wind coordinate frame and the aircraft body coordinate frame.


### 2.2. Dynamic Equations

The following quantities are studied in this paper:(1)
The forward, rightward, and downward airspeed componentsu, v, and w, respectively(2)
The rollϕ, pitch θ, and yaw ψ of the Euler attitude angles(3)
The northward, eastward, and downward wind velocity componentswx, wy, and wz(4)
The latitudeL, the longitude λ, and the altitude H of the aircraft(5)
The pitot tube velocityVpitot, angle of attack α, and sideslip angle β derived from the Air Data System(6)
The northward, eastward, and downward velocity (ground speed)VxmVymVzm and positions LmλmHm provided by the GNSS(7)
The accelerationsaxayaz and angular velocity pqr derived from the Inertial Measurement Unit (IMU).Firstly, based on the flight dynamics and the assumption that the plane is a rigid body, the motion equations of the aircraft are established as follows.The conversion from the aircraft body coordinate frame to the geographic coordinate system is defined by the Direction Cosine Matrix (DCM), which is derived from the continuous rotation of rollϕ, pitch θ, and yaw ψ of the aircraft.
(1)Cbn=DCMϕ,θ,ψ=cosψ−sinψ0sinψcosψ0001cosθ0sinθ010−sinθ0cosθ1000cosϕ−sinϕ0sinϕcosϕ.The following force equations can be obtained by decomposing the total aerodynamic forces and the thrustsFx,Fy,Fz into the aircraft body coordinate system and then decomposing the gravity into the aircraft body coordinate system, too:
(2)u̇=vr−wq−gsinθ+Fxm,v̇=−ur+wp+gcosθsinϕ+Fym,ẇ=uq−vp+gcosθcosϕ+Fzm,where m is the total mass of the aircraft which is assumed to be a constant. The relationship between the attitude rate ϕ̇θ̇ψ̇ and the three angular velocity components pqr in the body coordinate system can be described by the dynamic equations as follows [27]:
(3)ϕ̇=p+rcosϕ+qsinϕtanθ,θ̇=qcosϕ−rsinϕ,ψ̇=rcosϕ+qsinϕsecθ.According to the characteristic of randomness, the state model of wind speed is set as a random walk process [19].The dynamics for the position states are then given by groundspeeds of the aircraft:(4)L̇=VxRm+H,λ̇=VyRn+HsecL,Ḣ=Vz,where Rm=Re1−2e+3esin2L, Rn=Re1+esin2L, and Re=6378393 m, and e=1/297.3 is the oblation of rotating ellipsoid. VxVyVz are groundspeeds of the aircraft in navigation coordinates. They can be correlated with the wind speed in geographic coordinates and the airspeed in the aircraft body coordinate through the Direction Cosine Matrix as follows:
(5)VxVyVz=Cbnuvw+wxwywz.Secondly, the mathematical models related with pitot speed, angle of attack, and sideslip angle are presented. The pitot speedVpitot is calculated using the dynamic pressure measured from a sensor connected to a pitot tube. The calculation formula is based on Bernoulli equation, and it yields
(6)Vpitot2=K2ΔPρ=u2+v2+w2,where ΔP is the dynamic pressure, ρ is the air density, and K is the correction factor [28].Attack angleα and sideslip angle β can be derived from the airspeed components in the body coordinate system [29]:
(7)α=arctanwu,(8)β=arcsinvu2+v2+w2.Through the angle of attack, the sideslip angle, andVpitot, airspeed Va in navigation coordinates can be expressed as
(9)Va=CbncosαcosβsinβsinαcosβVpitot.Therefore, the wind speedwx, wy, and wz can be computed directly from the wind speed triangle as [30]
(10)wxwywz=VxmVymVzm−CbncosαcosβsinβsinαcosβVpitot.This method is traditionally used for wind information, and it is compared in the Section4 for the better accuracy of our method.Besides of speed, the information of the wind field around the aircraft also includes the wind direction, namely, the wind azimuth and the wind pitch. The angle between the wind speed and the north is defined as the wind azimuth:(11)waz=arctanwzwx.The angle between the wind speed and the horizontal plane is defined as the wind pitch:(12)wan=arctanwywx2+wz2.Their accuracy is also compared in the simulations.


## 2.1. Introduction of Coordinates

First, the coordinates that are used in the paper are defined.On−XnYnZn denotes the navigation frame which coincides with the geographic frame (Og−XgYgZg), north east and down (NED).Ob−XbYbZb denotes the aircraft body coordinate frame. Ob is at the center of gravity, Xb is in the plane of symmetry of the aircraft and parallel to the design axis of the aircraft and points to the nose, Yb is perpendicular to the plane of symmetry of the aircraft and points to the right of the fuselage, and Zb is in the plane of symmetry of the aircraft, perpendicular to the Xb and pointing below the fuselage. The rotation from NED to the aircraft body coordinate system is defined by Euler angles, as shown in Figure 1.Figure 1
Relationship between the aircraft body coordinate frame and the geographic system.Ow−XwYwZw denotes the wind coordinate frame. Ow is at the center of gravity and fixed connection with aircraft body; Xw is pointing at the direction of airspeed Va; Zw is in the plane of symmetry of the aircraft, perpendicular to the Xw and pointing below the fuselage; and Yb is defined by the right-hand rule. Rotation from the aircraft body coordinate frame to the wind coordinate frame is defined by the angle of attack α and sideslip β, as shown in Figure 2.Figure 2
Relationship between the wind coordinate frame and the aircraft body coordinate frame.


## 2.2. Dynamic Equations

The following quantities are studied in this paper:(1)
The forward, rightward, and downward airspeed componentsu, v, and w, respectively(2)
The rollϕ, pitch θ, and yaw ψ of the Euler attitude angles(3)
The northward, eastward, and downward wind velocity componentswx, wy, and wz(4)
The latitudeL, the longitude λ, and the altitude H of the aircraft(5)
The pitot tube velocityVpitot, angle of attack α, and sideslip angle β derived from the Air Data System(6)
The northward, eastward, and downward velocity (ground speed)VxmVymVzm and positions LmλmHm provided by the GNSS(7)
The accelerationsaxayaz and angular velocity pqr derived from the Inertial Measurement Unit (IMU).Firstly, based on the flight dynamics and the assumption that the plane is a rigid body, the motion equations of the aircraft are established as follows.The conversion from the aircraft body coordinate frame to the geographic coordinate system is defined by the Direction Cosine Matrix (DCM), which is derived from the continuous rotation of rollϕ, pitch θ, and yaw ψ of the aircraft.
(1)Cbn=DCMϕ,θ,ψ=cosψ−sinψ0sinψcosψ0001cosθ0sinθ010−sinθ0cosθ1000cosϕ−sinϕ0sinϕcosϕ.The following force equations can be obtained by decomposing the total aerodynamic forces and the thrustsFx,Fy,Fz into the aircraft body coordinate system and then decomposing the gravity into the aircraft body coordinate system, too:
(2)u̇=vr−wq−gsinθ+Fxm,v̇=−ur+wp+gcosθsinϕ+Fym,ẇ=uq−vp+gcosθcosϕ+Fzm,where m is the total mass of the aircraft which is assumed to be a constant. The relationship between the attitude rate ϕ̇θ̇ψ̇ and the three angular velocity components pqr in the body coordinate system can be described by the dynamic equations as follows [27]:
(3)ϕ̇=p+rcosϕ+qsinϕtanθ,θ̇=qcosϕ−rsinϕ,ψ̇=rcosϕ+qsinϕsecθ.According to the characteristic of randomness, the state model of wind speed is set as a random walk process [19].The dynamics for the position states are then given by groundspeeds of the aircraft:(4)L̇=VxRm+H,λ̇=VyRn+HsecL,Ḣ=Vz,where Rm=Re1−2e+3esin2L, Rn=Re1+esin2L, and Re=6378393 m, and e=1/297.3 is the oblation of rotating ellipsoid. VxVyVz are groundspeeds of the aircraft in navigation coordinates. They can be correlated with the wind speed in geographic coordinates and the airspeed in the aircraft body coordinate through the Direction Cosine Matrix as follows:
(5)VxVyVz=Cbnuvw+wxwywz.Secondly, the mathematical models related with pitot speed, angle of attack, and sideslip angle are presented. The pitot speedVpitot is calculated using the dynamic pressure measured from a sensor connected to a pitot tube. The calculation formula is based on Bernoulli equation, and it yields
(6)Vpitot2=K2ΔPρ=u2+v2+w2,where ΔP is the dynamic pressure, ρ is the air density, and K is the correction factor [28].Attack angleα and sideslip angle β can be derived from the airspeed components in the body coordinate system [29]:
(7)α=arctanwu,(8)β=arcsinvu2+v2+w2.Through the angle of attack, the sideslip angle, andVpitot, airspeed Va in navigation coordinates can be expressed as
(9)Va=CbncosαcosβsinβsinαcosβVpitot.Therefore, the wind speedwx, wy, and wz can be computed directly from the wind speed triangle as [30]
(10)wxwywz=VxmVymVzm−CbncosαcosβsinβsinαcosβVpitot.This method is traditionally used for wind information, and it is compared in the Section4 for the better accuracy of our method.Besides of speed, the information of the wind field around the aircraft also includes the wind direction, namely, the wind azimuth and the wind pitch. The angle between the wind speed and the north is defined as the wind azimuth:(11)waz=arctanwzwx.The angle between the wind speed and the horizontal plane is defined as the wind pitch:(12)wan=arctanwywx2+wz2.Their accuracy is also compared in the simulations.


## 3. Estimation Method of Wind Speed and Direction

Since the GNSS/INS/ADS integrated navigation system possesses nonlinear characteristics, nonlinear estimation algorithms such as EKF or UKF are needed. Because EKF is a first-order approximation of nonlinear system, while UKF has at least second-order approximation accuracy and Jacobian matrices are unnecessary, which reduce the computation complexity, UKF is adopted to estimate wind speed and direction online in this study. The research of Liu et al. [24] and Ding et al. [26] also trues this point.According to the dynamic equations detailed in Section2, the state vector X consists the following elements: the airspeeds uvwT in the body coordinate frame, Euler attitude angles ϕθψT, three-axis wind speed wxwywzT, and three-dimensional position LλHT. The input vector U consists of accelerations and angular rates obtained through IMU. The groundspeed and position measured by GNSS, airspeed, angle of attack, and sideslip angle provided by the ADS are adopted as the measurements Z:
(13)X=uvwϕθψwxwywzLλHT,U=axayazpqrT,Z=VxmVymVzmVpitotαβLmλmHmT.The state equation of the UKF is established as follows:(14)Xk=fXk−1,uk−1+Wk−1,where f is the nonlinear vector function detailed in equations (2), (3), and (4) and Xk and Xk−1 are the state vectors at k and k−1 moment, respectively. Wk‐1 is the zero-mean Gaussian process noise vector.The measurement equation is established as follows:(15)Zk=hXk+Vk,where h is a nonlinear observation function detailed in equations (5), (6), (7), and (8). Vk is the zero-mean Gaussian measurement noise vector.The UKF algorithms are as follows.First, select the initial filtering values:(16)X^0=EX0,P0=EX0−X^0X0−X∧0T,where E is the mathematical expectation, X0 is the initial value of the state vector, X^0 is the mean value of the initial value of the state vector, and P0 is the covariance matrix of the initial value of the state vector.Second, calculate the sigma sample points atk−1 moment:
(17)χ~k−10=X^k−1,χ~k−1i=X^k−1+n+λPk−1i,i=1,2,⋯,n,χ~k−1i=X^k−1‐n+λPk−1i‐n,i=n+1,n+2,⋯2n,where n is the state dimension and λ=α2n+κ−n, α∈10−4,1, κ=3‐n.Third, determine the weight:(18)W0m=λn+λ,W0c=λn+λ+1−α2+β,Wim=Wic=12n+λ,i=1,2,⋯,2n,where W0m and Wim are the weights of the state vector of the ith sigma point, W0c and Wic are the weights of the covariance matrix of the state vector of the ith sigma point, and β is the state distribution parameter. For Gaussian white noise, β=2.Fourth, calculate the one-step prediction model value ofk moment:
(19)χk/k−1i=fχ~k−1i,uk−1,i=0,1,2,⋯,2n,X^k/k−1=∑i=02nWimχk/k−1i,Pk/k−1=∑i=02nWicχk/k−1i−X^k/k−1χk/k−1i−X∧k/k−1T+Qk−1.X^k/k−1 is the prediction state vector, Pk/k−1 is the covariance matrix of the prediction state vector, and Qk−1 is the state vector noise matrix in the last moment.Fifth, calculate the one-step prediction sample point ofk moment:
(20)χk/k−10=X^k/k−1,χk/k−1i=X^k/k−1+n+λPk/k−1i,i=1,2,⋯,n,χk/k−1i=X^k/k−1‐n+λPk/k−1i‐n,i=n+1,n+2,⋯2n.Sixth is measurement update:(21)Zk/k−1i=hχk/k−1i,i=0,1,2,⋯,n,Z^k/k−1=∑i=02nWimZk/k−1i,PZZk/k−1=∑i=02nWiiZk/k−1i−Z^k/k−1Zk/k−1i−Z∧k/k−1T+Rk,PXZk/k−1=∑i=02nWicχk/k−1i−X^k/k−1Zk/k−1i−Z∧k/k−1T,where Z^k/k−1 is the prediction observation vector, Pzzk/k−1 is the covariance matrix of prediction vector and observation vector, PXZk/k−1 is the covariance matrix of prediction observation and state vector, and Rk is the observation noise matrix.Seventh is filter update:(22)Kk=PXZk/k−1PZZk/k−1−1,X^k=X^k/k−1+KkZk−Z^k/k−1,Pk=Pk/k−1−KkPZZk/k−1KkT,where Kk is the gain matrix, X^k is the estimated state vector, and Pk is the covariance matrix of the estimated state vector.


## 4. Simulation Verifications

In order to evaluate the feasibility and effectiveness of the estimation method for wind speed and direction, various simulations were carried out. The simulated trajectory lasts 900 s, and the three-dimension flight path is shown in Figure3.Figure 3
Flight path.In the simulation experiment, the sampling rate of the IMU is 100 Hz and the frequency of integrated navigation is 1 Hz. The measurement noise is assumed to be Gaussian white noise, and the standard deviations are shown in Table1.Table 1
Standard deviations of measurement noise.
ParameterVariableNoise standard deviationPosition (GNSS)Lm,λm,Hm1 mGround speed (GNSS)Vx,Vy,Vz0.2 m/sAcceleration (IMU)ax,ay,az0.001 gAngular rate (IMU)p,q,r0.01°/hAirspeedVa0.5 m/sAngle of attackα0.1°Sideslip angleβ0.1°The position error obtained by the integrated navigation is shown in Figure4. The standard deviations of latitude, longitude, and height errors are 0.014 m, 0.011 m, and 0.65 m, respectively.Figure 4
Position error obtained by integrated navigation.The speed error obtained by integrated navigation is shown in Figure5. The standard deviations of the north, east, and down speed errors are 0.070 m/s, 0.065 m/s, and 0.064 m/s, respectively.Figure 5
Speed error obtained by integrated navigation.In the first simulation, wind speed and direction are set to constants. The results of the proposed estimation method based on UKF estimation are depicted in Figure6. In the figure, the Measurements refer to raw wind speed information, which is calculated according to equation (10). The Estimates refer to the wind speed estimated by UKF, and Truth is the reference wind speed. As shown in Figure 6, the accuracy of wind speed obtained by the proposed method in this paper is much higher than that of the raw wind speed. The standard deviations of the wind speed error of the proposed method are 0.12 m/s, 0.14 m/s, and 0.13 m/s, respectively, while they are 0.33 m/s, 0.34 m/s, and 0.38 m/s for raw wind speed scenario.Figure 6
Wind speed measurements, estimates, and true values.After the wind speeds are obtained, the corresponding wind azimuth (as shown in Figure7) and pitch (as shown in Figure 8) can be acquired. The standard deviations of the wind direction error of the proposed method are 1.73° and 1.32°, respectively, while they are 4.70° and 3.28° for raw wind speed scenario.Figure 7
Wind azimuth measurements, estimates, and true values.Figure 8
Wind pitch measurements, estimates, and true values.In order to verify the reliability of the proposed method, we have carried out more than 1000 simulations and the statistical results are shown in Table2.Table 2
Standard deviation of estimated and measured wind speed for constant wind.
Speed: m/sAngle: °Group 1Group 2Group 3Standard deviation of measurementsStandard deviation of estimationsStandard deviation of measurementsStandard deviation of estimationsStandard deviation of measurementsStandard deviation of estimationsNorth0.330.120.430.120.530.11East0.330.130.310.130.310.12Down0.380.130.420.160.440.15Resultant0.330.120.430.130.540.12Azimuth4.781.694.701.834.701.73Pitch3.661.573.411.343.281.32Obviously, compared with the directly measured values, the accuracy of the estimated values is improved by nearly two times.In order to further demonstrate the effectiveness of the proposed method, we also carried out the simulations with randomly changing wind speed. The true wind gust fields are modeled using a frozen Dryden turbulence mode [31]. Air mass movement is expressed as a sum of sinusoids:
(23)w·=w·,0+∑n=1Na·,nsinΩ·,ns+ϕ·,n,where · represents the gust component and s is motion along the flight path. Random values of ϕ·,n simulate the random process, and the choice of coefficients a·,n defines the power spectral density.The simulation results are shown in Figure9. The standard deviations of the wind speed error of the proposed method are 0.50 m/s, 0.36 m/s, and 0.27 m/s, respectively, while they are 1.42 m/s, 0.72 m/s, and 0.53 m/s for raw wind speed scenario.Figure 9
Wind speed measurements, estimates, and true values.Then, the wind azimuth (as shown in Figure10) and pitch (as shown in Figure 11) are compared. The standard deviations of the wind direction error of the proposed method are 0.73° and 0.85°, respectively, while they are 1.41° and 1.23° for raw wind speed scenario.Figure 10
Wind azimuth measurements, estimates, and true values.Figure 11
Wind pitch measurements, estimates, and true values.Similarly, more than 1000 simulations are carried out and the statistical results are shown in Table3.Table 3
Standard deviation of estimated and measured value of random wind speed.
Speed: m/sAngle: °Group 1Group 2Group 3Standard deviation of measurementsStandard deviation of estimationsStandard deviation of measurementsStandard deviation of estimationsStandard deviation of measurementsStandard deviation of estimationsNorth4.620.511.420.501.710.50East1.350.370.720.360.580.37Down2.150.270.530.270.850.26Resultant5.230.601.460.591.930.59Azimuth1.350.731.410.731.320.75Pitch1.090.881.230.850.950.87In order to realize the establishment of the Dryden model of random wind field from low to high frequency and considering the amount of calculation,N=5 in equation (23) is selected. The resultant velocity of the multifrequency random wind field generated by the simulation is shown in Figure 12.Figure 12
Multifrequency random wind speed.The corresponding simulation results are shown in Figure13. The standard deviations of the wind speed error of the proposed method are 0.25 m/s, 0.35 m/s, and 0.27 m/s, respectively, while they are 0.69 m/s, 0.55 m/s, and 0.61 m/s for raw wind speed scenario.Figure 13
Multifrequency random wind speed measurements, estimates, and true values.The corresponding wind azimuth and pitch are shown in Figures14 and 15. The standard deviations of the wind direction error of the proposed method are 0.83° and 0.88°, respectively, while they are 1.65° and 1.43° for raw wind speed scenario.Figure 14
Multifrequency random wind azimuth measurements, estimates, and true values.Figure 15
Multifrequency random wind pitch measurements, estimates, and true values.Similarly, more than 1000 simulations are carried out and the statistical results are shown in Table4.Table 4
Standard deviation of estimated and measured value of multifrequency random wind speed.
Speed: m/sAngle: °Group 1Group 2Group 3Standard deviation of measurementsStandard deviation of estimationsStandard deviation of measurementsStandard deviation of estimationsStandard deviation of measurementsStandard deviation of estimationsNorth0.680.230.820.250.640.24East0.630.440.570.390.520.40Down0.760.330.630.340.590.34Resultant0.700.280.820.290.650.27Azimuth1.860.911.740.921.590.92Pitch1.551.081.650.951.320.94For all the set of simulations, the proposed wind speed estimation method provides wind speed with accuracy nearly two times better than raw measurements.


## 5. Conclusions

This paper proposes a novel online estimation method for the wind speed and direction around the fuselage, which can be used to improve the performance of on-time arrival in 4DT applications. The estimation method fuses INS, GNSS, and ADS in an UKF and provides great accuracy and robustness in nonlinearity conditions. The dynamic models of wind are established, and implementations of the UKF are detailed. Through flight simulations, the effectiveness of the method is verified by comparing with the traditional direct measurement method. Results demonstrate that the accuracy of wind speed and direction obtained by our method is nearly two times higher than the traditional direct measurement method.



---

*Source: 6844509-2021-12-20.xml*