hidelinks

# Improved Model and Analysis for RIS-Assisted Indoor Terahertz Wireless
Networks

Zhi Chai, Jiajie Xu, , Mohamed-Slim Alouini, , Justin P. Coon Z. Chai, and J.
P. Coon are with the Department of Engineering Science, University of Oxford,
Parks Road, Oxford, OX1 3PJ, United Kingdom. (e-mail<EMAIL_ADDRESS>justin.coon@eng.ox.ac.uk).
J. Xu and M.-S. Alouini are with the Computer, Electrical and Mathematical
Sciences and Engineering Division, King Abdullah University of Science and
Technology, Thuwal, 23955, Kingdom of Saudi Arabia. (e-mail:
<EMAIL_ADDRESS>slim.alouini@kaust.edu.sa).

###### Abstract

In this paper, we propose a new model for indoor THz communication assisted by
RIS. We conduct a realistic modeling of indoor obstacles and analyze their
impact on performance. Order statistics are applied to calculate the
cumulative distribution functions (CDFs) of distances from the transmitter to
the selected RIS, i.e., the nearest RIS in the bounded indoor environment to
the transmitter, and from the selected RIS to the receiver. We calculate the
coverage probability (CP) as a function of the RIS number, obstacle density,
room size, and the transmitter’s location. By comparing the numerical results
obtained from the analytical expressions with Monte Carlo simulations, we
verify the accuracy of our analysis. Through the numerical results, it is
observed that room size and obstacle density affect the CP in a significant
way. However, by optimizing the transmitter’s location and increasing the RIS
number deployed in the room, the CP can be significantly improved (e.g., an
increase of around 15% by optimizing the transmitter’s location, and an
increase of around 30% by increasing the RIS number deployed in the room).

###### Index Terms:

RIS, THz, Poisson point process, performance analysis, stochastic geometry,
order statistics

##  I Introduction

###  I-A Motivation

With the depletion of the current commercial spectrum, the soaring demand from
consumers, and the emerging techniques such as virtual reality (VR) [1, 2, 3]
and autonomous driving [4, 5, 6], the global mobile traffic is estimated to
increase by 55% each year, and reach 5016 exabytes, with the transmission rate
growing up to 1Tbps at the end of 2030 [7, 8]. Hence, more attention is paid
to the terahertz (THz) spectrum due to the substantial bandwidth, the ultra-
fast transmission rate [9, 10], and the ultra-low latency[11, 12]. THz
communication also has a wide range of applications such as wireless local
area networks (WLAN), wireless data centers, and space communication networks
[13]. However, a side effect of using the THz spectrum is that the signal will
experience severe path loss and be hugely attenuated by scatterers, which
leads to limited coverage. Reconfigurable intelligent surface (RIS)
technology, as widely studied, is regarded as one of the most potential
solutions to enhance THz communication’s coverage. RIS is defined as a
synthesized surface that has many antenna elements with tunable devices (e.g.,
positive-intrinsic-negative diodes, varactors) mounted on top to manipulate
the phase of the impinging waves to realize a constructive superposition after
reflection [14]. The reflected signal with an enhanced amplitude can improve
the SNR to achieve broader coverage. Because of the severe path loss issue,
the majority of THz applications focus on indoor environments, such as kiosk
downloads [15] (communication range: $\sim 1\textrm{m}$), WLAN [16]
(communication range: $\sim 10~{}\textrm{m}$), and data center networks [17]
(communication range: $<20~{}\textrm{m}$). To achieve a high coverage
probability (CP), high transmission speed, and low latency, multiple RISs can
be deployed in the indoor environment to facilitate THz communications. In
this paper, focusing on indoor THz communication, with the assistance of RIS,
we analyze the CP variation with respect to the room size, obstacle density,
RIS number, and the transmitter’s location. We also maximize CP by optimizing
the transmitter’s location and the deployed RIS number.

###  I-B Related Work

In this subsection, we first give an overview of the CP analysis in THz
communications and then introduce the state-of-the-art RIS-assisted network
performance analysis.

####  I-B1 CP Analysis in THz Communications

In [18], the outage probability is explicitly derived for a transmitter-relay-
receiver system model, and beam misalignment is taken into consideration at
THz frequencies. In [19], an unmanned aerial vehicle (UAV) is used as a
decode-and-forward (DF) relay in the THz communication system. In the proposed
system, the transmit power that minimizes the outage power is derived.
Compared with the fixed power allocation scheme, the outage probability
obtained through transmit power optimization is 20% lower. Both [18] and [19]
focus on the link-level outage probability of a THz wireless communication
system, which cannot capture the outage/coverage probability of the THz
wireless communication network, and cannot provide insights into THz
applications from a statistical perspective. In [20], a relay selection scheme
is proposed in radio frequency (RF) and THz hybrid wireless communication
networks. It can be concluded that the scheme would choose nodes supporting
THz communication when the demand for transmission rate is high or the
transmitter-receiver distance is short. Otherwise, it will choose nodes
supporting RF communication. In [21], a CP analysis of the downlink is derived
in an RF/THz hybrid wireless communication network in the presence of
interference using stochastic geometry. In [22], a coverage performance
analysis towards user-centric THz networks is carried out. It is concluded
that static or dynamic clustering in user-centric networks has advantages over
the conventional single base station (BS) association policy. It is worth
pointing out that [21] and [22] do not consider relays in the system, which
may cause an overestimation of the system as the transmission distance is
limited in THz communication systems. Meanwhile, active relays, for example,
DF-based relays considered in [18, 19, 20] and [23, 24, 25, 26, 27], require
higher energy consumption, and may not be an optimal candidate in THz
communication systems.

####  I-B2 Stochastic Geometry-Based Analysis of the RIS-Assisted Wireless
Communication Network

As a passive relay, RIS is widely researched to improve communication
performance. In [28], base stations (BSs) and RISs are assumed to follow two
independent homogeneous Poisson point processes (HPPPs), respectively. It is
also assumed that the density of the BS is lower than the density of the RIS,
and there is no LOS path. The association policy proposed in [28] follows the
shortest distance principle111The shortest distance principle means the user
will associate the nearest RIS in the indirect path.. Through a comparison
with conventional networks, it is concluded that the RIS-assisted network can
achieve higher capacity and energy efficiency under the condition that the RIS
density is higher than the BS density. It is shown in [28] that when the BS
density is higher than the RIS density, it is unnecessary to deploy RIS
devices. In [29], an application scenario with multiple users, BSs, and RISs
in the presence of interference from BSs and RISs is considered, the BSs are
assumed to follow an HPPP while the RIS is assumed to follow a binomial point
process (BPP), and the Laplace Transform of the received signal power,
approximated aggregate interference, CP, ergodic capacity, and the energy
efficiency are derived based on the assumptions. The analytical results are
validated by Monte-Carlo simulations. In [30], stochastic geometry is used to
analyze the performance of a RIS-aided multi-cell non-orthogonal multiple
access (NOMA) network. It is shown in [30] that the achievable data rate
reaches an upper limit with the increase in the RIS size and the path loss
intercept can be enhanced to improve the network performance without
influencing the bandwidth. In [31], the authors propose to coat RIS on a
subset of scatterers in the environment, and the light-of-sight (LOS)
probability is taken into consideration. The authors model the locations of
BSs and users following two independent HPPPs. The scatterers follow the line
boolean model proposed in [32]. It is concluded that the deployment of RIS
highly improves the coverage regions of the BSs, and it is shown that 6
RIS/$\textrm{km}^{2}$ is enough to maintain the ratio of blind spots over the
total area below $10^{-5}$ when the scatterer density is 300
scatterer/$\textrm{km}^{2}$. In [33], the spatial correlation between
transmitters and RISs is described by using the Gaussian point process (GPP),
and the CP expressions under a fixed association strategy and a nearest
association strategy are derived in closed form in the presence of limited
interference. It is also shown that the system performance is independent of
the density of transmitters with the nearest association strategy when the
scenario is interference-limited.

####  I-B3 RIS-Assisted indoor THz wireless network CP Analysis

In [34], a ray tracing simulation was conducted for the RIS-assisted THz
indoor system. The coverage probability variation with respect to the signal-
to-noise ratio threshold was given for up to four RISs. Yet, the deployment
strategy and the association policy were not considered. In [35], stochastic
geometry was applied to an indoor THz wireless network, and the optimal RIS
density that can maximize the coverage of users in THz wireless networks was
obtained in the presence of interference. It is worth pointing out that the
RISs were assumed to attach some of the obstacles in the indoor environment.
Besides the THz spectrum, similar investigations have been conducted in the
millimeter wave (mmWave) spectrum. [36] discussed small-cell densification and
RIS deployment for coverage improvement. For the RIS solution, the deployment
of RISs was assumed on the walls of a room with the locations following the
1-dimensional HPPP, but the blockage effect was not considered. [37]
considered a humanity mobility model to calculate the outage probability in an
indoor environment and the RIS locations (assumed to be on the wall) were
optimized with respect to coverage. However, the blockage and the RIS amount
showed in the simulation section was up to 2 (e.g., 1 and 2 humans as the
blockages, and 1 and 2 RISs being deployed on the wall). Such results cannot
provide insight when the network scales up. In addition, [37] assumes multi-
hop existed between RISs222The signal can be reflected twice by the RISs
before it reaches the receiver., which is unclear currently whether RIS is
feasible for this functionality. Compare [36, 37] with [35], [36, 37] is more
practical in terms of deployment because the RISs were assumed to be deployed
on the walls of a room rather than on obstacles333The obstacles may have
irregular shapes, which is hard for RIS to be deployed on.. What’s more, both
[35] and [36] did not consider the boundary effect (See Appendix A for
details) of an indoor environment when conducting the analysis.

###  I-C Contributions

In this paper, we investigate the THz wireless communication network coverage
performance in an indoor environment with boundaries considered. The main
contributions are given below:

  * •

We propose a new RIS-assisted indoor THz communication model. We explicitly
consider the effects of boundaries. In addition, a new blockage model
considering the size and density of the obstacle is proposed, which can
describe the indoor environment more realistically.

  * •

We analytically characterize the coverage probability for the RIS-assisted
indoor THz communication network. The communication link is classified as the
direct link (transmitter-receiver) and the indirect link (transmitter- RIS-
receiver). In the indirect link, the association policy is that the nearest
RIS will be selected as the desired RIS to reflect the signal. The combined
coverage performance of the direct and indirect link is analyzed with the room
size, RIS number, obstacle density, and transmitter’s location as variables.

  * •

We maximize the CP by optimizing the transmitter’s location and the deployed
RIS number. Based on the analytical results and simulations, it is shown that
the room size and the obstacle density are two factors that dominate the CP.
However, by optimizing the transmitter’s location and increasing the RIS
number, a significant improvement (15% to 30%) to the CP can be observed.
These relationships will serve as guidance in the future RIS-assisted indoor
THz wireless communication network deployment.

Figure 1: RIS-assisted indoor THz wireless communication network illustration.

##  II System Model

The RIS-assisted indoor THz wireless communication network is illustrated in
Fig. 1. We consider a square-bounded region. The length and width of the
square room are denoted as $a$. The deployment of RISs is restricted to the
wall. It is assumed that each wall has the same number of RISs and the
location of each RIS is an independent and identically distributed (i.i.d.)
random variable following the uniform distribution $\mathcal{U}(0,a)$. The
yellow circles with a constant diameter $d_{b}$ in Fig. 1 represent obstacles.
The transmitter’s location is given, and receivers and the centers of
obstacles are assumed to follow two independent HPPPs in the room. We define a
transmitter as $T_{0}$, the nearest RIS to $T_{0}$ as $S_{0}$, and a receiver
as $R_{0}$. We apply a two-step communication strategy. $T_{0}$ will first try
to connect to $R_{0}$ via the direct link, and the communication distance is
denoted as $D_{tr}$. $D_{tr}$ can be described as the distance between a given
location (the transmitter) and a point subject to an HPPP. If the direct link
is blocked, or the SNR is lower than the threshold, $T_{0}$ will choose
$S_{0}$, and then $S_{0}$ will reflect the signal to the desired receiver on
the condition that there are no obstacles between $T_{0}$ and $S_{0}$, as well
as $S_{0}$ and $R_{0}$. For better tractability, as shown in Fig. 1, we define
$D_{i}$ to be the random variable that represents the shortest distance from
$T_{0}$ to the $i$th wall. Define $R_{i}$ as the distance between the
projection of $T_{0}$ on the $i$th wall and the nearest RIS on that wall.
$h_{i}$ is a constant describing the distance from $T_{0}$ to the $i$th wall.
We define $D_{ts}$ as the minimum among $D_{i}$, $i=1,2,3,4$, and $D_{sr}$ as
the distance between $S_{0}$ and $R_{0}$. In addition, there are three
highlighted rectangular areas between any two nodes in Fig. 1. These
rectangular areas are related to the calculation of the LOS/NLOS probabilities
and will be explained in detail later.

##  III Distances Distributions

###  III-A CDF expression of $D_{ts}$

According to the association policy, $D_{ts}$ is the shortest distance between
the transmitter and the RIS on the walls. We define
$\mathcal{D}=\\{D_{1},D_{2},D_{3},D_{4}\\}$,
$\mathcal{R}=\\{R_{1},R_{2},R_{3},R_{4}\\}$, and
$\mathcal{H}=\\{h_{1},h_{2},h_{3},h_{4}\\}$ as three sets. With the help of
$\mathcal{D}$ and $\mathcal{R}$, $D_{ts}$ can be expressed as

$\displaystyle D_{ts}=\min\\{\mathcal{D}\\}$
$\displaystyle=\min\bigg{\\{}\sqrt{\mathcal{R}^{2}+\mathcal{H}^{2}}\bigg{\\}}.$
(1)

From (1), it is noted that acquiring the CDF of $R_{i}$ is necessary for
calculating the CDF of $D_{ts}$. Because the RIS’s location distribution is
the same on every wall, we only show the derivation process of $R_{1}$. The
CDF derivation processes of $R_{2}$, $R_{3}$, and $R_{4}$ are similar to
$R_{1}$. Define the lower left corner of the room as the origin of the room.
Define the location of the transmitter as $(x_{t},y_{t})$.
$X_{1},X_{2},\cdots,X_{n}$ are $n$ uniform random variables following
$\mathcal{U}(0,a)$ that represent the locations of the RISs on the first wall.
$Y_{j}$ represents the distance from the $j$th RIS to $(0,y_{t})$. If $y_{t}$
is smaller than $a/2$, the CDF of $Y_{j}$ is given by

$\displaystyle F_{Y_{j}}(y)=P(Y_{j}\leq y)$ $\displaystyle=P(|y_{t}-X_{j}|\leq
y)$ $\displaystyle=P(y_{t}-y\leq X_{j}\leq y_{t}+y)$
$\displaystyle=\frac{1}{a}\bigg{(}\min\\{a-y_{t},y_{t}+y\\}-\max\\{0,y_{t}-y\\}\bigg{)},$
(2)

if $y_{t}$ is larger than $a/2$, the CDF of $Y_{i}$ is given by

$F_{Y_{j}}(y)=\frac{1}{a}\bigg{(}\min\\{y_{t},y_{t}+y\\}-\max\\{0,y_{t}-y\\}\bigg{)}.$
(3)

Because $X_{1},X_{2},\cdots,X_{n}$ are i.i.d. uniform random variables,
$Y_{1},Y_{2},\cdots,Y_{n}$ are also i.i.d. random variables with the CDF of
(2) or (3). $R_{1}$ can be interpreted as the first order variable among
$Y_{1},Y_{2},\cdots,Y_{n}$. Based on order statistics, the CDF of the first
order variable is given by

$F_{R_{1}}(r)=1-(1-F_{Y_{j}}(r))^{n}.$ (4)

Note that, when deriving the CDFs of $R_{2}$, $R_{3}$, and $R_{4}$, the value
of $y_{t}$ should be replaced to the corresponding values (e.g., $x_{t}$ for
$R_{2}$, $y_{t}$ for $R_{3}$, and $x_{t}$ for $R_{4}$). With the CDF of
$R_{1}$, the CDF of $D_{1}$ is given by

$\displaystyle F_{D_{1}}(d)=P(D_{1}\leq d)$
$\displaystyle=P\big{(}\sqrt{R_{1}^{2}+h_{1}^{2}}\leq d\big{)}$
$\displaystyle=P\big{(}R_{1}\leq\sqrt{d^{2}-h_{1}^{2}}\big{)}$
$\displaystyle=F_{R_{1}}(\sqrt{d^{2}-h_{1}^{2}})\mathds{1}(d,h_{1}),$ (5)

where $\mathds{1}(a,b)$ is the indicator function. $\mathds{1}(a,b)$ is equal
to 1 when $a\geq b$, $\mathds{1}(a,b)$ is equal to 0 when $a<b$. As the CDF
derivation steps of $D_{2}$, $D_{3}$, and $D_{4}$ are similar to $D_{1}$, we
do not show the derivation details here. The CDF of $D_{i}$ is given by

$\displaystyle
F_{D_{i}}(d)=F_{R_{i}}(\sqrt{d^{2}-h_{i}^{2}})\mathds{1}(d,h_{i}).$ (6)

With the CDF of $D_{i}$, the CDF of $D_{ts}$ can be calculated as

$\displaystyle F_{D_{ts}}(d)=P(D_{ts}\leq d)$
$\displaystyle=P\bigg{(}\min\\{\mathcal{D}\\}\leq d\bigg{)}$
$\displaystyle=1-P\bigg{(}\min\\{\mathcal{D}\\}\geq d\bigg{)}$
$\displaystyle=1-\prod_{i=1}^{4}\bigg{(}1-F_{D_{i}}(d)\bigg{)}.$ (7)

A comparison between the analytical expression of the CDF of $D_{ts}$ and the
Monte Carlo simulations is given in Fig. 2. The dashed lines represent the CDF
of $D_{ts}$ obtained from the analytical expression, and the solid lines
represent the CDF of $D_{ts}$ obtained from the Monte Carlo simulations. It is
noted that there is a high match between the dashed lines and solid lines,
which indicates the accuracy of the derived analytical expression.

As we already obtain the CDF expression of $D_{ts}$ related to the first
communication step in the indirect link, we will forward to the second
communication step and analyze the related CDF expression in the following.

Figure 2: CDF of $D_{ts}$ comparison between the Monte Carlo simulations and
the analytical expression with $a=30~{}\textrm{m}$, $x_{t}=3~{}\textrm{m}$,
$y_{t}=27~{}\textrm{m}$ and the RIS number on one wall $n=1,3,9$.

###  III-B CDF Expression of $D_{sr}$

We first show the exact analytical expression of the CDF of $D_{sr}$. This is
given by

$\displaystyle F_{D_{sr}}(d)$
$\displaystyle=\sum_{i=1}^{4}F_{D_{sr}}\bigg{(}d\ |\
D_{i}=\min\\{\mathcal{D}\\}\bigg{)}$ $\displaystyle\times
P\bigg{(}D_{i}=\min\\{\mathcal{D}\\}\bigg{)}$
$\displaystyle=\sum_{i=1}^{4}F_{D_{sr}}\bigg{(}d\ |\
D_{i}=\min\\{\sqrt{\mathcal{R}^{2}+\mathcal{H}^{2}}\\}\bigg{)}$
$\displaystyle\times P\bigg{(}D_{i}=\min\\{\mathcal{D}\\}\bigg{)}.$ (8)

Then, we show how to calculate the fourth term in the summation in (8) (The
calculation steps for other terms are similar to the fourth term). For now, we
assume the selected RIS is on the fourth wall, and omit the probability term
in the fourth term. Define $(X_{r},Y_{r})$ as the coordinate of the receiver.
Because the receiver follows HPPP, $X_{r}$ and $Y_{r}$ are two i.i.d. random
variables following uniform distribution $\mathcal{U}(0,a)$. By observing the
geometry relationship between the transmitter, receiver, and selected RIS in
Fig. 3, we find that $D_{sr}$ can be represented by $X_{r}$, $Y_{r}$, $x_{t}$,
$y_{t}$, $h_{4}$ and $R_{4}$. The expression is given by

$D_{sr}=\sqrt{\big{(}R_{4}+(X_{r}-x_{t})\big{)}^{2}+\big{(}h_{4}+(Y_{r}-y_{t})\big{)}^{2}},$
(9)

and the expression of the fourth term (with the conditions omitted) is given
by

$\displaystyle F_{D_{sr}}(d\ |\ D_{4})=P(D_{sr}\leq d)$
$\displaystyle=P\bigg{(}\sqrt{\big{(}R_{4}+(X_{r}-x_{t})\big{)}^{2}+\big{(}h_{4}+(Y_{r}-y_{t})\big{)}^{2}}\leq
d\bigg{)}$ $\displaystyle=P\bigg{(}\big{(}R_{4}+(X_{r}-x_{t})\big{)}^{2}\leq
d^{2}-\big{(}h_{4}+(Y_{r}-y_{t})\big{)}^{2}\bigg{)}$
$\displaystyle\stackrel{{\scriptstyle(a)}}{{=}}\underbrace{P\bigg{(}x_{t}-X_{r}\leq
R_{4}\leq
x_{t}-X_{r}+\sqrt{d^{2}-\big{(}h_{4}+(Y_{r}-y_{t})\big{)}^{2}}\bigg{)}}_{P_{1}}$
$\displaystyle+\underbrace{P\bigg{(}x_{t}-X_{r}-\sqrt{d^{2}-\big{(}h_{4}+(Y_{r}-y_{t})\big{)}^{2}}\leq
R_{4}\leq x_{t}-X_{r}\bigg{)}}_{P_{2}},$ (10)

$\displaystyle
P_{1}=\bigintsss_{0}^{a}\frac{1}{a}\bigintsss_{0}^{a}\frac{1}{a}\bigintsss_{\max\big{\\{}0,x_{t}-x_{r}\big{\\}}}^{\min\bigg{\\{}a,x_{t}-x_{r}+\sqrt{\max\big{\\{}0,d^{2}-\big{(}h_{4}+(y_{r}-y_{t})\big{)}^{2}\big{\\}}}\bigg{\\}}}f_{R_{4}}(r)\textrm{d}r\textrm{d}x_{r}\textrm{d}y_{r}.$
(11) $\displaystyle
P_{2}=\bigintsss_{0}^{a}\frac{1}{a}\bigintsss_{0}^{a}\frac{1}{a}\bigintsss_{\max\bigg{\\{}0,x_{t}-x_{r}-\sqrt{\max\big{\\{}0,d^{2}-\big{(}h_{4}+(y_{r}-y_{t})\big{)}^{2}\big{\\}}}\bigg{\\}}}^{\min\big{\\{}a,x_{t}-x-r\big{\\}}}f_{R_{4}}(r)\textrm{d}r\textrm{d}x_{r}\textrm{d}y_{r}.$
(12)

where (a) in (10) comes from that the term $R_{4}+(X_{r}-x_{t})$ can be
negative or positive. The first term (defined as $P_{1}$) in (10) can be
expanded as (11) and the second term (defined as $P_{2}$) in (10) can be
expanded as (12) (shown on the next page). Under the assumption of the
selected RIS being on the fourth wall, three integrals are needed to calculate
the CDF of $D_{sr}$. When we incorporate the probability in the fourth term,
there will be six integrals, including the integrals with respect to $R_{1}$,
$R_{2}$, and $R_{3}$. This suggests the direct calculation of $F_{D_{sr}}(d)$
is not possible. We know most of the computation in the expression of
$F_{D_{sr}}(d)$ comes from deciding which wall the selected RIS is on. To
simplify this process, we can use $h_{1}$, $h_{2}$, $h_{3}$, and $h_{4}$ to
approximate which wall has the selected RIS. For instance, if $h_{1}$ is the
minimum among $h_{1}$, $h_{2}$, $h_{3}$, and $h_{4}$, we assume the selected
RIS is on the first wall444The reason this method is not 100% accurate is that
there may exist more than one minimum values. For instance,
$h_{1}=h_{2}<h_{3}<h_{4}$. Under this scenario, the selected RIS may show on
the first wall or the second wall.. We then show the comparisons between the
approximated CDF of $D_{sr}$ and the Monte Carlo simulations. In Fig. 4, the
dots represent the Monte Carlo simulations, and solid lines represent the
approximated analytical results. It is observed that there is a high match
between the approximation and simulations. In the zoom-in figure, it is shown
that the maximum deviation is around 0.05 at $n=3$. We give more comparisons
with different parameters in Appendix A to show the accuracy of the
approximation. With the CDF of $D_{ts}$ and $D_{sr}$, the analysis of the
indirect link is completed. We then forward to the direct link where the CDF
of $D_{tr}$ will be described in the following.

Figure 3: Geometric illustration for CDF of $D_{sr}$ calculation Figure 4: CDF
of $D_{sr}$ comparison between the Monte Carlo simulation and the analytical
expression with $a=30~{}\textrm{m}$,
$x_{t}=3~{}\textrm{m},y_{t}=27~{}\textrm{m}$, and the RIS number on a wall
$n=1,3,9$.

###  III-C CDF Expression of $D_{tr}$

The location of the receiver can be defined as $(X_{r},Y_{r})$, with $X_{r}$,
and $Y_{r}$ follow $\mathcal{U}(0,a)$, $D_{tr}$ can be expressed as

$D_{tr}=\sqrt{(X_{r}-x_{t})^{2}+(Y_{r}-y_{t})^{2}}.$ (13)

The CDF of $D_{tr}$ is given by

$\displaystyle F_{D_{tr}}=P(D_{tr}\leq d)$
$\displaystyle=P(\sqrt{(X_{r}-x_{t})^{2}+(Y_{r}-y_{t})^{2}}\leq d)$
$\displaystyle=P((X_{r}-x_{t})^{2}\leq d^{2}-(Y_{r}-y_{t})^{2}).$ (14)

Due to $X_{r}-x_{t}$ can be positive or negative, (14) can be expanded as

$\displaystyle P\big{(}(X_{r}-x_{t})^{2}\leq d^{2}-(Y_{r}-y_{t})^{2}\big{)}$
$\displaystyle=P\big{(}x_{t}\leq X_{r}\leq
x_{t}+\sqrt{d^{2}-(Y_{r}-y_{t})^{2}}\big{)}$
$\displaystyle+P\big{(}x_{t}-\sqrt{d^{2}-(Y_{r}-y_{t})^{2}}\leq X_{r}\leq
x_{t}\big{)}$ $\displaystyle=P\big{(}x_{t}\leq x_{r}\leq
x_{t}+\sqrt{d^{2}-(y_{r}-y_{t})^{2}}|x_{r},y_{r}\big{)}$
$\displaystyle+P\big{(}x_{t}-\sqrt{d^{2}-(y_{r}-y_{t})^{2}}\leq x_{r}\leq
x_{t}|x_{r},y_{r}\big{)}$
$\displaystyle=\int_{0}^{a}\frac{1}{a}\int_{x_{t}}^{\max\bigg{\\{}a,x_{t}+\sqrt{\max\big{\\{}0,d^{2}-(y_{r}-y_{t})^{2}\big{\\}}}\bigg{\\}}}\frac{1}{a}\textrm{d}x_{r}\textrm{d}y_{r}$
$\displaystyle+\int_{0}^{a}\frac{1}{a}\int_{\max\bigg{\\{}0,x_{t}-\sqrt{\max\big{\\{}0,d^{2}-(y_{r}-y_{t})^{2}\big{\\}}}\bigg{\\}}}^{a}\frac{1}{a}\textrm{d}x_{r}\textrm{d}y_{r}.$
(15)

######  Remark 1.

In [38], a similar CDF was derived, i.e., the distance between a fixed point
and a point generated by an HPPP. Interested readers can refer to page-224 to
page-225 [38] for more details.

###  III-D CP Expression

Before we derive the expression of CP, we briefly discuss the calculation of
LOS/NLOS probability. Previous work related to LOS/NLOS probability
calculation includes [32]. In [32], the LOS/NLOS probability of urban cellular
networks was analyzed using random shape theory, where the obstacles are
modeled as rectangles with the length and width being random variables, and
obstacles’ centers following PPP. It is worth pointing out that due to the
mathematical model of the obstacle, the CP calculation in [32] has a high
computational complexity. In this paper, we aim to simplify the obstacle
modeling process. Such a simplification can reduce the computation time when
the network scales up. We define the obstacle as a circle with a diameter of
$d_{b}$. If the communication link, which is assumed to be a line, is away
from the obstacle’s center more than the obstacle’s radius, it will not be
blocked. Otherwise, it will be blocked. The LOS probability then can be
interpreted as the center of the obstacle not appearing in the corresponding
rectangular area555The rectangular area between any two nodes can be drawn by
setting the node to be the middle point of the width, the obstacle diameter to
be the width, and the transmission distance to be the length. For a clearer
view, one can refer to Fig. 1, where the rectangular areas for any two nodes
are highlighted.. The NLOS probability can be interpreted as the center of the
obstacle appearing in the corresponding rectangular area. To be clearer, one
can refer to Fig. 1, where the rectangular area for any two nodes is
highlighted. Then, the LOS probability for any two nodes $P_{\textrm{LOS}}$,
whose expression can be given by

$P_{\textrm{LOS-PQ}}^{\textrm{RV}}=e^{-\lambda_{o}D_{\textrm{PQ}}d_{b}},$ (16)

where $D_{\textrm{PQ}},\ \textrm{PQ}=\\{\textrm{TS, SR, TR}\\}$ and
$\lambda_{o}$ represent the distance between two nodes and the obstacle
density. As the distances between any two nodes are characterized as random
variables in this context, the LOS probability in (16) is also a random
variable. We use the expected value of (16) to describe the LOS probability,
which is defined as

$P_{\textrm{LOS-PQ}}=\mathbb{E}[e^{-\lambda_{o}D_{\textrm{PQ}}d_{b}}].$ (17)

The NLOS probability of any two nodes is given by

$P_{\textrm{NLOS-PQ}}=1-\mathbb{E}[e^{-\lambda_{o}D_{\textrm{PQ}}d_{b}}].$
(18)

We then give the CP of the indirect link666The link reflected from the RIS.
When considering the indirect link, the direct link is assumed to be blocked,
or the SNR of the direct link is lower than the threshold, or both. and direct
link777The link from the transmitter to the receiver, without the assistance
of the RIS. When considering the direct link, there is no assumption of the
indirect link. below:

$\displaystyle P_{\textrm{cov}}^{i}=P_{\textrm{LOS-TS}}P_{\textrm{LOS-
SR}}P_{\textrm{NLOS-
TR}}P(\textrm{SNR}_{i}\geq\tau)P(\textrm{SNR}_{d}\geq\tau)$
$\displaystyle+P_{\textrm{LOS-TS}}P_{\textrm{LOS-SR}}P_{\textrm{LOS-
TR}}P(\textrm{SNR}_{i}\geq\tau)P(\textrm{SNR}_{d}\leq\tau)$
$\displaystyle+P_{\textrm{LOS-TS}}P_{\textrm{LOS-SR}}P_{\textrm{NLOS-
TR}}P(\textrm{SNR}_{i}\geq\tau)P(\textrm{SNR}_{d}\leq\tau),$ (19)
$\displaystyle P_{\textrm{cov}}^{d}=P_{\textrm{LOS-
TR}}P(\textrm{SNR}_{d}\geq\tau),$ (20)

where $\textrm{SNR}_{i}$ and $\textrm{SNR}_{d}$ denote the SNR for the
indirect link and the SNR for the direct link. In the indirect link, we assume
the RIS operates in the electrically small regime that the SNR is proportional
to the square of the surface area and inversely proportional to the square of
the multiplication of indirect distances [39]. Due to the beam being highly
directive and narrow in the THz communication system and the angular
selectivity [40] of the RIS, We assume that there is no interference from
other RISs to an intended receiver. Hence, the term
$P(\textrm{SNR}_{i}\geq\tau)$ in (19) can be expressed as

$\displaystyle
P(\textrm{SNR}_{i}\geq\tau)=P\bigg{(}\frac{P_{\textrm{Tx}}G_{t}G_{r}M^{4}}{\sigma^{2}}\bigg{(}\frac{\lambda}{4\pi
D_{ts}D_{sr}}\bigg{)}^{\alpha}\geq\tau\bigg{)}$
$\displaystyle=P\bigg{(}D_{sr}\leq\frac{\lambda}{4\pi
D_{ts}}\bigg{(}\frac{P_{\textrm{Tx}}G_{t}G_{r}M^{4}}{\tau\sigma^{2}}\bigg{)}^{-\alpha}\bigg{)}$
$\displaystyle=\int_{0}^{a}f_{D_{ts}}(d_{ts})\int_{0}^{\frac{\lambda}{4\pi
d_{ts}}\bigg{(}\frac{P_{\textrm{Tx}}G_{t}G_{r}M^{4}}{\tau\sigma^{2}}\bigg{)}^{-\alpha}}$
$\displaystyle\times f_{D_{sr}}(d_{sr})\textrm{d}d_{sr}\textrm{d}d_{ts},$ (21)

and the term $P(\textrm{SNR}_{d}\geq\tau)$ can be expanded as

$\displaystyle
P(\textrm{SNR}_{d}\geq\tau)=P\bigg{(}\frac{P_{\textrm{Tx}}G_{t}G_{r}}{\sigma^{2}}\bigg{(}\frac{\lambda}{4\pi
D_{tr}}\bigg{)}^{\alpha}\geq\tau\bigg{)}$
$\displaystyle=P\bigg{(}D_{tr}\leq\frac{\lambda}{4\pi}\bigg{(}\frac{P_{\textrm{Tx}}G_{t}G_{r}}{\tau\sigma^{2}}\bigg{)}^{\frac{1}{\alpha}}\bigg{)}$
$\displaystyle=\int_{0}^{\frac{\lambda}{4\pi}\big{(}\frac{P_{\textrm{Tx}}G_{t}G_{r}}{\tau\sigma^{2}}\big{)}^{\frac{1}{\alpha}}}f_{D_{tr}}(d_{tr})\textrm{d}d_{tr},$
(22)

where $P_{\textrm{Tx}}$, $G_{t}$, $G_{r}$, $M$, $\sigma^{2}$, $\lambda$,
$\tau$, and $\alpha$ denote the transmitted power, the antenna gain of the
transmitter, the antenna gain of the receiver, the RIS element number, noise
power, the wavelength of the transmitted wave, the SNR threshold, and the
attenuation factor. $f_{D_{ts}(d_{ts})}$, $f_{D_{sr}(d_{sr})}$, and
$f_{D_{tr}(d_{tr})}$ are probability density functions (PDFs) of $D_{ts}$,
$D_{sr}$, and $D_{tr}$ that can be obtained from the CDFs of $D_{ts}$,
$D_{sr}$, and $D_{tr}$, which have been given in the previous section. The
total CP is calculated as

$P_{\textrm{cov}}=P_{\textrm{cov}}^{i}+P_{\textrm{cov}}^{d}.$ (23)

Table I: Parameter settings for analytical expression validation. Parameter | Value
---|---
$a$ | [10, 20, 30, 40, 50, 60, 70, 80, 90, 100] m
$x_{t}$ | a/10
$y_{t}$ | 9a/10
RIS amount on four walls |  [20, 40, 60, 80, 100, 120, 140, 160, 180, 200]
$\lambda_{o}$ | [0.01, 0.02, 0.05] points/$\textrm{m}^{2}$
M | 16
$d_{b}$ | 1 m
$P_{\textrm{Tx}}$ |  $1\times 10^{-3}$ W
$\sigma^{2}$ |  $1\times 10^{-12}$ W

##  IV Simulations

###  IV-A Analytical Expression Validation

We first give the parameter settings for the validation of the analytical
expression in Table. I.
The comparison between the analytical results and the Monte Carlo simulations
for a square room with the transmitter at the upper left corner of the room is
shown in Fig. 5. Firstly, it is observed that there is a high match between
the analytical curves and the Monte Carlo simulations for different obstacle
densities. Secondly, it is noted that the CP will decrease with the increase
in the room size. It can be noticed that compared with $\lambda_{o}=0.01$ and
0.02, when $\lambda=0.05$, the slope of the CP curve is more steep. This
interesting phenomenon is expected as the result of the different LOS
probabilities in direct and indirect links in different scenarios. To
understand the trend of the three curves in Fig. 5 in a clearer way,
$P_{\textrm{LOS-TR}}$ and $P_{\textrm{LOS-TS}}P_{\textrm{LOS-SR}}$ are shown
in Fig. 5 and Fig. 5. As we can see, both $P_{\textrm{LOS-TR}}$ and
$P_{\textrm{LOS-TS}}P_{\textrm{LOS-SR}}$ decrease fast along with the increase
of the room size, which results in the steep slope of the curve for
$\lambda_{o}=0.05$. The relatively small slopes for $\lambda_{o}=0.01$ and
0.02 in Fig. 5 and Fig. 5 leads to the relatively gentle descent process on
the CP. Fig. 5 indicates that the room size will have a more critical effect
on the CP when the $\lambda_{o}$ is bigger. To explore more insights, CP
performance over different $\lambda_{o}$ and room size is shown in Fig. 6.
Observing Fig. 6 horizontally, CP decreases slowly over room size when
$\lambda_{o}$ is small, and fast when $\lambda_{o}$ is large. Observing Fig. 6
vertically, $\lambda_{o}$ dominates the CP performance when the room size is
small. However, with the increase of the room size, the variation of
$\lambda_{o}$ has less and less effect on the CP performance.

Figure 5: (a) Performance comparison between the analytical results and the
Monte Carlo simulations in a room with the transmitter being in the center of
the room. (b) Simulated direct link LOS probability. (c) Simulated indirect
link probability. Figure 6: Performance evaluation with respect to the room
size and the obstacle density with the RIS amount, transmitted power, and
transmitter location outlined in Table I.

###  IV-B Investigation on the Transmitter Location

In the above subsection, we investigate the CP performance with the Tx placed
at the upper left corner of the room. It is easy to understand that in the
small size indoor environment for the THz wireless communication networks, the
transmitter’s location hugely affects the CP performance. In this subsection,
we investigate four typical Tx locations as templates: (1) upper left corner
of the room with the coordinate $(x_{t},y_{t})=(a/10,9a/10)$, (2) center of
the room with the coordinate $(x_{t},y_{t})=(a/2,a/2)$, (3) upper-middle of
the room with the coordinate $(x_{t},y_{t})=(a/2,9a/10)$, (4) middle point
between location 1 and location 2 with the coordinate
$(x_{t},y_{t})=(3a/10,7a/10)$. We define these locations as the first, second,
third, and fourth locations. We show the CP performance of the network in Fig.
7 for the last three different locations of the transmitter (the performance
of the first location is shown in Fig. 5). Comparing Fig. 5 with three sub-
figures in Fig. 7, it is noted that the performance of Fig.7 is the best
compared to others. This can be explained by examining the expected values of
$D_{ts}$ and $D_{sr}$. Because the third location is closer to the wall than
the second and fourth locations (same compared to the first location), the
expected value of $D_{ts}$ of the third location is smaller than the second
and fourth locations and equal to the first location. In addition, the
expected value of $D_{sr}$ of the third location is smaller than the first and
fourth locations888The expected value of $D_{sr}$ of the third location is
larger than the second location. However, due to a large expected value of
$D_{ts}$, the performance of the second is worse than the third location..
Hence, Fig. 7 performs best. It is also noted that the performance in Fig. 7
is the worst. This is because the expected value of $D_{ts}$ in Fig. 7 is
larger than others, which leads to a worse CP performance. However, with the
increase of the obstacle density, the performance improvement introduced by
the transmitter location will diminish. To view the relationship between the
transmitter’s location and the performance of the network in a clearer way, we
also provide a three-dimensional diagram. In Fig. 8, the horizontal axis
represents the location of the transmitter on the $x$-axis, and the vertical
axis represents the location of the transmitter on the $y$-axis. The colors
indicate the CP in the corresponding locations. It is observed that the CP is
higher when the transmitter is closer to the wall, except in four corners of
the room. A certain degree of location symmetry can also be observed in Fig.
8. Hence, it is concluded that the optimal location for a transmitter is near
the wall except for the corners.

Figure 7: (a) Performance of the RIS-assisted indoor network with the first
transmitter location. (b) Performance of the RIS-assisted indoor network with
the second transmitter location. (c) Performance of the RIS-assisted indoor
network with the third transmitter location. Figure 8: Performance evaluation
with respect to the transmitter’s location in a 30-by-30 meters square room
with RIS number on one wall $n=15$, 0.02 $\textrm{points/m}^{2}$ obstacle
density, and $1\times 10^{-3}~{}\textrm{W}$ transmitted power.

###  IV-C CP variation with respect to $\lambda_{s}$ and $\lambda_{o}$

In this section, we investigate the CP performance variation over the RIS
number and obstacle density. In Fig. 9, the CP performance is evaluated at
four different Tx locations with respect to the RIS number at different room
sizes when the obstacle density is fixed. Firstly, it is noted that the CP
variation trend is similar (increase and then saturate) at different room
sizes for all sub-figures. Secondly, it is noted that the CP performance has a
significant increase in Fig. 9 and Fig. 9 (30% in Fig. 9 at $a=80$) with the
increase of the total RIS number, while Fig. 9, and Fig. 9 do not show similar
increase (maximum CP increase amount below 0.05). For Fig. 9 and Fig. 9, the
trend of the curves can be explained by examining the variation of the
expected values of $D_{ts}$ with the increase of the total RIS number. With
the increasing RIS number, the nearest RIS is approaching the projection point
of the transmitter on the wall, which will lead to the decrease of the
expected value of $D_{ts}$. This explains the CP increase in Fig. 9 and Fig.
9. Observing these two sub-figures carefully, we can also find that the
increased amount of Fig. 9 is higher than Fig. 9 for every room size. This is
because when optimizing the RIS location999By saying optimizing the RIS
location, we mean increasing the RIS number. Because the more RIS we have, the
closer the nearest RIS to the projection point. in Fig. 9, the expected value
of $D_{sr}$ will increase, which will counteract the benefit brought by RIS
location optimization. However, such phenomena do not exist when the
transmitter is placed in the third location. For Fig. 9 and Fig. 9, it can be
observed that the CP increased amount is not significant (below 0.05 for every
room size) with the increase of the RIS number. This is due to the large
expected value of $D_{ts}$ even after the RIS location optimization. To
investigate the improvement over CP made by the obstacle density and the RIS
number, more simulation results are shown in Fig.10, where the CP performance
variation is shown with respect to RIS number and obstacle density. Observing
vertically, we can see that the CP performance will decrease with the increase
in the obstacle density. However, the speed of decrease is smaller when the
RIS number is high (e.g., 40) and bigger when the RIS number is low (e.g.,
10). Observing horizontally, it is found that the obstacle density dominates
the CP performance when the obstacle density exceeds
0.05$\textrm{points}/\textrm{m}^{2}$. Below
0.05$\textrm{points}/\textrm{m}^{2}$, the increase of the RIS number can
improve the CP.

Figure 9: (a)-(d) Performance evaluation with respect to the RIS number at
different room sizes, with obstacle density being
$\lambda_{o}=0.02~{}\textrm{points/m}^{2}$, transmitter location being the
first, second, third, and fourth locations, and transmitted power being
$1\times 10^{-3}~{}\textrm{W}$. Figure 10: Performance evaluation with respect
to the RIS number and the obstacle density in a 30-by-30 meters square room
with the transmitter located at the first location. The transmitted power is
$1\times 10^{-3}~{}\textrm{W}$.

##  V Conclusion

In this paper, we described an order statistic-based method for indoor THz
wireless communication network performance analysis. We derived the analytical
expression of the CP. Through comparison with the Monte Carlo simulation, it
was verified that the analytical expression of the CP is accurate. We then
showed the performance variation with respect to the room size, transmitter’s
location, obstacle density, and RIS number. The impact on the CP performance
can be sorted as obstacle density, room size, total RIS number, and the
transmitter’s location. Although the effect of obstacles and room size is
dominant to the CP, a significant CP improvement (around 15% of CP increase by
optimizing the transmitter’s location, and around 30% by increasing the RIS
number) can also be observed with the transmitter location optimization and
the increase of the total RIS number.

##  Appendix A More $D_{sr}$ CDF comparisons with different parameters

Figure 11: (a) CDF of $D_{sr}$ comparison between the Monte Carlo simulation
and the analytical expression with $a=30$, $x_{t}=15\textrm{m}$,
$y_{t}=15\textrm{m}$, and the RIS number on a wall $n=1,3,9$. (b) CDF of
$D_{sr}$ comparison between the Monte Carlo simulation and the analytical
expression with $a=40$, $x_{t}=15\textrm{m}$, $y_{t}=15\textrm{m}$, and the
RIS number on a wall $n=1,3,9$.

In Fig. 11, the transmitter is placed at the center of a square room that has
a length of $30\textrm{m}$. We can observe that the largest deviation happens
when $n=1$, which is around 0.05. In Fig. 11, the transmitter is placed at the
center of a square room that has a length of $40\textrm{m}$, similarly, the
largest deviation happens when $n=1$, and it is around 0.05. With Fig. 4 and
Fig. 11, we can assure that the approximation is accurate under different
parameters.

## References

  * [1] K. Danyluk, T. Ulusoy, W. Wei, and W. Willett, “Touch and beyond: Comparing physical and virtual reality visualizations,” IEEE Transactions on Visualization and Computer Graphics, vol. 28, no. 4, pp. 1930–1940, 2020.
  * [2] N. T. Banerjee, A. J. Baughman, S.-Y. Lin, Z. A. Witte, D. M. Klaus, and A. P. Anderson, “Side-by-side comparison of human perception and performance using augmented, hybrid, and virtual reality,” IEEE Transactions on Visualization and Computer Graphics, vol. 28, no. 12, pp. 4787–4796, 2021.
  * [3] S. Ghosh, L. Winston, N. Panchal, P. Kimura-Thollander, J. Hotnog, D. Cheong, G. Reyes, and G. D. Abowd, “Notifivr: Exploring interruptions and notifications in virtual reality,” IEEE transactions on visualization and computer graphics, vol. 24, no. 4, pp. 1447–1456, 2018.
  * [4] J. Guo, U. Kurup, and M. Shah, “Is it safe to drive? An overview of factors, metrics, and datasets for driveability assessment in autonomous driving,” IEEE Transactions on Intelligent Transportation Systems, vol. 21, no. 8, pp. 3135–3151, 2019.
  * [5] L. Claussmann, M. Revilloud, D. Gruyer, and S. Glaser, “A review of motion planning for highway autonomous driving,” IEEE Transactions on Intelligent Transportation Systems, vol. 21, no. 5, pp. 1826–1848, 2019.
  * [6] D. Feng, C. Haase-Schütz, L. Rosenbaum, H. Hertlein, C. Glaeser, F. Timm, W. Wiesbeck, and K. Dietmayer, “Deep multi-modal object detection and semantic segmentation for autonomous driving: Datasets, methods, and challenges,” IEEE Transactions on Intelligent Transportation Systems, vol. 22, no. 3, pp. 1341–1360, 2020.
  * [7] I. Union, “IMT traffic estimates for the years 2020 to 2030,” Report ITU, vol. 2370, 2015.
  * [8] K. B. Letaief, W. Chen, Y. Shi, J. Zhang, and Y.-J. A. Zhang, “The roadmap to 6G: AI empowered wireless network,” IEEE communications magazine, vol. 57, no. 8, pp. 84–90, 2019.
  * [9] K. Liu, S. Jia, S. Wang, X. Pang, W. Li, S. Zheng, H. Chi, X. Jin, X. Zhang, and X. Yu, “100 Gbit/s THz photonic wireless transmission in the 350-GHz band with extended reach,” IEEE Photonics Technology Letters, vol. 30, no. 11, pp. 1064–1067, 2018.
  * [10] H.-H. Lu, C.-Y. Li, W.-S. Tsai, R.-D. Lin, Y.-S. Tang, Y.-X. Chen, Y.-S. Lin, and W.-C. Fan, “An Integrated Fiber-FSO-5G NR Sub-THz Link With 86.112 Gbps High Aggregate Data Rates,” Journal of Lightwave Technology, vol. 40, no. 24, pp. 7790–7798, 2022.
  * [11] C. Chaccour, M. N. Soorki, W. Saad, M. Bennis, and P. Popovski, “Can terahertz provide high-rate reliable low-latency communications for wireless VR?,” IEEE Internet of Things Journal, vol. 9, no. 12, pp. 9712–9729, 2022.
  * [12] L. Fan, H. Wang, Q. Yang, Y. Wang, B. Deng, and H. Xiao, “High frame-rate and low-latency video SAR based on robust Doppler parameters estimation in the terahertz regime,” IEEE Transactions on Geoscience and Remote Sensing, 2023.
  * [13] Z. Chen, X. Ma, B. Zhang, Y. Zhang, Z. Niu, N. Kuang, W. Chen, L. Li, and S. Li, “A survey on terahertz communications,” China Communications, vol. 16, no. 2, pp. 1–35, 2019.
  * [14] M. Di Renzo, A. Zappone, M. Debbah, M.-S. Alouini, C. Yuen, J. De Rosny, and S. Tretyakov, “Smart radio environments empowered by reconfigurable intelligent surfaces: How it works, state of research, and the road ahead,” IEEE journal on selected areas in communications, vol. 38, no. 11, pp. 2450–2525, 2020.
  * [15] H.-J. Song, H. Hamada, and M. Yaita, “Prototype of KIOSK data downloading system at 300 GHz: Design, technical feasibility, and results,” IEEE Communications Magazine, vol. 56, no. 6, pp. 130–136, 2018.
  * [16] H. Zhao, L. Wei, M. Jarrahi, and G. J. Pottie, “Extending spatial and temporal characterization of indoor wireless channels from 350 to 650 GHz,” IEEE Transactions on Terahertz Science and Technology, vol. 9, no. 3, pp. 243–252, 2019.
  * [17] J. M. Eckhardt, T. Doeker, S. Rey, and T. Kürner, “Measurements in a real data centre at 300 GHz and recent results,” in 2019 13th European Conference on Antennas and Propagation (EuCAP), pp. 1–5, IEEE, 2019.
  * [18] A.-A. A. Boulogeorgos and A. Alexiou, “Outage probability analysis of THz relaying systems,” in 2020 IEEE 31st Annual International Symposium on Personal, Indoor and Mobile Radio Communications, pp. 1–7, IEEE, 2020.
  * [19] S. Farrag, E. Maher, A. El-Mahdy, and F. Dressler, “Outage probability analysis of UAV assisted mobile communications in THz channel,” in 2021 16th Annual Conference on Wireless On-demand Network Systems and Services Conference (WONS), pp. 1–8, IEEE, 2021.
  * [20] Z. Lou, B. E. Y. Belmekki, and M.-S. Alouini, “Coverage analysis of hybrid RF/THz networks with best relay selection,” IEEE Communications Letters, 2023.
  * [21] J. Sayehvand and H. Tabassum, “Interference and coverage analysis in coexisting RF and dense terahertz wireless networks,” IEEE Wireless Communications Letters, vol. 9, no. 10, pp. 1738–1742, 2020.
  * [22] K. Humadi, I. Trigui, W.-P. Zhu, and W. Ajib, “Coverage analysis of user-centric dense terahertz networks,” IEEE Communications Letters, vol. 25, no. 9, pp. 2864–2868, 2021.
  * [23] A.-A. A. Boulogeorgos and A. Alexiou, “Error analysis of mixed THz-RF wireless systems,” IEEE Communications Letters, vol. 24, no. 2, pp. 277–281, 2019.
  * [24] P. Bhardwaj and S. M. Zafaruddin, “Performance of Dual-Hop Relaying for THz-RF Wireless Link Over Asymmetrical $\alpha$-$\mu$ Fading,” IEEE Transactions on Vehicular Technology, vol. 70, no. 10, pp. 10031–10047, 2021.
  * [25] P. Bhardwaj and S. M. Zafaruddin, “Fixed-gain AF relaying for RF-THz wireless system over $\alpha$-$\kappa$-$\mu$ shadowed and $\alpha$-$\mu$ channels,” IEEE Communications Letters, vol. 26, no. 5, pp. 999–1003, 2021.
  * [26] S. Li and L. Yang, “Performance analysis of dual-hop THz transmission systems over $\alpha$-$\mu$ fading channels with pointing errors,” IEEE Internet of Things Journal, vol. 9, no. 14, pp. 11772–11783, 2021.
  * [27] G. Stratidakis, E. N. Papasotiriou, H. Konstantinis, A.-A. A. Boulogeorgos, and A. Alexiou, “Relay-based blockage and antenna misalignment mitigation in THz wireless communications,” in 2020 2nd 6G Wireless Summit (6G SUMMIT), pp. 1–4, IEEE, 2020.
  * [28] Y. Zhu, G. Zheng, and K.-K. Wong, “Stochastic geometry analysis of large intelligent surface-assisted millimeter wave networks,” IEEE Journal on Selected Areas in Communications, vol. 38, no. 8, pp. 1749–1762, 2020.
  * [29] T. Shafique, H. Tabassum, and E. Hossain, “Stochastic geometry analysis of IRS-assisted downlink cellular networks,” IEEE Transactions on Communications, vol. 70, no. 2, pp. 1442–1456, 2022.
  * [30] C. Zhang, W. Yi, Y. Liu, K. Yang, and Z. Ding, “Reconfigurable intelligent surfaces aided multi-cell NOMA networks: A stochastic geometry model,” IEEE Transactions on Communications, vol. 70, no. 2, pp. 951–966, 2021.
  * [31] M. A. Kishk and M.-S. Alouini, “Exploiting randomly located blockages for large-scale deployment of intelligent surfaces,” IEEE Journal on Selected Areas in Communications, vol. 39, no. 4, pp. 1043–1056, 2020.
  * [32] T. Bai, R. Vaze, and R. W. Heath, “Analysis of blockage effects on urban cellular networks,” IEEE Transactions on Wireless Communications, vol. 13, no. 9, pp. 5070–5083, 2014.
  * [33] T. Wang, G. Chen, M.-A. Badiu, and J. P. Coon, “Performance analysis of RIS-assisted large-scale wireless networks using stochastic geometry,” IEEE Transactions on Wireless Communications, 2023.
  * [34] H. T. P. Da Silva, M. S. De Alencar, W. J. L. Queiroz, and H. S. Silva, “Coverage Probability Analysis of a RIS-Assisted THz Indoor System Based on Ray Tracing Simulations,” in 2023 IEEE Latin-American Conference on Communications (LATINCOM), pp. 1–6, IEEE, 2023.
  * [35] O. Abbas, N. Kouzayha, M. Kishk, H. Sarieddeen, M.-S. Alouini, and T. Y. Al-Naffouri, “Performance Analysis of Indoor THz Networks with Intelligent Reflective Surfaces,” in ICC 2023-IEEE International Conference on Communications, pp. 2816–2821, IEEE, 2023.
  * [36] Z. Li, H. Hu, J. Zhang, and J. Zhang, “Enhancing indoor mmWave wireless coverage: Small-cell densification or reconfigurable intelligent surfaces deployment?,” IEEE Wireless Communications Letters, vol. 10, no. 11, pp. 2547–2551, 2021.
  * [37] H. Qin, Z. Liu, and C. Yang, “Indoor mm-wave coverage enhancement: Reconfigurable intelligent surface deployment strategy based on human mobility model,” IEEE Communications Letters, vol. 26, no. 10, pp. 2475–2479, 2022.
  * [38] A. M. Mathai, An introduction to geometrical probability: distributional aspects with applications, vol. 1.  CRC Press, 1999.
  * [39] Ö. Özdogan, E. Björnson, and E. G. Larsson, “Intelligent reflecting surfaces: Physics, propagation, and pathloss modeling,” IEEE Wireless Communications Letters, vol. 9, no. 5, pp. 581–585, 2019.
  * [40] W. Tang, M. Z. Chen, X. Chen, J. Y. Dai, Y. Han, M. Di Renzo, Y. Zeng, S. Jin, Q. Cheng, and T. J. Cui, “Wireless communications with reconfigurable intelligent surface: Path loss modeling and experimental measurement,” IEEE Transactions on Wireless Communications, vol. 20, no. 1, pp. 421–439, 2020.

|  Zhi Chai received his B.Sc. and M.Sc from Xidian University and University
College London in 2020 and 2021. He is a Ph.D. student in the Engineering
Science Department at the University of Oxford. His current research interests
include liquid crystal-based reconfigurable intelligent surface modeling,
stochastic geometry-based reconfigurable intelligent surface-assisted network
performance analysis, and terahertz communications.
---|---
|  Jiajie Xu [S’19, M’24] received his M.Sc. and Ph.D. degrees from Yanshan
University and King Abdullah University of Science and Technology (KAUST) in
2019 and 2023. He is a postdoctoral research fellow in the communication
theory lab at KAUST. His current research interests include underwater
wireless acoustic communication, underwater target detection, underwater
cooperative networks, maritime communication, space-air-ground-sea integrated
communication systems and joint sensing and communication systems, stochastic
geometry, and energy harvesting wireless networks.
---|---
|  Mohamed-Slim Alouini [S’94-M’98-SM’03-F’09] was born in Tunis, Tunisia. He
received the Ph.D. degree in Electrical Engineering from the California
Institute of Technology (Caltech) in 1998. He served as a faculty member at
the University of Minnesota then in the Texas A&M University at Qatar before
joining in 2009 the King Abdullah University of Science and Technology (KAUST)
where he is now a Distinguished Professor of Electrical and Computer
Engineering. Prof. Alouini is a Fellow of the IEEE and of the OSA. He is
currently particularly interested in addressing the technical challenges
associated with the uneven distribution, access to, and use of information and
communication technologies in far-flung, rural, low-density populations, low-
income, and/or hard-to-reach areas.
---|---
|  Justin P. Coon Justin P. Coon received a BSc degree (with distinction) in
electrical engineering from the Calhoun Honours College, Clemson University,
USA and a PhD in communications from the University of Bristol, UK in 2000 and
2005, respectively. He held various research positions at Toshiba Research
Europe Ltd. (TREL) from 2004 until 2013, including the position of Research
Manager from 2010 to 2013, during which time he led all research on physical
layer communications and network science at TREL. Dr. Coon was a Visiting
Fellow with the School of Mathematics at the University of Bristol, UK from
2010 until 2012, and he held a position as Reader in the Department of
Electrical and Electronic Engineering at the same university from 2012 to
2013. He joined the University of Oxford in 2013 where he is currently a
Professor of Engineering Science as well as the Emmott Fellow in Engineering
Science at Oriel College. Dr. Coon is the recipient of Toshiba’s Distinguished
Research Award for his work on 4G systems and has received three “best paper”
awards. He has published more than 200 papers in IEEE and APS journals and
conferences and is a named inventor on more than 30 patents. He has served as
an Editor for several IEEE journals and has chaired or co-chaired various IEEE
conferences.
---|---