Open AccessArticle

Correction of Error of Airborne Anemometers Caused by Self-Excited Air Turbulence

Jianqiang Liu

^1,2,

Zhan Zhao

^1,2,

Zhen Fang

^1,2,

Yong Li

³ and

Lidong Du

^1,2,*

State Key Laboratory of Transducer Technology, Aerospace Information Research Institute, Chinese Academy of Sciences, Beijing 100194, China

University of Chinese Academy of Sciences, Beijing 100049, China

Key Laboratory of Transportation Meteorology, China Meteorological Administration, Nanjing 210008, China

Author to whom correspondence should be addressed.

Sensors 2023, 23(9), 4288; https://doi.org/10.3390/s23094288

Submission received: 8 March 2023 / Revised: 15 April 2023 / Accepted: 24 April 2023 / Published: 26 April 2023

(This article belongs to the Special Issue Advanced Sensing and Safety Control for Connected and Automated Vehicles: Volume Ⅱ)

Download

Browse Figures

Figure 1
The working status of the six rotors. "> Figure 2
MRF model diagram. "> Figure 3
The internal structure of the differential pressure anemometer. "> Figure 4
Schematic diagram of <math display="inline"><semantics> <mrow> <msub> <mrow> <mi>P</mi> </mrow> <mrow> <mi>D</mi> <mn>1</mn> </mrow> </msub> </mrow> </semantics></math> and <math display="inline"><semantics> <mrow> <msub> <mrow> <mi>P</mi> </mrow> <mrow> <mi>D</mi> <mn>2</mn> </mrow> </msub> </mrow> </semantics></math>. "> Figure 5
Display of AOA and <math display="inline"><semantics> <mrow> <mi mathvariant="sans-serif">α</mi> </mrow> </semantics></math>. "> Figure 6
Grid-independence test. "> Figure 7
Computing domain and the encryption area. "> Figure 8
The mesh diagram of the UAV. "> Figure 9
The mesh diagram of the dynamic area of the UAV rotors. "> Figure 10
The UAV velocity flow field diagram at ascent speed of 0 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s. "> Figure 11
The UAV velocity flow field diagram at ascent speed of 3 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s. "> Figure 11 Cont.
The UAV velocity flow field diagram at ascent speed of 3 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s. "> Figure 12
The UAV velocity flow field diagram at ascent speed of 5 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s. "> Figure 13
Curves between simulation crosswind speed and standard crosswind speed under different ascent speeds. "> Figure 14
Wind speed measurement results between the meteorological tower and airborne anemometer. "> Figure 15
Wind speed measurement errors of the airborne anemometer. ">

Versions Notes

Abstract

An airborne anemometer, which monitors wind on the basis of Meteorological Multi-rotor UAVs (Unmanned Aerial Vehicles), is important for the prevention of catastrophe. However, its performance will be affected by the self-excited air turbulence generated by UAV rotors. In this paper, for the purpose of the correction of an error, we developed a method for the elimination of the influence of air turbulence on wind speed measurement. The corresponding correction model is obtained according to the CFD (Computational Fluid Dynamics) simulation of a six-rotor UAV which is carried out with the sliding grid method and the S-A turbulence model. Then, the model is applied to the developed prototype by adding the angle of attack compensation model of the airborne anemometer. It is shown by the actual application that the airborne anemometer can maintain the original measurement accuracy at different ascent speeds.

Keywords:

air turbulence error; CFD simulation; multi-rotor UAVs; meteorological observation

1. Introduction

Multi-rotor UAVs have prevailed in many fields such as chemical [1], agricultural [2,3] and meteorological monitoring. By integrating miniaturized instruments, they have greatly promoted the development of scientific, industrial, and regulatory fields, especially in meteorological environment monitoring. It has great advantages over traditional automatic weather stations (AWS), satellites, remote sensing, and other measurement methods. As a platform for meteorological monitoring, multi-rotor UAVs can collect sensor data more sensitively and timely, and can obtain data with high spatial and temporal resolution [4]. A lot of research has been initiated in recent years. The US and Europe have begun to use UAVs as important instruments for disaster and environmental monitoring [5,6]. Brooke Potter et al. [7] made use of a UAV to collect data from a remote stream site. Zhewen Xing [8] and Ruisheng Ma [9] used multi-rotor UAVs to monitor meteorological disasters. Daniel Leuenberger et al. [10] used drones to improve the accuracy of weather forecasts.

Although multi-rotor UAVs have advantages in various measurement tasks, there is an urgent demand to resolve the effect of air turbulence generated by rotors. Many researchers have done a lot of meaningful work. Seokkwan Yoon et al. [11] calculated and simulated the airflow of rotors to study the best separation distance between the fuselage and the wings. Neal [12] solved the time-dependent Navier–Stokes equations for isolated rotors in hover and forward flight using detached eddy simulation and adaptive mesh refinement. Scott E. [13] used a fixed LBM grid and an adaptive refinement method to establish a simulation model for the four rotors of the drone. Qiwei Guo et al. [14] studied the formation process and flow distribution of the downwash airflow of the quadrotor UAV, and established the calculation model of the downwash airflow of the quadrotor agricultural UAV using CFD simulation. Hao Zhang et al. [15] studied the downwash airfield distribution of a six-rotor UAV when hovering at different flight speeds and altitudes, and performed numerical simulations on the airflow field. Most of the papers are about the simulation and modeling of the downwash airflow for multi-rotor UAVs. However, for many meteorological monitoring UAVs, the sensors are established on top of the multi-rotor UAVs. Upwash airflow excited by the multi-rotor UAVs will disturb the sensor even more.

Therefore, the research on the influence of the upwash airflow on the multi-rotor UAVs is more significant, especially the anemometer-involved application. When multi-rotor drones are used as airborne anemometers, the impact of rotor airflow should be compensated. The angle of attack (AOA) of the multi-rotor UAVs will also affect the performance of the anemometer as well. It is more urgent to resolve the union effects which come from these two weak points. Taro Nakai et al. [16,17] made a very prominent contribution to the correction of the AOA. They improved the accuracy of the correction method for ultrasonic wind sensors. In this paper, the differential pressure anemometer developed by Cheng Liu and Yichen Pan [18,19], is used. Although it can maintain its original measurement accuracy in the AOA range of 0–45 degrees, the union effects still need to be corrected when it is used in multi-rotor UAVs.

In this paper, for the purpose of the correction of error, we developed a method for the elimination of the influence of air turbulence on wind speed measurement. The corresponding correction model is obtained according to the CFD (Computational Fluid Dynamics) simulation of a six-rotor UAV which is carried out with the sliding grid method and the S-A turbulence model. Then, the model is applied in the developed prototype by adding the angle of attack compensation model of the airborne anemometer. The model has been verified in actual measurement, and it can make the airborne anemometer maintain the original wind speed measurement accuracy in the angle of attack range of 0–45° at various ascent speeds.

2. Methods

The UAV used in this work is a common six-rotor UAV, which has six propellers, and all its attitude and position control are achieved by adjusting the speed of the six driving motors. When the UAV is working normally, the three propellers are separated by 120 degrees rotate clockwise, and the other three propellers rotate counterclockwise, as shown in Figure 1. In general, the motion state of a six-rotor UAV is mainly divided into five types: hovering, vertical motion, rolling motion, pitching motion, and yaw motion. Only the hovering and vertical motions are simulated in this paper to study the impact on the anemometer because the two types often occur in measurement scenarios.

2.1. Basic Control Equation

In the process of UAV flight, it is difficult to study the complex flow field and phenomenon generated by the rotation of the rotor using traditional aerodynamics. With the continuous development of computer technology and numerical methods, the use of computational fluid dynamics to calculate and simulate the rotor flow field has become one of the important methods for studying the characteristics of the rotor flow field.

The flying speed of meteorological UAVs is low, and the ascent speed is within 5 m/s normally. Therefore, the air medium in the external flow field can be regarded as incompressible. Navier–Stokes (NS) equations are the most suitable differential equation to express incompressible fluid. The NS equation reflects the basic laws of viscous fluids, and it relies on differential equations to describe fluid motion. The three-dimensional incompressible N-S equation is expressed as follows:

\{\begin{matrix} ρ \frac{D u}{D t} = ρ f_{x} - \frac{\partial p}{\partial x} + μ \nabla^{2} u \\ ρ \frac{D v}{D t} = ρ f_{y} - \frac{\partial p}{\partial y} + μ \nabla^{2} v \\ ρ \frac{D w}{D t} = ρ f_{z} - \frac{\partial p}{\partial z} + μ \nabla^{2} w \end{matrix}

(1)

where u, v, and w are the components of the dimensionless velocity along the x, y, and z directions, p and t are the dimensionless pressure and time,

f_{x}

f_{y},

and

f_{z}

denote the components of the external force per unit volume of fluid in the x, y, and z directions, respectively. Multiply the above equations by the unit vectors

i

j

, and

k

in the three directions and add them to obtain the simpler vector form of the N-S equation for incompressible viscous fluid:

\frac{D \vec{V}}{D t} = \vec{f} - \frac{1}{ρ} \nabla p + \frac{μ}{ρ} \nabla \vec{V}

(2)

where

\vec{V}

is the velocity vector,

\nabla

is the Hamiltonian, and

\vec{f}

is the total external force per unit volume of fluid.

2.2. Calculation Method

When using the CFD method to simulate the rotor flow field, there are two main methods. The first method is to use the Actuator Disk theory [20] to equate the rotating blade with an actuator disk. The momentum source method [21] is a kind of actuator disk method. Its basic idea is that the action of the blade on the airflow is added to the governing equations (Euler or N-S) equivalent to the time-averaged momentum source term. In this way, the effect of the blade on the airflow is characterized by the change of the airflow. The second method is the sliding grid method, which generates a body-fitted grid around each blade, and uses the entire rotor grid system as a motion-nested grid. In this grid, the rotor flow field is simulated by solving the Euler equation or N-S equation. Essentially, the rotor rotation of a multi-rotor UAV belongs to the mechanical rotation, so a simple and adaptable sliding grid can be used to complete the calculation of various states with a multi-reference (MRF) system model.

The MRF model is one of the multi-region calculation methods, which uses a steady-state approximation. Different rotation or movement speeds can be assumed in each region. The equations of the motion reference system are used to solve the flow problem in each motion area grid. On the interface of the computational domain, a local reference system is used to calculate the flux of the flow variables in one area and convert them to adjacent areas. The schematic diagram of a typical MRF system model is shown in Figure 2. It is a coordinate system that rotates at a stable angular velocity

\vec{w}

for a stationary reference system. The origin of the rotating system is positioned by the position vector

\vec{r}

The position of any point in the calculation domain of the rotation system can be determined by the position vector

\vec{r}

and the origin of the rotation system. The implicated velocity can be expressed as follows:

\vec{u_{r}} = \vec{w} \times \vec{r}

(3)

The velocity

\vec{u_{r}}

can be converted from a stationary system to a rotating system by the following equation:

\vec{v_{r}} = \vec{v} - \vec{u_{r}}

(4)

where

\vec{v_{r}}

is the relative velocity and

\vec{v}

is the absolute velocity. When solving the problem of multiple moving individuals in a rotating coordinate system, the additional term in the momentum equation will cause the fluid acceleration to increase. The fluid governing equations in the form of relative velocity are shown as follows:

\frac{\partial ρ}{\partial t} + \nabla 🞄 ρ \vec{v_{r}} = 0

(5)

\frac{\partial}{\partial t} (ρ \vec{v_{r}}) + \nabla \cdot (ρ \vec{v_{r}} \vec{v_{r}}) + ρ (2 \vec{w} \times \vec{v_{r}} + \vec{w} \times \vec{w} \times \vec{v_{r}}) = - \nabla p + \nabla \cdot \vec{τ} + \vec{F}

(6)

\frac{\partial}{\partial t} (ρ E_{r}) + \nabla \cdot (ρ \vec{v_{r}} H_{r}) = \nabla \cdot (k \nabla T + \overset{̿}{τ_{r}} \cdot v_{r}) + s_{h}

(7)

where Equation (5) is the continuity equation, Equation (6) is the momentum equation, and Equation (7) is the energy equation. The momentum equation contains two additional acceleration terms: Coriolis acceleration

2 \vec{w} \times \vec{v_{r}}

and centripetal acceleration

\vec{w} \times \vec{w} \times \vec{v_{r}}

. Compared with the original equation, the viscous stress

\overset{̿}{τ_{r}}

uses the relative velocity derivative term. The energy equation uses relative internal energy

E_{r}

and relative total enthalpy

H_{r}

, and these variables are defined as:

E_{r} = h - \frac{p}{ρ} + \frac{1}{2} (v_{r}^{2} - u_{r}^{2})

(8)

H_{r} = E_{r} + \frac{p}{ρ}

(9)

In a sliding grid, the relative motion between the stationary and rotating parts causes transient interaction effects, which is a strong unsteady phenomenon, but these effects are ignored in the MRF system. The sliding grid technology uses two or more calculation areas, each area can generate a grid independently, which is extremely convenient for complex models. There is at least one interface between each area and adjacent areas. The interface of adjacent computing areas forms a “grid boundary”, and the dynamic domain will move along the interface. The grid on the interface does not need to be aligned, and the flux is calculated by the information interpolation between the grid nodes. A virtual grid layer is generated on both sides of the slip surface, which overlaps the computational domain grid on both sides of the sliding surface. During calculation, the nodes on the virtual grid layer are interpolated to realize the flux transfer on the computational domains on both sides of the interface.

When using the sliding grid technology for numerical simulation, the model needs to be divided into two parts: the rotor part and the stator part, and these two parts have meshed separately. In this paper, the rotor part is the cylindrical area where the propeller rotates, and the stator part is the entire computational domain minus other areas of the rotor part. In the modeling, the connecting parts of the rotor part and the stator part are paired to form multiple interfaces.

2.3. Turbulence Model

In this paper, the method of numerical simulation calculation of the UAV flow field is the S-A turbulence model which is widely used in aviation. Compared with the k-ε turbulence model, the S-A turbulence model is more robust in simulating and calculating complex flows and consumes fewer computing resources. The S-A turbulence model is based on a transport equation of eddy viscosity containing the convection term, diffusion term, and source term. This application was proposed by Spalart and Allmaras [22]. Ashford and Powell [23] improved this to avoid negative values in the generated term. The fluctuating amount

\tilde{v}

of turbulent kinetic energy can be obtained from the transport equation:

\frac{\partial v}{\partial t} + \vec{V} \cdot \nabla \tilde{v} = \frac{1}{σ} \{\nabla \cdot [v + (1 + c_{b 2}) \tilde{v} \nabla \tilde{v}] - c_{b 2} \tilde{v} \nabla \tilde{v}\} + Q

(10)

where

\vec{V}

is the mean velocity,

Q

is the source term,

σ

and

c_{b 2}

are constant. Source term

Q

contains the generating term and dissipative term as follows:

Q = \tilde{v} P (\tilde{v}) - \tilde{v} D (\tilde{v})

(11)

\tilde{v} P (\tilde{v}) = c_{b 1} S \tilde{v}

(12)

\tilde{v} D (\tilde{v}) = c_{w 1} f_{2} {(\frac{\tilde{v}}{d})}^{2}

(13)

The generating term can be obtained by Equations (14)–(16) in the following:

\tilde{S} = S f_{v 3} + \frac{\tilde{v}}{k^{2} d^{2}} f_{v 2}

(14)

f_{v 2} = \frac{1}{{(1 + \frac{χ}{c_{v 2}})}^{3}}

(15)

f_{v 3} = \frac{(1 + χ f_{v 1}) (1 - χ f_{v 2})}{χ}

(16)

where d is the minimum distance to the wall surface, S is the vorticity.

f_{m}

can be obtained by Equations (17)–(19) as follows:

f_{v 2} = g {(\frac{1 + c_{w}^{6}}{g^{6} + c_{w}^{6}})}^{6}

(17)

g = r + c_{w 2} (r^{6} - r)

(18)

r = \frac{\tilde{v}}{\tilde{S} k^{2} d^{2}}

(19)

The constant value in the S-A turbulence model is:

c_{w 1} = 3.239, c_{w 2} = 0.3, c_{w 3} = 2, c_{v 1} = 7.1, c_{v 2} = 5, c_{b 1} = 0.1355, c_{b 2} = 0.622, k = 0.41, σ = 0.667

(20)

2.4. Correction Model of the Angle of Attack

In this paper, a solid-state differential pressure anemometer is mounted on the multi-rotor UAV, and its structure is shown in Figure 3. The principle of the differential pressure anemometer in this work is that the differential pressure between the two ends of the cylinder varies with the wind speed. According to the variation in the differential pressure and distribution, the corresponding wind speed and wind direction can be calculated. Figure 4 is a schematic diagram of the anemometer measurement.

The relationship between differential pressure and wind speed and direction can be expressed as follows [24]:

U_{\infty} = 2 \sqrt{\frac{P_{D 2}}{ρ (a \sin 2 θ + a + 2 b)}}

(21)

θ = \frac{1}{2} [a r c c o s \frac{(R_{D} - 1) (1 + \frac{2 b}{a})}{\sqrt{{(R_{D})}^{2} + 1}} - a r c t a n R_{D}]

(22)

where

U_{\infty}

is the wind speed,

ρ

is the air density, and a and b are the correction coefficients obtained by fitting the measured data.

R_{D}

is the ratio of the two largest differential pressures(

P_{D 1} a n d P_{D 2}

), which is expressed by:

R_{D} = \frac{P_{D 1}}{P_{D 2}}

(23)

When the angle of attack is greater than 15°, the measurement result is affected [25]. In this paper, the angle of attack is the angle between the wind speed vector and the anemometer measurement plane where the eight holes are located, as shown in Figure 5. The tilt angle α between the anemometer and the vertical axis is used to replace the angle of attack because it is equivalent to it and can easily be obtained by the accelerometer inside the anemometer in practical applications.

According to previous work, the angle of attack error of the anemometer can be corrected and compensated by the model as shown in Equations (24)–(27), so that the anemometer can maintain the original measurement accuracy and range [19].

U_{T} = 2 \sqrt{\frac{P_{T D 2}}{ρ (a \sin 2 θ + a + 2 b)}}

(24)

P_{T D 2} = \frac{P_{D 2}}{T (g (α, θ))}

(25)

T (α_{r}) = a_{0} + a_{1} \cos α_{r} + a_{2} \cos {α_{r}}^{2}

(26)

g (α, θ) = α_{r} = \arcsin (\sin α \cdot \cos θ)

(27)

where

U_{T}

is the corrected wind speed under

α

P_{T D 2}

is the second-largest differential pressure under the tilt angle

α

and

T (g (α, θ))

represents the influence of the angle of attack on the pressure distribution. For Equations (24)–(27), there is a detailed derivation process and explanation in reference [19], which will not be introduced here.

3. Simulation and Modeling

3.1. Mesh and Boundary Conditions

In this paper, the geometric model is very complicated, there are small gaps between the rotating area and the static area. To better express as many detailed areas as possible, an unstructured grid method is used for numerical simulation calculations. For the calculation of the external flow domain of CFD, the larger the flow domain, the smaller the interference of the external flow field boundary on the flow field calculation. This requires the flow field to be set as large as possible during the calculation. However, a large computing domain needs to consume too many computing resources. When the calculated flow domain size reaches a certain range, the calculation accuracy remains stable. Before the formal simulation, we conducted a grid independence test to ensure the optimal grid size and distribution while maintaining calculation accuracy. The mesh size of the area where the airflow changes drastically and the area close to the surface is set to be smaller, and the mesh size of the area where the airflow is stable to be larger to keep the accuracy of calculation and save the calculation resources. We initially divided 5,376,248 mesh cells roughly according to the above rules, calculated the maximum wind speed error, and then refined the entire mesh four times. The final grid-independence test result is shown in Figure 6. The red triangle in the figure represents the number of grids divided in our five simulations and the corresponding maximum error.

According to the results of the grid-independence analysis, we ultimately used 10,803,973 grids for subsequent simulations. The dynamic calculation area is selected to cover the adjacent area of the propeller blades, as this area has the most significant impact on the air motion related to the rotation of the propeller. After further increasing the calculation area, the simulation results do not show significant changes, but the calculation time greatly increases. Therefore, after considering the calculation amount and simulation accuracy, we choose the adjacent area of the propeller blades as the dynamic calculation area. The external flow domain selected in this paper is shown in Figure 7. The “encryption area” in Figure 7 refers to the outer region of the “computing domain”, which is usually used to avoid the influence of boundary effects on the calculation results and to improve the computational efficiency. Since the “computing domain” can be divided into multiple regions for parallel computing, adding an external “encryption area” can expand the range of the computing domain, thereby improving the accuracy and efficiency of numerical simulation. The size of the flow domain is about eight times the size of the UAV simulation model. Due to the symmetry of the UAV model, symmetrical boundary conditions are used in the calculation, and only half of the UAV model is calculated for the flow field, which saves computing resources without sacrificing calculation accuracy. The mesh diagram of the UAV is shown in Figure 8. The copter has a length of 1 m and a width of 1 m, and each blade has a length of 30 cm and a width of 5 cm. The anemometer is located at a height of 5 cm above the rotor. The number of mesh cells here is 3,689,216.

For the sliding grid of the UAV, the rotors of the UAV should be wrapped in the rotation area. The static domain and each dynamic domain use interfaces to transfer data. The mesh diagram of the dynamic area of the UAV rotors is shown in Figure 9. The number of mesh cells here is 8,795,339.

According to the calculation requirements in this work, the entire calculation domain is divided into two parts, the dynamic domain, and the static domain. The relevant boundary conditions include the object boundary conditions, the far-field boundary conditions, and the interface boundary conditions. The surface of the aircraft model is set with no slippage and no penetration. The contact surface between the flow domain of the UAV and the outer flow domain is set as interfaces. Similarly, the contact surfaces between the flow domain of the rotors and the overall flow domain of the UAV are set as interfaces, which allows the two-flow domain to exchange data during the calculation process. Except for the symmetry plane, the surface of the flow domain of the UAV is all set as a velocity inlet to simulate the realistic flow field of the UAV during flight.

3.2. Simulation and Results

The measurement accuracy of the wind sensor on the top of the UAV will be affected by the airflow generated by the rotors. It is necessary to compare the simulation value and the standard value of the wind speed of the UAV under different ascent speeds and different crosswind conditions through CFD simulation. The wind speed measurement error of the UAV will be corrected by comparing the two values. In the paper, the UAV velocity flow field diagram is obtained through CFD simulation at different ascent speeds of 0 m/s, 3 m/s, and 5 m/s and different crosswind speeds of 0 m/s, 3 m/s, 5 m/s, 7 m/s, 10 m/s, 13 m/s, 15 m/s, 17 m/s, and 20 m/s. Figure 10, Figure 11 and Figure 12 are the velocity flow field diagrams of the UAV under different ascent speeds and crosswind speeds. The drone is hovering at the ascent speed of 0 m/s. When the crosswind speed is 0 m/s, the theoretical value of the wind speed measured by the UAV wind sensor should be 0 m/s. However, the flow field diagram shows a flow velocity exists at the wind sensor position, which indicates that the airflow driven by the rotation of the UAV rotors will affect the measurement results. When the UAV is rising at a constant speed, the UAV rotors have different effects on the flow velocity at the wind sensor position under different crosswind speeds.

To clearly express the influence of the UAV rotors on wind speed measurement, specific simulation crosswind speed values under different standard crosswind speeds and ascent speeds are listed in Table 1. As shown in Table 1, when the crosswind speed of the drone is the same, the higher the ascent speed, the closer the simulation speed is to the standard crosswind speed. As the speed in the UAV flow field increases, the influence of the motion of the UAV rotors has a smaller effect on the flow field near the wind sensor.

3.3. Modeling

The wind speed measurement error of the airborne anemometer comes from the angle of attack of the UAV and the air turbulence of the rotors. Combining the two correction models can well correct the wind speed measurement error of the airborne anemometer. Figure 13 shows curves between simulation crosswind speed and standard crosswind speed under different ascent speeds based on the data in Table 1. It shows that the curves under the three ascent speeds are almost the same. In other words, although the airflow of the UAV rotors has an influence on the wind speed measurement, the speed of the rotor within 5 m/s is not related to it.

The three curves in Figure 13 are fitted by the least-squares method to obtain Equation (28) as expressed:

V_{r} = c_{0} V_{m} + c_{1}

(28)

where

V_{r}

is the real crosswind speed and

V_{m}

is the measured crosswind speed.

c_{0}

and

c_{1}

are the fitting coefficients. The values of

c_{0}

and

c_{1}

are the result of averaging the coefficients of three formulas which is expressed in Equation (29):

c_{0} = 0.9456 c_{1} = - 0.1573

(29)

According to Equations (24) and (28), the wind speed measurement correction model of the airborne anemometer can be expressed as:

U_{R} = 2 c_{0} \sqrt{\frac{P_{T D 2}}{ρ (a \sin 2 θ + a + 2 b)}} + c_{1}

(30)

4. Test and Results

The wind speed correction model of the airborne anemometer is obtained by combining the angle of attack correction model and the air turbulence correction model, and the model was verified through a UAV flight test. The drone flies at different ascent speeds near the meteorological tower at a height of 70 m and performs the wind speed measurement. The measurement results of the cup anemometer in the meteorological tower are used as standard data for comparison with the measurement results from the airborne anemometer before and after correction, which are shown in Figure 14. The wind speed measured by the drone after correction is dynamically changing and consistent with the measurement result of the cup anemometer, while the measurement results from the airborne anemometer before correction had a larger error. It can be seen from the figure that there are some deviations between the test points and the standard value. This is because the cup anemometer and the airborne anemometer are in close positions but not absolutely the same. In the boundary layer, the uneven airflow causes this deviation, but these deviations are within a reasonable range.

To clearly verify the compensation model, the measurement error curve of the airborne anemometer is drawn, as shown in Figure 15. In this paper, the wind speed measurement error of the anemometer is ±(0.5 + 0.03 V) m/s (V is the standard wind speed). The error bar in Figure 15 is obtained according to the standard value, which is the reason for its dynamic change. It can be seen from the figure that the wind speed measurement errors of the airborne anemometer are all within the error bar, which verifies that the model has a good correction effect.

5. Conclusions

In this paper, the influence of the air turbulence generated by the rotors of UAVs on the measurement of the airborne anemometer is studied. The CFD simulation of the UAV is carried out using the sliding grid method and the S-A turbulence model. The relationship between the measured wind speed and the standard wind speed was obtained, and an air turbulence correction model was established. The angle of attack compensation model of the differential pressure anemometer is added to the air turbulence correction model to make it more practical.

The model is verified in the actual measurement, and the result shows the model has a good correction effect. The airborne anemometer maintains the original measurement accuracy at different ascent speeds. This study proves that for a six-rotor UAV, the air turbulence generated by the rotors has an impact on the measurement, but it is not related to the speed of rotors within 5 m/s. Since the ascent speed of meteorological UAVs is low, this paper does not study speeds above 5 m/s. The effect of the high-speed rotating rotors on the airflow needs to be further explored. Whether this model has the same corrective effect on UAVs with other rotor numbers or UAVs with other structures requires further research and verification.

Author Contributions

J.L., Z.Z., Z.F., Y.L. and L.D. collaborated on the research for various aspects of the paper. L.D. and Z.Z. designed the mathematical model; J.L. and L.D. designed and performed the experiments; Z.F. and Y.L. analyzed the data; J.L. wrote the paper; Y.L. polished the paper. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported by the National Key Research and Development Project (2021YFC3000304).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

No new data was created or analyzed in this study. Data sharing is not applicable to this article.

Conflicts of Interest

The authors declare no conflict of interest.

References

Javier, B.; Marco, S. Environmental chemical sensing using small drones: A review. Sci. Total Environ. 2020, 748, 141172. [Google Scholar] [CrossRef]
Ahirwar, S.; Swarnkar, R.; Bhukya, S.; Namwade, G. Application of drone in agriculture. Int. J. Curr. Microbiol. Appl. Sci. 2019, 8, 2500–2505. [Google Scholar] [CrossRef]
Mogili, U.M.R.; Deepak, B. Review on application of drone systems in precision agriculture. Procedia Comput. Sci. 2018, 133, 502–509. [Google Scholar] [CrossRef]
Li, C.; Sun, X. A novel meteorological sensor data acquisition approach based on unmanned aerial vehicle. Int. J. Sens. Netw. 2018, 28, 80–88. [Google Scholar] [CrossRef]
Malaver Rojas, J.; Gonzalez, F.; Motta, N.; Depari, A.; Corke, P. Towards the development of a gas sensor system for monitoring pollutant gases in the low troposphere using small unmanned aerial vehicles. In Proceedings of the Workshop on Robotics for Environmental Monitoring, 11 July 2012; pp. 1–3. Available online: http://eprints.qut.edu.au/52042/ (accessed on 5 December 2022).
Zhang, J.; Tang, J.; Wang, T.; Chen, F. Energy-efficient data-gathering rendezvous algorithms with mobile sinks for wireless sensor networks. Int. J. Sens. Netw. 2017, 23, 248–257. [Google Scholar] [CrossRef]
Potter, B.; Valentino, G.; Yates, L.; Benzing, T.; Salman, A. Environmental Monitoring Using a Drone-Enabled Wireless Sensor Network. In Proceedings of the Systems and Information Engineering Design Symposium (SIEDS), Charlottesville, VA, USA, 26 April 2019. [Google Scholar] [CrossRef]
Xing, Z.; Zhang, Y.; Su, C.Y.; Qu, Y.; Yu, Z. Kalman filter-based wind estimation for forest fire monitoring with a quadrotor UAV. In Proceedings of the 2019 IEEE Conference on Control Technology and Applications (CCTA), Hong Kong, China, 19–21 August 2019; IEEE: Piscataway, NJ, USA, 2019; pp. 783–788. [Google Scholar] [CrossRef]
Ma, R.; Li, X.; Sun, M.; Kuang, Z. Experiment of Meteorological Disaster Monitoring on Unmanned Aerial Vehicle. In Proceedings of the 2018 7th International Conference on Agro-Geoinformatics (Agro-Geoinformatics), Hangzhou, China, 6–9 August 2018; IEEE: Piscataway, NJ, USA, 2018; pp. 1–6. [Google Scholar] [CrossRef]
Leuenberger, D.; Haefele, A.; Omanovic, N.; Fengler, M.; Martucci, G.; Calpini, B.; Fuhrer, O.; Rossa, A. Improving high-impact numerical weather prediction with lidar and drone observations. Bull. Am. Meteorol. Soc. 2020, 101, E1036–E1051. [Google Scholar] [CrossRef]
Yoon, S.; Lee, H.C.; Pulliam, T.H. Computational analysis of multi-rotor flows. In Proceedings of the 54th AIAA Aerospace Sciences Meeting, San Diego, CA, USA, 4–8 January 2016; p. 0812. [Google Scholar] [CrossRef]
Chaderjian, N.M. Advances in rotor performance and turbulent wake simulation using DES and adaptive mesh refinement. In Proceedings of the 7th International Conference on Computational Fluid Dynamics, Big Island, HI, USA, 9–13 July 2012. [Google Scholar]
Thibault, S.; Holman, D.; Garcia, S.; Trapani, G. CFD Simulation of a quad-rotor UAV with rotors in motion explicitly modeled using an LBM approach with adaptive refinement. In Proceedings of the 55th AIAA Aerospace Sciences Meeting, Grapevine, TX, USA, 9–13 January 2017; p. 0583. [Google Scholar] [CrossRef]
Guo, Q.; Zhu, Y.; Tang, Y.; Hou, C.; He, Y.; Zhuang, J.; Zheng, Y.; Luo, S. CFD simulation and experimental verification of the spatial and temporal distributions of the downwash airflow of a quad-rotor agricultural UAV in hover. Comput. Electron. Agric. 2020, 172, 105343. [Google Scholar] [CrossRef]
Zhang, H.; Qi, L.; Wu, Y.; Musiu, E.M.; Cheng, Z.; Wang, P. Numerical simulation of airflow field from a six–rotor plant protection drone using lattice Boltzmann method. Biosyst. Eng. 2020, 197, 336–351. [Google Scholar] [CrossRef]
Nakai, T.; Van Der Molen, M.K.; Gash, J.H.C.; Kodama, Y. Correction of sonic anemometer angle of attack errors. Agric. For. Meteorol. 2006, 136, 19–30. [Google Scholar] [CrossRef]
Nakai, T.; Shimoyama, K. Ultrasonic anemometer angle of attack errors under turbulent conditions. Agric. For. Meteorol. 2012, 162, 14–26. [Google Scholar] [CrossRef]
Liu, C.; Du, L.; Zhao, Z. A directional cylindrical anemometer with four sets of differential pressure sensors. Rev. Sci. Instrum. 2016, 87, 035105. [Google Scholar] [CrossRef] [PubMed]
Yichen, P.; Zhenyu, L.; Zhen, F.; Lidong, D.; Zhan, Z. An angle-of-attack correction method for a 2-D miniature differential pressure solid-state anemometer. Measurement 2021, 171, 108833. [Google Scholar] [CrossRef]
Horlock, J.H. Actuator Disk Theory-Discontinuities in Thermo-Fluid Dynamics; McGraw-Hill, Inc.: New York, NY, USA, 1978. [Google Scholar]
Kim, Y.H.; Park, S.O. Navier-Stokes simulation of unsteady rotor-airframe interaction with momentum source method. Int. J. Aeronaut. Space Sci. 2009, 10, 125–133. [Google Scholar] [CrossRef]
Spalart, P.; Allmaras, S. A one-equation turbulence model for aerodynamic flows. In Proceedings of the 30th Aerospace Sciences Meeting and Exhibit, Reno, NV, USA, 6–9 January 1992; p. 439. [Google Scholar] [CrossRef]
Ashford, G.; Powell, K. Adaptive unstructured triangular mesh generation and flow solvers for the Navier-Stokes equations at high Reynolds number. In Proceedings of the 12th Computational Fluid Dynamics Conference, San Diego, CA, USA, 19–22 June 1995; p. 1724. [Google Scholar] [CrossRef]
Pan, Y.; Zhao, Z.; Zhao, R.; Fang, Z.; Wu, H.; Niu, X.; Du, L. High accuracy and miniature 2-D wind sensor for boundary layer meteorological observation. Sensors 2019, 19, 1194. [Google Scholar] [CrossRef] [PubMed]
Shao, W.Y.; Zhang, Y.P.; Zhu, D.Z.; Zhang, T.Q. Drag force on a free surface-piercing yawed circular cylinder in steady flow. J. Fluids Struct. 2013, 43, 145–163. [Google Scholar] [CrossRef]

Figure 1. The working status of the six rotors.

Figure 2. MRF model diagram.

Figure 3. The internal structure of the differential pressure anemometer.

Figure 4. Schematic diagram of

P_{D 1}

and

P_{D 2}

Figure 4. Schematic diagram of

P_{D 1}

and

P_{D 2}

Figure 5. Display of AOA and

α

Figure 5. Display of AOA and

α

Figure 6. Grid-independence test.

Figure 7. Computing domain and the encryption area.

Figure 8. The mesh diagram of the UAV.

Figure 9. The mesh diagram of the dynamic area of the UAV rotors.

Figure 10. The UAV velocity flow field diagram at ascent speed of 0 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s.

Figure 11. The UAV velocity flow field diagram at ascent speed of 3 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s.

Figure 12. The UAV velocity flow field diagram at ascent speed of 5 m/s. when the crosswind speed is: (a) 0 m/s, (b) 5 m/s, (c) 10 m/s, (d) 13 m/s.

Figure 13. Curves between simulation crosswind speed and standard crosswind speed under different ascent speeds.

Figure 14. Wind speed measurement results between the meteorological tower and airborne anemometer.

Figure 15. Wind speed measurement errors of the airborne anemometer.

Table 1. Simulation crosswind speed table under different standard crosswind speeds and ascent speeds.

Standard Crosswind Speed (m/s)	Simulation Crosswind Speed (m/s) at Ascent Speed of 0 m/s	Simulation Crosswind Speed (m/s) at Ascent Speed of 3 m/s	Simulation Crosswind Speed (m/s) at Ascent Speed of 5 m/s
0	0.168	0.015	0.009
3	3.344	3.333	3.222
5	5.928	5.758	5.419
7	7.612	7.552	7.438
10	10.970	10.670	10.661
13	13.807	13.794	13.682
15	16.148	15.989	15.871
17	18.282	18.099	17.979
20	21.483	21.263	21.142

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

© 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Liu, J.; Zhao, Z.; Fang, Z.; Li, Y.; Du, L. Correction of Error of Airborne Anemometers Caused by Self-Excited Air Turbulence. Sensors 2023, 23, 4288. https://doi.org/10.3390/s23094288

AMA Style

Liu J, Zhao Z, Fang Z, Li Y, Du L. Correction of Error of Airborne Anemometers Caused by Self-Excited Air Turbulence. Sensors. 2023; 23(9):4288. https://doi.org/10.3390/s23094288

Chicago/Turabian Style

Liu, Jianqiang, Zhan Zhao, Zhen Fang, Yong Li, and Lidong Du. 2023. "Correction of Error of Airborne Anemometers Caused by Self-Excited Air Turbulence" Sensors 23, no. 9: 4288. https://doi.org/10.3390/s23094288

APA Style

Liu, J., Zhao, Z., Fang, Z., Li, Y., & Du, L. (2023). Correction of Error of Airborne Anemometers Caused by Self-Excited Air Turbulence. Sensors, 23(9), 4288. https://doi.org/10.3390/s23094288

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

Correction of Error of Airborne Anemometers Caused by Self-Excited Air Turbulence

Abstract

1. Introduction

2. Methods

2.1. Basic Control Equation

2.2. Calculation Method

2.3. Turbulence Model

2.4. Correction Model of the Angle of Attack

3. Simulation and Modeling

3.1. Mesh and Boundary Conditions

3.2. Simulation and Results

3.3. Modeling

4. Test and Results

5. Conclusions

Author Contributions

Funding

Institutional Review Board Statement

Informed Consent Statement

Data Availability Statement

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI