Meteorology-Aware Path Planning for the UAV Based on the Improved Intelligent Water Drops Algorithm

A wide range of applications of the unmanned aerial vehicle (UAV) have been observed in the past few years, and path planning is one of the most critical issues that require to be resolved. UAVs are still prone to meteorological impediments such as thunderstorms, ice accumulation, and severe convective weather for the safety of flights. This paper proposes a meteorology-aware path planning method based on the improved intelligent water drops (IIWD) algorithm. The algorithm consists of both static and dynamic path planning. In the static path planning phase, ice accumulation and the Richardson number are considered when determining the trajectory with low risk. In the dynamic path planning phase, the latest forecast products are adopted to modify the planned path, and the virtual potential force is applied to adjust the flight direction of the UAV when encountered with severe convective weather. The validation and efficiency of the proposed algorithm are verified via simulations in comparison with the ant colony optimization, genetic algorithm, and Q-learning algorithm, where the quality of our algorithm’s performance in terms of flight time and risk degree is determined. Meanwhile, the risk degree of the UAV flight path at different altitudes is analyzed. The simulation results show that the average flight speed decreases, and the risk degree increases along with the descending of the flight altitude, respectively, which is found to be in consensus with the theory of meteorology.


I. INTRODUCTION
With the development of satellite communication, microelectronics, and composite materials, we have observed the widespread applications of unmanned aerial vehicles (UAVs) in various fields in the past few years [1]. As human operators are not required to operate on board, UAVs are more suitable for hazardous or mechanized tasks. For example, in military terms, UAVs are often used to perform reconnaissance, fire strikes, decoys, and other tasks [2]. In the field of agriculture, UAVs can significantly improve the efficiency of farmland irrigation [3]. In the scenario of environmental detection, UAVs are capable of obtaining high-resolution data with a low flight altitude. Moreover, UAVs have great potential to become the carrier of the logistics distribution in the future, thus reducing the cost of manual distribution.
The associate editor coordinating the review of this manuscript and approving it for publication was Arturo Conde .
According to the definition given by the Department of Defense (DoD) of the United States in the roadmap of the unmanned system [4], the UAV system is composed of necessary equipment, network, and UAV operators. In some cases, the UAV system also includes launch devices. It can be seen that the UAV plays the role of mission executor in the entire unmanned system [5]. Before carrying out the mission, it is necessary to plan the UAV path according to the mission requirements, which is considered one of the most important parts of UAV mission execution. UAV path planning aims to find the optimal path between the starting point and the destination point. The path is required to achieve the shortest length, fastest time and solve the problems of collision-free and obstacle avoidance.
Generally, there are four steps in UAV path planning. Firstly, the target flight area is divided into dispersed grid points, i.e., the search graph. After designing the search graph, the information related to the grid points can be determined to transform the path planning problem into a graph search problem. Secondly, the grid points' characteristics are updated in the search graph corresponding to the decision conditions of the path search at the beginning of each time slot. Thirdly, the optimal or the near-optimal path is calculated according to the path planning algorithm. Finally, the trajectory smoothing algorithm is applied for smoothening the obtained path.
Meteorological factors play a vital role in the navigation of the UAV. An accident survey of 2062 aviation accidents from 1999 to 2010 showed that weather-related accidents occupy 7% among all cases [6], which means that meteorological factors have a significant impact on the safe flight of aircraft. Moreover, reasonable use of wind can significantly reduce flight time and achieve the economic benefits of reducing fuel consumption. With the enrichment of various forecast products such as satellite cloud images, radar echoes, and ground/high-altitude observation data, the relationship between flight activities and meteorological conditions is changing from the meteorological conditions determine whether the flight mission is able to carry out to how to fly under complex meteorological conditions [7]. Although there are a few reports on the impact of meteorological conditions on the UAV flight, UAV, like manned aircraft, can be affected by severe convective weather, thus causing safety risks. However, in previous studies, most researchers have not considered the impact of meteorological factors on UAV path planning.
In this paper, the UAV path planning model based on meteorological factors is constructed, and an improved intelligent water drops (IIWD) algorithm is proposed to solve this NP-hard problem. The intelligent water drops (IWD) algorithm is inspired by the behavior of rivers and the interaction of water drops and soils in the river bed. This algorithm has been applied to many fields of natural science and engineering science, demonstrating great advantages and potential [8], [9].
To the best of the authors' knowledge, the basic IWD algorithm has not been adopted for meteorology-aware path planning. This fact motivated us to develop this algorithm to find a near-optimal solution to the problem. To solve this problem, an improved intelligent water drops algorithm for UAV path planning is proposed. This algorithm consists of both static and dynamic path planning. In the static path planning phase, ice accumulation and the Richardson number are taken into consideration when constructing the cost function. In the dynamic path planning phase, the latest forecast products are adopted to modify the planned path, and the virtual potential force is applied to adjust the flight direction of the UAV when encountered with severe convective weather. The main contributions of this paper are as follows: 1) Presenting improved intelligent water drops (IIWD) algorithm as a new swarm-based nature-inspired optimization method, which has not been adopted in the literature so far to tackle the meteorology-aware path planning.
2) Taking into account of meteorological factors such  as wind speed, temperature, and relative humidity,  the two indices, namely, ice accumulation and Richardson number, are selected to define the risk degree of the  UAV path planning. The cost function is formulated by  considering the path risk factors, wind direction, and  wind speed jointly. 3) The path planning framework of the UAV in the dynamic environment is constructed. Firstly, when the UAV receives new forecast products, it replans the path. Secondly, when encountered with sudden dangerous weather, the virtual potential force algorithm is adopted to dynamically adjust the flight direction of the UAV, thus ensuring flight safety.
The rest of the paper is organized as follows. The literature review of the state-of-the-art is presented in Section II. The problem description and modeling are presented in Section III. The improved intelligent water drops algorithm is described in Section IV. The effectiveness of the proposed algorithm is verified via simulation in Section V. Finally, the conclusion of this paper is drawn in Section VI.

II. LITERATURE REVIEW
The research on UAV path planning has achieved fruitful results. Scholars have studied this problem from different perspectives using various methods. This section reviews the literature from the following two aspects: the path planning algorithms and the UAV path planning under meteorological conditions.

A. PATH PLANNING ALGORITHMS
Path planning algorithms search for a feasible trajectory between the starting point and the destination point after determining the environmental conditions. The existing path planning algorithms are primarily divided into the following four categories: the search algorithm based on (1) geometric models, (2) probability sampling, (3) artificial potential field, and (4) intelligent algorithm [10]. The search algorithms based on geometric models are classical path search algorithms, which belong to discrete optimal planning. Chen et al. [11] compared the A* algorithm with the genetic algorithm (GA) and ant colony optimization (ACO) algorithm and considered reducing the fuel consumption of the UAV in path planning. Song and Hu [12] proposed an improved A* algorithm, where they constructed a weighted graph according to the cost incurred in the path. Simulations show that the planned path distance can be effectively reduced. Koenig and Likhachev [13] investigated the D* algorithm, which is a combination of the dynamic A* algorithm and the incomplete replanning and combines global planning with local information to avoid dynamic obstacles. The search algorithms based on the probability sampling show certain theoretical advantages according to the probability integrity or the asymptotic optimality. However, this method requires to know the global environment information and considers the environment as a set sampling, and finally VOLUME 9, 2021 randomly searches to find the path [14]. Huang et al. [15] presented the probabilistic route map method for the automatic path search, and this method has been applied to the deployment of submarine optical cables. When constructing the route map, since the randomness is included in the route search, the planned route is found to be near-optimal. The method based on the artificial potential field shows good performance in online path planning since the algorithm only accepts the local information of the UAV, and the calculation cost is low. However, it may be difficult to generate a feasible trajectory in a complex environment. Guang et al. [16] combined the artificial potential field with the virtual force and applied it to the two-dimensional unknown environment for path planning. The simulation results show that the algorithm can avoid local optimum. Intelligent algorithms can effectively plan a suitable track in a complex dynamic environment. Nevertheless, they are still encountered with the same challenging issues, such as calculation time and stability. Additionally, intelligent algorithms are prone to be trapped into local optimum. Zhang et al. [17] proposed the particle swarm optimization algorithm to plan the UAV path in the dynamic environment and verified the effectiveness of the algorithm in Monte Carlo tests. Zhang et al. [18] introduced the artificial neural network to plan the flight route of the UAV so that the UAV can handle emergencies in advance, and they can also use the path planning algorithm to find a near-optimal path. Challita et al. [19] proposed a deep learning reinforcement learning approach for interference management for cellular-connected UAVs. The proposed algorithm could achieve a balance between energy efficiency and interference caused on the ground network along the UAV's path.

B. UAV PATH PLANNING UNDER THE INFLUENCE OF METEOROLOGICAL FACTORS
The influence of meteorological factors on UAV flight is considered in some studies. Rubio et al. [20] proposed the basic framework of the oceanic search mission by UAVs, considering the impact of the aircraft ice accumulation and the wind on the UAV flight, and solved the problem by using the evolutionary algorithm (EA). Wirth et al. [21] applied the dynamic programming (DP) method to solve the problem of the long endurance solar-powered UAV path planning and comprehensively considered the influence of nine weather factors, such as temperature and relative humidity. However, when considering various meteorological factors, the impact of dangerous weather phenomena on UAV flight safety is not considered. Lee et al. [22] also considered the path planning problem of a solar-powered longendurance unmanned aerial vehicle. They obtained the maximum energy by adjusting the flight altitude of the UAV. Oettershagen et al. [23] modeled the multi-objective longendurance UAV path planning problem with the influence of terrain and adopted the improved A* algorithm to solve the problem. The results demonstrate that a reasonable flight takeoff time could effectively reduce the whole flight time from 106 hours to 52 hours. Hu et al. [24] studied the influence of the flight altitude on the aircraft fuel efficiency and achieved the purpose of fuel-saving by controlling the flight altitude. Coombes et al. [25] and Kim et al. [26] considered the influence of the wind on the flight of rotorcraft and solved the path planning of the UAV by using a geometric method. Sun et al. [27] investigated the resource allocation algorithm designed for multicarrier solar-powered UAV communication systems. The results show that the UAV first climbs up to a high altitude to harvest solar energy and then descents to a lower altitude to reduce the path loss of communication links could achieve a good performance.
In previous studies, few specific dangerous weather phenomena have been considered for the UAV flight's safety. This paper designs a UAV flight model and comprehensively considers various meteorological factors. Focusing on aircraft ice accumulation, Richardson number, and UAV performance constraints, the improved intelligent water drops algorithm is proposed to solve the UAV path planning problem.

III. PROBLEM DESCRIPTION AND MODELING A. PROBLEM DESCRIPTION
The general description of the UAV path planning problem is as follows. A feasible trajectory is calculated for the UAV moving from the starting point to the destination point in a given situation restricted by environmental factors and risk areas. In order to guarantee the safety of the UAV flight, the trajectory is not only limited to obstacle avoidance, but it also requires to fulfill the performance constraints of the UAV. In addition, the trajectory is optimized under certain evaluation indexes in order to ensure the minimum cost. The static path planning is performed based on the task demands and meteorological forecast products of the target area before taking off. However, with the changes in atmospheric circulation, which result in the flight area's changes in meteorological conditions, dangerous weather phenomena may be formed along the flight path. Hence, it is necessary to plan the flight path according to the dynamically changing environmental conditions. The overall framework is illustrated in Fig. 1.
For the UAV flight's safety, the flight altitude remains consistent after taking off, which means the UAV maintains the flight height after ascending to the operating altitude. The UAV path planning can be mapped to a two-dimensional space, except for the taking off and landing. In this study, a certain altitude of the UAV flight is mapped to a twodimensional grid space z, m×n grid points are available in the space, and the planned path is a collection of the grid points, which is expressed by formula (1).
where p 1 = source and p k = destination, respectively represent the starting point and the destination point, and the coordinate of each grid p i is (x i , y i ). The flight time t i of the UAV moving from position p i to p i+1 can be calculated by the flight speed v i , i.e., is determined by the length of the flight path and the flight speed. However, the UAV flight speed can be affected by factors such as the wind speed. Therefore, the shortest path does not necessarily guarantee the shortest flight time. The path planning problem can be expressed by formula (2).
In this study, the objective of the UAV path planning is to find a path with the minimum flight time.

B. UAV PERFORMANCE CONSTRAINTS
UAVs' mechanical properties determine that performance constraints need to be satisfied in path planning. The UAV performance constraints primarily include the minimum step size constraint, maximum steering angle constraint, and maximum flight distance constraint.

1) MINIMUM STEP SIZE CONSTRAINT
The maneuvering ability of the aircraft determines that the UAV must fly a short straight line before changing to the next flight state, which means the distance between p i and p i+1 should be greater than a specific value, which is expressed by formula (3).
where d(p i , p i+1 ) represents the distance between p i and p i+1 , and p min represents the minimum step size.

2) MAXIMUM STEERING ANGLE CONSTRAINT
Due to the effect of inertia, the heading of the UAV is constrained by the time and the steering angle. That is, the UAV has the maximum steering angle. The steering angle of the UAV moving from p i to p i+1 is expressed by formula (4).
where ϕ(p i , p i+1 ) represents the steering angle of UAV moving from p i to p i+1 , p i = (x i , y i ) is the position of p i , and the direction vector is . The steering angle constraint of the UAV is expressed in formula (5), where ϕ(max) represents the maximum steering angle.

3) MAXIMUM FLIGHT DISTANCE CONSTRAINT
The fuel that the UAV carries is limited to the size of the UAV. In order to return safely, the UAV must fly with the fuel allowed, which is expressed by formula (6).
where L max represents the limitation of the flight distance.

C. COMPLEX ENVIRONMENTAL IMPACT
The meteorological environment is complex and changeable, and many factors affect the flight safety of the UAV. The meteorological factors that affect the UAV are ice accumulation, vertical wind shear, thunderstorms, intense rainfall, and other severe convective weather.

1) AIRCRAFT ICE ACCUMULATION
Aircraft ice accumulation is primarily caused by supercooled water droplets hitting the aircraft, and it can also be formed by the condensation of water vapor on the surface of the aircraft. It has been proved that when the temperature is −14 ∼ 0 • C, and the relative humidity is above 50%, the aircraft is found most susceptible to ice accumulation when the fuselage hits large supercooled water droplets. The World Meteorological Organization developed the Ic index to predict the ice accumulation conditions based on the relative humidity, and the relative temperature [28], which is expressed by formula (7).
where RH is the relative humidity and T is the relative temperature. The first half of the formula represents the growth process of the content and the size of the supercooled water droplets in the cloud. The second half represents the growth rate of the water droplets. When RH increases from 50 to 100, the Ic index also increases from 0 to 100. The value of the second half goes up to 1 when T is −7 • C, which means the ice accumulation mostly occurs at this temperature. When either part of the two parts is negative, it is considered that ice accumulation does not occur.

2) RICHARDSON NUMBER
In the atmosphere's static equilibrium state, some air masses are disturbed by dynamic factors or thermal factors and produce upward or downward vertical movement. Whether VOLUME 9, 2021 the vertical movement deviates from its equilibrium position can continue to develop depends on the vertical distribution of atmospheric temperature and humidity. The Richardson number represents the atmospheric static stability ratio to the vertical wind shear [29]. Suppose an air cluster rises from the height of z − l to the height of z, the energy of the average unit volume of air in the unit time resists the gravitational field is expressed by formula (8).
In the turbulent motion, the turbulent friction consumption is converted from the Reynolds stress to the pulsating kinetic energy. The magnitude of this part of the energy is expressed by formula (9).
When W 1 < W 2 , the pulsating kinetic energy converted from the Reynolds stress is found to be greater than the pulsating kinetic energy consumed by the energy of resisting the gravitational field when the air cluster moves and the turbulent pulsation continues to grow. On the other hand, when W 1 > W 2 , the turbulent pulsation weakens.
Richardson number is defined as a dimensionless parameter that judges whether the turbulence can be enhanced or not, which is defined in formula (10).
The Richardson number represents the relative magnitude between the pulsating kinetic energy consumed by overcoming the gravitational field and the pulsating kinetic energy converted from the Reynolds stress. It has been observed that when R(p i ) < 0.5, it indicates that the turbulence increases. In particular, when R(p i ) < −2, the rain accumulation is prone to occur. Meanwhile, when R(p i ) < −1, clouds are prone to thunderstorms, and the systemic convection is prone to occur when −1 ≤ R(p i ) ≤ −1/4 [7].
Based on the above description, we define the degree of the Richardson number DoR in formula (11).

3) WIND SPEED
In the case of a fixed flight speed, the only way to minimize the flight time between two points is to fully use the influence of the wind to increase the ground speed.
Assume that the flight speed of the UAV is V , the flight altitude of the UAV is gpm, the angle between the UAV flight direction and the due east direction is θ, the radial wind velocity component at p i is u(p i ), and the zonal velocity component is v(p i ). Then, the influence of the wind on the UAV ground speed can be expressed in formula (12).
The time required for the UAV to fly from p i to p i+1 is: Compared with manned aircraft, UAV is prone to stall due to vertical wind shear when adjusting flight altitude. For the UAV flight's safety, once the level of flight altitude is selected, the flight altitude cannot be changed unless exceptional circumstances occur. Therefore, the difference between different altitudes should be considered when planning the path. The meteorological data we applied in this paper is the high-resolution numerical prediction product supported by the European Centre for Medium-Range Weather Forecasts (ECMWF). This data will be received twice a day, containing the meteorological data every three hours in the next seven days. When the missing data area is large, these points will not be chosen. Otherwise, an interpolation algorithm can be used to fill the data of the area. Meteorological factors and symbols used in this study are shown in Table 1. Therefore, once the level flight height is determined, the entire UAV path planning problem can be transformed into an optimization problem expressed in formula (14).
subject to the following constraints: where C1 represents that the distance between p i and p i+1 should be greater than the minimum step size p min ; C2 represents that the steering angle should be smaller than the maximum value; C3 represents that the flight distance should be smaller than the flight distance limitation; C4 and C5 represent are the constraints of meteorological factor, which guarantee the safety of UAV flight. The above expression is the mathematical description of the meteorology-aware path planning problem. The improved intelligent water drops algorithm is applied to solve the problem in the next section.

IV. IMPROVED INTELLIGENT WATER DROPS ALGORITHM A. STATIC UAV PATH PLANNING BASED ON THE IMPROVED INTELLIGENT WATER DROPS ALGORITHM
The design idea of the intelligent water drops (IWD) algorithm originates from the flow of water drops in nature. It is a swarm-based intelligent computing method that imitates the process of water drops interacting with the surrounding environment to form a river channel, which is proposed by Hamed Shah-Hosseini in 2007 [30].
In nature, due to the influence of the environment, the flow path of water drops encounters different types of obstacles. If there are no other resistance forces, water drops reach the destination in a straight line, which is considered the ideal path from the starting point to the destination point. The motion of water drops follows the below-listed rules: 1. The water drops with higher velocity carry more soil than slower ones. 2. The water drops gain more velocity in the path with less soil. 3. The water drops tend to choose the path with less soil.
The intelligent water drops algorithm is primarily divided into the following four steps: initial phase, path selection phase, soil update phase, and the optimal path selection phase.

1) INITIAL PHASE
The number of IWD N IWD , the initial speed of the IWD v ini , the initial soil of the water drop soil IWD , path soil soil(p i , p j ) from grid point p i to p j , and the maximum number of iterations I max are determined for the follow-up process.

2) PATH SELECTION PHASE
We imitate each of the water drops as the UAV to solve the UAV path planning problem. Every water drop selected on the grid point from the starting point to the destination point forms different tracks and gets optimized in continuous iterations. When the water drop is located at p i , due to the constraints such as the minimum step size and the maximum steering angle, the next selectable set of grid point is δ(i). The IWD algorithm selects the next grid point p j with a certain probability. The probability has an impact on the path's soil between p i and all elements in δ(i), which is expressed in formula (15).
where ε is a small non-negative constant, which prevents the function value from being less than 0. h(p i , p j ) represents the difficulty to move from p i to p j , in the conventional IWD algorithm, h(p i , p j ) = d(p i , p j ), which represents the distance between p i to p j . We consider the meteorological impact in path choosing, and if the probability of ice accumulation or the Richardson number is high, the water drop will have a low probability of choosing this path. We define h(p i , p j ) in formula (17) and (18).
The influence of meteorological factors is considered in formula (17), and D(p j ) is defined as the risk degree of the grid point, which is expressed in formula (18). If the risk degree of the next point p j is high, h(p i , p j ) gets a smaller value compared to other optional points, which indicates that the probability of this path being selected is significantly low.

3) SOIL UPDATE PHASE
Assume that at moment t the water drop is located at p i , thereafter it moves to p j at t +1 as shown in Fig. 3. The change in the water drop velocity is determined by the amount of soil in the path, which is expressed in formula (19).
where, a v , b v , and c v are all non-negative constants. The velocity of the water drop is affected by the amount of soil in the path. Meanwhile, the amount of soil in the path changes due to the flow of water drops. When the water drop moves from p i to p j , the amount of soil in the path is updated by formula (20). (20) where ρ n is a local soil updated parameter, which is a constant between 0 and 1, and soil(p i , p j ) represents the amount of VOLUME 9, 2021 soil change, which is determined by the time that the water drop moves from p i to p j and it is expressed in formula (21).
where, a s , b s , and c s are all non-negative constants. p j ) i.e., expressed in formula (17) and v g (p i ) represents the ground speed, which is expressed in formula (12).
Simultaneously, the amount of soil carried by water drop is updated by formula (23).

4) OPTIMAL PATH SELECTION PHASE
For all water drops in one iteration, when they have completed the path from the starting point to the destination point, according to our optimization objective, we select the optimal path and update it to the global path soil, which is expressed in formula (24).
The optimal path is generated in each iteration. Suppose the optimal path generated in the next iteration is better than before. In that case, it is saved until the maximum number of iterations so that the near-optimal solution of the trajectory planning can be obtained.
The computational complexity of the IIWD is determined by the iterations of the algorithm, the number of the water drops, and the maximum step size of the water drop. Therefore, the computational complexity is O(N iter N IWD Step max ), since there is a linear relationship among them, the computational complexity is rewritten as O(n 3 ). The pseudocode of the IIWD algorithm is expressed in Algorithm 1.

B. DYNAMIC PATH PLANNING OF THE UAV UNDER SUDDEN DANGEROUS WEATHER CONDITIONS
The UAV can receive instructions from the ground station by the satellite link during the flight. Due to the timeliness of forecast products, dangerous and severe convective weather Algorithm 1 Improved Intelligent Water Drops Algorithm Input: Starting position, Destination position, Environmental parameters, UAV performance parameters. Output: Trajectory.
---------Initialization---------1: p 1 ← Starting position, p k ← Destination position 2: Initialize environmental parameters of each grid, including T , RH , u and v. 3: Initialize IIWD parameters, including N IWD , v ini , soil IWD , soil(p i , p j ) and I max . 4: while k ≤ I max do 5: while i ≤ Step max do --------Path Selection--------- 6: for each IWD do 7: current_position ← p i 8: calculate p i (j) 9: select next_position by p i (j) ------Local parameters update------- 10: update v IWD 11: update soil(p i , p j ) 12: update soil IWD 13: current_position ← next_position 14: Track IWD ← Track IWD + next_position 15: find the best Track best ← Track IWD ------Global parameters update------ 16: update global path of soil 17: update the best solution Track global 18: return Track global phenomena may occur on the initial statically planned path after the UAV takes off and sails for a period of time. The ground station can receive forecast products including temperature, humidity, pressure wind, and other meteorological factors at different altitudes of the target area every 3 hours. According to the UAV's current location and flight direction, the improved intelligent water drops algorithm is used to adjust the trajectory dynamically. As shown in Fig. 4, the solid black line indicates the UAV's initial planned path. The ground station receives new forecast products at time t i , replans the path, and sends control commands to the UAV by the satellite link at t i + t, and the dotted line shows the UAV path after the first replanning. Similarly, the UAV receives the adjustment instruction again at t i+1 + t, and the segment line shows the UAV path after the second replanning.
However, some severe convective weather phenomena develop rapidly and often disappear before the relevant weather factors, such as thunderstorms and downbursts, are captured. These weather phenomena are difficult to predict, but they pose a major threat to aircraft flight safety since the UAV's wings carry sensors that can detect changes in temperature, humidity, etc. However, due to the limited computing power of the UAV, when the temperature and the wind speed change significantly, it can only roughly determine the location of dangerous weather phenomena. Due to the tremendous destructive power of dangerous weather phenomena, an accident can be easily caused by the crash of the aircraft once it enters the dangerous zone, so it is necessary to use a simple dynamic path adjustment method with low computational cost to adjust the flight direction of the UAV. We use the virtual potential force method to quickly avoid dangerous weather phenomena and achieve online UAV path planning.
When sudden dangerous weather phenomenon occurs, assuming that the UAV is located at the grid point p i and the dangerous weather is located at τ j , and its influence radius is r j . The repulsive force of the UAV by the dangerous weather is defined in formula (25).
The impact of dangerous weather on the UAV is only when the distance between the two is less than the radius of influence of the dangerous weather. The relative position of the two determines the direction of the repulsion. For multiple dangerous weather phenomena, the combined force of the impact on the UAV is : Meanwhile, in order to reach the destination, the attractive force is determined by formula (27).
This resultant force determines the direction of the UAV's flight. After moving out of the range of the dangerous weather, the ground station replans the flight trajectory. Figure 5 shows an example of a UAV that encounters sudden dangerous weather phenomena, and the red line represents the attractive force between the UAV and the destination, and the blue line represents the virtual repulsive force between the UAV and the risk area, and the black line represents the resultant force that indicates the moving direction of the UAV. When the UAV moves to position A, a dangerous weather phenomenon is detected. Under the joint action of repulsive force and attractive force, the UAV changes its moving direction to avoid the obstacle. When the UAV moves to position B, since the angle between the virtual repulsion force and gravity is less than 90, the UAV will report the current position and request to replan the path. When the UAV moves to position C, it detects two dangerous weather phenomena. Similarly, the UAV determines the flight direction according to the resultant direction of the resultant force, and finally reaches the destination.
The entire dynamic path planning process is expressed as follows: STEP 1: Determine the starting point and the destination point, and initialize the meteorological parameters and IIWD algorithm parameters. STEP 2: Estimate whether the UAV has reached the destination point; if not, proceed with the static path planning, otherwise end the process. STEP 3: Determine whether there is a new forecast product or replanning request; if so, return to STEP 1; otherwise, go to STEP 4. STEP 4: Estimate whether there is a sudden dangerous weather phenomenon; if so, perform the dynamic path planning; otherwise, return to STEP 1.
The moving direction of the UAV is determined by the repulsive force of the dangerous phenomena and the attractive force of the destination position. Therefore, the computational complexity is O(N τ + 1), that is O(n).
Once the path planning is completed, the resultant path requires to be smoothed, and the third B-spline curve is used to smooth the obtained path. B-spline is a curve generated by approximating a set of control points. The corresponding calculation is expressed in formula (28). where, P i is the characteristic point of the curve, and F i,k (t) is the order B-spline basis function. For cubic B-spline curve, the basis function in the equation is: A sample of the B-spline is shown in Fig. 6, where the blue dot is the original point, and the green line is the polygonal line, and the red line is the smoothed curve.

V. SIMULATION AND RESULTS
In order to verify the effectiveness of the proposed algorithm, we conduct simulation on the static path planning and the dynamic path planning. The performance constraints of the UAV are p min = 10km, ϕ min = 45 • , and L max = 4000km. The parameters of the improved intelligent water drops algorithm are expressed as follows: N IWD = 50, v ini = 120km/h, soil(p i , p j ) = 1000, I max = 500, a v = 1, b v = 0.01, c v = 1, ρ n = 0.8, a s = 1, b s = 0.01, and c s = 1.  [31], genetic algorithm(GA) [32] and Q-learning algorithm [33] are applied for comparative simulation. Figures 7 and 8 illustrate the resultant path obtained by using the four algorithms. After calculation, the Richardson  number of the grid points at the altitude of 6000 m in the target area is much greater than 0, which means the possibility of the convective weather phenomenon is close to 0. Under this circumstance, the main threat to the UAV flight is ice accumulation. The color bar at the bottom of the figure indicates the possibility of ice accumulation, which is consistent with Ic/100. The blue arrows represent the direction of the wind, and their lengths represent the wind speed. A number of performance metrics have been investigated to evaluate the optimality of the proposed algorithm. These metrics involve the flight distance, flight time, average flight speed, and risk degree. The improved intelligent water drops algorithm shows a better performance than the other three algorithms in simulation 1 and simulation 2. The simulation results are presented in Table. 2.
From the simulation results in Table. 2, it can be seen that the flight time and risk degree obtained by IIWD are smaller than the other algorithms. The mechanism of the IIWD shares a similar characteristic with ACO, and the latter chooses the path according to the pheromone in the path. However, the reason that IIWD shows a better performance than ACO in this scenario is that the velocity of water drop will have feedback to the soil in the path, which means the water drop with higher velocity will take more soil from the path. This path is likely to be chosen by another water drop, and this positive feedback mechanism can make the algorithm find a better solution. The genetic algorithm does not get a good solution probably because, in the process of crossover and mutation, GA needs to consider the constraints of path continuity and steering angle, which leads to the need to generate feasible solutions repeatedly randomly. When we use the Q-learning algorithm to solve the problem, if the wind speed, wind direction, and meteorological factors are brought into the cost function, the difference in the value of the grid points will lead to the failure of generating a feasible solution. Therefore, like the traditional path planning problem, we assign a positive value to all available points and a negative value to the infeasible nodes (the probability of ice accumulation is greater than 50%). The results show that the obtained path is basically a straight line between the starting point and the destination point. Although the total path length is short, the effect of risk degree, flight time, and other evaluation indicators are not acceptable.

B. PATH PLANNING AT DIFFERENT ALTITUDES
In the previous section, we only mentioned ice accumulation, the wind direction, and the wind speed on the path planning. This is because the atmosphere is found to be very stable at the altitude of 6000 m, and the obtained Richardson number is found to be much greater than 0. However, UAVs may operate at different altitudes according to the requirements of exploration missions. The UAV path planning at different altitudes (5000 m, 4000 m, 3000 m, and 2000 m) is considered in this section. Figures 9, 10, 11, and 12 show the flight distance, flight time, flight speed, and risk degree of the UAV path by using the improved intelligent water drops algorithm at different altitudes. In terms of time scale, the path planned every day is not constant, and the choice of the UAV takeoff time exhibits a greater impact on the path. Simultaneously, the flight time and the flight speed of the planned path show certain similarities. For example, on June 2, the flight time at four altitudes is found less than that at other times, and the corresponding average flight speed is found larger. Similarly, on June 3, the average flight speed decrease also led to a longer flight time.
As shown in Fig. 9, the total length of the trajectory planned at each altitude ranges from 2400 km to 2900 km. The path  planned at the altitude of 2000 m shows the lowest total path length. The flight time of each altitude ranges from 16 hours to 22 hours in Fig. 10. In Fig. 11, the average flight speed at the altitude of 5000 m is higher than that of other altitudes, while the flight speed at 2000 m is the lowest. The total risk degree experienced by the aircraft at 3000 m and 5000 m are smaller than that at the other two altitudes, and the total risk degree at 2000 m is found to be larger, expressed in Fig. 12. The horizontal wind speed increases with the ascending of flight altitude, which is the reason that the average flight speed at a high altitude is generally greater than that at a low altitude. IIWD algorithm selects the grid point with a higher velocity in order to increase  the ground speed when planning the path. In a low-altitude environment, the wind speed at each grid point is generally low, so the average wind speed remains at a low level. From the aspect of meteorological factors in IIWD, we define the risk degree DoR, which is determined by the two factors, namely, ice accumulation and Richardson number. At the altitude of 2000 m, due to the influence of a lower cushion surface, convective weather phenomenon occurs frequently. Therefore, we have to select the path with a high-risk level. Besides, due to the proximity to the ground, except at high latitudes, the temperature is much greater than 0, so the possibility of ice accumulation is close to 0. In the middle troposphere above 5000 m, severe convective weather phenomenon still occurs at this altitude. However, compared to 2000 m, the occurrence of the convective weather phenomenon significantly reduces, and the probability of aircraft ice accumulation increases. The height affected by convective weather phenomenon is generally 8000 m, above which the temperature is about −20 • C, and ice accumulation does not occur easily. Therefore, ordinary civil aircraft generally fly at an altitude of 10000 m. However, for the environment mapping or the reconnaissance UAV, the lower flight altitude can obtain higher-resolution data, so it is necessary to select the appropriate flight altitude taking into consideration the mission objectives and specific requirements.

C. DYNAMIC PATH PLANNING
As the weather changes with time, we carry out the static path planning for the UAV according to the weather forecast products before takeoff in the process of a mission. When the new forecast products are received, the path gets dynamically adjusted. Fig. 13 shows the path adjustment process during a voyage.
As shown in Fig. 13, the green points are the replanning points which means that the ground station sends the new path to the UAV by the satellite link via the calculation of IIWD after it receives the new forecast products. The blue line shows the original planned path, and the red line shows the adjusted path. As can be seen from the figure, five path adjustments are made after the UAV takes off. The new path is basically found consistent with the original path since the whole environment does not change much. Figure 14 shows the schematic diagram of the path planning when the UAV encounters sudden dangerous weather phenomena, which are detected at the position of green points. The blue line shows the original planned path. If the  path is not changed, the flight risk of the UAV significantly increases. The red line shows the adjustment by using the artificial potential field method. When the UAV detects the risk area, the flight direction is adjusted by the resultant force of both the repulsive force of the risk area and the attractive force of the destination. It can be seen from Fig. 14 that the UAV could effectively avoid dangerous areas and achieve the purpose of safe flight.

VI. CONCLUSION
In this paper, a meteorology-aware path planning algorithm of the UAV based on the improved intelligent water drops algorithm is proposed. The path planning model of the UAV in a dynamic environment is established. The influence of the wind direction, wind speed, ice accumulation, Richardson number, and other meteorological factors on the UAV flight are comprehensively considered. The intelligent water drops algorithm is improved by introducing parameters of risk degree in order to adjust the path selection probability to choose the trajectory with low risks. The simulation results show that the improved water drop algorithm exhibits better performance than the ACO, GA, and Q-learning algorithm in solving the problem of path planning. The results show that the higher the altitude level is, the more the flight time can be shortened by using wind speed. Moreover, the risk degree of the UAV decreases with the increase of the flight altitude. Finally, the flight effect of the virtual potential force algorithm in a dynamic environment is verified. The algorithm can effectively avoid the dangerous area and ensure the UAV flight's safety.
In the future study, the influence of wind speed, temperature, altitude, and relative humidity on fuel consumption and the efficiency of the engine will be considered. The minimum energy consumption will be considered the optimization objective in the path planning to improve the economic benefits while ensuring safety. Simultaneously, the sudden convective weather phenomenon considered at present is static, which does not change with the development of the atmospheric circulation. Therefore, we will take the dynamic of the convective weather phenomenon into account in the future.