A Multi-Frequency Investigation of Air-to-Ground Urban Propagation Using a GPU-Based Ray Launching Algorithm

Unmanned Aerial Vehicles (UAV), also known as “drones”, are attracting increasing attention as enablers for many technical applications and services, and this trend is likely to continue in the next future. When compared to conventional terrestrial communications, those making use of UAVs as base- or relay-stations can definitely be more useful and flexible in reaction to specific events, like natural disasters and terrorist attacks. Among the many and different fields, UAV enabled communications emerge as one of the most promising solutions for next-generation mobile networks, with a special focus on the extension of coverage and capacity of mobile radio networks. Motivated by the air-to-ground (A2G) propagation conditions which are likely to be different than those experienced by traditional ground communication systems, this paper aims at investigating the narrowband properties of the air-to-ground channel for 5G communications and beyond by means of GPU accelerated ray launching simulations. Line of sight probability as well as path loss exponent and shadowing standard deviations are analysed for different UAV flight levels, frequencies and dense urban scenarios, and for different types of on board antennas. Thanks to the flexibility of the ray approach, the role played by the different electromagnetic interactions, namely reflection, diffraction and diffuse scattering, in the air-to-ground propagation process is also investigated. Computation time is reported as well to show that designing UAV communication networks and optimising their performances in a fast and reliable manner, might avoid exhausting – multiple - measurement campaigns.


I. INTRODUCTION
Unmanned Aerial Vehicles (UAVs) have been historically used in military applications for more than twenty years but nowadays the excitement for extending their use to civilian and commercial applications is growing worldwide, thanks to the advancement of manufacturing technologies in batteries, electronics and lightweight materials together with their cost reductions, making UAVs more easily accessible to the commercial market. Endless developments in artificial intelligence technology and big data analytics have paved the The associate editor coordinating the review of this manuscript and approving it for publication was Yiming Huo .
way to a multitude of application fields, where the UAVs can do their jobs autonomously [1], [2] without any need for ground pilots to be in sight of the drone.
Among the many and different applications -ranging from public safety surveillance and law enforcement actions to precision agriculture or search and rescue operations [3], [4], many of them also quite challenging to set-up [5] -the concept of one or more UAVs operating as airborne Base Stations (BSs) or flying relays has recently received great attention in the framework of 5G systems [6]- [11] to expand wireless connectivity to unserved areas and/or to increase capacity in busy areas, or even to back up a damaged ground network in emergency situations. VOLUME 9, 2021 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ The quite short flight time, less than a hour in most cases, is possibly the major drawback of this idea, in spite of the on-going advances in battery technology and charging techniques [12]. Besides, concerns with public safety still exist and thus the deployment of flying UAV swarms for whatever commercial applications must often come to terms with strict limitations enforced by local regulations [13].
In order to tackle these manifold issues, UAV-aided communication networks require careful and thorough design, in order to optimally set the number and the routes of the simultaneously flying BSs. To this aim, some awareness about the air-to-ground (A2G) propagation channel is necessary, as it may contribute to plan the flying paths / hovering locations where the UAVs can provide the most helpful wireless access to users at street level and/or indoor.
In urban environments, as the UAVs may well stand out against the building layer, the A2G channel may benefit from a better line of sight (LOS) occurrence compared to standard cellular networks, where the antenna visibility is instead more easily prevented even in case the terrestrial BSs are placed on tall buildings or high masts. At the same time, LOS conditions cannot be counted on always and everywhere, and multipath propagation might therefore still play an important role. In particular, scatterers at ground levels, like buildings, bridges, overpass roads and other manmade and natural clutters can heavily affect electromagnetic signals propagation.
Although the different papers on the subject, as referenced in the next paragraphs, the reader should bear in mind that measurements campaigns for the characterization of the A2G channel represent a rather challenging task, as they undergo the aforesaid limitations about the payload UAVs can rise, the battery lifetime and the permission to fly, in addition to the common shortcomings of the experimental approach, i.e. expensive equipment and long execution time.
For these reasons, dense urban scenarios have not been primarily targeted in previous studies, and A2G propagation measurements have been mainly carried out in rural/suburban environments [14], and often limited to few and quite short flying paths in clear visibility conditions to the pilot [15]- [17].
Because of such hindrances to extensive measurement campaigns, established path loss empirical models, seem also lacking, at least for the dense urban case and for a wide range of frequencies.
Preliminary and basic studies have been done in [18], where the authors studied different path loss models with potential applications in the urban scenario -especially in emergency communication systems -and compared their performances. In [19], the authors extracted the main parameters of path loss model in urban environment and propose a simplified Saleh-Valenzuela model to describe the stochastic properties of the arrival delays and amplitudes of resolvable multipath components (MPCs) in wireless transmission systems. In [20], the researchers studied the propagation characteristics of air-to-air channels in a virtual urban scenario generated with the ITU-R model.
Generally, for practical cases there is still a trend to rely on very simple propagation models, like the simple free-space or ''two-rays'' formulas [3], [21]. Nevertheless, they cannot be reliable when multipath components (besides the ground reflection) are significant. Similarly, empirical/statistical propagation models for terrestrial cellular networks (like the Okumura-Hata model [22]- [24] and the many ''Hata-like'' models) may fit the air-to-ground scenario to a limited extent, as their reliability is likely to be doubtful when the drone flight level is remarkably higher than the buildings level [25]. In conclusion, empirical/statistical models specifically conceived for the A2G channel can be hardly found in the open literature [26] and this paper aims at bringing insights to a broader extent.
Deterministic propagation tools, either Ray Launching (RL) or Ray Tracing (RT), may represent an effective solution for A2G channel assessment for two main reasons. Firstly, they are inherently site specific and yet much faster than actual measurements: therefore, they can replace measurements to some extent and might be run over the huge spatial domain -including a broad range of distances and heights -necessary to characterize A2G propagation over large urban areas, thus providing a suitable data set for statistical analysis and/or the development of empirical A2G path loss models. Secondly, the ray approach can provide a significant insight into the A2G propagation process, e.g., about LOS/NLOS probability or the composition of the major multipath contributions at different flight levels and/or different communication frequencies. Unfortunately, running RT/RL simulations often involves a fairly high computational effort, which turns out particularly heavy in the A2G urban scenario, as the drone might be in clear visibility with a multitude of buildings underneath, and therefore a very large number of propagation paths might be present. In order to limit the computational burden, ray-based investigation of A2G propagation has been often limited to restricted areas with few ground locations, thus downgrading the actual validity of the propagation analysis [27], [28].
In this paper, narrowband characterization of the A2G channel in urban environment is carried out by means of a GPU-based, discrete RL model (also called DED-RL) specifically conceived for fast coverage assessment over large areas and already validated vs. measurements in urban environment [29]. Wideband characterization -which focuses on key parameters like power delay profiles or root-mean-square (RMS) delay spreads -is also necessary to fully characterize the A2G wireless channel. As existing studies are mainly limited to rural/suburban scenarios [30], wideband description of A2G urban links is by far not consolidated yet and requires further specific investigations. This is not in the current scope of this work and is going to be addressed in a different study. With respect to previously published studies [31]- [34], the narrowband analysis is here extended over three different urban scenarios, considering a broad range of 54408 VOLUME 9, 2021 altitudes and with communication frequencies spanning up to the millimetre-wave band. In particular, the considered frequencies range from the 700 MHz band to the 3.5 GHz and 26 GHz bands (bands assigned worldwide to 5G technology) up to the 70 GHz band (already existing unlicensed band worldwide). Besides path loss and shadowing, attention is here devoted to the role played by the different electromagnetic interactions, namely reflection, diffraction and diffuse scattering, when multipath occurs in the A2G channel.
After a brief background and motivation in Section I, Section II and Section III respectively provide some general remark about the A2G channel and a short description of the ray-launching software together with the simulation setup. Section IV deals with LOS probability while Section V shares the outcomes of the narrowband air-to-ground channel study for three typical dense urban environments with omnidirectional and directional antennas. Section VI is dedicated to the role played by the different mechanisms in the propagation process. Section VII deals with computational time assessments related to A2G channel simulation, while some conclusions are finally drawn in Section VIII.

II. PRELIMINARY CONSIDERATIONS ABOUT A2G PROPAGATION MODELING
Compared to terrestrial cellular networks, UAV-assisted communication systems are characterized by some peculiar aspects that should be accounted for in the propagation modelling procedure. They are briefly discussed in the following sub-sections.

A. EFFECT OF THE ANTENNA RADIATION PROPERTIES
Propagation measurement/simulation activities should be always planned with the transmitting-receiving locations mutually included inside the main radiation lobe of the antennas. This is done to minimise the possible side propagation effects like the multiple direction of departure / arrival of the different multipath components. Therefore, the simulation/measurement area shrinks down for increasing directivity of the on-board antenna and/or for lower flight levels of the UAV (from red to green colour in Fig.1), since the UAV is expected to be equipped with an antenna commonly pointing towards the ground. This difference in the area extension -which depends thus on the on-board antenna radiation properties -may correspond to different best-fitting lines matching the collected attenuation values as shown in Fig. 2, again moving from red to green colour. As a consequence, this means that the on-board antenna directivity can significantly affect propagation parameters like the Path Loss Exponent (PLE) α and the shadowing standard deviation coefficient σ of a Hata-like path-loss model [24], respectively representing the slope of the best-fit line and the standard deviation of the attenuation samples with respect to the fitting line.
It is worth pointing out that the same impact does not automatically occur in terrestrial cellular networks, where the main lobe of the base station antennas is fundamentally parallel to the ground (or slightly tilted downwards), and therefore increasing the directivity primarily restricts the size of the measurement/simulation area mainly in the azimuth domain rather than reducing the distance range.

B. NOTES ON THE PATH LOSS FITTING FORMULAS
According to the Hata-like approach, the PLE can be derived by fitting with a straight line ( Fig. 2) a large set of pointspecific attenuation values displayed in a log-log graph. The attenuations samples often come from measurement investigation, but they can be provided by deterministic propagation model as well, as achieved in the study herein.
Two different equations are usually considered for the best-fit line [35]: PL dB (d) = 10 · α 2 · log 10 (d) + β the former is a ''fixed intercept'' model (also called in literature ''close in''), where d is the distance between the drone and the ground tile, PL free (d 0 ) is the free space loss at a reference distance d 0 (set equal to 1 m for simplicity) and α 1 represents the PLE; the latter as in (2) consists of a ''floating intercept model'', where a second fitting coefficient (β) is introduced in addition to the PLE (α 2 ). As long as the distance range is similar to or even larger than the corresponding calculated attenuation range (in relative terms), Equations (1) and (2) turn out to be quite similar, with α 1 ≈ α 2 . This may be not the case in UAV to ground communications, when the flight altitude is low and/or the UAV antenna directivity is large, corresponding to a quite narrow beam and thus a small antenna footprint on the ground. Ground locations falling within the main antenna lobe would basically share the same link distance, whereas attenuation can still undergo large fluctuations, especially in case both LOS and NLOS conditions occur. In Equation (2), the corresponding optimal value of α 2 might become awkwardly large, with questionable physical meaning although mathematically correct. Conversely, in Equation (1), α 1 would be much less affected under the same conditions, i.e., it keeps physically sounded whatever the UAV height and its antenna directivity are. This is shown as an example by the continuous and dashed lines in Fig. 2. Because of this greater robustness, the close-in model only is considered in this work.

C. MAXIMUM TOLERABLE LOSS
Wireless communications are always affected by signal attenuation and a maximum acceptable propagation loss (PL max ) can be easily set based on the transmitter power level, the receiver sensitivity and the noise threshold at the receiver side. Depending on the attenuation law, the maximum communication range can be easily estimated. The same limitations do not automatically apply when the propagation analysis is carried out by means of a software simulation tool, which can virtually compute any received signal intensity value, whatever the attenuation experienced by the propagating signal.
In order for the channel simulations to effectively replace channel measurement, i.e., returning similar data set with a high level of confidence, a PL max value must be properly set: the analysis of the simulations data should be restricted to the distance range where the probability to compute attenuation values greater than PL max is basically negligible. In accordance with [36], the statistical analyses carried out in this work are limited to a maximum threshold distance d max .
Besides, due to the broadband nature of the analysis in this paper, it has been further extended to a dual PL max definition: PL max = 160 dB for the lower frequencies − 700 MHz and 3.5 GHz -and PL max = 180 dB for the higher ones − 26 and 70 GHz. This strategy stems from the need to better address such a broad range of frequencies, where a single PL max value, as the one used in [36] and [37] for a single frequency, would not have been sufficient to fit them all.

III. SIMULATION SET-UP A. RAY LAUNCHING TOOL
The Discrete Environment-Driven RL model (DED-RL), has been introduced and extensively described for the first time in [29]. The simulation tool relies on a digitalised urban model of the city (Fig. 3) where each building is modelled as a polygon with a defined shape, material, position, and height. The model is discrete, i.e., the building walls are properly divided into tiles with a predetermined size. In addition, other advanced features are implemented in the model to achieve very high accuracy while drastically reducing computation time. The main advantage of discretization is that the tile centres can be assumed as fixed points, therefore all the visibility relations among the tiles can be pre-computed and properly stored into a visibility matrix.
This visibility pre-processing is carried out only once, for a given simulation scenario. Of paramount relevance is the fact that both the visibility pre-processing and the bouncing of the ray tubes in RL, are suitable to be parallelized into Graphics Processing Units (GPUs) and can be thus run by means of NVIDIA cards using the CUDA framework. Typical computation times for complete predictions range from a few seconds to tens of minutes, depending on the size of the urban scenario, the characteristics of propagation and other elements on board of the drone, like the antenna or the transmitting power, as it will be shown in Section VII.
In order to perform propagation studies with general validity and not influenced by specific antenna characteristics, in this paper we consider UAVs equipped either with an 54410 VOLUME 9, 2021 isotropic (omnidirectional) antenna or with a directional one, installed below the drone body and pointing perpendicular towards the ground, with a fixed aperture angle θ. This angle θ is simply the aperture of an ideal and symmetrical radiation cone within which, we assumed an ideal constant gain and consequently negligible side lobes to keep the RL simulations affordable. This can be also seen as a best-case scenario where rays are launched only within a specific cone angle θ, thus drastically reducing the computation effort of calculating multiple interactions. In addition, this choice allows us to possibly extend our results to a wider range of directional antennas as a function of the simple aperture angle θ and not linked to any specific antenna radiation diagram. Consequently, for the spot area covered on the ground (i.e. the drone antenna footprint) we assume a circular shape whose radius depends on the aperture angle and the UAV altitude.

B. SIMULATION SCENARIO AND ITS PARAMETERS
The DED-RL simulations take into account both antenna patterns (isotropic and directive) and different parameters of the transmitting system such as frequency and power, as well as the number and types of the allowed propagation interactions (i.e., reflections, diffractions, diffuse scatterings) to be considered. Specular reflection and diffraction are modelled according to the Geometrical Theory of Propagation (Geometrical Optics + Uniform Theory of Diffraction), while diffuse scattering is modelled through the Effective Roughness (ER) model [38], [39], whose main parameter is the inherited diffuse scattering coefficient S in simple way. The output of the simulations we are interested in is the received power level in dBm at ground, while path loss and other major propagation figures are retrieved by means of specific post-processing scripts.
The three urban environments have been picked up as good representatives of dense inhabited areas: one American and two European cities with towers, churches, lengthy porticoes and a specific narrow maze of streets that lead to a well-preserved historical centre or -conversely -wide streets full of skyscrapers that end up into a modern metropolitan centre. Concerning resolution, a discretization with 10m × 10m tiles has been considered; further details about the simulation scenarios can be found in Table 1.
Ray Launching simulations have been extensively carried out for 8 different hovering positions at 7 different altitudes and 4 communication frequencies. Table 2 shows the key simulation parameters used for the three dense urban scenarios.

IV. LOS PROBABILITY ANALYSIS
The availability of Line of Sight has a great impact on the properties of the wireless channel, and therefore on the performance of the wireless link. Propagation with a high Probability of LOS (PLOS) with a strong direct path intensity corresponds to free space like conditions (i.e. PLE ≈ 2) with limited σ values, whereas increase in both the PLE and the shadowing effect is expected when a LOS condition seldom occurs. Investigation of LOS probability is especially important in the mm-wave band, where multipath is often made of few, sparse components for the user to rely on in case the LOS path is obstructed.
Whether the LOS path exists or not depends on various factors, e.g., terrain features, density of buildings, their height and mutual distance or distribution. As ground users are inherently placed in the tiles centre by the RL algorithm, i.e., almost everywhere along the streets, the LOS probability discussed in this section basically depends on the UAV flying altitude and on the geometrical properties of the environment. It is worth mentioning that the geometrical LOS is independent of the system frequency.
Models for the computation of the LOS probability are discussed in [40]- [43], whereas the possibility of identifying LOS and NLOS conditions based on the narrow/wide band channel parameters is proposed in [44].
Most of the referred models stem from ITU [42], where the range dependence of the LOS probability is mainly considered, and no specific attention is paid to the antenna height. This limitation does not fit the A2G channel characterization, where investigating the LOS occurrence at the different possible height of the UAV is instead of primary interest. In this VOLUME 9, 2021 context, the formula proposed in [43] seems to be the only exception, as it accounts for both the link distance and the antenna height.
Therefore, LOS probability vs. UAV height for isotropic antennas is presented In Fig. 4 for the three cities of interest and compared with the PLOS trends calculated from the analytical expression in [43].
For the sake of comparison, i.e., to avoid the computation being affected by the different areas of the simulation maps (Table 1), the PLOS occurrence was calculated via RL over circular spots of 2 km 2 underneath each UAV position. As expected, it can be seen in Fig.4 that for an isotropic antenna, the higher the drone is flying, the higher the PLOS. Conversely, when the UAV is flying at very low altitudes, PLOS drops down, because of the increased shadowing level. In both cases, the trend is clearly monotonous, while the difference in LOS probability between the three urban RL samples may be attributed to the different properties of the building layers. San Francisco exhibits the greatest building density, and it includes some high-rise building areas together with residential neighbourhood with lower houses. Furthermore, San Francisco has a hilly terrain profile. These reasons explain the corresponding lower values of LOS probability in Fig. 4. The analytical model in Fig.4 also returns greater PLOS values at higher UAV altitudes -that is not actually surprising -as well as it keeps the same general relationship between the different scenarios, with San Francisco and Munich showing the heaviest and the lightest degree of obstruction, respectively, and Bologna in between. Nevertheless, the analytical and the RL assessment for the same urban layout are not always in fair agreement. It can be noted that the model in [43] is based on the description of the environment through few, rather simple parameters, like the percentage of built area and the building density. On the one hand, this leads to a very friendly, ready-to-use final expression, on the other hand it may represent an oversimplification, thus limiting in some cases the accuracy of the retrieved LOS probability.
In Fig. 5 the same LOS probability vs. UAV height is presented for a directive antenna, with an aperture θ = 80 • as an example -but similar trends have been observed also for other aperture values. As a directive antenna on the UAV is expected to provide wireless connectivity to the users primarily through its main radiation lobe, the evaluation of the LOS probability is now carried out on the tiles inside the antenna footprint at ground level. Compared to Fig. 4, LOS probability in Fig. 5 is always greater, irrespective of the drone height. When the UAV is very close to ground, the footprint underneath is rather small and most of the few ground tiles inside it are in line of sight, corresponding to a probability equal to 80-85% in Fig. 5. While the UAV is flying at a greater altitude, more tiles fall within the ground spot. At first, they are likely to be shadowed by buildings rather than in sight of the UAV, and the LOS probability therefore drops down; then, the situation is reversed, and the probability begins to increase. This trend reaches then a steady state for which the PLOS becomes nearly constant at higher altitudes. According to our investigations, the minimum in this case can be related to a new parameter called ''urbanization threshold'' -expressed in meters -that gives a feeling of the complexity of the 3D urban environment, in terms of number of buildings, heights, standard deviation, as per Table 2: while the UAV is flying higher -but still below this threshold -the PLOS drops down. Once the UAV overtakes the same threshold, the visibility of the UAV over the ground spot and consequently the PLOS starts to increase.
This trend reaches then a steady state when the height of the surrounding buildings becomes somehow negligible with respect to the flying altitude of the drone. This value is between 50 to 75 meters in Bologna and Munich, while it is definitely higher in San Francisco where it reaches a value close to 300 meters.
The achieved results -with takeaways in both figuresshould draw our attention to the key fact that having an A2G link does not always imply perfect LOS conditions. Since the UAVs might be expected to fly at altitudes lower than 100m in many applications, LOS turns out to be a rather uncommon condition in the isotropic case, whereas it does not occur in approximately 30% of the cases, when the directive antenna is considered.
The A2G channel behaviour is therefore less simple than it may appear, with multipath and shadowing effects that shall not be automatically negligible, and there is thus a need for further investigations about the A2G propagation mechanism.

V. NARROWBAND CHANNEL PARAMETERS
In this section, we extract and analyse the narrowband channel parameters. The choice of an isotropic antenna is made in the following sub-section to focus our attention on the channel itself, with no misleading information coming from the antenna on board of the drone. Furthermore, it is worth noting that the isotropic case can somehow represent real situations where the radiation lobe of the antenna is wide enough to illuminate a large area below the UAV frame, e.g., when UAVs are deployed to set-up temporary emergency communication infrastructure after natural disasters. As drones can be as well equipped with directive antennas, e.g., to improve wireless connectivity at the cell boundary in cellular networks, the narrowband analysis is also carried out taking into account different beamwidths of the on-board antennas. To collect a common value for a specific altitude and frequency, regardless of the specific UAV hovering position, PLE and σ were calculated by joining the different cloud of points for each of the 8 different drone positions ( Table 2) and then calculating the best fit line according to (2). This strategy makes particularly sense for low UAV flight levels, where values can differ significantly among the chosen positions due to specific local obstructions; the higher the UAV flight level, the less this effect is seen, as all the location experience similar mean propagation conditions and obstructions from the environment.

A. PATH LOSS EXPONENT AND SHADOWING STANDARD DEVIATION -ISOTROPIC
The first parameter we have focused on has been the path loss exponent α, which accounts for attenuation and for the ''mean degree of obstruction'' on the propagation of the paths. The main results about PLE are shown in Fig. 6 for different UAV altitudes and different frequencies, using an isotropic antenna. As seen in previous works [26], [32], PLE has a decreasing trend with height: the higher the UAV is flying, the lower the PLE. It is worth noting that for each environment the values of α achieved at the lowest UAV height and frequency -respectively equal to 30m and 700MHz -are in quite good agreement with those retrieved from the wellknown Hata model [24].
Nevertheless, the path loss factor reduction in Fig. 6 for increasing flight height is faster compared to the trend highlighted by the Hata formula. Although the number of receiving LOS locations increases as the UAV flies at higher altitudes, the PLE at 450m is still greater than 2, according to LOS probability values which were well lower than 1 (Fig. 4).
Concerning frequencies, PLE is faintly depending on them although we would have expected not, as already revealed in [24]: α might increase with frequency due to more challenging overall propagation conditions. This dependency is slightly higher for the lower bands, between 700 MHz and 3.5 GHz, than for the higher ones, between 26 GHz and 70GHz, where the gap is marginally noticeable.
The achieved σ values as a function of the UAV flying altitude are reported in Fig.7. It is observed that σ also has a decreasing trend with height: the higher the UAV is flying, the lower the σ value is, since the effects of propagation at local level become negligible in favour of a wider covered area. The same Fig. 7 shows the dependency of the shadowing component from frequency. This dependency seems to be more visible for the lower bands than for the higher ones, where the gap between 26GHz and 70GHz is again narrower.
Shadowing levels in Fig. 7 look quite large compared to what usually assumed for wireless cellular networks, where σ hardly exceeds 10 dB in urban scenarios at UHF frequencies [45]. This difference can be related to the use of an isotropic transmitting antenna in dense urban contexts, where VOLUME 9, 2021 both line-of-sight conditions and heavily obstructed ground locations are likely to be simultaneously present, thus increasing the signal spread. This is especially stressed at lower UAV altitudes, where pathloss at the locations in LOS is reduced by the limited distance, whereas it is instead particularly large at the heavily obstructed furthest tiles.
To achieve mathematical accuracy while keeping the necessary filtering on distance introduced in Section II, it is vital to ensure that the maximum path loss in the ray launching model exceeds the values expected in the actual radio system [37] and for this reason we had to account for PL max = 180 dB and 210 dB for the case of San Francisco.

B. PATH LOS EXPONENT AND SHADOWING STANDARD DEVIATION -DIRECTIVE
In this subsection, the UAV is then assumed to be equipped with a directional antenna of fixed aperture angle θ, placed under the UAV fuselage and pointing towards the ground, as introduced previously in Section III.
In Fig.8 and Fig.9 -which refer only to Munich city and to only one frequency at 3.5 GHz for simplicity -it is possible to see respectively the trend of the path loss exponent α and the shadowing standard deviation σ as a function of the antenna aperture angle θ, ranging from 40 • to 180 • , and for different flight altitudes.
As a general trend, the increase in the antenna directivity reduces the size of the service area at ground, resulting in a greater LOS probability ( Fig. 4 and Fig. 5). The weaker average obstruction on the propagating paths results in the lower values of PLE in Fig. 8, compared to Fig. 6(b). At the same time, the extent of shadowing effects is also reduced, corresponding to a decrease in σ for decreasing antenna aperture values. It is interesting to report that for θ < 80 • the simulation samples available in our model are not statistically significative at low UAV altitudes but -being already close to α = 2 for the same θ value at higher altitudes -we can reasonably extrapolate them as those in the free space LOS since all the tiles within the antenna footprint happen to be in clear line of sight to the drone. This constraint does not exist anymore for higher altitudes, where statistical samples are available also for θ < 80 • . The more the antenna aperture angle increases -and thus the directivity decreases -the more both α and σ parameters tend to the isotropic case (corresponding to antenna beamwidth equal to 180 • in Fig. 8 and Fig. 9). From these trends at different UAV altitudes, it is possible to see again that the higher the drone is flying, the lower both α and σ values tend to be. The same trend was observed for the other frequencies, also in Bologna and San Francisco, except for a different gradient which accounts in our study for the different urban environments and the different frequencies.

VI. WEIGHT OF MULTIPATH COMPONENTS
To the best of the authors' knowledge, no article has shed some light so far on the different components that contribute to the A2G propagation. It is known that the received power in a specific location (tile) is the result of different contributions: mainly one direct ray (when existing) and many indirect rays, which consist of multipath signal components generated by reflection, diffraction and diffuse scattering phenomena.
As a matter of fact, investigating the impact of the different electromagnetic interactions on the A2G channel through in field measurement is an overwhelming task, as UAVs could hardly carry the necessary heavy equipment like channel sounders and processing hardware. Conversely, ray-launching simulations can easily provide insight into the power brought to the receivers by the direct path (when present) or indirect ones or even combinations of them.
By running several ''single events only'' simulations, i.e., enabling only direct rays or only diffracted rays or only reflected rays and so on, it is possible to evaluate the incoherent power contribution of each event -or combination of events -to the total received power.
To further investigate the A2G link, the analysis has been focused on the isotropic and directive antenna cases, bearing in mind the same configuration parameters as in Table 2. For the sake of simplicity and clarity, figures have been focused only on 700MHz and 26 GHz, as significative representative frequencies for the low and high 5G bands.
At the lower frequency of 700 MHz (Fig. 10), diffraction is the major propagation mechanism overall in case of isotropic antenna: as the UAV flying altitudes are always higher than the average building heights (Table 1) while single diffractions on the rooftop boundaries can convey significant power to most of the ground locations. Anyway, the lower the UAV height, the further the distance of the diffracted rays from their shadow boundary; this may explain why A2G propagation better relies on combinations of different interaction at the lowest altitudes (Fig. 10). In particular, combined events challenge the leadership with diffraction in San Francisco (Fig. 10c), where hilly terrain and high-rise buildings can support propagation through mixed interactions along longer wireless paths. This seems especially true at the higher altitudes, where the UAV has better visibility of the whole urban area, and rays experiencing few bounces on the taller buildings are likely to provide radio coverage almost everywhere. Not surprisingly, the strength of the LOS path increases with the UAV altitude to the detriment of diffraction relative weight (Fig. 10 a-b). This is not as much true in San Francisco (Fig. 10c) where the high skyscrapers prevent LOS to dominate even at high flying altitudes and A2G link is mainly relying also on diffraction and combined events.
Scattering as well as reflection alone contribute to the total received power to a relative extent always lower than 10%.
At higher frequencies, namely 26 GHz (Fig. 12), power contributions from diffraction undergo a great reduction, that is of course physically sounded in case of isotropic antenna. Conversely, scattering rises to greater importance (20%-40% of the total), to the extent that it stands out as the major propagation interaction, at least at the lower flying levels. This seems in contrast with the general statement that the millimetre wave channel is sparse and dominated by few powerful contributions commonly classified as reflections.
In this respect, it should be noted that: i) A2G propagation is mostly happening in the longitudinal direction, perpendicular to the ground since the antenna is pointing downwards ii) sparsity of the mm-wave channel has been mainly investigated indoor and with the antenna pointing parallel to the floor iii) the scattering model embedded into the RL tool aims at modelling not only actual scattering from surface roughness, but also other multipath components -regardless of their specific nature -spread by the architectural and structural elements on the facades of buildings (i.e. balconies, pipes, windows frames, etc.) usually not included in the simulation input files for practical reasons. LOS power is still increasing with the UAV height and becomes remarkably uppermost at the higher altitudes. It is worth noting that up to few hundred meters in altitude, propagation seems never driven by the LOS contribution only, thus highlighting the need for multipath propagation models for the A2G channel.
The outcomes in case of a directive antenna are reported in Fig.11 and Fig.13, in good agreement with the research expectations. Fig. 11 and Fig. 13 show the same trends discussed before but for a directive antenna with an aperture angle θ = 100 • , taken as a representative example. Not surprisingly, the LOS path in both cases is the predominant contribution with an outstanding value ranging between 55% and 60% along the different UAV altitudes, with peaks between 70% and 80% at low UAV altitudes. This result matches the trends seen in Section IV, since most of the tiles within the antenna footprint are in clear line of sight. At the lower frequencies, diffraction supports propagation for a good 30% while it drops down to 10% at higher frequencies.

VII. NOTES ON COMPUTATION TIME
In addition to getting a huge number of data samples and to investigating different situations, the same runs gave us a better insight on the computation time and the speed-ups over GPU accelerated simulations, according to the parallelization concept introduced in Section III and which forms the basis of the DED-RL tool used in this paper. Simulations were run on an Intel R Xeon R CPU E5-2620 v3 @ 2.40GHz [6c/12t], 48 GB RAM, cooperating with a Nvidia Titan Xp GPU with 12 Gb and 3840 CUDA Cores. Further details of the NVIDA GPU card can be found in Table 3.   were observed also for the other two cities of Munich and San Francisco.
This clearly shows a dependency of the computation time over frequency and flying altitude of the UAV. As a matter of fact, the lower the frequency the longer the computation time, which is seen more than halved when moving from 700 MHz to 70 GHz. Conversely, the higher the flying altitude, the longer the simulation time although simulation time keeps basically steady from a certain height of the UAV onwards.
In addition to the isotropic case, the investigation has been carried out also in case of directive antennas. Replacing the isotropic antenna with a directive one, having namely an aperture angle θ = 40 • , 60 • and 80 • -to take advantage of the most significant runs already calculated in the previous sections -leads to more realistic study cases. Fig. 15 shows the same trend as in Fig.14, i.e., computation time as a function of drone altitude, but this time the results have been highlighted with regards to the aperture angle θ at 3.5 GHz instead of the whole frequency span.   This figure shows an interesting connection between directivity and computation times: the wider is the antenna aperture θ or the transmitting frequency, the longer the computation time.
Both situations described in Fig. 14 and 15 are in good agreement with the research activities introduced in [46] under the assumption that there are more rays that propagate at lower frequencies, for low directivity values or for higher altitudes of the drone: as reported in section II, the antenna footprint shrinks down for increasing directivity of the on-board antenna. The explanation of this trend relays on the ray-launching algorithm and its strong relationship between the number of rays that are launched -and kept ''alive'' along the propagation way -and the system parameters, configured under these specific circumstances.

VIII. CONCLUSION
In this paper, air-to-ground narrowband propagation is investigated over three different reference urban scenarios by means of a fast ray-launching propagation tool based on Nvidia GPU acceleration. Probability of LOS, path loss exponent, shadowing characteristics as well as the impact of different propagation mechanisms have been explored for several drone hovering altitudes, directive and non-directive antennas as well as transmitting frequencies, thanks to the computation agility of the DED-RL software tool.
Radiation properties of the on-board antenna turned out to greatly affect line of sight probability, which steadily increases with the flight altitude only in case of a large radiation lobe. Conversely, increasing the antenna directivity corresponds to a somehow convex relationship between lineof-sight probability and UAV height.
Narrowband analysis has proved that both path loss exponent and shadowing standard deviation decreases with both the UAV altitude and the on-board antenna aperture. Path loss exponents turn out to be remarkably greater than 2 to for flying heights up to 150m.
Diffraction appears as the leading propagation mechanism at lower frequencies, while scattering at higher frequencies for the larger aperture of the UAV antenna, whereas most of the power is likely to be conveyed through the direct path for greater directivity. Anyway, the intensity of the direct field hardly exceeds the 50%-60% of the total power, thus meaning multipath propagation is not negligible on the average.
In this regard, next steps in the research will be related to further simulations of additional dense urban environments, extending the ray launching tool to wideband simulations -to investigate the trend of the power delay profiles and RMS delay spreads -and updating the tool code to the latest Nvidia GPU architectures for possible enhancements. Field measurements, which are currently ongoing only for specific cases, will try to check simulation prediction with experimental data. MARINA BARBIROLI received the Laurea degree in electronic engineering and the Ph.D. degree in computer science and electronic engineering from the University of Bologna, in 1995 and 2000, respectively.
Since 2000, she has been an Associate Professor with the University of Bologna. Her research interests are on propagation models for mobile communications systems, with focus on wideband channel modeling for 5G systems. Her research interests include investigation of planning strategies for mobile systems, broadcast systems and broadband wireless access systems, analysis of exposure levels generates by all wireless systems, and for increasing spectrum efficiency. The research activity includes the participation to European research and cooperation programs (COST 259, COST 273 COST2100, COST IC004, and COST IRACON) and the European Networks of Excellence FP6-NEWCOM and FP7-NEWCOM++. In 1998, he joined the University of Bologna, where he is currently an Associate Professor in electromagnetic fields. He has authored more than 60 scientific publications on peer-reviewed international journals and more than 100 scientific publications on proceedings of international conferences. His research interests include nonlinear microwave circuit simulation and design, with emphasis on nonlinear/electromagnetic co-design of integrated radiating subsystems/systems for wireless power transfer and energy harvesting applications. He is a member of the Paper Review Board of the main Journals of the microwave sector. He serves on the Editorial Board of the International Journal of Antennas and Propagation, the Cambridge journal of Wireless Power Transfer, IEEE ACCESS, and Electronics Letters.
FRANCO FUSCHINI received the M.Sc. degree in telecommunication engineering and the Ph.D. degree in electronics and computer science from the University of Bologna, in March 1999 and in July 2003, respectively. He is currently an Associate Professor with the Department of Electrical, Electronic and Information Engineering ''G. Marconi,'' University of Bologna. He is the author or co-author of more than 30 journal articles on radio propagation and wireless system design. His main research interests include radio systems design and radio propagation channel theoretical modeling and experimental investigation. In April 1999, he got the 'Marconi Foundation Young Scientist Prize' in the context of the XXV Marconi International Fellowship Award. VOLUME 9, 2021