A Unified Spatial Framework for UAV-aided MmWave Networks

In this paper, we propose a unified three-dimensional (3D) spatial framework in order to evaluate the average performance of unmanned aerial vehicle (UAV) aided networks with millimeter wave (mmWave) communications. More specifically, the locations of transceivers in downlink and uplink are modeled through the homogeneous Poisson point processes and Poisson cluster processes (PCPs), respectively. For PCPs, Matern cluster and Thomas cluster processes, are analyzed. Furthermore, both 3D blockage processes and 3D antenna beamforming patterns are introduced for appraising the effect of altitudes. Based on this unified framework, several closed-form expressions for the coverage probability in the uplink as well as in the downlink, are derived. By investigating the entire communication process, which includes the two aforementioned phases and the cooperative transmission between them, tractable expressions of system coverage probabilities are derived. Next, three practical applications in UAV networks are provided as case studies of the proposed framework. The results reveal that the impact of thermal noise and non-line-of-sight mmWave transmissions is negligible. In the considered networks, mmWave outperforms sub-6~GHz in terms of the data rate, due to the sharp direction beamforming and large transmit bandwidth. Additionally, there exists an optimal altitude of UAVs, which maximizes the system coverage probability.

(LOS) connections to users by adjusting their locations and altitudes [4]. Under this condition, mmWave bands become the best choice for UAV-aided networks due to its larger available bandwidth and consequent higher data rate transmission in LOS than sub-6GHz [7]. Moreover, with the aid of large antenna scales deployed at mmWave devices, sharp directional beams can be generated to increase network capacity [8] and mitigate the Doppler spread [4].
In order to analyze the average performance of wireless networks, stochastic geometry has rekindled the strong interest of academia [9]. More specifically, the locations of transceivers are distributed according to different random spatial point processes, which are capable of statistically depicting the nature of wireless systems [10]. One popular point pattern is Poisson point process (PPP), which has been widely adopted in recent works for modeling network nodes such as BSs and multiple pieces of user equipment [8], [10], [11]. The main advantage of this pattern is that all nodes are uniformly distributed in a plane. However, in practical scenarios, most wireless systems have a clustered property. For example, in device-to-device (D2D) communications or cognitive networks, users with the same demands are frequently gathered to form a cluster and the PPP model fails to accurately describe this property. As a result, another point pattern named Poisson cluster process (PCP), which includes multiple clusters, becomes increasingly attractive in the literature [12], [13].

A. State-of-the-Art and Motivation
The evaluation of the performance of UAV networks plays a vital role in proposing new relative protocols and designing efficient network structures. Accordingly, Al-Hourani et al. [14] provided a theoretical approach to obtain the optimal altitude of UAVs, with an aim to maximize the coverage on the ground. However, this primary work was based on deterministic locations. Unfortunately, acquiring the accurate location information is a difficult and costly task. Therefore, researchers turned to analyze randomly deployed spatial models for UAV networks. Mozaffari et al. [15] investigated the average coverage probability and throughput in a UAVaided network with underlaid D2D communications, where the D2D users were modeled as a homogeneous PPP and the downlink users were uniformly distributed in a finite area forming a variant of PPP, known as binomial point process (BPP). However, this article only considered a single-UAV scenario. Chetlur and Dhillon [16] extended this work to a multi-UAV case, where the locations of UAVs obey a BPP.
With the aid of this framework, the coverage performance of downlink transmissions was analysed. Recently, Zhang and Zhang in [17] introduced a three-dimensional (3D) PPP incorporating the adjustable altitude, in order to model the locations of UAVs. Analytical expressions for the coverage probability of drone small-cells were derived.
In terms of the frequency band utilized for UAV networks, the aforementioned research contributions [14]- [17] focused on the sub-6 GHz scenarios. However, when considering mmWave scenarios, blockage environment, antenna patterns, and fading channel model should be carefully modified. Instead of using the ray tracing method as discussed in [18], Bai et al. [19] leveraged the random shape theory to propose a mathematical blockage model for high frequencies. This model was approximated by a fixed line-of-sight (LOS) disc in [8], which had acceptable accuracy with high calculation efficiency. By invoking the altitude information of devices, the conventional two-dimensional (2D) blockage model was extended to the 3D case with the aid of Rayleigh distributed buildings [20]. In addition to the widely used omnidirectional antenna pattern in sub-6 GHz, Zhu et al. [21] introduced a 3D sectorized beamforming pattern for mmWave antennas to depict the sharp directional beam. Furthermore, due to the huge difference between LOS and non-LOS (NLOS) mmWave transmissions, the Nakagami-m fading channel was preferable in mmWave-enabled networks. This small-scale fading model was considered in [8], [11], [13], [21], [22].
In UAV-aided networks, for the communications between terrestrial BSs and UAVs, the locations of transceivers can be modeled as two independent PPPs to capture the randomness of the networks. However, for communications between UAVs and users, in most of the cases one UAV is dedicated to serve a cluster of users with the same requirements. To this end, PCP becomes a popular pattern for modeling these scenarios. Furthermore, Ganti and Haenggi [12] first proposed PCP to model the node locations in clustered wireless ad hoc networks. In this work, two specific models, i.e., Matern cluster processes (MCPs) and Thomas cluster processes (TCPs), were provided. Then, several networks were studied with these two types of PCPs, e.g. cognitive networks [23], heterogeneous networks [24], and device-to-device (D2D) networks [13]. To the best of our knowledge, the research in the context of UAV networks by utilizing PCP is still in its infancy.
When considering the altitude of UAVs and the cluster property of users, a unified 3D spatial framework including both uplink and downlink phases is needed, which can be effortless adjusted to fit most practical applications. Additionally, the evaluation of average coverage performance in UAV networks, especially for mmWave scenarios, is still at the very early stage. These two factors are key motivations for this article.

B. Contributions and Organization
We propose a unified spatial framework for UAV-aided networks, where a UAV first collects messages from one user in an uplink phase and then it retransmits the desired information to the requiring BS in a downlink phase. The transceivers in the downlink phase and uplink phase are distributed with the aid of a PPP and a PCP, respectively. Particularly, two commonly used PCPs, namely TCP and MCP, are investigated, and mmWave communications are considered as well. In summary, the main contributions of this paper are as follow: • We propose a novel 3D analytical framework for UAVaided networks by invoking a tractable 3D blockage model and a 3D sectorized antenna pattern, in which the altitude of all transceivers and the cluster property of users are jointly considered. Additionally, both uplink and downlink phases are assessed in order to model the entire transmission process of UAV-aided communications. • We investigate the distributions of communication distances in the downlink phase. Based on these distributions, closed-form expressions for the coverage probability are derived to enhance the evaluation efficiency.
Regarding the uplink phase, we derive a general expression for coverage probabilities. Moreover, this expression is able to characterize TCP as well as MCP. Note that the number of users in each cluster can be a constant or to follow a Poisson distribution depending on different network requirements. We theoretically demonstrate that these deployment scenarios perform similarly when the number of interferers is large. • We provide a practical cooperation transmission to model the behaviours of UAVs from the uplink phase to the downlink phase. Three typical applications of UAVs are examined with the aid of the unified framework, and thus validating the flexibility of our work. In order to study the performance of the entire communication process, we also derive tractable expressions for the system coverage probability. • We show that, 1) when the geographical size of clusters and the interfering number of users are small, the inter-cluster interference can be ignored to simplify the analysis; 2) the effect of thermal noise is negligible in the proposed UAV networks. In terms of mmWave communications, NLOS transmissions can be ignored as well, especially in a low-density blockage environment; 3) there exists an optimal altitude of UAVs for achieving the maximum system coverage probability; and 4) a large antenna scale is able to enlarge its main beam gain and narrow the beamwidth for compensating path loss, thereby enhancing coverage performance. The rest of this paper is organized as follows: Section II proposes the network model including the PPP distributed downlink phase and the PCP distributed uplink phase. Section III first describes the coverage performance in the proposed UAVaided networks. Then, three practical applications are studied with the aid of the unified framework. Section IV provides the validating and numerical results. Section V presents the conclusions.

II. NETWORK MODEL
We introduce a unified spatial framework for UAV-aided networks by integrating mmWave communications, which describes a general case where a typical BS is located far from the  served users, and hence, UAVs are deployed to establish the desired communication link. Due to the mobility of UAVs, we assume that the serving UAV is able to fly from the user region to the BS region. Therefore, one corresponding UAV receives the message from a transmitting user and then it flies toward the typical BS to forward this message 1 .In addition to the flight assumption, more practical cooperation transmission scenarios are discussed at the end of this section. The entire transmission process can be divided into two phases: 1) Uplink Phase, where the transmitting user uploads the desired message to the corresponding UAV; and 2) Downlink Phase, where the typical BS downloads the required message from the corresponding UAV. These two phases are processed in different time slots. The details are presented in the following.

A. Spatial Distribution
Different spatial distributions are assumed for the two phases. For ease of understanding, we first introduce the downlink and then the uplink phase.
1) Downlink Phase: In this phase, we consider multiple transceivers. To capture the randomness of networks, we assume that macro BSs and UAVs are modeled as two independent PPPs with density λ b and λ down v , denoted by Φ b and Φ down v , respectively. Regarding the altitude of each device, UAVs are hovering at a height of h down v , while the receiving antennas of macro BSs are located at an altitude h b , as illustrated in Fig. 1.
2) Uplink Phase: For the uplink phase, since users with similar requests frequently stay together, we utilize PCPs to model this cluster property [12]. In one PCP, parent points are distributed following a homogeneous PPP Φ up v with density λ up v . Around each parent point at v ∈ Φ up v (v ∈ R 2 ), daughter points (denoted by U v ) are independently and identically distributed (i.i.d.) forming a cluster. The number of the daughter points can be constant, N , or a Poisson random variable with meann for different purposes. Based on this spatial model, we assume that parent and daughter points represent locations of UAVs and users, respectively.
Note that TCP and MCP are two widely used patterns of PCPs. The main difference between them is the distribution of daughter points in each cluster. For TCP, these points are modeled as a symmetrical normal distribution around the central parent point, with a standard deviation σ. For MCP, the daughter points are uniformly distributed in a disc with radius R. Therefore, the two density functions for the distance between a user at u ∈ U v (u ∈ R 2 ) to the above UAV at v can be expressed as where the superscript T ho represents TCP and M at represents MCP. The U(.) is the unit step function, which is given by Regarding the altitude in the uplink phase, we assume that all users are located at the ground with an altitude h u and each UAV hovers above one cluster center at an altitude h up v to serve the intra-cluster users.

B. Blockage Model
The authors in [20] have proposed a tractable blockage model. The extracted theoretical expressions are derived based on a three dimensional obstacle environment, where the average height and the density of obstacles are adjustable. Therefore, this model is suitable for various practical scenarios, e.g., the suburban with low-density and low-altitude obstacles, the urban with dense and high buildings, and so forth. In this model, the density of obstacles and the ratio of the obstacle area to the total area are represented by β b m −2 and β a , respectively. The height of each obstacle is modeled as a Rayleigh distribution with a scale parameter ε. When  the horizontal transmission distance is r, then the probability density function (PDF) of a LOS link is given by [25] p L (γ |h t , h r ) = max(0,γ) where γ = r √ β a β b and ⌊.⌋ is the floor function. The subscript t and r represent transmitters and receivers, respectively 2 . Intuitively, the PDF of a NLOS link is p N (.) = (1 − p L (.)). As a result, the path loss law in the UAV network can be expressed as where B(x) is a Bernoulli random variable with parameter of success probability x. The parameters C κ and α κ are the considered intercept and path loss exponent, respectively. Additionally, κ ∈ {L, N }, where L represents LOS links and N denotes NLOS links.

C. Directional Beamforming
We consider the typical sectorized antenna pattern, as mentioned in [26] at all transceivers to accomplish the sharp directional beamforming (see Fig. 2). The numbers of antenna elements at macro BSs, UAVs and users are N b , N v and N u , respectively. Four main characteristics of antennas are introduced in this pattern, which are the half-power beamwidth in the azimuth plane θ a c , the counterpart in the elevation plane θ e c , the main beam gain M c , the side lobe gain m c (c ∈ {t, r}). Based on the hierarchical beam search scheme [27], the location information of devices can be obtained at the transmitter. Then, it adjusts the antenna direction toward the corresponding receiver in order to achieve the maximum beamforming gain G 0 . Regarding one interfering transmission, the angles deviating from the boresight direction in the azimuth plane and the elevation plane are ψ c and ϕ c , respectively. We assume that ψ c is uniform in the range [0, 2π) and ϕ c is uniformly distributed over the range [0, π] for all interfering transmissions [21]. The directivity gain at one receiver with the interfering transmitter located at l can be expressed as follows: 2 In the downlink phase, the transmitters are UAVs and the receivers are BSs, while in the uplink phase, t and r represent users and UAVs, respectively. We omit this explanation in the rest of this paper. where G(θ a c , θ e c , M c , m c ) denotes the directional antenna gain. Therefore, G l has four patterns as shown in Table. I. Each pattern has the value o i with the probability p i , where i ∈ {1, 2, 3, 4} and G 0 = o 1 .

D. Signal Model 1) Downlink Phase:
We assume that the number of BSs is greater than that of UAVs, namely λ b > λ down v , such that all UAVs in the downlink phase are active at the considered time slot. The typical BS at b 0 (b 0 ∈ Φ b , b 0 ∈ R 2 ) is randomly selected and it is fixed at the origin of the downlink plane. To achieve the best quality of service, the corresponding UAV, which is located at , is the closest UAV to the typical BS. Since the interfering signals are offered by the remaining UAVs, the received SINR at the typical BS can be expressed as where P v is the transmit power for each UAV,ĥ is the small scale gain for Nakagami fading channel with parameter N κ such that |ĥ| 2 is a normalized Gamma variable, and n 0 is the thermal noise.
2) Uplink Phase: Since the downlink and uplink phases are analyzed in two independent planes, the corresponding UAV can belong to any of the clusters in the uplink plane. We assume that this UAV is located at the origin of the uplink plane and its location is denoted by v 0 (v 0 ∈ Φ up v ). In one time slot, we randomly select an intra-cluster user at u 0 (u 0 ∈ U v0 ) to be the transmitting user. In contrast to the downlink phase, the interference in the uplink phase is originated by both intracluster and inter-cluster users. Then, the received SINR at the corresponding UAV is given at the top of next page. In (8), P u is the transmit power for each user.

III. PERFORMANCE EVALUATION FOR UAV-AIDED NETWORKS
In this section, we first evaluate the performance of UAVaided networks with the aid of the proposed unified framework. Then, we provide several practical scenarios which can be modeled by the proposed framework after minor adjustments.

A. Performance of the Downlink Phase
In this part, we analyze the coverage performance of downlink transmissions based on the distribution of distances 3 and Laplace transform of interference.

1) Distance Distributions in Poisson Point Processes:
Assuming that the distance realizations of UAVs in the downlink phase form a set r n = v n n=1,2,...,N dl , where v n ∈ Φ down v and N dl is the number of UAVs. The subscript n is the parameter of ranked distances, namely r 1 < r 2 < ... < r N dl . Note that the corresponding UAV is the closest node to the reference BS. Therefore, the PDF of the nearest communication distance r 1 is given by [10] In addition to the corresponding UAV, the rest of UAVs which are located further than r 1 are interfering transmitters. Due to the Slivnyak's theorem, the density of interfering UAVs is still λ down v in the area O(r 1 , +∞), where O(a, b) represents an annulus with the inner radius a and outer radius b [28]. Moreover, each interferer is distributed independently such that the index n can be dropped from r n . Therefore, the PDF of the interfering distance r can be derived via the probability generating function of PPP [29].
2) Laplace Transform of Interference in Poisson Point Processes: For simplifying the notation, we define the Laplace transform of interference as follows: where E [.] is the expectation function, s is the transform parameter, and I represents the received power of interference, which is given by Next, we provide closed-form expressions for the Laplace transform of interference in the following lemma and corollary.
Lemma 1. Since the corresponding UAV at v 1 is the closest node to the typical BS with a communication distance r 1 = v 1 , all interfering UAVs are located in the area O(r 1 , +∞), namely r > r 1 . A closed-form expression for the Laplace transform of interference in the downlink phase conditional on r 1 is given by where and when α κ > 2, F κ ακ (.) can be expressed as with 2 F 1 (., .; .; .) being Gauss hypergeometric function. For most of mmWave frequencies, the path loss exponent of LOS transmissions equals to two [30]- [32]. Therefore, when Proof. See Appendix A.

3) Coverage Probability for the Downlink Phase:
In order to insure the quality of the downlink communicating service, a targeted rate R th down is pre-decided in most of the cases. Then, the corresponding desired SINR threshold is given by where B down is the bandwidth for each downlink resource block. Therefore, the coverage probability is defined as the proportion of received SINR that exceeds Υ th down , which can be expressed as where P[.] denotes probability. With the aid of the aforementioned distance distribution and Laplace transform of interference, we present the coverage probability in the following theorem.
Theorem 1. When choosing the nearest UAV at v 1 as the corresponding UAV, then the serving communication distance is Accordingly, a closed-form expression for the coverage probability at the typical BS is given by where and η κ = N κ (N κ !) −1/Nκ . The set of non-negative integers Z * = {0, 1, 2, ...}, and ζ is a Gauss-Chebyshev node, which equals cos 2k−1 2m π | k=1,2,...,m . When m → ∞, the equality holds. We are able to change the value of m to balance the complexity and efficiency [10].
Proof. See Appendix B.
Remark 1. Since we randomly select the typical BS, the coverage provability in the downlink phase is independent of the distribution of the macro BSs Φ b .
Remark 2. Although the range of γ is infinite due to γ ∈ Z * , we are able to choose the first three values, namely γ = 0, 1, 2, for enhancing the computation efficiency. The reason is that when the communication distance increases, the received power of both signals and interference decreases, which have a negligible impact on the received SINR.
Remark 3. Note that the point process of UAVs is same with that of macro BSs. When considering the uplink transmissions from the BSs to UAVs in the proposed downlink phase, the expression of coverage probability in Theorem 1 is also valid.

B. Performance of the Uplink Phase
Due to the cluster property of PCP, the analysis of the uplink phase is more challenging compared with the downlink phase. Similarly, we first present the distance distribution and Laplace transform of the interference in the uplink phase. Then, the coverage performance at the corresponding UAV is analyzed.

1) Distance Distribution in the Corresponding Cluster:
The communication distances in the same cluster with the corresponding UAV form a group w n = u v0 n n=1,2,3,...,Nuser and u v0 n ∈ U v0 , where N user is the number of users in the considered cluster. By utilizing the fair selection strategy, the subscript n can be removed from w n and u v0 n . Regarding TCP, the PDF of the distance (w ≥ 0) between one user in the corresponding cluster to the corresponding UAV is given by [13] f T ho Then, for MCP, the PDF of that distance can be expressed as [10] f Mat 2) Distance Distributions in Other Clusters: In other clusters, the distances between the users with a cluster center v and the corresponding UAV compose a set g n = u v n n=1,2,3,...,Nuser and u v n ∈ U v . As all users are i.i.d, we are able to drop the subscript n from g n and u v n . For TCP, the PDF of this distance (g ≥ 0) is conditional on the distance q = v ≥ 0, which is given by [13] f T ho where I 0 (.) is the first kind for the modified Bessel function with order zero. Then, for MCP, the PDF of this distance is [33] f Mat Q (g |q ) = 2g πR 2 arccos Remark 4. In fact, eqs. (22) and (23) are the PDF of the distance between a daughter point to the origin. The origin can be anywhere on the plane. For example, (22) and (23) are also valid for D2D communications, where one daughter point is fixed at the origin instead of one parent point as mentioned in the proposed system.

3) Laplace Transform of Interference in Poisson Cluster
Processes: There exist two kinds of interference in the uplink phase. One is intra-cluster interference I intra and the other is inter-cluster interference I inter . As a result, Laplace transform of interference can be defined as where In several special networks, e.g., D2D communications and non-orthogonal multiple access (NOMA) networks, the number of users in one cluster should be fixed as they are paired. For other networks, the number of users should be random across different clusters to enhance the generality. Therefore, we consider two cases in this paper: 1) Fixed Case: the number of users in each cluster is fixed as N ; and 2) Random Case: the number of users in each cluster follows Poisson distribution with the meann. Lemma 2. When the distance between an intra-cluster interferer and the corresponding UAV is w, then the Laplace transform of the interference in the fixed case is given by where and ∆h u = h u − h up v0 , γ 1 = w √ β a β b , and Ω ∈ {T ho, M at}.
On the other hand, if we consider the random case, the Laplace transform of the intra-cluster interference can be expressed as L a (s) = exp (− (n − 1) (1 − O a (s))) .
Proof. See Appendix C.
In terms of the inter-cluster interference, the Laplace transform of this interference can be expressed in the following lemma and corollary. Lemma 3. When the distance between an inter-cluster interferer and the corresponding UAV is q, for the fixed case, the Laplace transform of inter-cluster interference can be expressed as and for the random case, the corresponding Laplace transform of inter-cluster interference is changed to q))))qdq , where O e (s, q) Proof. See Appendix D.

Remark 5. Note that when
if the number of users in each cluster is large and N =n, the difference between Laplace transform of interference in the fixed case and the random case is negligible.

4) Coverage Probability for the Uplink Phase:
We predecide a targeted rate R th up for the uplink phase such that the corresponding SINR threshold is given by Υ th up = 2 R th up /Bup −1 and B up is the bandwidth for each uplink resource block. Then, the coverage probability for the uplink phase is defined as Note that in the corresponding cluster, the transmitting user is randomly selected. The communication distance between the corresponding UAV and the transmitting user is w 0 = u 0 . Based on eqs. (8) and (33), the analytical expression for the uplink coverage probability can be expressed in the following theorem and corollary.
Theorem 2. When the SINR threshold for the uplink phase is Υ th up , then the coverage probability at the corresponding UAV is given by where and Proof. In terms of coverage probability, there exist two differences between the downlink phase and uplink phase. One is the distribution of communication distance and the other is the received interference power. Based on Theorem 1, we replace (9) by (20) and (21). Additionally, we utilize L a (.)L e (.) instead of L down (.). After that, we obtain Theorem 2.
By analyzing Theorem 1 and Theorem 2, we are able to abstract a general expression of the coverage probability with an SINR threshold Υ th for all point processes, including PPP and PCP. Proposition 1. Assuming that Nakagami-m fading channel is utilized in the considered wireless networks, then the coverage probability for both PPP and PCP can be expressed as where the index g represents the general case, P noise (.) is the signal-to-noise-ratio (SNR) coverage probability part, L I (.) is the Laplace transform of interference, and f Rg (.) is the PDF of communication distance r g between the serving transmitter and the served receiver.

C. System Coverage Probability
Without loss of generality, we assume that the probability of the message being successfully transmitted from the uplink to the downlink phase is P s . In the proposed system, we assume that the corresponding UAV is able to successfully flies from the uplink plane to the downlink plane, such that the messages can be guaranteed to reach the downlink phase, namely P s ≡ 1. Therefore, we are capable to deduce the system coverage probability for the proposed UAV networks.
Proposition 2. Assuming that the SINR thresholds for the downlink phase and uplink phase are Υ th down and Υ th up , respectively, then the system coverage probability is given by

D. Practical Applications
Next, we provide three typical UAV applications: UAVaided ubiquitous coverage, UAV-aided information dissemination and data collection, and UAV-aided relaying [2]. The proposed unified framework is able to model them after minor modifications. As it can be seen from Fig. 3, for the ubiquitous coverage case, since the corresponding UAV aims to offload the data traffic of the central BS, all transceivers should be modeled in one plane. For the information dissemination and data collection case, the proposed framework can be utilized directly. For the relaying case, since two corresponding UAVs act as two relays, the success probability P s is mainly decided by their communication conditions. In order to make the analysis more complete, we introduce a simplified cooperative transmission between two relays in the following part, which can be extended to other complicated scenarios. Moreover, the time slot for this transmission is located between the uplink and downlink phases. Special modifications for modeling these three applications are summarized in Table II.

1) Simplified Cooperative Transmission:
In most of the cases, two relaying UAVs are connected through the air. Due to limited obstacles in the air, these links can be regarded as LOS transmissions. Assuming that the distance between the two relays is fixed, y 0 = v 1 − v 0 and for simplicity, the cochannel interferers can be ignored, then the received SNR at the receiver is given by In this cooperative transmission, P s represents the corresponding coverage probability, which can be calculated in the following lemma.
Lemma 4. When the considered cooperative transmission distance is fixed as y 0 , the closed-form coverage probability with an SNR threshold Υ th link can be expressed as follows: and Γ(., .) is the upper incomplete gamma function.
Proof. By using the cumulative distribution function (CDF) of the Gamma distribution and following similar proof as in Theorem 1, we obtain Lemma 4.

2) Coverage Probability for the Relaying Case:
The system coverage probability for the relaying case can be represented by Proposition 2 as well. The only difference is the probability P s . Therefore, the system coverage probability is expressed as follows: Proposition 3. Assuming that the SINR thresholds for the downlink phase, uplink phase, and cooperative transmission are Υ th down , Υ th up , and Υ th link respectively, the system coverage probability for the relaying case with a fixed cooperative transmission distance y 0 is given by Remark 7. Note that the LOS transmissions outperform NLOS transmissions in mmWave communications. When the altitude of UAVs increases, the PDF p L (.) in (4) enlarges, which enhances the system coverage probability P sc . However, the high altitude of UAVs also increases the communication distance such that the path loss L(.) in (5) raises. Therefore, P sc can be maximized by selecting optimal altitude of UAVs.

IV. NUMERICAL RESULTS
Since the used blockage model utilizes the bands from 3 to 60 GHz [20], the results can be extended to sub-6 GHz scenarios. The main difference between sub-6 GHz and mmWave is the antenna pattern and channel model. Omnidirectional antennas and Rayleigh fading are commonly used in sub-6 GHz communications [14]- [17]. More specifically, for the omnidirectional antenna, we assume that  and M c = m c = 0 dB. For Rayleigh fading, only NLOS transmissions are considered, namely p N (γ |h t , h r ) ≡ 1 and N N = 1. In this section, sub-6 GHz scenarios are regarded as the benchmark.

A. Simulations and Discussions
We use Monte Carlo (MC) simulations, which include noise, LOS and NLOS transmissions, to appraise the accuracy of derived expressions. The antenna pattern is modeled by a uniform planar square array (UPA) (see Table III) [26]. General network settings are illustrated in Table IV [8], [13]. Additionally, we assume that P u = P v = 1 W and P b = 10 W for simplicity. In terms of the intercept C κ , we consider the reference distance with one meter. The validation of the derived expressions are illustrated in Fig. 4, with N L = 2, N N = 1, β a = 0.2, and β b = 10 × 10 −6 . As shown in Fig. 4(a), the theoretical results for the downlink phase from Theorem 1 match MC simulations perfectly. Regarding the uplink phase in Fig. 4(b), expressions in Theorem 2 are the tight approximations of the exact MC simulations. Particularly, the deviation is negligible in high coverage probability regions.

B. Comparison and Analysis of MCP and TCP
Regarding the clustered property of the proposed framework, we compare two typical PCPs (MCP and TCP) in Fig. 5. More specifically, Fig. 5(a) demonstrates that the clusters in MCP are appropriate to model bounded regions, while the counterparts in TCP are suitable for open-boundary scenarios. In addition, Fig. 5(b) shows that the coverage probability in the uplink phase has a negative correlation with the number of users in one cluster. In terms of intra-cluster and intercluster interference, exact results can be approximated by the simplified scenario that only considers intra-cluster interference, which significantly improve the analysis efficiency. When σ increases from 50 to 150, the effect of inter-cluster interference is gradually enhanced. Moreover, when R = σ, coverage probabilities with MCP are higher than those with TCP in small N regions. Although the number of users in each cluster can be a constant N or a variable with the mean n, coverage probabilities for two cases are nearly the same under the condition N =n in high N regions, which indicates that two cases can be replaced by each other as discussed in Remark 5.

C. The Impact of Blockage Environment and Noise
In this part, we analyze the system coverage performance affected by the blockage environment and the thermal noise, with N = 30, h v = 100 m, Υ th up = Υ th down = Υ th link = −20 dB, and the bandwidth for all scenarios equals 100 MHz. Since MCP has similar trends with TCP and thus we only study TCP here. Regarding blockage effects, Fig. 6(a) shows that an optimal altitude of UAVs exists for maximizing the system coverage probability as discussed in Remark 7. It is interesting that all optimal values are larger than the height of traditional BSs, which is around 30 m (dash line in Fig. 6(a)). Therefore, UAV networks are capable of achieving better performance than terrestrial cellular networks by adjusting the serving altitude. When the density of UAVs increases from λ up v = λ down v = 1/(250 2 π) to 3/(250 2 π) m −2 , the optimal altitude rises. Moreover, the hight density of obstacles β a also enlarges the optimal height. In terms of signal-to-interferenceratio (SIR), Fig. 6(b) demonstrates that the thermal noise has a negligible effect on system coverage probabilities in both mmWave and sub-6 GHz scenarios. NLOS transmissions can

Cluster Centers Matern Processes Thomas Processes
(a) Comparing the spatial difference between MCP and TCP, with the fixed number of daughter points in each cluster N = 10. be ignored in mmWave communications in the considered framework. Additionally, both figures in Fig. 6 illustrate that mmWave outperforms sub-6 GHz.

D. The Impact of Antennas and Carrier Frequencies in MmWave Scenarios
Since sub-6 GHz can be regarded as one special case of mmWave scenarios, we study mmWave scenarios in this part to comprehensively evaluate the proposed framework. Fig. 7(a) concentrates on the impact of antenna scales. For UPA, the large number of antenna elements is capable of increasing the main lobe gain and narrowing the main lobe beamwidth. Therefore, the desired signal is enhance and the interference is weakened. The system coverage probability has a positive correlation with the number of antenna elements N c as shown in Fig. 7(a). Then, we focus on different carrier frequencies. It is worth noting that compared with 28 GHz, higher carrier frequencies allows more antenna elements to be deployed at transceivers. We provide the path loss laws and the estimated number of antennas of three typical carrier frequencies for mmWave communications in Table V [31], [32]. The comparison of those candidates are illustrated in Fig. 7(b). When the number of antennas is fixed as two, namely N c = 2, it is obvious that 28 GHz, with the best path loss law and the intercept, achieves the highest coverage probabilities in low SINR threshold region, while 38 GHz outperforms other two candidates in high SINR threshold area. With the increasing of the antenna scale, 60 GHz becomes the best among three candidates.

V. CONCLUSION
In this article, a unified 3D spatial framework for UAV networks has been provided, where the stochastic geometry has been utilized for modeling the locations of BSs, UAVs and users. Especially, in the uplink phase, two typical PCPs (MCP and TCP) have been analyzed to enhance the generality. Tractable expressions in terms of coverage probabilities have been deduced for mmWave communications, which can be extended to sub-6 GHz scenarios. For the proposed system, there exists an optimal altitude of UAVs to achieve the maximum coverage probability. Another remark is that the effects of the thermal noise and NLOS transmissions are negligible in mmWave-aided UAV networks. Moreover, large number of antenna elements have the capability of enhancing the coverage performance. Without considering the antenna scales, 28 GHz is the best choice for mmWave scenarios in low SINR regions.
where (f) depends on step property of (4). (g) calculates the expectation of the variable o i for antenna gain. Then, with the similar process, Laplace transform of interference via NLOS links Ξ N is given by For LOS links, the corresponding coverage probability (B.2) can be expressed at the top of next page. In (B.2), (a) follows the fact that the serving UAV is the closest node to the reference BS. (b) depends on a tight upper bound for the normalized gamma variable |ĥ v1 | 2 , which is P |ĥ v1 | 2 < χ ≤  [10]. With the similar deriving procedure, the coverage probability for NLOS transmissions is given by If this number follows a Poisson distribution with the mean n, the number of interfering devices, in this case, has the mean (n − 1). Therefore, (C.1) is changed to L a (s) = Regarding mmWave scenarios, the Laplace transform of inter-cluster interference can be derived with the aid of (24) and (27), which is given by Oe(s,q) , (D.1) where g = u + v and q = v . (a) follows the similar procedure of (a) and (b) in Appendix C. The difference is the distance distributions in this case obey (22) and (23). If the number of users in each cluster is fixed as N , the number of interfering devices in every other cluster is same with N . Therefore the Laplace transform of inter-cluster interference can be expressed as where (b) follows the probability generating function of PPP [13], [29].
In terms of the Poisson distributed number of users in one cluster, the average number of interfering users in each cluster isn. (D.1) can be further deduced as Oa(s) . (C.1)