Low-altitude atmospheric turbulence sounding on the basis of unmanned aerial vehicle

Based on the theory of turbulence, equations are derived for estimations of turbulent fluctuations of the longitudinal and transverse components of the wind velocity during ideal hovering of a quadcopter in a turbulent atmosphere. We present the results of experiments which were carried out on the territory of the Geophysical Observatory of Institute of Monitoring of Climatic and Ecological Systems, Siberian Branch, Russian Academy of Sciences, located in Tomsk Akademgorodok on the territory with complex orography, in a parkland zone with buildings of research institutes and motorways. Time series of turbulent fluctuations of the longitudinal and transverse components of wind velocity fluctuations were received with the use of an automated weather station, and time series of estimates of these components, from data of a DJI Phantom 4 Pro quadcopter during hovering. According to the automated weather station data, anisotropy was observed in one experiment during measurements in the atmosphere, but this phenomenon was not observed in the other experiment: the fluctuation spectra of all components of wind speed fluctuations coincide. The spectra of fluctuations of the longitudinal and transverse wind velocity components based on the automated weather station data and UAV telemetry are presented. The fluctuation spectra of these components for the automated weather station data and quadcopter generally coincide. The behavior of the spectra coincides with the spectrum which corresponds to Kolmogorov–Obukhov “–5/3” law within the inertial range. The turbulent spectra of the wind velocity fluctuations obtained with the use of the automatic weather station and with the unmanned aerial vehicle differ in the high-frequency spectral region.


Introduction
Low-altitude sounding of atmospheric turbulence is of great importance for many scientific and applied problems. Traffic control of unmanned aerial vehicles (UAV) of a "smart city", study of changes in the urban climate, monitoring of air pollution in urban environment and other regions with complex orography are among problems which could be solved in the presence of reliable high-resolution data on the state of atmospheric turbulence. The use of UAV, both fixed-wing and rotary-wing, is one of the main directions of development of the low-altitude sounding of turbulent flows in the atmosphere. The paper [1] presents the historic review and describes promises of a system of atmospheric observations. In addition, this paper discusses the widely used specialized observation systems, such as radars, lidars, and research aircraft, and presents some examples of their application for weather forecasting and prediction of climate.
The analysis given in [2][3][4] showed multi-rotor UAS to be more promising instrument for remote low-altitude sounding than observation balloons with sensors, since the former are cheap, easy to operate, longeval in a wide range of atmospheric conditions, and reusable. The advantages of multirotor UAS include vertical take-off and landing, as well as the capability of hovering at a fixed point. These advantages make this UAS type an ideal platform for measuring vertical profiles or The technique for measuring an atmospheric turbulence and temperature at a space point on the basis of integration of fast-response multihole pressure probe and a thermocouple fixed at a multirotor UAV-based platform is described in [20]. The technique allows accurate measurement of time series of the three components of the velocity vector and temperature at any space point where the UAV can fly within the limits of its autonomy.
A method of measuring wind speed using the data logged by the UAV is proposed in [21]. Experiments were conducted at a test site with laminar wind conditions with the UAV hovering next to a static 2D ultrasonic anemometer. During the flight tests, the wind velocity varied within 0-5 m/s. Horizontal wind estimates provide good results with root mean square error (RMSE) values between 0.26-0.29 m/s for wind speed, as well as between 4.1-4.9 for wind direction.
In [22], the air temperature was measured by attaching a temperature sensor 40 cm above the rotary-wing of the UAV and the UAV-based wind velocity was estimated using attitude data (pitch, roll, and yaw angles) recorded using the UAV's inertial measurement unit. The wind velocity was close to that observed by the automatic weather system within an error of approximately 10%.
The aim of this work is to carry out the comparative analysis of the spectra of fluctuations of the longitudinal and transverse components of wind velocity, measured with an AMK-03 automated weather station (AMK-03 is a 3D acoustic anemometer) [23,24]and a DJI Phantom 4 Pro quadcopter in hovering mode. The experiments were carried out in an urban area with a complex orography on the basis of the Geophysical Observatory of the IMCES SB RAS. The main equations of the turbulence theory are given in the next section, which are required for the correct statement of the problem of low-altitude sounding of atmospheric turbulence using quadcopter without wind sensor. The third and fourth section the equations of estimates of the wind velocity components and estimates of longitudinal and transverse components are provided for quadcopter without wind sensor. The general information about the experiment and the measurement results obtained with the DJI Phantom 4 Pro quadcopter and the AMK-03 automatic meteorological complex are presented in the 5th and 6th sections. The fluctuation spectra of the longitudinal and transverse components in a turbulent atmosphere using the AMK-03 data and quadcopter telemetry are compared.

Model of atmospheric turbulence
Airflow in the atmosphere is accompanied by random temporal and spatial fluctuations. Therefore, the turbulent air flow velocity field is represented as the sum of the mean and fluctuation components [25,26]: where ⟨ ( , )⟩ is the average wind velocity; ′ ( , ) are the wind velocity fluctuations; ⟨. . . ⟩ is the operator of statistical averaging.
In the coordinates commonly accepted in meteorology [26], the X axis is eastward directed (E), the Y axis is northward directed (N), and the Z axis is upward directed and is normal to the actual horizon surface. When mathematically describing the wind field fluctuations, the second-rank correlation tensor is introduced, where and enumerate the orthogonal components of the vector of velocity field fluctuations ′ ( , ) along the X, Y, and Z axes, i.e., , ∈ {X,Y,Z}. For isotropic turbulence, the correlation tensor is expressed via the longitudinal ( ) and transverse ( ) components: Where = � � = , = | |. Summation in Eqs (4) and (5) is carried over double subscripts.
As is well known, horizontal air mass transfer often prevails over vertical motion in the atmosphere; therefore, the behavior of the correlation tensor is studied with respect to the mean horizontal wind velocity direction � || ( , )�. In this case, the vector n in Eqs. (3)-(5) has the form where ⟨ ⟩ and ⟨ ⟩ are the mean wind velocity along the X and Y axes; || = �⟨ ⟩ 2 + ⟨ ⟩ 2 is the mean speed of the horizontal turbulent flow. Substitution of Eq. (6) in longitudinal and transverse correlation tensor components (4) and (5) results in the equations Here, is the longitudinal component of the wind velocity fluctuations and is its transverse components. Equations (9) and (10) where ( ) = 2 ( ) is the spectral function; ( ) is the spectral density of the average kinetic energy of a unit air mass [25]; are the coordinates of the vector ; is the Kronecker delta. For isotropic fluctuations, the Kolmogorov-Obukhov "-5/3" law [25]: is fulfilled in the inertial range of wavenumbers, where ε is the dissipation rate of turbulent kinetic energy.
The correlation between the spatiotemporal and spatial parameters of the wind field can be derived on the basis of the Taylor hypothesis about "frozen" turbulent fluctuations [25], the essence of which reduces to the fact that the spatial pattern ′ ( , ) moves in time with the mean wind velocity ⟨ ⟩, i.e., The equation which connects the frequency spectrum Φ( ) of a random field with its spatial spectrum ( ) was derived in [25]. Using this equation, one can write the equation for the frequency spectrum in the case of the Kolmogorov-Obukhov "-5/3" law: Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 25 September 2020 doi:10.20944/preprints202009.0595.v1 Equation (15) shows that the frequency spectrum of a random wind velocity damps proportional to the frequency in the power -5/3.

Estimates of the wind velocity components
Dynamic equations for the center of gravity of a quadcopter can be written in the inertial coordinates associated with the Earth as [27][28][29][30] where (•) = sin(•), (•) = cos(•); is the roll angle; is the pitch angle; is the yaw angle; is the aerodynamic force generated by propellers; , , and are the drag force components along the x, y, and z axes; is the quadcopter mass; is the gravity acceleration. The components of the drag force along the x, y, and z axes, which arises during a quadcopter flight, have the form in the linear case and in the square-law case. In Eqs. (20) and (21), and are the drag coefficient along the x, y, and z axes; is the subscript for enumeration of the orthogonal components of vectors, i.e., ∈ { , , }; are the quadcopter speed component; are the wind velocity components in the atmosphere; is the air density; are the projections of the quadcopter area on the corresponding axes; sgn(•) is the sign function.
As is well known, the ideal hover is possible when = 0 and all forces which act on a quadcopter are compensated. The roll, pitch, and yaw angles in a turbulent atmosphere are sums of the average and fluctuation components: (19) can be transformed to the case of ideal hover by means of their linearization. In the small-angle approximation , ≪ π, at ′ ≪ and under the conditions ̈=̈=̈= 0 and = 0, the equations for estimates of the horizontal components of the wind field � = and � = take the forms in the linear case and in the square-law case. Equations

Estimates of longitudinal and transverse components
Let us define estimates the longitudinal and transverse components of wind field fluctuations by analogy with Eqs. (9) and (10) as where � ′ ( , ) is the estimate of the longitudinal component, and � ′ ( , ) , of the transverse component in the XoY plane. The estimate of the vector direction � in Eqs. (26) and (27) can be written as where ⟨ � ⟩ and � � � are the estimates of the components of the horizontal wind velocity along the x and y axes; � || = �⟨ � ⟩ 2 + � � � 2 is the estimate of the horizontal wind velocity.
It follows from Eqs. (26) and (27) that the estimates of the longitudinal and transverse components could be used for calculation of fluctuation spectra of the its components using quadcopter data. AMK-03 is the 3D acoustic anemometer that measures the longitudinal and transverse components, so Eqs. (9), (10), (26) and (27) allows to calculate the fluctuation spectra and compare them in the turbulent atmosphere.

General information about experiments
Experiments were carried out in the Geophysical Observatory of IMCES SB RAS on February 20 and August 13, 2020. The Observatory is situated in Tomsk Akademgorodok, on the territory with complex orography, in a parkland zone with buildings of research institutes and motorways. Figure 1 shows the quadcopter flight trajectory during the experiments. A DJI Phantom 4 Pro drone was launched on February 20, 2020, in close proximity to the Observatory building. After takeoff, it climbed to an altitude of 28 m and flew up to the AMK-03 automated weather station, which was mounted at a mast on the Observatory roof. The quadcopter axis x was approximately northward oriented. The experiment started at 10:21 UTC and ended at 10:34 UTC. The quadcopter held the altitude near AMK-03 for about 10 min and then returned to the starting point.
On August 13, 2020, the drone started near the 30-m meteorological mast situated between the Observatory and a forest. After takeoff, it climbed to an altitude of 27 m and flew up to AMK-03 automated weather station mounted at the mast. The quadcopter axis x was also approximately northward oriented. The experiment was carried out from 07:01 to 07:22 UTC. The quadcopter held the altitude near AMK-03 for about 18 min and then returned to the starting point.
AMK-03 is designed to measure and record wind speed and direction, air temperature and relative humidity, and atmospheric pressure. We used two AMK-03 stations. The station on the Observatory roof measured wind speed and direction with a frequency of 10 Hz, and the weather station at the 30-m meteorological mast recoded the time series with a frequency of 80 Hz. The quadcopter flight data (time series of the roll, pitch, and yaw angles) were recorded with a frequency of 10 Hz. When calculating the spectra of wind velocity fluctuations, the AMK-03 data were reduced to a frequency of 10 Hz by means of partial averaging of the time series. According to data of the Tomsk International Airport, which is located at a distance of ~10 km from the building of the Geophysical Observatory of IMCES SB RAS, the weather conditions were difficult for quadcopter flight during the experiment on February 20: south-southwestern wind at a speed of 6.0 m/s, air temperature of -1.5°С, air humidity of 100%, and snowing. During the experiment on August 13, the southeastern wind at a speed of 4 m/s, air temperature of 26°C, and air humidity of 70%, and no precipitation were recorded. Figure 2 shows variations in the quadcopter velocity components along the x, y, and z axes during hovering. The quadcopter velocity components were generally zero during the measurements. Sometimes, the forces which acted on the quadcopter exceeded the capabilities of the control system for a short time, and high-precision positioning was perturbed. After control recovery, the quadcopter began moving to its original position and stopped upon reaching it. Thus, the time periods when the spatial positioning of the quadcopter was perturbed can be ignored due to their insignificance, and the hover could be considered ideal from 10:23:00 to 10:33:00 UTC on February 20 and 07:01 to 07:22 UTC on August 13, 2020. Let us consider the behavior of the estimates of the longitudinal and transverse wind velocity components in a turbulent atmosphere in the altitude holding mode. Figure 3 shows the time series of these components (red curve) measured by AMK-03 automated weather station and their estimates (black curve) found on the basis of the DJI Phantom 4 Pro quadcopter data. It follows from Fig. 3 that the time series of � and and the time series of � and generally coincide, and the differences are observed in the region of high-frequency fluctuations.

Measurement results
For objective assessment of the spectra of fluctuations of the wind field components and spectra of their estimates in a turbulent atmosphere, we used averaging over the hover period. The spectra were calculated with the use of standard FFT algorithms and a smoothing procedure.  Kolmogorov-Obukhov "-5/3" law; σ 2 is the normalization coefficient Figure 4 shows relative fluctuation spectra of the components calculated on the basis of AMK-03 data for two experiments. It is seen from figure that anisotropy of the wind velocity fluctuations was observed on February 20, but on August 13 this phenomenon was not observed: the fluctuation spectra of all components ′ ′ and ′ coincide. In the inertial interval, the spectra of the ′ , ′ , and ′ components behave similarly to the spectrum which corresponds to the Kolmogorov-Obukhov "-5/3" law. Figure 5 shows the results of comparison between the fluctuation spectra of the estimates of transverse and longitudinal components. One can see that the fluctuation spectra of the estimates of transverse and longitudinal components w � r ′ and w � t ′ coincide in the altitude holding mode. In the inertial range, the spectra of the w � r ′ and w � t ′ components coincide with the spectrum which corresponds to the Kolmogorov-Obukhov "-5/3" law.  Figure 5 shows the results of comparison between the fluctuation spectra of the estimates of transverse and longitudinal components. One can see that the fluctuation spectra of the estimates of transverse and longitudinal components w � r ′ and w � t ′ coincide in the altitude holding mode. In the inertial range, the spectra of the w � r ′ and w � t ′ components coincide with the spectrum which corresponds to the Kolmogorov-Obukhov "-5/3" law.  7 show the results of comparison between the fluctuation spectra of the transverse and longitudinal components and of their estimates. The fluctuation spectra of ′ and � ′ and ′ and � ′ generally coincide, and significant differences are seen in the high-frequency spectral region. Here, the frequency spectra of turbulent fluctuations of w � r ′ and w � t ′ damp with the frequency faster that according to the Kolmogorov-Obukhov "-5/3" law.

Discussion
Based on the study of the spectra of turbulent fluctuations of the longitudinal and transverse wind velocity components, measured with the automated weather station AMK-03 and DJI Phantom 4 Pro quadcopter, the following conclusions can be drawn.
We used DJI Phantom 4 Pro quadcopter in our experiments, for which we have shown that the fluctuations spectra of the transverse and longitudinal components measured with the 3D acoustic anemometer and an UAV generally coincide in a wide spectral range. In the inertial range, the behavior of the spectra coincides with the spectrum, which corresponds to the Kolmogorov-Obukhov "-5/3" law. The fluctuations spectra of these components measured with the 3D acoustic anemometer and an UAV coincide regardless of the state of the turbulent atmosphere: for isotropic and anisotropic turbulence. This means that the state of a turbulent atmosphere can be controlled with the use of a rotary-wing UAV instead of contact sensors, e.g., 3D acoustic anemometer.
The low-altitude sounding of a turbulent atmosphere with the use of a fixed-wing UAV [8][9][10][11][12][13][14][15] allow retrieving information from spatially averaged data. However, for certain scientific and applied problems, for example, traffic control of UAV of "smart city", detailed information is required on local processes which run in a turbulent atmosphere. For such problems, rotary-wing UAV capable of hovering at different altitudes for a long time have an undoubted advantage.
Differences in the spectra of turbulent fluctuations obtained with AMK-03 and DJI Phantom 4 Pro are observed in the high-frequency spectral region, which corresponds to a small-scale low-energy turbulence. One possible reason for the difference in the spectra is the presence of an induced turbulence: in the hovering mode, the quadcopter creates high-energy turbulent vortices around it, which destroy the small-scale structure of atmospheric turbulence. As a result, the measured spectrum of turbulent fluctuations of the wind field in the high-frequency spectral region may differ from the spectrum which exists in the atmosphere.
Another possible reason could be large mass and size of the quadcopter used. A small-scale turbulence has low kinetic energy; hence, the quadcopter is insensitive to weak high-frequency impacts due to high inertia and large size. It is possible that a lighter and smaller quadcopter will be more sensitive to small-scale turbulent fluctuations and will produce weaker local turbulent distortions, which is to provide a correct pattern in the high-frequency spectral region. Therefore, the future direction of our research is to study the effect of atmospheric turbulence on the operation of a quadcopter with lower mass and smaller size than the DJI Phantom 4 Pro.