5G-NR

    5G New Radio
3GPP

    3rd Generation Partnership Project
ABS

    aerial base station
AC

    address coding
ACF

    autocorrelation function
ACR

    autocorrelation receiver
ADC

    analog-to-digital converter
AIC

    Analog-to-Information Converter
AIC

    Akaike information criterion
ARIC

    asymmetric restricted isometry constant
ARIP

    asymmetric restricted isometry property
ARQ

    Automatic Repeat Request
AUB

    asymptotic union bound
AWGN

    Additive White Gaussian Noise
AWGN

    additive white Gaussian noise
PSK

    asymmetric PSK
AWRICs

    asymmetric weak restricted isometry constants
AWRIP

    asymmetric weak restricted isometry property
BCH

    Bose, Chaudhuri, and Hocquenghem
BCHSC

    BCH based source coding
BEP

    bit error probability
BFC

    block fading channel
BG

    Bernoulli-Gaussian
BGG

    Bernoulli-Generalized Gaussian
BPAM

    binary pulse amplitude modulation
BPDN

    Basis Pursuit Denoising
BPPM

    binary pulse position modulation
BPSK

    Binary Phase Shift Keying
BPZF

    bandpass zonal filter
BSC

    binary symmetric channels
BU

    Bernoulli-uniform
BER

    bit error rate
BS

    base station
BW

    BandWidth
BLLL

     binary log-linear learning
CP

    Cyclic Prefix
CDF

    cumulative distribution function
CDF

    Cumulative Distribution Function
CDF

    cumulative distribution function
CCDF

    complementary cumulative distribution function
CCDF

    complementary CDF
CCDF

    complementary cumulative distribution function
CD

    cooperative diversity
CDMA

    Code Division Multiple Access
ch.f.

    characteristic function
CIR

    channel impulse response
CoSaMP

    compressive sampling matching pursuit
CR

    cognitive radio
CS

    compressed sensing
CS

    Compressed sensing
CS

    compressed sensing
CSI

    channel state information
CCSDS

    consultative committee for space data systems
CC

    convolutional coding
COVID-19

    Coronavirus disease
DAA

    detect and avoid
DAB

    digital audio broadcasting
DCT

    discrete cosine transform
DFT

    discrete Fourier transform
DR

    distortion-rate
DS

    direct sequence
DS-SS

    direct-sequence spread-spectrum
DTR

    differential transmitted-reference
DVB-H

    digital video broadcasting – handheld
DVB-T

    digital video broadcasting – terrestrial
DL

    DownLink
DSSS

    Direct Sequence Spread Spectrum
DFT-s-OFDM

    Discrete Fourier Transform-spread-Orthogonal Frequency Division Multiplexing
DAS

    Distributed Antenna System
DNA

    DeoxyriboNucleic Acid
EC

    European Commission
EED

    exact eigenvalues distribution
EIRP

    Equivalent Isotropically Radiated Power
ELP

    equivalent low-pass
eMBB

    Enhanced Mobile Broadband
EMF

    ElectroMagnetic Field
EU

    European union
EI

    Exposure Index
eICIC

    enhanced Inter-Cell Interference Coordination
FC

    fusion center
FCC

    Federal Communications Commission
FEC

    forward error correction
FFT

    fast Fourier transform
FH

    frequency-hopping
FH-SS

    frequency-hopping spread-spectrum
FS

    Frame synchronization
FS

    frame synchronization
FDMA

    Frequency Division Multiple Access
GA

    Gaussian approximation
GF

    Galois field
GG

    Generalized-Gaussian
GIC

    generalized information criterion
GLRT

    generalized likelihood ratio test
GPS

    Global Positioning System
GMSK

    Gaussian Minimum Shift Keying
GSMA

    Global System for Mobile communications Association
GS

    ground station
GMG

     Grid-connected MicroGeneration
HAP

    high altitude platform
HetNet

    Heterogeneous network
IDR

    information distortion-rate
IFFT

    inverse fast Fourier transform
IHT

    iterative hard thresholding
i.i.d.

    independent, identically distributed
IoT

    Internet of Things
IR

    impulse radio
LRIC

    lower restricted isometry constant
LRICt

    lower restricted isometry constant threshold
ISI

    intersymbol interference
ITU

    International Telecommunication Union
ICNIRP

    International Commission on Non-Ionizing Radiation Protection
IEEE

    Institute of Electrical and Electronics Engineers
ICES

    IEEE international committee on electromagnetic safety
IEC

    International Electrotechnical Commission
IARC

    International Agency on Research on Cancer
IS-95

    Interim Standard 95
KPI

    Key Performance Indicator
LEO

    low earth orbit
LF

    likelihood function
LLF

    log-likelihood function
LLR

    log-likelihood ratio
LLRT

    log-likelihood ratio test
LoS

    Line-of-Sight
LRT

    likelihood ratio test
LWRIC

    lower weak restricted isometry constant
LWRICt

    LWRIC threshold
LPWAN

    Low Power Wide Area Network
LoRaWAN

    Low power long Range Wide Area Network
NLoS

    Non-Line-of-Sight
Li-Fi

    light-fidelity
LED

    light emitting diode
LABS

    LoS transmission with each ABS
NLABS

    NLoS transmission with each ABS
MB

    multiband
MC

    macro cell
MDS

    mixed distributed source
MF

    matched filter
m.g.f.

    moment generating function
MI

    mutual information
MIMO

    Multiple-Input Multiple-Output
MISO

    multiple-input single-output
MJSO

    maximum joint support cardinality
ML

    maximum likelihood
MMSE

    minimum mean-square error
MMV

    multiple measurement vectors
MOS

    model order selection
${M}$-PSK

     $M$-ary phase shift keying
${M}$-PSK

     $M$-ary asymmetric PSK
MP

     multi-period
MINLP

    mixed integer non-linear programming
$M$-QAM

     $M$-ary quadrature amplitude modulation
MRC

    maximal ratio combiner
MSO

    maximum sparsity order
M2M

    Machine-to-Machine
MUI

    multi-user interference
mMTC

    massive Machine Type Communications
mm-Wave

    millimeter-wave
MP

    mobile phone
MPE

    maximum permissible exposure
MAC

    media access control
NB

    narrowband
NBI

    narrowband interference
NLA

    nonlinear sparse approximation
NLOS

    Non-Line of Sight
NTIA

    National Telecommunications and Information Administration
NTP

    National Toxicology Program
NHS

    National Health Service
LOS

    Line of Sight
OC

    optimum combining
OC

    optimum combining
ODE

    operational distortion-energy
ODR

    operational distortion-rate
OFDM

    Orthogonal Frequency-Division Multiplexing
OMP

    orthogonal matching pursuit
OSMP

    orthogonal subspace matching pursuit
OQAM

    offset quadrature amplitude modulation
OQPSK

    offset QPSK
OFDMA

    Orthogonal Frequency-division Multiple Access
OPEX

    Operating Expenditures
OQPSK/PM

    OQPSK with phase modulation
PAM

    pulse amplitude modulation
PAR

    peak-to-average ratio
PDF

    probability density function
PDF

    probability density function
PDF

    probability distribution function
PDP

    power dispersion profile
PMF

    probability mass function
PMF

    probability mass function
PN

    pseudo-noise
PPM

    pulse position modulation
PRake

    Partial Rake
PSD

    power spectral density
PSEP

    pairwise synchronization error probability
PSK

    phase shift keying
PD

    power density
$8$-PSK

     $8$-phase shift keying
PPP

    Poisson point process
PCP

    Poisson cluster process
FSK

    Frequency Shift Keying
QAM

    Quadrature Amplitude Modulation
QPSK

    Quadrature Phase Shift Keying
OQPSK/PM

    OQPSK with phase modulator
RD

    raw data
RDL

    ”random data limit”
RIC

    restricted isometry constant
RICt

    restricted isometry constant threshold
RIP

    restricted isometry property
ROC

    receiver operating characteristic
RQ

    Raleigh quotient
RS

    Reed-Solomon
RSSC

    RS based source coding
r.v.

    random variable
R.V.

    random vector
RMS

    root mean square
RFR

    radiofrequency radiation
RIS

    Reconfigurable Intelligent Surface
RNA

    RiboNucleic Acid
RRM

    Radio Resource Management
RUE

    reference user equipments
RAT

    radio access technology
RB

    resource block
SA-Music

    subspace-augmented MUSIC with OSMP
SC

    small cell
SCBSES

    Source Compression Based Syndrome Encoding Scheme
SCM

    sample covariance matrix
SEP

    symbol error probability
SG

    sparse-land Gaussian model
SIMO

    single-input multiple-output
SINR

    signal-to-interference plus noise ratio
SIR

    signal-to-interference ratio
SISO

    Single-Input Single-Output
SMV

    single measurement vector
SNR

    signal-to-noise ratio
SP

    subspace pursuit
SS

    spread spectrum
SW

    sync word
SAR

    specific absorption rate
SSB

    synchronization signal block
SR

    shrink and realign
tUAV

    tethered Unmanned Aerial Vehicle
TBS

    terrestrial base station
uUAV

    untethered Unmanned Aerial Vehicle
PDF

    probability density functions
PL

    path-loss
TH

    time-hopping
ToA

    time-of-arrival
TR

    transmitted-reference
TW

    Tracy-Widom
TWDT

    TW Distribution Tail
TCM

    trellis coded modulation
TDD

    Time-Division Duplexing
TDMA

    Time Division Multiple Access
Tx

    average transmit
UAV

    Unmanned Aerial Vehicle
URIC

    upper restricted isometry constant
URICt

    upper restricted isometry constant threshold
UWB

    ultrawide band
UWB

    Ultrawide band
URLLC

    Ultra Reliable Low Latency Communications
UWRIC

    upper weak restricted isometry constant
UWRICt

    UWRIC threshold
UE

    User Equipment
UL

    UpLink
WiM

    weigh-in-motion
WLAN

    wireless local area network
WM

    Wishart matrix
WMAN

    wireless metropolitan area network
WPAN

    wireless personal area network
WRIC

    weak restricted isometry constant
WRICt

    weak restricted isometry constant thresholds
WRIP

    weak restricted isometry property
WSN

    wireless sensor network
WSS

    Wide-Sense Stationary
WHO

    World Health Organization
Wi-Fi

    Wireless Fidelity
SpaSoSEnc

    sparse source syndrome encoding
VLC

    Visible Light Communication
VPN

    Virtual Private Network
RF

    Radio Frequency
FSO

    Free Space Optics
IoST

    Internet of Space Things
GSM

    Global System for Mobile Communications
2G

    Second-generation cellular network
3G

    Third-generation cellular network
4G

    Fourth-generation cellular network
5G

    Fifth-generation cellular network
gNB

    next-generation Node-B Base Station
NR

    New Radio
UMTS

    Universal Mobile Telecommunications Service
LTE

    Long Term Evolution
QoS

    Quality of Service

#  End-to-End Uplink Performance Analysis of Satellite-Based IoT Networks: A
Stochastic Geometry Approach

Jiusi Zhou, Ruibo Wang, Basem Shihada, (Senior Member, IEEE), and Mohamed-Slim
Alouini, Fellow, IEEE The authors are with King Abdullah University of Science
and Technology (KAUST), CEMSE division, Thuwal 23955-6900, Saudi Arabia.
Corresponding author: Ruibo Wang. (e-mail<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>slim.alouini@kaust.edu.sa).

###### Abstract

With the deployment of satellite constellations, Internet-of-Things (IoT)
devices in remote areas have gained access to low-cost network connectivity.
In this paper, we investigate the performance of IoT devices connecting in up-
link through low Earth orbit (LEO) satellites to geosynchronous equatorial
orbit (GEO) links. We model the dynamic LEO satellite constellation using the
stochastic geometry method and provide an analysis of end-to-end availability
with low-complexity and coverage performance estimates for the mentioned link.
Based on the analytical expressions derived in this research, we make a sound
investigation on the impact of constellation configuration, transmission
power, and the relative positions of IoT devices and GEO satellites on end-to-
end performance.

###### Index Terms:

Satellite-based IoT network, coverage probability, availability probability,
stochastic geometry.

##  I Introduction

Recently, non-terrestrial networks have garnered significant attention due to
their seamless global coverage [1], enabling even remote Internet-of-Things
(IoT) devices to connect to the network [2, 3]. Remote areas usually lack
sufficient terrestrial communication infrastructure, such as fiber optic
networks or cellular towers. Therefore, IoT devices in many cases can only
access the network via satellite or other aerial links [4]. However, the power
sources for IoT devices are typically limited, relying on batteries, solar
energy, or wireless charging technology [5, 6]. This implies that the signal
transmission power of IoT devices is ultimately weaker [7]. Facing substantial
path loss in long-distance communications [8], for instance, in direct
transmission with geosynchronous equatorial orbit (GEO) satellites, IoT
devices need to rely on alternative means to maintain stable communication.
The research has demonstrated the feasibility of direct communication between
IoT devices and low Earth orbit (LEO) satellites [9]. Therefore, establishing
an IoT device to GEO satellite link with LEO satellites as relays is a
reasonable approach [10].

Figure 1: Diagram of satellite-based IoT networks.

For the uplink IoT device-LEO satellite-GEO satellite communication link,
analyzing its end-to-end performance is a critical and meaningful research
topic [11]. Performance analysis allows us to know the probability of
connection for a given transmission power from IoT devices. These insights aid
in better energy management and allocation for IoT devices. However, due to
the non-geostationary nature of LEO satellites, accurately modeling and
analyzing the aforementioned end-to-end communication link becomes challenging
[12].

So far, the vast majority of literature achieves the dynamic variation of LEO
satellite positions through numerical simulations [13, 14]. However, modeling
and simulating massive LEO satellite constellations are computationally
expensive. When system parameters, such as the transmission power of IoT
devices, are adjusted, simulations have to be re-executed. Therefore, we aim
to represent metrics such as coverage probability as functions of system
parameters like transmission power through desired analytical methods,
reducing computational complexity. In the literature providing analytical
results for air-to-ground communication systems, very few studies have
considered the randomness of the positions of GEO and LEO satellites in the
satellite-relayed communication system [15]. Although the authors in [15, 16,
17] considered the randomness in the distribution of satellites or other
aerial devices, they did not model the satellite constellation exactly. As a
result, they only accounted for the dynamic variation of satellite positions
within a given orbit, lacking a comprehensive view of the entire satellite
system’s performance.

###  I-A Related Works

Based on the above discussion, we introduce stochastic geometry (SG) as a
mathematical tool suitable for dynamic network performance analysis [18].
Under the SG framework, LEO satellites are modeled as a stochastic point
process and the spherical binomial point process (BPP) is one of the most
common models among them [19, 20]. The core idea of the SG analytical
framework is to achieve strong tractability by assuming that satellites follow
a specific distribution. Since the BPP model does not account for the modeling
of orbits, it differs from actual satellite constellation models. Fortunately,
both in terms of network topology and performance estimation outcomes, the
differences between BPP and deterministic constellations are neglected [21,
22]. Therefore, although the spatial distribution of the satellite
constellation is regular but not completely random, we can still model
satellite constellations as a BPP, accepting a small cost and reasonable
modeling differences in exchange for the strong tractability it provides.

Figure 2: Diagram of satellite-based IoT networks.

There are several SG-based studies analyzing the performance of satellite
networks under the stochastic geometry framework. The authors in [23]
considered the uplink performance analysis of IoT-over-satellite link. The
above satellite-based IoT network coverage analysis framework is further
expanded in [24, 25]. Specifically, Authors in [24] adopted a more practical
modeling method, and authors in [25] proposed an adaptive coverage enhancement
strategy. However, research on satellite relay systems remains limited. This
is due to the complex distance distribution in spherical multi-hop networks,
posing significant challenges in providing low-complexity analytical results.
To our best knowledge, within the SG framework, only our previous research
[26] has analyzed the distance distribution of satellite-relay systems. In
this study, we considered a dual-hop communication scenario where signals are
transmitted from a ground transmitter to a ground receiver via an LEO
satellite as the relay. We respectively derived the analytical expressions for
the distance distribution of each hop. However, our study did not involve
performance analysis nor encompass modeling of GEO satellites or channel
modeling, thus differing significantly from the scope of this article.

###  I-B Contribution

As the first study to analyze the performance of IoT-GEO satellite link’s
performance under the SG framework, the contributions of this article are as
follows.

  * •

We derive the distance distributions of the IoT-LEO satellite link, as well as
the LEO satellite-GEO satellite link (LGL). Due to the non-independence of the
distance distributions between the two links, the distance distribution
derivation in this paper differs significantly from existing studies based on
spherical SG.

  * •

For the first time, we incorporate the inter-satellite link channel model into
the modeling of SG-based satellite coverage issues. Taking pointing errors
into account as small-scale fading, we analyze the coverage performance of the
LGL.

  * •

Unlike most SG-based studies on satellite networks focusing on point-to-point
performance, we analyze the end-to-end performance of the dual-hop link.
Specifically, we first derive the probability of having at least one LEO
satellite being within the common reliable communication range of the IoT
device and the GEO satellite, called the availability probability. We derive
the analytical expression for coverage probability, which is the probability
of end-to-end successful communication given that the availability is ensured.

  * •

Through numerical simulations, we validate the accuracy performance of the
aforementioned analytical expressions and present the impact of constellation
configuration, transmission power, and the relative positioning of IoT devices
and GEO satellites, on these metrics.

##  II System Model

In this section, we establish the spatial model of the IoT device, GEO
satellite, and LEO satellite constellation. Then, channel models of the IoT
device-LEO satellite link (ILL) and LGL are provided.

###  II-A Spatial Model

We consider a system including an IoT device, a GEO Satellite, and
$N_{\mathrm{LEO}}$ LEO satellites. LEO satellites are distributed on a sphere
around Earth with a radius $R_{\mathrm{LEO}}=R_{\oplus}+h_{\mathrm{LEO}}$ form
a spherical homogeneous BPP, where $R_{\oplus}$ is the radius of the Earth and
$h_{\mathrm{LEO}}$ denotes the altitude of LEO satellites. According to
Slivnyak’s theorem [27], the rotation of the coordinate system does not affect
the distribution of homogeneous point processes. Without loss of generality,
we take the Earth’s center as the origin and establish a spherical coordinate
system. The coordinates of the IoT device and the GEO satellites are denoted
as $x_{\mathrm{IoT}}(R_{\oplus},0,0)$ and
$x_{\mathrm{GEO}}(R_{\mathrm{GEO}},\Theta,0)$, where
$R_{\mathrm{GEO}}=R_{\oplus}+h_{\mathrm{GEO}}$ is the radius of the sphere
that the GEO satellite is located. Here $\Theta$ represents the central angle
between the IoT device and the GEO satellite, and the definition of the
central angle is given as follows.

######  Definition 1 (Central Angle).

The central angle between two communication devices, A and B, refers to the
angle created by the lines extending from the Earth’s center to device A and
device B respectively [28].

###  II-B Channel Model

In this article, we consider the uplink channel models to follow the free
space propagation model experienced with large-scale fading and small-scale
fading, and the received power is given as [29]

$\rho_{Q}^{r}=\rho_{Q}^{t}G_{Q}\left(\frac{\lambda_{Q}}{4\pi
l_{Q}}\right)^{2}\zeta_{Q}W_{Q},\ l_{Q}\leq l_{Q}^{\max},$ (1)

where $\rho_{Q}^{t}$, $G_{Q}$, $\lambda_{Q}$, $\zeta_{Q}$, and $W_{Q}$ denote
the transmission power, antenna gain, wavelength, additional attenuation
during propagation, and the power of small-scale fading, respectively;
$Q\in\\{{\mathrm{IL},{\mathrm{LG}}}\\}$. Parameters with subscripts
${\mathrm{IL}}$ and ${\mathrm{LG}}$ represent labels for the two links, ILLs
and LGLs respectively. We consider $\zeta_{\mathrm{IL}}$ is mainly caused by
rain attenuation [19] and the attenuation of satellite-GEO satellite link is
negligible ($\zeta_{\mathrm{LG}}=0$ dB).

Considering the energy constraints of the IoT device, we assume that the IoT
device will select the LEO satellite closest to it as the relay.
$l_{\mathrm{IL}}$ and $l_{\mathrm{LG}}$ denote the distance between the IoT
device and the selected relay LEO satellite, and the distance between this LEO
satellite and the GEO satellite. Furthermore, $l_{\mathrm{IL}}^{\max}$ and
$l_{\mathrm{LG}}^{\max}$ are the maximum distances that can maintain stable
communication for ILLs and LGLs. To ensure these links are not blocked by the
Earth, $l_{\mathrm{IL}}^{\max}\leq\sqrt{R_{\mathrm{LEO}}^{2}-R_{\oplus}^{2}}$
and $l_{\mathrm{LG}}^{\max}\leq
2\sqrt{R_{\mathrm{GEO}}^{2}-R_{\mathrm{LEO}}^{2}}$ need to be satisfied.

The small-scale fading of ILL is assumed to follow the shadowed-Rician (SR)
fading. Considering the shadowing and multi-path effects, SR fading stands out
as one of the most precise models describing small-scale fading in space-
terrestrial communication links. The cumulative distribution function (CDF) of
the SR fading power $W_{\mathrm{IL}}$ is given as follows [30]:

$\begin{split}F_{W_{\mathrm{IL}}}(w)&=\left(\frac{2b_{0}m}{2b_{0}m+\Omega}\right)^{m}\sum_{z=0}^{\infty}\frac{(m)_{z}}{z!\Gamma(z+1)}\\\
&\times\left(\frac{\Omega}{2b_{0}m+\Omega}\right)^{z}\Gamma_{l}\left(z+1,\frac{w}{2b_{0}}\right),\end{split}$
(2)

where $(m)_{z}$ is the Pochhammer symbol, $m$, $b_{0}$ and $\Omega$ are
parameters of the SR fading. $\Gamma(\cdot)$ and $\Gamma_{l}(\cdot,\cdot)$
denote the gamma function and lower incomplete gamma function, respectively.

Furthermore, the shadowing and multi-path effects have minimal impact in the
space [31]. In this case, the pointing error becomes the primary factor in
small-scale fading due to the rapid movement of LEO satellites. Therefore, we
consider $W_{\mathrm{LG}}$ follows the pointing error model [32]. Given that
the deviation angle of the beam as $\theta_{d}$, the conditional probability
density function (PDF) of $W_{\mathrm{LG}}$ can be written as

$f_{W_{\mathrm{LG}}\,|\,\theta_{d}}\left(w\right)=\frac{\eta_{s}^{2}w^{\eta_{s}^{2}-1}\cos\left(\theta_{d}\right)}{A_{0}^{\eta_{s}^{2}}},\
0\leq w\leq A_{0},$ (3)

where $\eta_{s}$ and $A_{0}$ are parameters of the pointing error, and the
deviation angle $\theta_{d}$ is subject to Rayleigh distribution with variance
$\varsigma^{2}$[32],

$f_{\theta_{d}}\left(\theta_{d}\right)=\frac{\theta_{d}}{\varsigma^{2}}\exp\left(-\frac{\theta_{d}^{2}}{2\varsigma^{2}}\right),\
\theta_{d}\geq 0.$ (4)

##  III Performance Analysis

In this subsection, we give the definitions and analytical expressions of two
end-to-end performance metrics, namely satellite availability probability and
coverage probability in order.

###  III-A Satellite Availability

In the system model, both ILL and LGL have been set with upper bounds for
reliable communication distances. Therefore, the presence of a relay LEO
satellite within the common reliable communication range of IoT devices and
GEO satellites is a prerequisite for communication. If the required above
conditions are achieved, the LEO satellite is considered available.

######  Definition 2 (Satellite Availability Probability).

The satellite availability probability is the probability that the selected
relay LEO satellite’s distance to the IoT device is less than
$l_{\mathrm{IL}}^{\max}$ and its distance to the GEO satellite is less than
$l_{\mathrm{LG}}^{\max}$.

Because the closest LEO satellite is selected as a relay to the IoT device and
the availability of LGL depends on the relay satellite’s position, we present
the following lemma for the PDF of the central angle in our IoT-LEO model
regarding the distance between the IoT and the relay satellite. As we have
employed a spherical coordinate system for modeling, expressing distance in
terms of central angles proves to be more concise.

######  Lemma 1.

The PDF of the central angle between the IoT device and the LEO satellite is

$f_{\theta_{\mathrm{IL}}}(\theta)=\frac{N_{\mathrm{LEO}}\sin\theta}{2}\left(\frac{1+\cos\theta}{2}\right)^{N_{\mathrm{LEO}}-1},$
(5)

where $N_{\mathrm{LEO}}$ represents the number of LEO satellites,
$\theta\leq\theta_{\mathrm{IL}}^{\max}$, and $\theta_{\mathrm{IL}}^{\max}$ can
be expressed as,

$\theta_{\mathrm{IL}}^{\max}=\arccos\left(\frac{R_{\mathrm{LEO}}^{2}+R_{\oplus}^{2}-\left(l_{\mathrm{LEO}}^{\max}\right)^{2}}{2R_{\mathrm{LEO}}R_{\oplus}}\right).$
(6)

###### Proof:

See Appendix A ∎

Based on the above central angle distribution, also known as the contact angle
distribution [33], the satellite availability probability can be determined by
the following theorem.

######  Theorem 1.

The LEO satellite availability probability can be given by the following
expression

$P^{A}=\int_{0}^{\theta_{\mathrm{IL}}^{\max}}f_{\theta_{\mathrm{IL}}}(\theta)P_{\mathrm{LG}}^{A}(\theta){\mathrm{d}}\theta,$
(7)

where $f_{\theta_{\mathrm{IL}}}(\theta)$ and $\theta_{\mathrm{IL}}^{\max}$ are
defined in (5) and (6), respectively. $P_{\mathrm{LG}}^{A}(\theta)$ is called
the conditional availability probability of LGL, given that the contact angle
is $\theta$, and its expression is given as

$\displaystyle P_{\mathrm{LG}}^{A}(\theta)=\left\\{\begin{array}[]{lll}0,\ \ \
\
\sqrt{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-2R_{\mathrm{LEO}}R_{\mathrm{GEO}}(\sin\theta\sin\Theta+\cos\theta\cos\Theta)}>l_{\mathrm{LG}}^{\max},\\\
1,\ \ \ \
\sqrt{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-2R_{\mathrm{LEO}}R_{\mathrm{GEO}}(\cos\theta\cos\Theta-\sin\theta\sin\Theta)}<l_{\mathrm{LG}}^{\max},\\\
\frac{1}{\pi}\arccos\left(\frac{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-\left(l_{\mathrm{LG}}^{\max}\right)^{2}}{2R_{\mathrm{LEO}}R_{\mathrm{GEO}}\sin\theta\sin\Theta}-\frac{\cos\theta\cos\Theta}{\sin\theta\sin\Theta}\right),\
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ {\mathrm{otherwise}}.\end{array}\right.$
(11)

###### Proof:

See Appendix B. ∎

As the theorem demands specific requirements for the distances of both links,
thus it can be considered an end-to-end performance metric. Considering that
the availability probabilities of LGL and dual-hop have been analyzed, we can
derive the following corollary about the availability probability of ILL based
on the above theorem.

######  Corollary 1.

The availability probability of ILL, which is the probability that the
distance between the LEO satellite and the IoT device is less than
$\theta_{\mathrm{IL}}^{\max}$, is given by

$\begin{split}P_{\mathrm{IL}}^{A}=1-\left(\frac{1+\cos\theta_{\mathrm{IL}}^{\max}}{2}\right)^{N_{\mathrm{LEO}}},\end{split}$
(12)

where $\theta_{\mathrm{IL}}^{\max}$ is defined in (6).

###  III-B Coverage Probability

In this subsection, we first provide the definition and mathematical
expression of coverage probability, followed by the deductions for analytical
expression.

######  Definition 3 (Coverage Probability).

The coverage probability is the probability that the signal-to-noise ratio
(SNR) of ILL is greater than the coverage threshold $\gamma_{\mathrm{IL}}$ and
the SNR of LGL greater than the coverage threshold $\gamma_{\mathrm{LG}}$.

Mathematically, the coverage probability can be given by

$P^{C}=\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}},\frac{\rho_{\mathrm{LG}}^{r}}{\sigma_{\mathrm{GEO}}^{2}}>\gamma_{\mathrm{LG}}\right],$
(13)

where $\sigma_{\mathrm{LEO}}^{2}$ and $\sigma_{\mathrm{GEO}}^{2}$ denote the
noise power at the LEO satellite and the GEO satellite.  Coverage probability
is one of the most widely used performance metrics for evaluating network
performance. It represents the ability of a communication receiver to
successfully demodulate the signal. Before providing the analytical expression
for coverage probability, we need to derive the unconditional CDF of pointing
errors as a lemma.

######  Lemma 2.

The unconditional CDF of pointing errors is approximately given by

$\displaystyle F_{W_{\mathrm{LG}}}(w)=\left\\{\begin{array}[]{lll}0,&w<0,\\\
\frac{w^{\eta_{s}^{2}}}{A_{0}^{\eta_{s}^{2}}}\left(1-\varsigma^{2}\right)&0\leq
w\leq A_{0},\\\ 1,&w>A_{0}.\\\ \end{array}\right.$ (17)

###### Proof:

See Appendix C. ∎

Based on the above lemma, the end-to-end coverage probability can be deduced
by the following theorem.

######  Theorem 2.

The coverage probability is given in as

$\begin{split}&P^{C}=\int_{0}^{\theta_{\mathrm{IL}}^{\max}}\int_{0}^{P_{\mathrm{LG}}^{A}(\theta)}\left(1-F_{W_{\mathrm{IL}}}\left(\frac{\gamma_{\mathrm{IL}}\sigma_{\mathrm{LEO}}^{2}}{\rho_{\mathrm{IL}}^{t}G_{\mathrm{IL}}\zeta_{\mathrm{IL}}}\left(\frac{4\pi}{\lambda_{\mathrm{IL}}}\right)^{2}\left(R_{\mathrm{LEO}}^{2}+R_{\oplus}^{2}-2R_{\mathrm{LEO}}R_{\oplus}\cos\theta\right)\right)\right)\\\
&\times\left(1-F_{W_{\mathrm{LG}}}\left(\frac{\gamma_{\mathrm{LG}}\sigma_{\mathrm{GEO}}^{2}}{\rho_{\mathrm{LG}}^{t}G_{\mathrm{LG}}\zeta_{\mathrm{LG}}}\left(\frac{4\pi}{\lambda_{\mathrm{LG}}}\right)^{2}d^{2}\left(R_{\mathrm{GEO}},\Theta,0;R_{\mathrm{LEO}},\theta,\varphi\right)\right)\right)\frac{\sin\theta
f_{\theta_{\mathrm{IL}}}(\theta)}{4\pi}{\mathrm{d}}\varphi{\mathrm{d}}\theta.\end{split}$
(18)

where $\theta_{\mathrm{IL}}^{\max}$, $P_{\mathrm{LG}}^{A}$,
$F_{W_{\mathrm{IL}}}$ and $F_{W_{\mathrm{LG}}}$ are defined in (6), (11), (2)
and (17), respectively.
$d\left(R_{1},\theta_{1},\varphi_{1};R_{2},\theta_{2},\varphi_{2}\right)$
calculates the Euclidean distance between position
$(R_{1},\theta_{1},\varphi_{1})$ and position $(R_{2},\theta_{2},\varphi_{2})$
in the spherical coordinate system,

$\small\begin{split}&d\left(R_{1},\theta_{1},\varphi_{1};R_{2},\theta_{2},\varphi_{2}\right)\\\
&=\sqrt{R_{1}^{2}+R_{2}^{2}-2R_{1}R_{2}(\sin\theta_{1}\sin\theta_{2}\cos(\varphi_{1}-\varphi_{2})+\cos\theta_{1}\cos\theta_{2})}.\end{split}$
(19)

###### Proof:

See Appendix D. ∎

It can be shown that coverage probability is determined by the SNR of two
links, thus it is also an end-to-end performance metric.

TABLE I: Simulation Parameters [19, 32].

Notation | Meaning | Default Value
---|---|---
$N_{\mathrm{LEO}}$ | Number of LEO satellites | 1000
$R_{\oplus},R_{\mathrm{LEO}},R_{\mathrm{GEO}}$ | Radius of the Earth, LEO satellites, and GEO satellites |  $6371,7371,35860$ km
$\sigma_{\mathrm{LEO}}^{2},\sigma_{\mathrm{GEO}}^{2}$ | Noise power |  $5\times 10^{-10}$ mW
$\Omega,b_{0},m$ | Parameters of the SR fading | $1.29,0.158,19.4$
$\eta_{s}$, $A_{0}$ | Parameters of the pointing error | $1.00526,3.2120$
$\varsigma$ | Variance of Rayleigh distribution |  $15$ mrad
$\zeta_{\mathrm{IL}}$, $\zeta_{\mathrm{LG}}$ | Additional attenuation during propagation |  $-2$ dB, $0$ dB
$G_{\mathrm{ST}},G_{\mathrm{SS}}$ | Antenna gain |  $41.7$ dBi
$\lambda_{\mathrm{IL}},\lambda_{\mathrm{LG}}$ | Wavelength |  $1550$ nm
$\Theta$ | Central angle between the IoT device and the GEO satellite | $\pi/4$
$l_{\mathrm{IL}}^{\max},l_{\mathrm{LG}}^{\max}$ | Maximum distance of communication |  $3000$ km, $35000$ km
$\rho_{\mathrm{IL}}^{t},\rho_{\mathrm{LG}}^{t}$ | Transmission power |  $15$ dBW, $50$ dBW

##  IV Numerical Results

This section demonstrates the numerical results of availability probability
and coverage probability. The alignment of the results obtained from Monte
Carlo simulations (lines) with the derived analytical results (marks) in the
figures confirms the accuracy of the analytical expressions presented in this
article. Unless otherwise specified, the parameters in this section will be
set to their default values given in Table I.

Figure 3: Satellite availability probability with different constellation
configurations. Figure 4: Satellite availability probability with different
central angles.

In Fig. 3, we analyze the impact of the number of LEO satellites and
constellation altitude on the availability. For the same number of satellites,
constellations at lower altitudes have a larger availability probability. In
the subgraph, we show the number of satellites required to achieve an
availability probability of $1$ at different altitudes. In Fig. 4, the number
of satellites is fixed at $100$. When the central angle between IoT devices
and GEO satellites is less than $60$ degrees, the influence of the central
angle on the availability can be neglected. At this point, availability
probability is mainly determined by the ILL. When the central angle is larger
than $60$ degrees, the impact of central angles on satellite availability
probability is much larger compared to the constellation altitude. The
availability probability decreases rapidly as the central angle further
increases.

Figure 5: Satellite coverage probability with different constellation
configurations.

Fig. 5 shows the results of satellite coverage probability with different
constellation configurations. The lower altitude of LEO satellites and the
greater number of satellites result in shorter distances for ILL, thereby
improving coverage performance. Under default parameters, compared to LGL, the
quality of ILL plays a more decisive role in coverage probability.

Figure 6: Satellite coverage probability with different transmission power of
IoT device. Figure 7: Satellite coverage probability with different
transmission power of LEO satellite.

Fig. 6 and Fig. 7 present end-to-end satellite coverage probability with
different transmission power of IoT devices and LEO satellites. In general,
constellations with lower altitudes have better coverage performance. As shown
in Fig. 7, with the increase in LEO transmission power, the end-to-end
coverage probability first increases and then converges. The final converged
value depends on the coverage probability of the ILL.

##  V Additional Insights

In this section, we analyze the influence of transmission policy/association
strategy on coverage probability, compare different network architectures, and
describe the methods of counteracting Doppler shift.

###  V-A Transmission Policy

In this subsection, we discuss the impact of the following three different
transmission policies on coverage probability:

  * •

Policy 1 (proposed): Select the nearest LEO satellite to the IoT device for
relaying.

  * •

Policy 2: Select the nearest LEO satellite to the GEO satellite for relaying.

  * •

Policy 3: Select the LEO satellite with the shortest distance to the direct
link between the IoT device and the GEO satellite for relaying.

Fig. 8 compares the coverage probabilities of the three policies when the
central angle between the IoT device and the GEO satellite is fixed at
$\Theta=\frac{\pi}{12}$. Regardless of the transmission power of the IoT
device, our proposed transmission policy (policy 1) consistently outperforms
the other two policies. In terms of coverage probability, the performance of
Policy 1 and Policy 3 is similar, with both significantly outperforming Policy
2.

Figure 8: Satellite coverage probability with different transmission policies
and transmission powers.

In Fig. 9, we fix the transmission power and vary the central angle to observe
the impact of different transmission policies on coverage probability. The
change in central angle has no impact on the position of the relay LEO
satellite selected by Policy 1, so the central angle has little impact on the
coverage probability for Policy 1. In contrast, the relay positions selected
by Policy 2 and Policy 3 are influenced by the central angle, leading to a
significant decline in their coverage performance as the central angle
increases. Especially for Policy 2, when $\Theta>20$ degrees, the chosen relay
LEO satellite becomes unavailable to the IoT device, resulting in a coverage
probability of $0$.

Figure 9: Satellite coverage probability with different transmission policies
and central angles.

###  V-B Network Architecture

Here, we compare the advantages and disadvantages of the network
architectures: $(a_{1})$ IoT-LEO-GEO (proposed), $(a_{2})$ IoT-GEO direct
link, and $(a_{3})$ IoT-LEO direct link.

Compare $(a_{1})$ and $(a_{2})$: Direct transmission from IoT devices to GEO
satellites can simplify network configuration by eliminating the need for LGL
synchronization. Introducing LEO satellites as relays requires additional
capital and operational expenditures. Moreover, the Doppler effect is
introduced with LEO satellites serving as relays, which is discussed in detail
in the next subsection, is required to be overcome. However, the distance from
IoT devices to GEO satellites is several dozen times greater than to LEO
satellites, meaning that achieving the same coverage probability or data rate
would require an increase in transmission power by several thousand times.
Under the default parameters, with an IoT device transmission power of $18$
dB, nearly $100\%$ coverage probability can be achieved with LEO satellites as
relays, whereas without LEO satellites, the coverage probability is
approximately $0\%$. Therefore, direct communication between IoT devices and
GEO satellites requires stable energy sources and powerful hardware support,
which is often challenging to achieve in remote areas where IoT devices are
energy-constrained.

Compare $(a_{1})$ and $(a_{3})$: The path loss for a direct link between IoT
devices and LEO satellites is undoubtedly smaller compared to communication
with GEO satellites. Single-hop communication between IoT devices and LEO
satellites improves both availability and coverage probability compared to
dual-hop communication. Although the IoT-LEO direct link offers better
communication performance, it is typically not a complete communication. LEO
satellites have limited computational and storage capabilities, and thus they
are often not the final destination for the data. Consequently, after
receiving data, LEO satellites need to relay it back to ground gateways or
transmit it through other satellites before it reaches the ground.

###  V-C Doppler Effect

The rapid relative motion between satellites necessitates consideration of the
Doppler shift’s impact. Satellite communications typically operate at lower
transmission rates, which simplifies the process of bit synchronization [34].
Consequently, bit synchronization adjustments can be used to mitigate the
Doppler effect in satellite communications. Additionally, because bandwidth
constraints are less critical in satellite systems, broadening the range of
receiving frequencies does not lead to significant co-frequency interference.
Therefore, spectrum spreading emerges as an effective technique to counteract
Doppler shift [35].

##  VI Conclusion

In this paper, we provide analytical expressions for the end-to-end
availability and coverage performance of satellite-based IoT networks. Based
on these analytical expressions, the main conclusions are as follows. Lowering
the altitude of the relay, and increasing the number of satellites and
transmission power can enhance network performance. The central angle between
IoT devices and GEO satellites has little impact on network performance when
it is small. However, once the central angle exceeds a threshold, network
performance deteriorates rapidly with the central angle.

##  Appendix A Proof of Lemma 1

Considering the rotational invariance of the BPP, the PDF of the distance
between the IoT device and the LEO satellite should be independent of the
azimuthal angle of the LEO satellite, relying solely on the central angle.
Denote $\mathcal{S}(\theta)$ as the spherical cap with a central angle as
$2\theta$, the CDF of the central angle is given by

$\begin{split}&F_{\theta_{\mathrm{LEO}}}(\theta)=\mathbb{P}\left[\theta_{\mathrm{LEO}}\leq\theta\right]\\\
&=1-\mathbb{P}\left[{\mathcal{N}\left({{\mathcal{S}(\theta)}}\right)=0}\right]\\\
&=1-\left(1-\frac{\mathcal{A}\left(\mathcal{S}(\theta)\right)}{\mathcal{A}\left(\mathcal{S}(\pi)\right)}\right)^{N_{\mathrm{LEO}}}\\\
&=1-\left(1-\frac{2\pi R_{\mathrm{LEO}}^{2}(1-\cos\theta)}{4\pi
R_{\mathrm{LEO}}^{2}}\right)^{N_{\mathrm{LEO}}}\\\
&=1-\left(\frac{1+\cos\theta}{2}\right)^{N_{\mathrm{LEO}}},\end{split}$ (20)

where $\mathcal{N}\left(\mathcal{S}(\psi)\right)$ counts the number of
satellites in the spherical cap $\mathcal{S}(\psi)$, and
$\mathcal{A}\left(\mathcal{S}(\psi)\right)$ denotes the area of
$\mathcal{S}(\psi)$. Then, the PDF of the central angle can be derived by

$\begin{split}&f_{\theta_{\mathrm{LEO}}}(\theta)=\frac{\mathrm{d}}{\mathrm{d}\theta}F_{\theta_{\mathrm{LEO}}}(\theta)\\\
&=\frac{N_{\mathrm{LEO}}\sin\theta}{2}\left(\frac{1+\cos\theta}{2}\right)^{N_{\mathrm{LEO}}-1}.\end{split}$
(21)

Finally, when the distance between the LEO relay satellite and the IoT device
is $l_{\mathrm{LEO}}^{\max}$, the central angle between them reaches its upper
bound. Denote the upper bound as $\theta_{\mathrm{IL}}^{\max}$, which can be
expressed by the cosine rule,

$\cos\theta_{\mathrm{IL}}^{\max}=\frac{R_{\mathrm{LEO}}^{2}+R_{\oplus}^{2}-\left(l_{\mathrm{LEO}}^{\max}\right)^{2}}{2R_{\mathrm{LEO}}R_{\oplus}}.$
(22)

##  Appendix B Proof of Theorem 1

Because the distance distribution of LGL is not independent of the distance
distribution of ILL, we denote the probability that the distance from the LEO
satellite to the GEO satellite is less than $l_{\mathrm{LG}}^{\max}$ as
$P^{A}_{\mathrm{LG}}(\theta)$, given that the central angle between the IoT
device and LEO satellite is $\theta$. From the definition, the satellite
availability probability can be written as,

$P^{A}=\mathbbm{E}_{\theta_{\mathrm{IL}}}\left[P^{A}_{\mathrm{LG}}(\theta_{\mathrm{IL}})\right]=\int_{0}^{\theta_{\mathrm{IL}}^{\max}}f_{\theta_{\mathrm{IL}}}(\theta)P^{A}_{\mathrm{LG}}(\theta){\mathrm{d}}\theta.$
(23)

Therefore, the next step is to derive the analytical expression for
$P^{A}_{\mathrm{LG}}(\theta)$ to complete the proof. As the polar angle of the
LEO satellite has been fixed as $\theta$, we need to know the range of the
azimuth angle $\varphi$ of the LEO satellite when the LGL is available. When
$\varphi=0$, the LEO satellite is closest to the GEO satellite, and the
distance between them is

$\small
d_{1}=\sqrt{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-2R_{\mathrm{LEO}}R_{\mathrm{GEO}}(\sin\theta\sin\Theta+\cos\theta\cos\Theta)}.$
(24)

When $d_{1}>l_{\mathrm{LG}}^{\max}$, $P_{\mathrm{LG}}^{A}=0$. When
$\varphi=\pi$, the LEO satellite is farthest to the GEO satellite, and the
distance between them is

$\small
d_{2}=\sqrt{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-2R_{\mathrm{LEO}}R_{\mathrm{GEO}}(\cos\theta\cos\Theta-\sin\theta\sin\Theta)}.$
(25)

When $d_{2}\leq l_{\mathrm{LG}}^{\max}$, $P_{\mathrm{LG}}^{A}=1$. Otherwise,
when $d_{1}<l_{\mathrm{LG}}^{\max}<d_{2}$, the distance between the LEO
satellite located at $\left(R_{\mathrm{LEO}},\theta,\varphi\right)$ and the
GEO satellite is

$\begin{split}d_{3}&=\bigg{(}R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-2R_{\mathrm{LEO}}R_{\mathrm{GEO}}\\\
&\times\left(\sin\theta\sin\Theta\cos\varphi+\cos\theta\cos\Theta\right)\bigg{)}^{\frac{1}{2}}.\end{split}$
(26)

When $d_{3}=l_{\mathrm{LG}}^{\max}$, $\varphi$ reaches its critical values

$\cos\varphi=\frac{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-\left(l_{\mathrm{LG}}^{\max}\right)^{2}}{2R_{\mathrm{LEO}}R_{\mathrm{GEO}}\sin\theta\sin\Theta}-\frac{\cos\theta\cos\Theta}{\sin\theta\sin\Theta}.$
(27)

Since the value of $\varphi$ follows a uniform distribution, therefore

$\small
P_{\mathrm{LG}}^{A}(\theta)=\frac{1}{\pi}\arccos\left(\frac{R_{\mathrm{LEO}}^{2}+R_{\mathrm{GEO}}^{2}-\left(l_{\mathrm{LG}}^{\max}\right)^{2}}{2R_{\mathrm{LEO}}R_{\mathrm{GEO}}\sin\theta\sin\Theta}-\frac{\cos\theta\cos\Theta}{\sin\theta\sin\Theta}\right).$
(28)

Finally, we present a corollary of the above result. When
$P_{\mathrm{LG}}^{A}(\theta)$ is replaced by $1$, LGL is always available. At
this point, $P^{A}$ degenerates into the availability probability of ILL,

$\begin{split}P_{\mathrm{IL}}^{A}&=\int_{0}^{\theta_{\mathrm{IL}}^{\max}}f_{\theta_{\mathrm{IL}}}(\theta){\mathrm{d}}\theta=F_{\theta_{\mathrm{IL}}}\left(\theta_{\mathrm{IL}}^{\max}\right)\\\
&=1-\left(\frac{1+\cos\theta_{\mathrm{IL}}^{\max}}{2}\right)^{N_{\mathrm{LEO}}}.\end{split}$
(29)

##  Appendix C Proof of Lemma 2

The unconditional PDF of the pointing error can be derived based on the
definition:

$\begin{split}&f_{W_{\mathrm{LG}}}(w)=\int_{0}^{\infty}f_{W_{\mathrm{LG}}|\theta_{d}}(w)f_{\theta_{d}}(\theta_{d}){\mathrm{d}}\theta_{d}\\\
&=\frac{\eta_{s}^{2}w^{\eta_{s}^{2}-1}}{A_{0}^{\eta_{s}^{2}}}\int_{0}^{\infty}\cos\theta_{d}\frac{\theta_{d}}{\varsigma^{2}}\exp\left(-\frac{\theta_{d}^{2}}{2\varsigma^{2}}\right)\mathrm{d}\theta_{d}\\\
&\overset{(a)}{\approx}\frac{\eta_{s}^{2}w^{\eta_{s}^{2}-1}}{A_{0}^{\eta_{s}^{2}}}\left(1-\frac{1}{2}\int_{0}^{\infty}\frac{\theta_{d}^{3}}{\varsigma^{2}}\exp\left(-\frac{\theta_{d}^{2}}{2\varsigma^{2}}\right)\mathrm{d}\theta_{d}\right)\\\
&\overset{(b)}{=}\frac{\eta_{s}^{2}w^{\eta_{s}^{2}-1}}{A_{0}^{\eta_{s}^{2}}}\left(1-\int_{0}^{\infty}z\varsigma^{2}\exp\left(-z\right)\mathrm{d}z\right)\\\
&\overset{(c)}{=}\frac{\eta_{s}^{2}w^{\eta_{s}^{2}-1}}{A_{0}^{\eta_{s}^{2}}}\left(1-\varsigma^{2}\right),\end{split}$
(30)

where step $(a)$ follows the second-order Taylor expansion of
$\cos\theta_{d}\approx 1-\frac{\theta_{d}^{2}}{2}$, since $\theta_{d}$ is
generally a small value. Step $(b)$ is derived by the substitution of
$z={\theta_{d}^{2}}/{2\varsigma^{2}}$, and step $(c)$ is obtained by taking
the expectation of the exponential distribution.

Then, the unconditional CDF obtained by integrating the unconditional PDF,

$\begin{split}F_{W_{\mathrm{LG}}}(w)=\int_{0}^{w}\frac{\eta_{s}^{2}z^{\eta_{s}^{2}-1}}{A_{0}^{\eta_{s}^{2}}}\left(1-\varsigma^{2}\right){\mathrm{d}}z=\frac{w^{\eta_{s}^{2}}}{A_{0}^{\eta_{s}^{2}}}\left(1-\varsigma^{2}\right).\end{split}$
(31)

##  Appendix D Proof of Theorem 2

Due to the correlation between the distances of ILL and LGL,

$\begin{split}&P^{C}=\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}},\frac{\rho_{\mathrm{LG}}^{r}}{\sigma_{\mathrm{GEO}}^{2}}>\gamma_{\mathrm{LG}}\right]\\\
&\neq\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}}\right]\times\mathbbm{P}\left[\frac{\rho_{\mathrm{LG}}^{r}}{\sigma_{\mathrm{GEO}}^{2}}>\gamma_{\mathrm{LG}}\right].\end{split}$
(32)

Considering that, given the position of relay LEO satellite
$(R_{\mathrm{LEO}},\theta_{\mathrm{IL}},\varphi)$, the conditional coverage
probabilities of the two links are independent of each other. Therefore, we
calculate the end-to-end coverage probability by sequentially obtaining the
conditional coverage probabilities and then taking the expectation over
$\theta_{\mathrm{IL}}$. From the definition, the coverage probability can be
expressed as,

$\begin{split}&P^{C}=\mathbbm{E}_{\theta_{\mathrm{IL}},\varphi}\left[\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}},\frac{\rho_{\mathrm{LG}}^{r}}{\sigma_{\mathrm{GEO}}^{2}}>\gamma_{\mathrm{LG}}\,\bigg{|}\,\theta_{\mathrm{IL}},\varphi\right]\right]\\\
&=\int_{0}^{\theta_{\mathrm{IL}}^{\max}}\int_{0}^{P_{\mathrm{LG}}^{A}(\theta)}f_{\theta_{\mathrm{IL}}}(\theta)\frac{R_{\mathrm{LEO}}^{2}\sin\theta}{4\pi
R_{\mathrm{LEO}}^{2}}\\\
&\times\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}},\frac{\rho_{\mathrm{LG}}^{r}}{\sigma_{\mathrm{GEO}}^{2}}>\gamma_{\mathrm{LG}}\,\bigg{|}\,\theta,\varphi\right]{\mathrm{d}}\varphi{\mathrm{d}}\theta\\\
&=\int_{0}^{\theta_{\mathrm{IL}}^{\max}}\int_{0}^{P_{\mathrm{LG}}^{A}(\theta)}f_{\theta_{\mathrm{IL}}}(\theta)\frac{\sin\theta}{4\pi}\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}}\,\bigg{|}\,\theta\right]\\\
&\times\mathbbm{P}\left[\frac{\rho_{\mathrm{LG}}^{r}}{\sigma_{\mathrm{GEO}}^{2}}>\gamma_{\mathrm{LG}}\,\bigg{|}\,\theta,\varphi\right]{\mathrm{d}}\varphi{\mathrm{d}}\theta.\end{split}$
(33)

where $\theta_{\mathrm{IL}}^{\max}$ and $\pi P_{\mathrm{LG}}^{A}(\theta)$ are
the upper bound of the LEO satellite’s polar angle and azimuth angle. These
two upper bounds ensure the availability of the LEO satellite.

Next, we explicitly represent the first conditional probability as

$\begin{split}&\mathbbm{P}\left[\frac{\rho_{\mathrm{IL}}^{r}}{\sigma_{\mathrm{LEO}}^{2}}>\gamma_{\mathrm{IL}}\,\bigg{|}\,\theta\right]\\\
&=\mathbbm{P}\left[\frac{\rho_{\mathrm{LG}}G_{\mathrm{LG}}\zeta_{\mathrm{LG}}W_{\mathrm{LG}}}{\sigma_{\mathrm{LEO}}^{2}}\left(\frac{4\pi}{\lambda_{\mathrm{LG}}}l_{\mathrm{LG}}\right)^{2}>\gamma_{\mathrm{LG}}\,\bigg{|}\,l_{\mathrm{LG}}^{2}=R_{\mathrm{LEO}}^{2}+R_{\oplus}^{2}-2R_{\mathrm{LEO}}R_{\oplus}\cos\theta\right]\\\
&=\mathbbm{P}\left[W_{\mathrm{LG}}>\frac{\gamma_{\mathrm{IL}}\sigma_{\mathrm{LEO}}^{2}}{\rho_{\mathrm{IL}}^{t}G_{\mathrm{IL}}\zeta_{\mathrm{IL}}}\left(\frac{4\pi}{\lambda_{\mathrm{IL}}}\right)^{2}\left(R_{\mathrm{LEO}}^{2}+R_{\oplus}^{2}-2R_{\mathrm{LEO}}R_{\oplus}\cos\theta\right)\right]\\\
&=1-F_{W_{\mathrm{LG}}}\left(\frac{\gamma_{\mathrm{IL}}\sigma_{\mathrm{LEO}}^{2}}{\rho_{\mathrm{IL}}^{t}G_{\mathrm{IL}}\zeta_{\mathrm{IL}}}\left(\frac{4\pi}{\lambda_{\mathrm{IL}}}\right)^{2}\left(R_{\mathrm{LEO}}^{2}+R_{\oplus}^{2}-2R_{\mathrm{LEO}}R_{\oplus}\cos\theta\right)\right).\end{split}$
(34)

As for the derivation of the second conditional probability, given the
position of LEO satellite is $(R_{\mathrm{LEO}},\theta,\varphi)$, the distance
between the LEO satellite and GEO satellite is
$d\left(R_{\mathrm{GEO}},\Theta,0;R_{\mathrm{LEO}},\theta,\varphi\right)$. The
specific explicit expression of the operator
$d\left(\cdot,\cdot,\cdot;\cdot,\cdot,\cdot\right)$ is provided in (19). The
remaining parts follow a similar procedure and, therefore omitted here.

## References

  * [1] K. An, Y. Sun, Z. Lin, Y. Zhu, W. Ni, N. Al-Dhahir, K.-K. Wong, and D. Niyato, “Exploiting multi-layer refracting RIS-assisted receiver for HAP-SWIPT networks,” _IEEE Transactions on Wireless Communications_ , May 2024, early access.
  * [2] Z. Lin, M. Lin, T. De Cola, J.-B. Wang, W.-P. Zhu, and J. Cheng, “Supporting IoT with rate-splitting multiple access in satellite and aerial-integrated networks,” _IEEE Internet of Things Journal_ , vol. 8, no. 14, pp. 11 123–11 134, Jan. 2021.
  * [3] Z. Lin, M. Lin, B. Champagne, W.-P. Zhu, and N. Al-Dhahir, “Secrecy-energy efficient hybrid beamforming for satellite-terrestrial integrated networks,” _IEEE Transactions on Communications_ , vol. 69, no. 9, pp. 6345–6360, June 2021.
  * [4] J. A. Fraire, O. Iova, and F. Valois, “Space-terrestrial integrated internet of things: Challenges and opportunities,” _IEEE Communications Magazine_ , vol. 60, no. 12, pp. 64–70, Dec. 2022.
  * [5] X. Liu and N. Ansari, “Toward green IoT: Energy solutions and key challenges,” _IEEE Communications Magazine_ , vol. 57, no. 3, pp. 104–110, Mar. 2019.
  * [6] M. A. Kishk and H. S. Dhillon, “Downlink performance analysis of cellular-based IoT network with energy harvesting receivers,” in _IEEE Global Communications Conference (GLOBECOM)_. IEEE, Dec. 2016, pp. 1–6.
  * [7] R. Ma, W. Yang, X. Guan, X. Lu, Y. Song, and D. Chen, “Covert mmwave communications with finite blocklength against spatially random wardens,” _IEEE Internet of Things Journal_ , vol. 11, no. 2, pp. 3402–3416, Jan. 2024.
  * [8] O. Kodheli, N. Maturo, S. Chatzinotas, S. Andrenacci, and F. Zimmer, “PNB-IoT via LEO satellites: An efficient resource allocation strategy for uplink data transmission,” _IEEE Internet of Things Journal_ , vol. 9, no. 7, pp. 5094–5107, Sep. 2021.
  * [9] Z. Qu, G. Zhang, H. Cao, and J. Xie, “LEO satellite constellation for internet of things,” _IEEE access_ , vol. 5, pp. 18 391–18 401, Aug. 2017.
  * [10] H. Huang, S. Guo, W. Liang, and K. Wang, “Online green data gathering from GEO-distributed IoT networks via LEO satellites,” in _IEEE International Conference on Communications (ICC)_. IEEE, July 2018, pp. 1–6.
  * [11] A. K. Dwivedi, S. P. Chokkarapu, S. Chaudhari, and N. Varshney, “Performance analysis of novel direct access schemes for LEO satellites based IoT network,” in _IEEE 31st Annual International Symposium on Personal, Indoor and Mobile Radio Communications_. IEEE, Oct. 2020, pp. 1–6.
  * [12] R. Wang, M. A. Kishk, and M.-S. Alouini, “Reliability analysis of multi-hop routing in multi-tier LEO satellite networks,” _IEEE Transactions on Wireless Communications_ , Mar. 2024, early Access.
  * [13] H.-C. Lim, J. U. Park, M. Choi, C.-S. Choi, J.-D. Choi, and J. Kim, “Performance analysis of DPSK optical communication for LEO-to-ground relay link via a GEO satellite,” _Journal of Astronomy and Space Sciences_ , vol. 37, no. 1, pp. 11–18, Jan. 2020.
  * [14] J. Chen, X. Di, R. Xu, H. Qi, L. Cong, K. Zhang, Z. Xing, X. He, W. Lei, and S. Zhang, “A remote sensing data transmission strategy based on the combination of satellite-ground link and GEO relay under dynamic topology,” _Future Generation Computer Systems_ , vol. 145, pp. 337–353, Aug. 2023.
  * [15] X. Lin, H. Zhang, G. Pan, S. Wang, and J. An, “LEO relay-aided GEO satellite-terrestrial transmissions,” _IEEE Transactions on Vehicular Technology_ , pp. 1–7, Jul. 2023.
  * [16] Y. Ma, T. Lv, G. Pan, Y. Chen, and M.-S. Alouini, “On secure uplink transmission in hybrid RF-FSO cooperative satellite-aerial-terrestrial networks,” _IEEE Transactions on Communications_ , vol. 70, no. 12, pp. 8244–8257, Oct. 2022.
  * [17] H. Zhang, G. Pan, S. Ke, S. Wang, and J. An, “Outage analysis of cooperative satellite-aerial-terrestrial networks with spatially random terminals,” _IEEE Transactions on Communications_ , vol. 70, no. 7, pp. 4972–4987, Jun. 2022.
  * [18] Z. Lou, B. E. Y. Belmekki, and M.-S. Alouini, “Coverage analysis of hybrid RF/THz networks with best relay selection,” _IEEE Communications Letters_ , vol. 27, no. 6, pp. 1634–1638, 2023.
  * [19] A. Talgat, M. A. Kishk, and M.-S. Alouini, “Stochastic geometry-based analysis of LEO satellite communication systems,” _IEEE Communications Letters_ , vol. 25, no. 8, pp. 2458–2462, Aug. 2021.
  * [20] Z. Lou, B. E. Y. Belmekki, and M.-S. Alouini, “HAPS in the non-terrestrial network nexus: Prospective architectures and performance insights,” _IEEE Wireless Communications_ , vol. 30, no. 6, pp. 52–58, 2023.
  * [21] R. Wang, M. A. Kishk, and M.-S. Alouini, “Evaluating the accuracy of stochastic geometry based models for LEO satellite networks analysis,” _IEEE Communications Letters_ , vol. 26, no. 10, pp. 2440–2444, Oct. 2022.
  * [22] N. Okati, T. Riihonen, D. Korpi, I. Angervuori, and R. Wichman, “Downlink coverage and rate analysis of low Earth orbit satellite constellations using stochastic geometry,” _IEEE Transactions on Communications_ , vol. 68, no. 8, pp. 5120–5134, Aug. 2020.
  * [23] C. C. Chan, B. Al Homssi, and A. Al-Hourani, “A stochastic geometry approach for analyzing uplink performance for IoT-over-satellite,” in _ICC 2022-IEEE International Conference on Communications_. IEEE, May. 2022, pp. 2363–2368.
  * [24] A. Talgat, M. A. Kishk, and M.-S. Alouini, “Stochastic geometry-based uplink performance analysis of IoT over LEO satellite communication,” _IEEE Transactions on Aerospace and Electronic Systems_ , Oct. 2024, early access.
  * [25] T. Hong, X. Yu, Z. Liu, X. Ding, and G. Zhang, “Narrowband internet of things via low earth orbit satellite networks: An efficient coverage enhancement mechanism based on stochastic geometry approach,” _Sensors_ , vol. 24, no. 6, p. 2004, Mar. 2024.
  * [26] R. Wang, A. Talgat, M. A. Kishk, and M.-S. Alouini, “Conditional contact angle distribution in LEO satellite-relayed transmission,” _IEEE Communications Letters_ , vol. 26, no. 11, pp. 2735–2739, Nov. 2022.
  * [27] W. Feller, _An introduction to probability theory and its applications, Volume 2_. John Wiley & Sons, 1991, vol. 81.
  * [28] R. Wang, M. A. Kishk, and M.-S. Alouini, “Ultra-dense LEO satellite-based communication systems: A novel modeling technique,” _Communications Magazine_ , vol. 60, no. 4, pp. 25–31, Apr. 2022.
  * [29] D.-H. Na, K.-H. Park, Y.-C. Ko, and M.-S. Alouini, “Performance analysis of satellite communication systems with randomly located ground users,” _IEEE Transactions on Wireless Communications_ , vol. 21, no. 1, pp. 621–634, Jul. 2021.
  * [30] A. Abdi, W. C. Lau, M.-S. Alouini, and M. Kaveh, “A new simple model for land mobile satellite channels: First-and second-order statistics,” _IEEE Transactions on Wireless Communications_ , vol. 2, no. 3, pp. 519–528, May. 2003.
  * [31] P. Gopal, V. Jain, and S. Kar, “Modulation techniques used in Earth-to-satellite and inter-satellite free space optical links,” in _Unmanned/Unattended Sensors and Sensor Networks X_ , vol. 9248. SPIE, Oct. 2014, pp. 210–219.
  * [32] Y. Ata and M.-S. Alouini, “Performance of integrated ground-air-space FSO links over various turbulent environments,” _IEEE Photonics Journal_ , vol. 14, no. 6, pp. 1–16, Dec. 2022.
  * [33] R. Wang, M. A. Kishk, and M.-S. Alouini, “Stochastic geometry-based low latency routing in massive LEO satellite networks,” _IEEE Transactions on Aerospace and Electronic Systems_ , vol. 58, no. 5, pp. 3881–3894, Oct. 2022.
  * [34] G. Kaddoum and F. Gagnon, “Design of a high-data-rate differential chaos-shift keying system,” _IEEE Transactions on Circuits and Systems II: Express Briefs_ , vol. 59, no. 7, pp. 448–452, June 2012.
  * [35] N. I. Ziedan and J. L. Garrison, “Bit synchronization and doppler frequency removal at very low carrier to noise ratio using a combination of the Viterbi algorithm with an extended Kalman filter,” in _Proceedings of the 16th International Technical Meeting of the Satellite Division of The Institute of Navigation (ION GPS/GNSS 2003)_ , Sep. 2003, pp. 616–627.


# Probing barrier transmission in ballistic graphene

Daniel Gunlycke<EMAIL_ADDRESS>Naval Research Laboratory,
Washington, D.C. 20375, USA  Carter T. White  Naval Research Laboratory,
Washington, D.C. 20375, USA

###### Abstract

We derive the local density of states from itinerant and boundary states
around transport barriers and edges in graphene and show that the itinerant
states lead to mesoscale undulations that could be used to probe their
scattering properties in equilibrium without the need for lateral transport
measurements. This finding will facilitate vetting of extended structural
defects such as grain boundaries or line defects as transport barriers for
switchable graphene resonant tunneling transistors. We also show that barriers
could exhibit double minima and that the charge density away from highly
reflective barriers and edges scales as $x^{-2}$.

###### pacs:

72.80.Vp, 73.20.At, 73.23.Ad, 73.63.Bd

One challenge currently preventing widespread use of graphene in
nanoelectronic devices is the absence of a band gap at the Fermi level.
Without a practical band gap, other ways to switch on and off electron and
hole currents are needed. A promising possibility is to use graphene transport
barriers formed by extended structural defects such as grain boundaries Malola
_et al._ (2010); Yazyev and Louie (2010a, b); Huang  _et al._ (2011); Yu  _et
al._ (2011); Kim  _et al._ (2011); Tapasztó  _et al._ (2012); Ahmad  _et al._
(2012); Koepke  _et al._ (2013); Clark  _et al._ (2013); Tison  _et al._
(2014) or line defects.Appelhans  _et al._ (2010); Lahiri  _et al._ (2010);
Gunlycke and White (2011); Chen  _et al._ (2014) It has been shown that two
such barriers in a parallel configuration produces a graphene resonant
tunneling transistor with an appreciable transport gap and perfect valley
filtering.Gunlycke  _et al._ (2013) This approach, however, requires transport
barriers that are both penetrable and fairly reflective. In vetting potential
candidates, it would be advantageous to be able to probe the barrier
transmissivity without having to perform lateral transport measurements.

While it might seem impossible to probe the conduction through a barrier with
no current, we show herein that the quantum nature of the charge carriers let
us do exactly that. We combine two properties enabled by the wave-particle
duality: quantum tunneling and quantum interference. Quantum tunneling allows
carriers to transmit across narrow barriers as evanescent waves. Line defects
and grain boundaries with a limited transmission probability are examples of
such barriers in graphene.Gunlycke  _et al._ (2013); Gargiulo and Yazyev
(2014) As illustrated in Fig. 1(a), these barriers in equilibrium are
surrounded by density undulations. These mesoscale undulations are related to
Friedel oscillations from isolated impurity sites Friedel (1952); Cheianov and
Fal’ko (2006); Cheianov (2007); Hwang and Das Sarma (2008); Ádám Bácsi and
Virosztek (2010) and arise because the limited number of wave vectors allowed
by the band structure is unable to describe sharp real-space features. We show
that the undulations result from quantum interference between incoming and
outgoing waves, but not all outgoing waves—and this is the key—only the
reflected waves. The undulations on the left side of the barrier in Fig. 1(a)
are therefore the same as those in Fig. 1(b), which illustrates non-
equilibrium with carriers originating from the left side only. This is the
connection that allows us to probe the transmission probability through the
barrier, even in equilibrium.

Figure 1: Schematic illustrations of the relationship between the barrier
transmissivity and the undulations in the LDOS itinerant from (a) both sides
and (b) the left side, representing equilibrium and non-equilibrium,
respectively. Because carriers originating from the left and right do not
interfere, the undulations are the same in (a) and (b), except for the
transmitted side in (b), where they are absent. The transmission probability
through the barrier is related to the undulations, which could be probed in
equilibrium, even through no net current flows through the barrier.

To understand the relationship between the transmissivity and the undulations
in the local density of states (LDOS) in equilibrium, we express the LDOS at
energy $E$ and coordinate $x$ (cf. Fig. 1), centered at the barrier, as

$\rho_{E}(x)=\rho_{E}^{\rightarrow}(x)+\rho_{E}^{\leftarrow}(x)+\rho_{E}^{\mathrm{(b)}}(x),$
(1)

where the terms, respectively, represent the local densities of states
itinerant from the left and right sides of the barrier, and boundary states at
the barrier.

It can be shown that the LDOS per unit area originating from the left side is

$\rho_{E}^{\rightarrow}(x)=\frac{\rho_{E}}{4\pi}\sum_{\tau}\int_{-\pi/2}^{\pi/2}\operatorname{Tr}\,\langle\vec{r}|\Psi_{\eta\tau\vec{q}}\rangle\langle\Psi_{\eta\tau\vec{q}}|\vec{r}\rangle\mathrm{d}\alpha,$
(2)

where $\rho_{E}\equiv\frac{2|E|}{\pi(\hbar v_{F})^{2}}$ is the graphene
density of states per unit area with the graphene Fermi velocity $v_{F}\approx
8.5\times 10^{5}$ ms-1, and $\langle\vec{r}|\Psi_{\eta\tau\vec{q}}\rangle$ is
the wave function in the presence of the barrier. The wave function indices
are carrier type $\eta=\pm 1$ for electrons and holes, respectively, valley
index $\tau=\pm 1$ representing the two inequivalent graphene symmetry points
$K$ and $K^{\prime}$, and wave vector $\vec{q}=q\hat{q}$ centered at the
considered symmetry point. For elastic scattering, we have
$E=\varepsilon_{\eta\tau\vec{q}}$, where
$\varepsilon_{\eta\tau\vec{q}}=\eta\hbar v_{F}q$ assuming a linear graphene
dispersion. Because the carrier group velocity
$\vec{v}\equiv(\eta/\hbar)\nabla_{\vec{q}}\varepsilon_{\eta\tau
q}=v_{F}\hat{q}$, the wave vector is locked to the propagation angle $\alpha$,
defined in Fig. 1, so that $\vec{q}=q(\hat{x}\cos\alpha+\hat{y}\sin\alpha)$.
The scattering is specular near the Dirac point, owing to energy and wave
vector conservation along the barrier, and for the majority of extended
structural defects, these conservation laws also prohibit intervalley
scattering.Yazyev and Louie (2010b); Gunlycke and White (2014) We can express
the wave function for a state originating from the left as

$\displaystyle\langle\vec{r}|\Psi_{\eta\tau\vec{q}}\rangle$
$\displaystyle=\frac{1}{\sqrt{2}}\left(\begin{array}[]{c}1\\\ i\eta
e^{i\tau\alpha}\end{array}\right)e^{iq_{x}x}e^{iq_{y}y}$ (3c)
$\displaystyle+\frac{r_{\eta\tau\alpha}}{\sqrt{2}}\left(\begin{array}[]{c}1\\\
-i\eta e^{-i\tau\alpha}\end{array}\right)e^{-iq_{x}x}e^{iq_{y}y},$ (3f)
$\displaystyle\langle\vec{r}|\Psi_{\eta\tau\vec{q}}\rangle$
$\displaystyle=\frac{t_{\eta\tau\alpha}}{\sqrt{2}}\left(\begin{array}[]{c}1\\\
i\eta e^{i\tau\alpha}\end{array}\right)e^{iq_{x}x}e^{iq_{y}y},$ (3i)

on the left and right sides, respectively, where $r_{\eta\tau\alpha}$ and
$t_{\eta\tau\alpha}$ are the reflection and transmission amplitudes, and
phases associated with the microscopic structure has been dropped, for
clarity. Inserting Eq. (3) into Eq. (2) yields

$\rho_{E}^{\rightarrow}(x)=\left\\{\begin{array}[]{lr}\Big{.}\left(\rho_{E}/2\right)\left(1+\bar{R}\right)+\Delta\rho_{E}(x),&\quad(x<0),\\\
\Big{.}\left(\rho_{E}/2\right)\bar{T},&\quad(x>0),\end{array}\right.$ (4)

where $\bar{R}$ and $\bar{T}$ are the reflection and transmission
probabilities, $R=|r|^{2}$ and $T=|t|^{2}$, respectively, averaged over all
angles of incidences, and

$\Delta\rho_{E}(x)=\frac{\rho_{E}}{4\pi}\sum_{\tau}\int_{-\pi/2}^{\pi/2}\operatorname{Re}\Delta_{\eta\tau\alpha}\mathrm{d}\alpha,$
(5)

describes the undulations resulting from the quantum interference term
$\operatorname{Re}\Delta_{\eta\tau\alpha}$ on the incoming and reflected side,
where $\Delta_{\eta\tau\alpha}\equiv
r_{\eta\tau\alpha}\left(1-e^{-2i\tau\alpha}\right)e^{2iq|x|\cos\tau\alpha}$.
No such undulations are present on the transmitted side, where
$\rho_{E}^{\rightarrow}(x)$ is proportional to the average transmission
probability.

Let us now take advantage of present symmetry. First, we apply the parity
operator $\mathcal{P}_{x}$ to $\rho_{E}^{\rightarrow}(x)$, which leads to

$\rho_{E}^{\leftarrow}(x)=\rho_{E}^{\rightarrow}(-x).$ (6)

Next, we note that both the valley index $\tau$ and the propagation angle
$\alpha$ are odd under the parity operator $\mathcal{P}_{y}$. The reflection
amplitude, on the other hand, is even and can be shown to only depend on
$\tau$ and $\alpha$ through the product $\tau\alpha$. Because
$\Delta_{\eta\tau\alpha}$ is also even, we conclude, as one might expect, that
the two valleys contribute the same amount to the undulations. Applying the
time-reversal operator $\mathcal{T}$ to the wave function
$\langle\vec{r}|\Psi_{\eta\tau\vec{q}}\rangle$ yields
$r_{\eta,-\tau,\pi+\alpha}=r^{*}_{\eta\tau\alpha}$, and concomitantly
$\Delta_{\eta,-\tau,\pi+\alpha}=\Delta^{*}_{\eta\tau\alpha}$. This allows the
integral in Eq. (5) to be formulated as a contour integral using $z\equiv
ie^{i\tau\alpha}$. Identifying the generating function for Bessel functions of
the first kind and expressing it as a Laurent series lets us express the
undulation term as

$\Delta\rho_{E}(x)=\frac{\rho_{E}}{2}\sum_{n=-\infty}^{\infty}I_{n}[r]\,J_{n}\left(2q|x|\right),$
(7)

where the Bessel function coefficients

$I_{n}[r]\equiv\frac{1}{2\pi
i}\oint_{\mathrm{u.c.}}r_{\eta}(z)\left(z^{n-1}+z^{n-3}\right)\mathrm{d}z$ (8)

are functionals of the reflection amplitude. This functional dependence on the
reflection amplitude establishes the connection between the scattering
properties and the LDOS that enables the transmission probability to be probed
in equilibrium. Poles within the reflection amplitude could also lead to LDOS
contributions from boundary states. These contributions could be expressed as

$\rho_{E}^{\mathrm{(b)}}(x)=\frac{\rho_{E}}{4\pi
i}\int_{\mathcal{C}}r_{\eta}(z)\left(z^{-1}+z^{-3}\right)e^{q|x|(z-z^{-1})}\mathrm{d}z,$
(9)

where the contour $\mathcal{C}$ is shown in Fig. 2(a), and have a localization
length $\xi\equiv\max_{z_{p}}[q(z_{p}-z_{p}^{-1})]^{-1}$ with $z_{p}$ being
poles.

Figure 2: Contours to determine LDOS around a graphene barrier. (a) The unit
circle and real axis contours account for the itinerant and boundary states in
graphene, respectively. (b) These contributions are evaluated for a barrier
with an effective coupling $\lambda$ and effective potential $\varepsilon$.

Summing everything up, we obtain the equilibrium LDOS

$\rho_{E}(x)=\rho_{E}+\Delta\rho_{E}(x)+\rho_{E}^{\mathrm{(b)}}(x)$ (10)

from Eq. (1), where we have used $\bar{R}+\bar{T}=1$, required by carrier
conservation. The charge density has the same form

$n_{\mu}(x)=n_{\mu}+\Delta n_{\mu}(x)+n_{\mu}^{\mathrm{(b)}}(x),$ (11)

where each term is $-e$ times the integral over the corresponding term in Eq.
(10) from the Dirac point $E=0$ to the electrochemical potential $E=\mu$.

Figure 3: Local densities of states around neutral barriers ($\varepsilon=0$)
at the carrier energy $E=50$ meV that include (a) all contributions and (b)
contributions from states originating from the left. Contributions from
transmitted carriers in (b) are proportional to the average transmission
probability and do not exhibit undulations.

The expressions above follow from from the properties of graphene and the
barrier geometry. They can be used with different barrier models.Gunlycke and
White (2011); Rodrigues  _et al._ (2012); Ebert  _et al._ (2014); Gunlycke and
White (2014) For illustration, we adopt below a model Gunlycke and White
(2014) that consists of an effective coupling parameter across the barrier
$\lambda$ and an effective barrier potential parameter $\varepsilon$ [cf. Fig.
2(b)], both in units of the effective graphene nearest-neighbor hopping
parameter $\gamma=-2\hbar v_{F}/\sqrt{3}a$ with $a\approx 0.246$ nm being the
graphene lattice constant. The reflection and transmission amplitudes are then
given by

$\displaystyle r_{\eta}(z)$
$\displaystyle\equiv-1+\frac{1}{2}\left[\frac{1-z^{2}}{1+(\eta\varepsilon+\lambda)z}+\frac{1-z^{2}}{1+(\eta\varepsilon-\lambda)z}\right],$
(12a) $\displaystyle t_{\eta}(z)$
$\displaystyle\equiv\frac{1}{2z}\left[\frac{1-z^{2}}{1+(\eta\varepsilon+\lambda)z}-\frac{1-z^{2}}{1+(\eta\varepsilon-\lambda)z}\right],$
(12b)

respectively. From the reflection amplitude, we find that the undulation term
for this barrier model is

$\Delta\rho_{E}(x)=-\frac{\rho_{E}}{2}\left[J_{0}+J_{2}-\sum_{n=1}^{\infty}\frac{z_{+}^{n}+z_{-}^{n}}{2}(J_{n-2}-J_{n+2})\right],$
(13)

where the Bessel function argument $2q|x|$ is implicit and

$z_{\pm}=\left\\{\begin{array}[]{cc}\Big{.}\big{(}\eta\varepsilon\pm\lambda\big{)},&\qquad|\eta\varepsilon\pm\lambda|\leq
1,\\\
\Big{.}-\big{(}\eta\varepsilon\pm\lambda\big{)}^{-1},&\qquad|\eta\varepsilon\pm\lambda|>1.\\\
\end{array}\right.$ (14)

The LDOS contribution from the boundary states is

$\rho_{E}^{\mathrm{(b)}}(x)=\frac{\rho_{E}}{4}\sum_{\forall
z_{p}\in\\{z_{\pm}>0\\}}\left(z_{p}^{-2}-z_{p}^{2}\right)e^{q|x|(z_{p}-z_{p}^{-1})}.$
(15)

Further insight into the density undulations can be obtained from the two
barrier limits: neutral barriers ($\varepsilon=0$) and decoupled barriers
($\lambda=0$).

Figure 4: Local densities of states around decoupled barriers ($\lambda=0$) at
the carrier energy $E=50$ meV. The effective potential shifts the undulation
peaks in (a), and if positive, leads to states localized around the barrier.
The boundary states take their measure from the itinerant states, as shown in
(b) for $\varepsilon=0.8$.

Neutral barriers ($\varepsilon=0$): These barriers have smooth undulations,
owing to the cancellation of all odd order Bessel functions. This smoothness,
however, is concealed by the presence of a boundary state, except in the weak
coupling limit $|\lambda|\rightarrow 0$. In this limit, the boundary state
becomes restricted to $E=0$ and
$\Delta\rho_{E}(x)\rightarrow-(\rho_{E}/2)(J_{0}+J_{2})$ as
$|\lambda|\rightarrow 0$, as previously found for a neutral graphene
edge.Fujita  _et al._ (1996); Jaskólski  _et al._ (2011); Sasaki  _et al._
(2010) As $|\lambda|$ increases, the undulations fade, as shown in Fig. 3, and
vanish, $\Delta\rho_{E}(x)\rightarrow 0$, in the limit $|\lambda|\rightarrow
1$, where there is no barrier, and hence no quantum interference.

Decoupled barriers ($\lambda=0$): This case also describes isolated graphene
edges. The undulations are generally not smooth at the barrier and the
undulation peaks shift as a function of the potential, as can be seen in Fig.
4(a). As $\eta\varepsilon\rightarrow 0$, we approach the case of neutral
decoupled barriers or edges mentioned above. We additionally find that the
charge density undulations for these barriers or edges are given by

$\Delta
n_{\mu}(x)=-\frac{n_{\mu}}{2q_{\mu}^{2}x^{2}}\left[1-J_{0}(2q_{\mu}|x|)\right].$
(16)

The predicted quadratic decay of the charge density undulations away from the
edge is different from the cubic decay of the Friedel oscillations away from
isolated impurity sites.Cheianov and Fal’ko (2006)

Figure 5: Properties in parameter space. The derivative of the density
undulations (a) at the barrier $\Delta\rho_{E}^{\prime}(0^{+})$ correlates
with the number of boundary states (b), a result of the boundary states taking
measure from the itinerant states. The boundary states are maximally localized
at the bold cross. The average transmission probability (c) ranges from $0$
for decoupled barriers (or edges) to $1$ for nonexistent barriers.

Any quantum interference at the barrier must be destructive, which follows
from $\Delta\rho_{E}(0)=-(\rho_{E}/4)(2-z_{+}^{2}-z_{-}^{2})\leq 0$. Because
the undulations are generally not smooth at the barrier, the derivative
$\Delta\rho_{E}^{\prime}(0^{+})=(q\rho_{E}/4)(z_{+}^{3}+z_{-}^{3}-z_{+}-z_{-})$,
shown in Fig. 5(a), is not necessarily positive. In fact, as
$\Delta\rho_{E}^{\prime}(0^{+})$ depends only on odd powers of $z_{\pm}$, it
is odd under the electron-hole operator $\mathcal{C}$: $\eta\mapsto-\eta$, and
thus the signs of the derivative for electrons and holes are opposite.
Therefore, the undulations must exhibit a double minimum in half the parameter
space [cf. the itinerant LDOS contribution in Fig. 4(b)]. Furthermore, these
barriers could appear as two closely situated barriers. It is an open question
if this is related to the double-barrier features observed in graphene plasmon
experiments.Fei  _et al._ (2013)

The undulations near the barrier could be masked by boundary states. In this
barrier model, there are $0$, $1$, or $2$ boundary states present. Figure 5
shows that the number of boundary states correlates with
$\Delta\rho_{E}^{\prime}(0^{+})$ so that barriers with undulations exhibiting
double-minima, also have boundary states. This is a result of measure being
transferred from the itinerant states to the boundary states, which can be
shown in Fig. 4(b), where the sum of the LDOS contributions from itinerant and
boundary states hovers around the graphene density of states.

The LDOS in Eq. (10) and its undulations in Eq. (13) and boundary
contributions in Eq. (15) can be fitted to experiment to extract the barrier
parameters $\lambda$ and $\varepsilon$. These parameters could then be
inserted into Eq. (12) to obtain the scattering properties, including the
average transmission probabilityGunlycke and White (2014)

$\bar{T}=\left\\{\begin{array}[]{cc}\frac{\big{.}2\lambda^{2}}{\big{.}1-\varepsilon^{2}+\lambda^{2}},&\quad|\eta\varepsilon\pm\lambda|<1,\\\
\frac{\big{.}\lambda}{\big{.}\eta\varepsilon+\lambda},&\quad|\eta\varepsilon-\lambda|<1,~{}|\eta\varepsilon+\lambda|>1,\\\
-\frac{\big{.}\lambda}{\big{.}\eta\varepsilon-\lambda},&\quad|\eta\varepsilon-\lambda|>1,~{}|\eta\varepsilon+\lambda|<1,\\\
\frac{\big{.}2\lambda^{2}}{\big{(}\varepsilon^{2}-\lambda^{2}\big{)}\big{(}\varepsilon^{2}-\lambda^{2}-1\big{)}},&\quad|\eta\varepsilon\pm\lambda|>1,\\\
\end{array}\right.$ (17)

plotted in Fig. 5(c).

For highly reflective barriers with weak coupling and weak potential, we could
also use approximate expressions. From the second order Born approximation, we
obtain

$\Delta\rho_{E}(x)\approx-\frac{\rho_{E}}{2q|x|}\left[J_{1}+2\eta\varepsilon
J_{2}-\big{(}\varepsilon^{2}+\lambda^{2}\big{)}\big{(}J_{1}-3J_{3})\right].$
(18)

The corresponding average transmission is $\bar{T}=2\lambda^{2}$.

Figure 6: LDOS comparison between the graphene barrier model and the nearest-
neighbor tight-binding (TB) model at the carrier energy $E=50$ meV. There are
three sets of TB curves for the three barrier structures shown as insets, each
with a different barrier orientation with respect to the graphene lattice. The
highlighted barrier sites have a potential $-1.2$ eV. The barrier model is
computed for $(\lambda,\varepsilon)=(-1,-0.5)$.

To test our analytical approach above, we have also calculated the LDOS
numerically using an exact formalism within the tight-binding
approximation.Lee and Joannopoulos (1981) Figure 6 shows the results for three
barriers with a constant potential along a chain of sites. Although the
graphene lattice orientations for the three barriers are quite different, the
LDOS contributions are consistent and in excellent agreement with the barrier
model.

In summary, we have derived the local density of states around a generic
transport barrier described by an effective coupling parameter and an
effective potential parameter. We showed that the mesoscale undulations are
related to the average transmission probability through the barrier, therefore
making it possible to probe the transport properties through the barrier
without the need for lateral transport measurements. Rather, scanning probe
techniques could be used to estimate the transmission probability.Clark  _et
al._ (2013) This could aid the search for suitable transport barriers in
graphene, and down the road lead to new switchable graphene nanoelectronics.

###### Acknowledgements.

The authors acknowledge support from the U.S. Office of Naval Research,
directly and through the U.S. Naval Research Laboratory.

## References

  * Malola  _et al._ (2010) S. Malola, H. Häkkinen, and P. Koskinen, Phys. Rev. B 81, 165447 (2010).
  * Yazyev and Louie (2010a) O. V. Yazyev and S. G. Louie, Phys. Rev. B 81, 195420 (2010a).
  * Yazyev and Louie (2010b) O. V. Yazyev and S. G. Louie, Nat. Mat. 9, 806 (2010b).
  * Huang  _et al._ (2011) P. Y. Huang, C. S. Ruiz-Vargas, A. M. van der Zande, W. S. Whitney, M. P. Levendorf, J. W. Kevek, S. Garg, J. S. Alden, C. J. Hustedt, Y. Zhu, J. Park, P. L. McEuen, and D. A. Muller, Nature 469, 389 (2011).
  * Yu  _et al._ (2011) Q. Yu, L. A. Jauregui, W. Wu, R. Colby, J. Tian, Z. Su, H. Cao, Z. Liu, D. Pandey, D. Wei, T. F. Chung, P. Peng, N. P. Guisinger, E. A. Stach, J. Bao, S.-S. Pei, and Y. P. Chen, Nat. Mat. 10, 443 (2011).
  * Kim  _et al._ (2011) K. Kim, Z. Lee, W. Regan, C. Kisielowski, M. F. Crommie, and A. Zettl, ACS Nano 5, 2142 (2011).
  * Tapasztó  _et al._ (2012) L. Tapasztó, P. Nemes-Incze, G. Dobrik, K. J. Yoo, C. Hwang, and L. P. Biró, Appl. Phys. Lett. 100, 053114 (2012).
  * Ahmad  _et al._ (2012) M. Ahmad, H. An, Y. S. Kim, J. H. Lee, J. Jung, S.-H. Chun, and Y. Seo, Nanotechnology 23, 285705 (2012).
  * Koepke  _et al._ (2013) J. C. Koepke, J. D. Wood, D. Estrada, Z.-Y. Ong, K. T. He, E. Pop, and J. W. Lyding, ACS Nano 7, 75 (2013).
  * Clark  _et al._ (2013) K. W. Clark, X.-G. Zhang, I. V. Vlassiouk, G. He, R. M. Feenstra, and A.-P. Li, ACS Nano 7, 7956 (2013).
  * Tison  _et al._ (2014) Y. Tison, J. Lagoute, V. Repain, C. Chacon, Y. Girard, F. Joucker, R. Sporken, F. Gargiulo, O. V. Yazyev, and S. Rousset, Nano Lett. 14, 6382 (2014).
  * Appelhans  _et al._ (2010) D. J. Appelhans, L. D. Carr, and M. T. Lusk, New Journal of Physics 12, 125006 (2010).
  * Lahiri  _et al._ (2010) J. Lahiri, Y. Lin, P. Bozkurt, I. I. Oleynik, and M. Batzill, Nat. Nanotech. 5, 326 (2010).
  * Gunlycke and White (2011) D. Gunlycke and C. T. White, Phys. Rev. Lett. 106, 136806 (2011).
  * Chen  _et al._ (2014) J.-H. Chen, G. Autès, N. Alem, F. Gargiulo, A. Gautam, M. Linck, C. Kisielowski, O. V. Yazyev, S. G. Louie, and A. Zettl, Phys. Rev. B 89, 121407 (2014).
  * Gunlycke  _et al._ (2013) D. Gunlycke, S. Vasudevan, and C. T. White, Nano Letters 13, 259 (2013).
  * Gargiulo and Yazyev (2014) F. Gargiulo and O. V. Yazyev, Nano Lett. 14, 250 (2014).
  * Friedel (1952) J. Friedel, Phil. Mag. 43, 153 (1952).
  * Cheianov and Fal’ko (2006) V. V. Cheianov and V. I. Fal’ko, Phys. Rev. Lett. 97, 226801 (2006).
  * Cheianov (2007) V. V. Cheianov, Eur. Phys. J. Special Topics 148, 55 (2007).
  * Hwang and Das Sarma (2008) E. H. Hwang and S. Das Sarma, Phys. Rev. Lett. 101, 156802 (2008).
  * Ádám Bácsi and Virosztek (2010) Ádám Bácsi and A. Virosztek, Phys. Rev. B 82, 193405 (2010).
  * Gunlycke and White (2014) D. Gunlycke and C. T. White, Phys. Rev. B 90, 035452 (2014).
  * Rodrigues  _et al._ (2012) J. N. B. Rodrigues, N. M. R. Peres, and J. M. B. Lopes dos Santos, Phys. Rev. B 86, 214206 (2012).
  * Ebert  _et al._ (2014) D. Ebert, V. C. Zhukovsky, and E. A. Stepanov, J. Phys. Cond. Mat. 26, 125502 (2014).
  * Fujita  _et al._ (1996) M. Fujita, K. Wakabayashi, K. Nakada, and K. Kusakabe, J. Phys. Soc. Jpn. 65, 1920 (1996).
  * Jaskólski  _et al._ (2011) W. Jaskólski, A. Ayuela, M. Pelc, H. Santos, and L. Chico, Phys. Rev. B 83, 235424 (2011).
  * Sasaki  _et al._ (2010) K. Sasaki, K. Wakabayashi, and T. Enoki, New J. Phys. 12, 083023 (2010).
  * Fei  _et al._ (2013) Z. Fei, A. S. Rodin, W. Garnett, S. Dai, W. Regan, M. Wagner, M. K. Liu, A. S. McLeod, G. Dominguez, M. Thiemens, A. H. C. Neto, F. Keilmann, A. Zettl, R. Hillenbrand, M. M. Fogler, and D. N. Basov, Nat. Nanotech. 8, 821 (2013).
  * Lee and Joannopoulos (1981) D. H. Lee and J. D. Joannopoulos, Phys. Rev. B 23, 4988 (1981).