* [AK91] A. Albano and S. Katz, _Lines on the Fermat quintic threefold and the infinitesimal generalized Hodge conjecture_ , Trans. Amer. Math. Soc. 324(1), 353-368 (1991).
  * [Art69a] M. Artin, _Algebraic approximation of structures over complete local rings_ , Publ. Math-Paris. 36, 23-58 (1969).
  * [Art69b] M. Artin, _Algebraization of formal moduli: I_ , Global Analysis (Papers in Honor of K. Kodaira), University of Tokyo Press, 21-71 1969, ISBN:978-1-4008-7123-0.
  * [BO01] A. Bondal and D. Orlov, _Reconstruction of a variety from the derived category and groups of autoequivalences_ , Compos. Math. J. 125, 327-344 (2001).
  * [BV03] A. Bondal and M. Van den Bergh, _Generators and representability of functors in commutative and noncommutative geometry_ , Mosc. Math. J. 3(1), 1-36 (2003).
  * [BFN10] D. Ben-zvi, J. Francis, and D. Nadler, _Integral transforms and Drinfeld centers in derived algebraic geometry_ , J. Am. Math. Soc. 23(4), 909-966 (2010).
  * [Coh] L. Cohn, _Differential graded categories are k-linear stable infinity categories_ , arXiv:1308.2587
  * [Căl00] A. Căldăraru, _Derived categories of twisted sheaves on Calabi–Yau manifolds_ , Ph.D. thesis, Cornell University (2000).
  * [CNS] A. Canonaco, A. Neeman and P. Stellari, _Uniqueness of enhancements for derived and geometric categories_ , arXiv:2101.04404
  * [DL13] O. D. Deken and W. Lowen _On deformations of triangulated models_ , Adv. Math 243, 330-374 (2013).
  * [DLL17] H. Dinh Van, L. Liu and W. Lowen, _Non-commutative deformations and quasi-coherent modules_ , Sel. Math. New Ser. 23, 1061–1119 (2017).
  * [Dri04] V. Drinfeld, _DG quotients of DG categories_ , J. Algebra 272(2), 643-691 (2004).
  * [GD61] A. Grothendieck and J. Dieudonné, _Elements de géométrie algebriqué: III. Étude cohomologique des faisceaux cohérent, première partie._ , Publ. Math-Paris. 11, 5-167 (1961).
  * [GJ] E. Getzler and J.D.S. Jones, _Operads, homotopy algebra and iterated integrals for double loop spaces_ , arXiv:hep-th/9403055
  * [GS87] M. Gerstenhaber and S. D. Schack, _A Hodge-type decomposition for commutative algebra cohomology_ , J. Pure. Appl. Algebra. 48(3) 229-247 (1987).
  * [GS88] M. Gerstenhaber and S. D. Schack, _Algebraic cohomology and Deformation theory_ , Deformation Theory of Algebras and Structures and Applications. NATO ASI Series (Series C: Mathematical and Physical Sciences), 247, Springer, https://doi.org/10.1007/978-94-009-3057-5$\\_2$
  * [Har77] R. Hartshorne, _Algebraic geometry_ , Graduate Texts in Mathematics, 52, Springer-Verlag, 1977, ISBN: 0-387-90244-9.
  * [Har10] R. Hartshorne, _Deformation theory_ , Graduate Texts in Mathematics, 257, Springer-Verlag, 2010, ISBN: 978-1-4419-1596-2.
  * [HMS11] D. Huybrechts, E. Macrì, and P. Stellari, _Formal deformations and their categorical general fibre_ , Comment. Math. Helv. 86(1), 41-71 (2011).
  * [Kel] B. Keller, _Derived invariance of higher structures on the Hochschild complexes_ , http://www.math.jussieu.fr/ $\tilde{}$ keller/publ/dih.dvi
  * [KL09] B. Keller and W. Lowen, _On hochschild cohomology and morita deformations_ , Int. Math. Res. Not. 17, 3221–3235 (2009)
  * [Lie06] M. Lieblich, _Moduli of complexes on a proper morphism_ , J. Algebraic Geom. 15, 175-206 (2006).
  * [Low08] W. Lowen, _Hochschild cohomology, the characteristic morphism and derived deformations_ , Compos. Math. 144(6), 1557-1580 (2008).
  * [LV06a] W. Lowen and M. Van den Bergh, _Hochschild cohomology of abelian categories and ringed spaces_ , Adv. Math. 198(1), 172-221 (2006).
  * [LV06b] W. Lowen and M. Van den Bergh, _Deformation theory of abelian categories_ , Trans. Amer. Math. Soc. 358(12), 5441-5483 (2006).
  * [LST13] A. C. López Martín, D. Sánchez Gómez, and C. Tejero Prieto, _Relative Fourier–Mukai transforms for Weierstrass fibrations, abelian schemes and Fano fibrations_ , Math. Proc. Cambridge. 155(1), 129-153 (2013).
  * [Lur09] J. Lurie, _Higher topos theory_ , Annals of Mathematics Studies, 170, Princeton University Press, 2009, ISBN: 978-0-6911-4049-0
  * [Man09] M. Manetti, _Differential graded Lie algebras and formal deformation theory_ , Proceedings of Symposia in Pure Mathematics, 80(2), American Mathematical Society, Providence, RI, 2009. ISBN: 978-0-8218-4703-9
  * [Mora] H. Morimura, _Algebraic deformations and Fourier–Mukai transforms for Calabi–Yau manifolds_ , to appear in Proc. Amer. Math. Soc.
  * [Morb] H. Morimura, _Categorical generic fiber_ , arXiv:2111.00239
  * [Por10] M. Porta, _The Popescu–Gabriel theorem for triangulated categories_ , Adv. Math. 225(3), 1669-1715 (2010).
  * [Qui70] D. Quillen, _On the (co-)homology of commutative rings_ , Applications of categorical algebra (Proc. Sympos. Pure. Math., Vol. XVII, New York, 1968), Amer. Math. Soc., Providence, R.I., 65-87 (1970)
  * [Shu61] U. Shukla, _Cohomologie des algèbres associatives_ , Ann. Sci. École Norm. Sup. 78(3), 163-209 (1961).
  * [STZ14] N. Sibilla, D. Treumann, and E. Zaslow _Ribbon Graphs and Mirror Symmetry I_ , Selecta Mathematica (N.S.) 20(4), 979-1002 (2014).
  * [Tab05a] G. Tabuada, _Une structure de catégorie de modèles de Quillen sur la catégorie des dg-catégories_ , Comptes Rendus de lAcadadémie de Sciences de Paris. 340, 15–19 (2005).
  * [Tab05b] G. Tabuada, _Invariants additifs de DG-catégories_ , Int. Math. Res. Not. 53, 3309-3339 (2005).
  * [Tod09] Y. Toda, _Deformations and Fourier–Mukai transforms_ , J. Differ. Geom. 81(1), 197-224 (2009).
  * [Töe] B. Töen, _Lectures on DG-categories_ , http://www.math.univ-toulouse.fr/ $\tilde{}$ toen/


# CubeSat-Enabled Free-Space Optics:
Joint Data Communication and Fine Beam Tracking

Hossein Safi, Mohammad Taghi Dabiri, Julian Cheng, Fellow, IEEE,
Iman Tavakkolnia, Member, IEEE, and Harald Haas, Fellow, IEEE H. Safi, I.
Tavakkolnia, and H. Haas are with the Electrical Division, Department of
Engineering, University of Cambridge, Cambridge, UK, e-mails: {hs905, it360,
<EMAIL_ADDRESS>M. T. Dabiri is with the Department of Electrical
Engineering, Qatar University, Doha, Qatar, e-mail: (m.dabiri@qu.edu.qa). J.
Cheng is with the School of Engineering, The University of British Columbia,
Kelowna, BC, Canada (e-mail: julian.cheng@ubc.ca).
This work is a contribution by Project REASON, a UK Government funded project
under the Future Open Networks Research Challenge (FONRC) sponsored by the
Department of Science Innovation and Technology (DSIT).

###### Abstract

The integration of CubeSats with Free Space Optical (FSO) links accelerates a
major advancement in high-throughput, low-Earth orbit communication systems.
However, CubeSats face challenges such as size, weight, and power (SWaP)
limitations, as well as vibrations that cause fluctuations in the angle-of-
arrival (AoA) of the optical beam at the receiver. These practical challenges
make establishing CubeSat-assisted FSO links complicated. To mitigate AoA
fluctuations, we expand the receiver’s field of view and track the location of
the focused beam spot using an array of avalanche photodiodes at the receiver.
Initially, we model the optical channel between the transmitter and the
detector array. Furthermore, to reduce the computational load of maximum
likelihood sequence detection, which is infeasible for CubeSats due to SWaP
constraints, we propose a sub-optimal blind sequence data detection approach
that relies on the generalized likelihood ratio test (GLRT) criterion. We also
utilize combining methods such as equal gain combining (EGC) and maximal ratio
combining (MRC) for data detection, benchmarking their performance against the
GLRT-based method. Numerical results demonstrate that the proposed low-
complexity GLRT-based method outperforms the combining methods, achieving
performance close to that of the ideal receiver.

###### Index Terms:

Angle of arrival fluctuations, blind data detection, CubeSats, detectors,
free-space optics.

##  I Introduction

In the sixth-generation (6G) era, the race to offer and commercialize high
data rate services faces challenges due to limited radio frequency (RF)
spectrum availability. As a result, there is a growing interest in exploiting
large blocks of spectrum in higher frequency bands to accommodate the need for
increased peak rates and capacity in wireless communication links [1]. Free-
space optical (FSO) communication has emerged as a promising solution to
address these challenges and enhance modern communication systems [2]. Unlike
RF communication, FSO uses optical wavelengths to transmit data through free
space, offering advantages such as low power consumption, significantly higher
unlicensed bandwidth, and immunity to electromagnetic interference. In
satellite communication, where traditional commercial bandwidth allocations
(i.e., S, X, Ku, and Ka bands) struggle to meet the demands of data-intensive
6G and beyond applications [3], FSO technology offers the potential for higher
data rates, enhanced security, and reduced size, weight and power (SWaP) due
to smaller antenna requirements [4]. These characteristics make FSO-equipped
satellite systems well-suited for a wide range of applications [5].

Interestingly, the space industry has recently witnessed a remarkable shift
from traditional large and costly spacecraft to smaller, more accessible
missions [6]. This transformation was catalyzed by advancements in commercial-
off-the-shelf (COTS) technology, which enabled the miniaturization of
spacecraft components [7]. Accordingly, a class of nanosatellites, also known
as CubeSats, has emerged as a cost-effective platform. Their compact size, low
cost, and rapid development timelines have led to their adoption for a variety
of applications. CubeSats typically operate in low earth orbit (LEO) to offer
low latency communication, reduced lunch cost, and to support rapid response
missions. They also adhere to standardized dimensions, measured in “units”
(U), where 1U corresponds to a cube measuring 10 $\times$ 10 $\times$ 10
centimetres. Notably, the global CubeSat market, valued at USD 461 million in
2023, is projected to reach USD 1445.04 million by 2031, with a compound
annual growth rate of 15.34% [8]. Consequently, major companies like Google,
SpaceX, OneWeb, and Facebook have shown interest in CubeSats for various
applications, including earth monitoring, disaster prevention, and providing
connectivity to internet of things devices in remote areas [7].

The capabilities of CubeSat are, however, limited by SWaP limitations,
particularly for high-bandwidth communications. Furthermore, the vast
potential of many CubeSat missions, is constrained not only by their
communication capabilities but also by regulatory burdens in obtaining
licenses for RF spectrum usage [9]. The integration of CubeSats with FSO links
represents a significant advancement in space-based communication systems by
offering access to huge unlicensed bandwidth [10, 11]. Indeed, FSO technology
finds a natural synergy with the compact size and agility of CubeSats. By
leveraging FSO links, CubeSats can establish high-speed, line-of-sight
communication links with ground stations or other satellites, enabling real-
time high-throughput data transmission without interfering with RF-based
terrestrial networks. This capability is particularly valuable for
applications requiring high-bandwidth communication, such as Earth
observation, disaster monitoring, and remote sensing. Early studies by the Jet
Propulsion Laboratory (JPL) have identified CubeSat laser communications as
pivotal for a diverse range of exploratory missions [12]. Significant efforts
are currently underway in numerous prominent research organizations around the
globe to miniaturize next-generation space-based laser communication systems,
while simultaneously achieving greater data transmission rates and enhancing
CubeSat capabilities [11, 13].

Although the integration of FSO links into CubeSats opens up new opportunities
for advanced capabilities, FSO-based CubeSat communication encounters several
challenges that need to be carefully considered when designing link parameters
and evaluating the system performance [7, 10, 4, 14]. In both uplink and
downlink scenarios, the transmitted optical beam propagates through Earth’s
atmosphere to reach the receiver (Rx), where multiple attenuations impact the
amount of received power. Furthermore, the random nature of atmospheric
turbulence induces signal fading, resulting in random fluctuations in received
optical power. Additionally, background noise from the Sun and other stars
makes it challenging to accurately detect transmitted data. Importantly, due
to the narrow beamwidth of lasers, establishing and maintaining optical links
becomes challenging for long-range FSO systems with moving terminals. Here,
the optical Rx subsystem mounted on the CubeSat must initially acquire the
pointed narrow beam laser signal from the optical transmitter (also referred
to as acquisition), before establishing the link. Once successful acquisition
occurs, the Rx must continuously track the angle-of-arrival (AoA) of the
incoming beam to maintain link alignment (also referred to as fine beam
tracking). Consequently, the successful implementation of FSO-based CubeSat
communication heavily relies on the performance of pointing, acquisition, and
tracking (PAT) systems. Furthermore, the power and payload constraints of
CubeSats limit the size and capabilities of FSO systems that can be
integrated, requiring innovative solutions to achieve efficient and reliable
communication. Overcoming these challenges will be crucial for realizing the
full potential of FSO-based CubeSat communications and harnessing its benefits
for a wide range of space-based applications.

###  I-A Major Contributions and Novelty

In this paper, we study the problem of joint fine beam tracking and optical
signal detection for a ground-to-CubeSat FSO link. To this end, we first
derive a novel analytical channel model for the considered link when the
CubeSat is equipped with a converging lens and a photodetector array with an
arbitrary size. Because the CubeSat is a mobile platform and moves from one
point to another during the communication time interval, the coherence time of
the underlying channel is reduced. Therefore, compared to fixed FSO
communications, we need to perform channel estimation at shorter intervals to
ensure a reliable communication. Under the circumstances, using pilot symbols
for channel estimation results in a signaling overhead and spectral efficiency
loss. Thus, it seems reasonable to suggest the use of blind signal detection
for CubeSat-based FSO links. Consequently, we consider a practical scenario in
which the Rx has no information about the instantaneous channel coefficient
and blindly detect data at the Rx. In particular, data detection is performed
based on the assumption that avalanche photodiodes (APDs) are employed at the
array. APDs are known for their high sensitivity and ability to detect weak
optical signals, which is crucial for maintaining reliable communication in
challenging conditions. We take a shot-noise limited scenario into account
[15] (i.e., the distribution of the outputs of the APD varies as the channel
states change), and thus the performance analysis of the considered FSO system
is more challenging than those with positive-intrinsic-negative (PIN) photo-
diodes and thermal or background noise limited scenarios.

In summary, our key contributions encompass the following aspects.

  * •

We first develop a model for the optical channel between the optical source at
the Tx and the Rx detector array. This model comprehensively considers various
factors, including atmospheric turbulence and attenuation, geometrical loss,
CubeSat vibrations, the number and size of photodetectors in the array, and
the junction width between two adjacent detectors in the array.

  * •

Utilizing the derived channel model, we proceed to assess the link performance
through comprehensive mathematical analysis. Our findings demonstrate that the
angular instability of the CubeSat significantly influences system
performance. Furthermore, we illustrate that employing an array of detectors
can alleviate the adverse effects of AoA fluctuations and enhance system
performance. These analytical results serve as a benchmark for determining the
optimal number of APDs in the array to attain desired performance levels under
varying degrees of instability without resorting to time-consuming
simulations.

  * •

We propose a sub-optimal blind sequence data detection method based on the
generalized likelihood ratio test (GLRT) criterion, which detects data over a
sequence of transmitted optical signals with length $L$ bits. Our analysis
demonstrates that when $L$ is sufficiently large, the GLRT-based method can
achieve performance close to that of the ideal Rx. Furthermore, although a
larger value of $L$ is required to achieve performance comparable to the ideal
receiver, our proposed method significantly reduces the computational load
compared to the maximum likelihood (ML)-based method, which is on the order of
$2^{L}$. This reduction in computational complexity translates to lower energy
consumption and allows for the use of simpler processing units in the CubeSat
payload. We also employ two combining techniques at the Rx, i.e., equal gain
combining (EGC) and maximal ratio combining (MRC). However, our findings
indicate that neither combining method can achieve performance comparable to
the GLRT-based method.

  * •

Finally, to fully leverage the benefits of utilizing a photodetector array at
the Rx, we estimate the position of the beam hot spot by comparing the output
of APDs. Subsequently, the instantaneous orientation of the mounted Rx with
respect to the arrival direction of the received optical beam can be
determined. This information can be fed back through a control message to the
mechanical subsystem of the CubeSat (e.g., fast steering mirror) to correct
any instantaneous orientation errors, allowing for fine-tracking of the beam
while simultaneously performing data communication.

The rest of this paper is structured as follows. Section II presents the prior
art. In Section III, we detail the system model of the ground-to-CubeSat FSO
link. Section IV outlines our proposed signal detection and beam tracking
schemes when utilizing an APD array at the Rx. Following that, Section V shows
the numerical results of the proposed schemes and system performances.
Finally, Section VI provides the conclusion of this paper.

##  II Prior Art

Exploring the integration of FSO communication with CubeSats opens up a new
frontier in research, particularly in system design, managing power efficiency
and reliable data communication, that needs to be further investigated. Prior
research in this context can be categorized into two main directions: i)
design and prototyping of different optical components mounted on CubeSats, as
demonstrated in [13, 10, 11, 16, 17, 18, 19, 20], and ii) data communication,
as evidenced by [21, 22, 23, 24, 25, 26, 27].

###  II-A Advancements in CubeSat Laser Technology

The latest advancements in technology development and the system design of
CubeSat laser communication terminals (LCTs) were demonstrated in [13, 10,
11]. For instance, NASA’s optical communication and sensor demonstration
(OCSD) was the first CubeSat laser payload mission in which laser terminals
were mounted on a 1.5U CubeSat weighing 2.5 kg [16]. Also, the CubeSat laser
infrared crosslink (CLICK) mission showcased terminals capable of conducting
full-duplex, high-data-rate crosslinks and enabling precise ranging on 3U
CubeSats in LEO. The CLICK-A payload had a mass of less than 1.2 kg and fit
within a 1.2U envelope [17]. Most of these LCTs were designed with the aim of
reducing the cost and weight of the CubeSat payload. Consequently, it becomes
feasible to utilize off-the-shelf components in ground-based optical networks,
such as detector arrays, which are commercially available [18, 19]. However,
challenges such as AoA fluctuations arise due to CubeSat vibrations resulting
from imperfectly stabilized commercial products, which often exhibit minor
mechanical instabilities and less precise attitude control compared to custom-
designed systems [7, 4]. Additionally, CubeSats rely on a combination of
gyroscopes, magnetometers, accelerometers, sun position sensors, Earth horizon
sensors, or star trackers to determine their orientation with limited
accuracy, primarily due to SWaP constraints, which induce pointing errors. In
addressing these challenges, detector arrays, typically composed of APDs for
long-range satellite links, play a crucial role. They enable the acceptance of
the received optical beam across a wider field of view (FoV) to compensate for
AoA fluctuations. Moreover, they facilitate received optical beam position
sensing (AoA estimation) for fine beam tracking [19, 28, 20].

###  II-B Data Communication

In the literature concerning CubeSat FSO data communication, researchers have
explored various aspects to enhance the capabilities of these miniature
satellite systems across different link configurations [21, 22, 23, 24, 25,
26, 27]. For instance, the study in [21] examined the capacity and error rate
performance of an integrated transmission RF/FSO system over the uplink
scenario where a high-altitude platform (HAP) is deployed as a relay between a
ground station and the LEO satellite. Furthermore, the authors extended their
findings from [21] by incorporating various channel parameters such as beam
divergence loss, free space loss, and atmospheric attenuation. They analyzed
the end-to-end outage probability and average symbol error probability
performance in [22]. Moreover, the tradeoff between satellite transmission
power and network latency for different laser inter-satellite links was
explored in [23]. Addressing physical layer security concerns, the authors in
[24] focused on FSO-based satellite communication in LEO. Meanwhile, high-
altitude ground stations as a transformative element towards an all-optical
satellite megaconstellation was proposed in [26]. The authors developed a
comprehensive analytical model-based cross-layer approach for evaluating
transmission control protocol (TCP) performance in the scenario of FSO-based
satellite-assisted internet of vehicles [25]. More recently, the optimal
design of a modulated retroreflector (MRR) laser link to establish a high-
speed downlink for CubeSat FSO communication was studied in [27].

###  II-C Rationale for Our Research

Although the idea of deploying a detector array for beam tracking has been
proposed in the literature, most of them are limited to 2$\times$2 arrays,
also known as quad detectors, except for a few recent works [29]. However,
these works only consider the problem of fine beam tracking with the ideal
assumption of neglecting the effect of the junction width of the detector in
the array (dead space), thus failing to fully capture the actual detector
characteristics in the system model. Furthermore, the aforementioned works
primarily focus on the channel between the Tx and receiver Rx lens. Given that
CubeSat vibrations cause shifted diffraction patterns at the photodetector
array [4], which can attenuate the received optical power, it is essential to
consider an additional channel parameter to address the channel coefficient
between the receiver lens and the detector array. Additionally, these prior
works mainly concentrate on the problem of optical beam tracking using
detector arrays, while optimal signal detection in this context when using
detector arrays at the Rx has yet to be investigated. Maximizing the CubeSat’s
battery life is crucial for keeping replacement costs low, as a satellite’s
mission duration depends heavily on its battery performance. Therefore, it is
vital for detection algorithms to achieve desirable performance with reduced
complexity, thereby minimizing energy consumption and addressing SWaP
limitations.

##  III System Model

We consider a ground-to-CubeSat FSO link with the propagation length
$Z=\frac{H_{s}-H_{g}}{\cos(\xi)}$, where $H_{s}$ is the sattelite altitude,
$H_{g}$ is the ground node altitude, and $\xi$ is the satellite’s zenith
angle. The Tx and Rx nodes are positioned at coordinates ($0,0,0$) and
($0,0,Z$), respectively, in a Cartesian coordinate system ($x,y,z$)
$\in\mathbb{R}^{1\times 3}$. The optical beam is propagated over the $z$-axis.
Here, we assume that coarse pointing and acquisition has been achieved, and
the coarse control loop is closed, allowing the ground station to locate and
track the satellite in space [30, 31, 32]. However, the mechanical vibrations
of the satellite cause fluctuations in the Rx aperture away from the boresight
direction along the $z$-axis. This results in performance degradation due to
fluctuations in the AoA of the received optical beam at the photodetector
plane. Hence, a fine tracking stage is required to supplement the coarse stage
in order to maintain the stability of the link. We employ an array of APDs at
the Rx to increase its FoV to combat the adverse effects of AoA fluctuations.
We also simultaneously perform fine beam tracking and data detection. The
details of the proposed schemes are provided in Section IV. In the sequel, we
first model the FSO channel and then present the signal model of the
considered FSO setup.

###  III-A Channel Model

During transmitting the optical signal by a ground node toward the CubeSat,
the received signal must be collected by the CubeSat’s aperture lens. The
collected signal is then focused on the detector plane. Finally, the array of
APDs converts the optical signal into the electrical signal. Hence, as
depicted in Fig. 1, the channel coefficient for the considered aerial FSO link
can be divided in two different parts as

$\displaystyle h_{c}=h_{1}h_{2}$ (1)

where $h_{1}$ represents the channel coefficient between the optical
transmitter at the ground station and the receiver lens, while $h_{2}$ denotes
the channel coefficient between the receiver lens and the array of APDs.

Figure 1: Graphical illustration of a ground-to-CubeSat FSO link for 6G and
beyond wireless networks.

####  III-A1 Channel Modelling Between the ground station and CubeSat Rx Lens

Let $P_{r}^{\prime}=h_{1}P_{t}$ be the collected optical power by the receiver
lens where $P_{t}$ denotes the transmit optical power. Let
$\theta_{ex}\sim\mathcal{N}(0,\sigma_{\theta_{e}}^{2})$ and
$\theta_{ey}\sim\mathcal{N}(0,\sigma_{\theta_{e}}^{2})$ denote pointing angle
errors in the directions of $x$ and $y$ axis, respectively, with
$\sigma_{\theta_{e}}$ being the standard deviation of the angle errors [27].
These errors can be induced by the limited accuracy of the CubeSat’s position
and orientation sensors, as well as beam wander resulting from inhomogeneities
in air temperature and pressure [33]. The radial displacement between the
received beam center and the receiver aperture center can be represented as
$d_{d}=\sqrt{d_{dx}^{2}+d_{dy}^{2}}$, where $d_{dx}$ and $d_{dy}$ are the
displacement in the directions of $x$ and $y$ axes, respectively, and they are
obtained as

$\displaystyle d_{dx}=Z\sin\left(\theta_{ex}\right),\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \&\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
d_{dy}=Z\sin\left(\theta_{ey}\right).$ (2)

We consider a Gaussian beam at the ground station, for which the normalized
spatial distribution of the received intensity at distance $Z$ is given by

$\displaystyle I_{r}(d,Z)=\frac{2}{\pi
w_{z}^{2}}\exp\left(-\frac{2(x^{2}+y^{2})}{w_{z}^{2}}\right),$ (3)

where $d=[x,y]$ is the radial distance vector from the beam center. Also,
$w_{z}$ is the beamwidth at distance $Z$ and can be approximated as [34, eq.
5]

$\displaystyle w_{z}\simeq w_{0}\sqrt{1+\left(\frac{\lambda Z}{\pi
w_{0}^{2}}\right)},$ (4)

where $\lambda$ is the optical wavelength and $w_{0}$ is the beamwidth at the
output of optical transmitter. For lower values of $w_{0}$, (4) can be
simplified as $w_{z}\simeq\theta_{div}Z$, where $\theta_{\textrm{div}}$ is the
transmit optical beam divergence angle.

For a receiver aperture with the effective area of $A_{r}$, the effective
channel coefficient due to the geometric spread with pointing error is
obtained as

$\displaystyle h_{p_{u}}$ $\displaystyle=\frac{2}{\pi
w_{z}^{2}}\int\int_{p_{A_{r}}(x,y)}$ $\displaystyle\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
\times\exp\left(-\frac{2((x-d_{dx})^{2}+(y-d_{dy})^{2})}{w_{z}^{2}}\right)\textrm{d}x\textrm{d}y,$
(5)

where $p_{A_{r}}(x,y)$ is the position of the effective aperture area in the
$x-y$ plane. Unlike terrestrial links, the laser beamwidth for satellite
communications is on the order of several tens of meters. In contrast, the
dimensions of the aperture in the CubeSat payload are on the order of several
centimeters. Therefore, the active area of the receiver lens is usually much
smaller than the received beamwidth $w_{z}$ in practical satellite FSO links.
From this perspective, one can conclude that the optical beam maintains its
plane wave nature locally at the receiver lens [33], and thus, (III-A1) can be
approximated as

$\displaystyle h_{p_{u}}$ $\displaystyle\approx\frac{2A_{r}}{\pi
w_{z}^{2}}\exp\left(-\frac{2(d_{dx}^{2}+d_{dy}^{2})}{w_{z}^{2}}\right).$ (6)

In addition to the aforementioned factors, the transmitted optical power from
the ground station is affected by atmospheric turbulence-induced fading and
atmospheric loss. The atmospheric attenuation is typically modeled by the
Beer-Lambert law as

$\displaystyle h_{l}=\exp(-Z\zeta),$ (7)

where $\zeta$ is the scattering coefficient and is a function of visibility.
Under moderate to strong turbulence conditions, we use the Gamma-Gamma (GG)
distribution to model the random variable $h_{a}$ representing turbulence
conditions as

$\displaystyle
f_{G}\left(h_{a}\right)=\frac{2(\alpha\beta)^{\frac{\alpha+\beta}{2}}}{\Gamma(\alpha)\Gamma(\beta)}h_{a}^{{}^{\frac{\alpha+\beta}{2}-1}}k_{\alpha-\beta}\left(2\sqrt{\alpha\beta
h_{a}}\right),$ (8)

where $\Gamma(\cdot)$ is the Gamma function and $k_{m}(\cdot)$ is the modified
Bessel function of the second kind of order $m$. Also, $\alpha$ and $\beta$
are respectively the effective number of large-scale and small-scale eddies,
which depend on Rytov variance $\sigma_{R}^{2}$ [15]. The parameter
$\sigma_{R}^{2}$ is also obtained for two nodes with different altitude as
[15, p. 509]

$\displaystyle\sigma_{R}^{2}=2.25\left({2\pi}/{\lambda}\right)^{7/6}\left(H_{s}-H_{g}\right)^{5/6}\sec(\xi)^{11/6}$
(9) $\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\
\times\int_{H_{g}}^{H_{s}}C_{n}^{2}(h)\left(1-\frac{h-H_{g}}{H_{s}-H_{g}}\right)^{5/6}\left(\frac{h-H_{g}}{H_{s}-H_{g}}\right)^{5/6}\textrm{d}h.$

In (9), $C_{n}^{2}(h)$ denotes the variation of the refractive index structure
parameter used to characterize the atmospheric turbulence described by the
Hufnagel-Valley model as

$\displaystyle
C_{n}^{2}(h)=0.00594(\mathcal{V}/27)^{2}\left(10^{-5}h\right)^{10}\exp\left(-\frac{h}{1000}\right)$
(10) $\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ +2.7\times
10^{-16}\exp\left(-\frac{h}{1500}\right)+C_{n}^{2}(0)\exp\left(-\frac{h}{100}\right),$

where $\mathcal{V}$ (in m/s) is the speed of strong wind, and $C_{n}^{2}(0)$
(in $\textrm{m}^{-2/3}$) is a strong nominal ground turbulence level. Finally,
$h_{1}$ can be obtained as

$\displaystyle h_{1}=h_{p_{u}}h_{l}h_{a}.$ (11)

Figure 2: A schematic of the deviated received beam due to the AoA
fluctuations on a $2\times 2$ APD detector.

####  III-A2 Channel Modelling Between the Rx Lens and the Photodetector Array

Here, we develop a mathematical model for the channel coefficient between the
receiver lens and the APDs. At the Rx side, a circular lens focuses both the
collected optical signal and undesired background noise, which fall within the
Rx FoV, onto the photodetector array. As discussed, vibrations in the CubeSat
cause the center of the Rx lens to deviate from the center of the received
optical beam, resulting in AoA fluctuations, as graphically illustrated in
Fig. 2. The AoA fluctuations increase the probability that the received
optical beam lies outside of the Rx FoV, significantly degrading system
performance. Even when the optical signal is collected by the aperture lens,
AoA fluctuations may cause the diffracted patterns to shift out of the
photodetector array. This can attenuate the amount of optical power received
at the Rx.

To compensate for the effect of AoA fluctuations on the link performance and
to achieve a wider FoV, we utilize an array of APDs at the Rx. Specifically,
we consider an array of size $N=N_{a}\times N_{a}$, where APDs are arranged
next to each other in a rectangular shape, as depicted in Fig. 3. Each APD has
an active width of $w_{a}$, determined by a pair of indices ($i,j$).
Additionally, there is a dead space between adjacent APDs, represented by the
junction width $w_{f}$.

Figure 3: The Graphical illustration of a $4\times 4$ photodetector array
where each photodetector is specified by $i$ and $j$ indices. Moreover,
$w_{f}$ denotes the junction width (dead space) and $w_{a}$ denotes the active
width of each photodetector.

We assume that both the Rx lens and the detector array are situated on the
$x-y$ plane, while the optical beam propagates along the $z$-axis. Let’s
define the AoA of the signal as the incidence angle relative to the Rx axis,
denoted by $\theta_{r}$. Consequently, $\theta_{r}$ can be closely
approximated as (see Fig. 2)

$\displaystyle\theta_{r}\simeq\sqrt{\theta_{x}^{2}+\theta_{y}^{2}}$ (12)

where $\theta_{x}$ and $\theta_{y}$ denote the deviations of received laser
beam in $x-z$ and $y-z$ planes, respectively. At the Rx, the converging lens
focuses the collected optical signals onto the area of the APD array placed at
the focal plane. However, due to AoA fluctuations, the center of the optical
beam at the APD area deviates, represented by
$r_{d}=\sqrt{x_{d}^{2}+y_{d}^{2}}$, where $x_{d}=f_{c}\tan(\theta_{x})\simeq
f_{c}\theta_{x}$, $y_{d}=f_{c}\tan(\theta_{y})\simeq f_{c}\theta_{y}$, and
$f_{c}$ denotes the focal length. To illustrate these parameters related to
the optical Rx, refer to the schematic of a $2\times 2$ APD array shown in
Fig. 2. The random variables (RVs) $\theta_{x}$ and $\theta_{y}$ are modeled
by zero-mean Gaussian distributions with variances $\sigma_{x}^{2}$ and
$\sigma_{y}^{2}$, respectively, and their joint probability density function
(PDF) is obtained as [35]

$\displaystyle
p_{\theta}\left(\theta_{x},\theta_{y}\right)=\frac{1}{2\pi\sigma_{x}\sigma_{y}}\exp\left(-\frac{\theta_{x}^{2}}{2\sigma_{x}^{2}}-\frac{\theta_{y}^{2}}{2\sigma_{y}^{2}}\right).$
(13)

We denote the Rx FoV in $x-z$ and $y-z$ planes by $\theta_{x\textrm{FoV}}$ and
$\theta_{y\textrm{FoV}}$, respectively. In this setup, we have
$\theta_{x\textrm{FoV}}=\theta_{y\textrm{FoV}}={\rm
arctan}\left(\frac{N_{a}w_{a}}{f_{c}}\right)$. Therefore, the Rx FoV in the
spherical coordinate system can be represented by [36]

$\displaystyle\Omega_{FoV}$
$\displaystyle=4\arccos\left(\sqrt{\dfrac{1+u^{2}+v^{2}}{(1+u^{2})(1+v^{2})}}\right)$
(14)

where $u=\frac{\tan\left(\theta_{x\textrm{FoV}}\right)}{2}$, and
$v=\frac{\tan\left(\theta_{y\textrm{FoV}}\right)}{2}$.

The intensity of the incident optical beam on the detector array can be
approximated as a two-dimensional Gaussian-shaped function [37, 38]

$\displaystyle
I_{p}(x,y)=\frac{h_{1}P_{t}}{2\pi\sigma_{I}^{2}}\exp\left(-\frac{(x-f_{c}\theta_{x})^{2}+(y-f_{c}\theta_{y})^{2}}{2\sigma_{I}^{2}}\right)$
(15)

where $\sigma_{I}^{2}$ is the variance of the intensity of optical beam on the
detector array. According to (1) and (15), the channel coefficient
corresponding to the $(i,j)$th detector can be obtained as

$\displaystyle h_{ij}=h_{1}h_{2,ij}$ (16)

where

$\displaystyle h_{2,ij}$
$\displaystyle=\int_{(i-1-\frac{N_{a}}{2})w_{a}+w_{f}/2}^{(i-\frac{N_{a}}{2})w_{a}-w_{f}/2}\int_{(j-1-\frac{N_{a}}{2})w_{a}+w_{f}/2}^{(j-\frac{N_{a}}{2})w_{a}-w_{f}/2}\frac{1}{2\pi\sigma_{I}^{2}}$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\times\exp\left(-\frac{(x-f_{c}\theta_{x})^{2}+(y-f_{c}\theta_{y})^{2}}{2\sigma_{I}^{2}}\right)\textrm{d}x\,\textrm{d}y$
$\displaystyle=\Bigg{[}Q\left(\frac{(i-1-N_{a}/2)w_{a}+w_{f}/2-f_{c}\theta_{x}}{\sigma_{I}}\right)$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
-Q\left(\frac{(i-N_{a}/2)w_{a}-w_{f}/2-f_{c}\theta_{x}}{\sigma_{I}}\right)\Bigg{]}$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\times\Bigg{[}Q\left(\frac{(j-1-N_{a}/2)w_{a}+w_{f}/2-f_{c}\theta_{y}}{\sigma_{I}}\right)$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
-Q\left(\frac{(j-N_{a}/2)w_{a}-w_{f}/2-f_{c}\theta_{x}}{\sigma_{I}}\right)\Bigg{]}.$
(17)

###  III-B Received Signal Model

Due to its relatively low implementation complexity, most current commercial
FSO systems utilize intensity modulation with direct detection (IM/DD) based
on on-off keying (OOK). In practical FSO links, the mean number of absorbed
photons is typically sufficiently large, allowing the distribution of the
number of APD output electrons to be well approximated by a Gaussian
distribution [39]. Accordingly, the photo-current corresponding to the $k$th
bit interval on the $(i,j)$th detector can be expressed as

$r_{ij}[k]=h_{ij}\mu s[k]+n_{ij}[k]$ (18)

where $i$ and $j\in\\{1,...,N_{a}\\}$, $\mu=\frac{eG\eta}{h_{p}\nu}$, $e$
denotes the charge of electron, $G$ is the average APD gain, $\eta$ denotes
the APD quantum efficiency, $h_{p}$ denotes the Planck's constant, and $\nu$
is the optical frequency. Moreover in (18), $s[k]$ and $n_{ij}[k]$,
respectively, denote the transmitted symbol with optical power $P_{t}$, and
the photo-current noise of the $(i,j)$th detector that is an additive Gaussian
noise having zero-mean and variance $\sigma_{ij,k}^{2}$ given by

$\displaystyle\sigma_{ij,k}^{2}=\sigma_{s}^{2}h_{ij}s[k]+\sigma_{0}^{2},$ (19)

where $\sigma_{s}^{2}$ is the variance of the shot noise111This noise is
attributed to the quantum nature of light, where the number of photons emitted
by a coherent optical source in a given time is never constant. Known as
quantum noise or photon noise, this phenomenon arises from the random arrival
rate of photons from the data-carrying optical source, acting as a shot noise
present in all photodetectors [40]. due to transmitted signal, and
$\sigma_{0}^{2}=\sigma_{b}^{2}+\sigma_{th}^{2}$ is the total noise variance
due to the variance of background radiation, $\sigma_{b}^{2}$, and receiver
thermal noise, $\sigma_{th}^{2}$. We have $\sigma_{b}^{2}=2eGF\mu BP_{b}$,
where $F$ denotes the APD excess noise factor, and $B$ is the bandwidth of the
receiver low-pass filter (in Hz). Furthermore,
$\sigma_{th}^{2}=\frac{4K_{B}T_{r}B}{R_{l}}$, where $K_{B}$ is Boltzmann
constant, $T_{r}$ is the receiver equivalent temperature in Kelvin, and
$R_{l}$ is the load resistance. The background power $P_{b}$ is a function of
the photodetector area and can be obtained as [35, 4]

$\displaystyle P_{b}=N_{b}(\lambda)\,B_{o}\,\Omega_{FoV}\,A_{a}$ (20)

where $N_{b}(\lambda)$ is the spectral radiance of the background radiations
at wavelength $\lambda$ (in Watts/${\rm cm}^{2}$-$\micro$m-srad), $B_{o}$ is
the bandwidth of the optical filter at the Rx (in $\micro$m), and $A_{a}=\pi
r_{a}^{2}$ is the lens area (in ${\rm cm}^{2}$).

##  IV Data Detection and Spatial Beam Tracking

We assume that the transmitted optical signal is captured during an
observation window consisting of $L$ bits. Therefore, the received signal
vector $\underline{r}_{ij}=\left[{r_{ij}[1],\ldots,r_{ij}[L]}\right]$ at the
$(i,j)$th APD is related to the $L$ transmitted signal vector
$\underline{s}=\left[{s[1],s[2],\ldots,s[L]}\right]$. We also assume slow
fading channel, i.e., the channel remains constant during the observation
window [4].

In this paper, we address a practical scenario where channel state information
(CSI) and the beam center position (or equally beam AoA) are unknown at the Rx
side. Our objective is to detect OOK signals without the need for pilot
symbols and to determine the position of the incident optical beam on the
array of APDs. To achieve this, we propose an efficient data-aided channel
estimation method without requiring pilot symbols. We then detect the sequence
of OOK signals and perform spatial beam tracking. We evaluate the performance
of our sub-optimal data detection method in terms of bit error rate BER,
considering the ideal scenario where perfect CSI and exact beam position are
available at the Rx as a benchmark. Further details of these methods are
presented in the following sections.

###  IV-A Ideal Data Detection and Spatial Beam Tracking

In this subsection, we first study the ideal data detection under the
assumption of perfect knowledge of the channel coefficient $h_{1}$, and the
beam angular deviations $\theta_{x}$ and $\theta_{y}$, as a benchmark.

####  IV-A1 Ideal Data Detection

For the ideal detection method, we assume that the Rx perfectly knows the
instantaneous channel coefficient $h_{1}$. We also assume that the values of
$\theta_{x}$ and $\theta_{y}$ are available with high accuracy (i.e., no need
to spatial beam tracking). Hence, the position of the center of optical beam
spot is available, which leads to the perfect knowledge of the channel
coefficient $h_{2,ij}$.

Under this ideal scenario, the ideal detection method based on the ML
criterion can be performed symbol-by-symbol as

$\displaystyle\hat{s}[k]=\underset{s[k]\in\\{0,1\\}}{\operatorname{\text{\leavevmode\nobreak\
arg\leavevmode\nobreak\ max\leavevmode\nobreak\ }}}\leavevmode\nobreak\
\prod_{i=1}^{N_{a}}\prod_{j=1}^{N_{a}}p(r_{ij}[k]|h_{ij};s[k])$ (21)

where $p(\cdot)$ denotes the probability of an event. Accordingly, from (18),
(19), and (21), the ideal data detection method can be obtained as

$\displaystyle\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\frac{-2\sigma_{0}^{2}\mu
h_{ij}r_{ij}[k]+\sigma_{0}^{2}\mu^{2}h_{ij}^{2}-\sigma_{s}^{2}h_{ij}r_{ij}[k]^{2}}{\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4}}\leavevmode\nobreak\
\begin{subarray}{c}\hat{s}[k]=0\\\ >\\\ <\\\
\hat{s}[k]=1\end{subarray}\leavevmode\nobreak\
\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right).$
(22)

Moreover, the BER of the ideal detection method is obtained as

$\displaystyle\mathbb{P}^{\textrm{ML}}_{e}=\int_{0}^{\infty}\int_{0}^{\infty}\int_{0}^{\infty}\mathbb{P}^{\textrm{ML}}_{e|h_{1},\theta_{x},\theta_{y}}\times
f_{h_{1}}(h_{1})f_{\theta_{x}}(\theta_{x})f_{\theta_{y}}(\theta_{y})\textrm{d}h_{1}\textrm{d}\theta_{x}\textrm{d}\theta_{y}$
(23)

where

$\displaystyle\mathbb{P}^{\textrm{ML}}_{e|h_{1},\theta_{x},\theta_{y}}=p_{0}\mathbb{P}^{\textrm{ML}}_{e|01}+p_{1}\mathbb{P}^{\textrm{ML}}_{e|10}.$
(24)

In (24), $p_{0}$ and $p_{1}$ denote the probability of transmitting bit “0”
and “1”, respectively, and it is assumed that
$p_{0}=p_{1}=\nicefrac{{1}}{{2}}$. Also, $\mathbb{P}_{e|01}$ and
$\mathbb{P}_{e|10}$ denote the BER conditioned on $s[k]=0$ and $s[k]=1$,
respectively, and they can be obtained as

$\displaystyle\mathbb{P}^{\textrm{ML}}_{e|10}$
$\displaystyle=\textrm{Prob}\Bigg{\\{}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\bigg{[}-2\sigma_{0}^{2}\mu
h_{ij}(\mu h_{ij}+n_{ij}[k])$ (25) $\displaystyle\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
+\sigma_{0}^{2}\mu^{2}h_{ij}^{2}-\sigma_{s}^{2}h_{ij}(\mu^{2}h_{ij}^{2}+2\mu
h_{ij}n_{ij}[k]+n_{ij}^{2}[k])$ $\displaystyle\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
-(\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4})\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)\bigg{]}>0\Bigg{\\}}$
$\displaystyle\simeq\textrm{Prob}\Bigg{\\{}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\bigg{[}n_{ij}[k]+\frac{\mu
h_{ij}}{2}+\frac{(\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4})}{2\mu
h_{ij}(\sigma_{0}^{2}+\sigma_{s}^{2}h_{ij})}$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\times\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)\bigg{]}<0\Bigg{\\}}$

and

$\displaystyle\mathbb{P}^{\textrm{ML}}_{e|01}$
$\displaystyle=\textrm{Prob}\Bigg{\\{}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\bigg{[}\frac{-2\sigma_{0}^{2}\mu
h_{ij}n_{ij}[k]+\sigma_{0}^{2}\mu^{2}h_{ij}^{2}}{\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4}}$
(26) $\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\
-\frac{\sigma_{s}^{2}h_{ij}n_{ij}[k]^{2}}{\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4}}-\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)\bigg{]}<0\Bigg{\\}}$
$\displaystyle\simeq\textrm{Prob}\Bigg{\\{}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\bigg{[}{n_{ij}[k]-\frac{\mu
h_{ij}}{2}}+\frac{\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4}}{2\sigma_{0}^{2}\mu
h_{ij}}$ $\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\
\times\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)\bigg{]}>0\Bigg{\\}}.$

The approximations utilized in (25) and (26) stem from neglecting the impact
of $n^{2}_{ij}[k]$ in comparison to $n_{ij}[k]$. As demonstrated in the
numerical results, this approximation proves to be reasonable, particularly in
high signal-to-noise ratio (SNR) regimes. Furthermore, considering that the
error probability expression in (25) is formulated for the scenario of
transmitting bit “1,” it becomes evident from (19) that the term
$\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}n_{ij}[k]$ constitutes a zero-mean
Gaussian noise with a variance of
$N_{a}^{2}\sigma_{0}^{2}+\sigma_{s}^{2}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}h_{ij}$.
Likewise, the term $\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}n_{ij}[k]$ in (26)
corresponds to Gaussian noise with a mean of zero and a variance of
$N_{a}^{2}\sigma_{0}^{2}$, as it is derived for the transmission of bit “0.”
Therefore, we have

$\displaystyle\mathbb{P}^{\textrm{ML}}_{e|10}\simeq
Q\left(\frac{\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\left[\frac{\mu
h_{ij}}{2}+\frac{(\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4})}{2\mu
h_{ij}(\sigma_{0}^{2}+\sigma_{s}^{2}h_{ij})}\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)\right]}{N_{a}^{2}\sigma_{0}^{2}+\sigma_{s}^{2}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}h_{ij}}\right)$
(27)

and

$\displaystyle\mathbb{P}^{\textrm{ML}}_{e|01}\simeq
Q\left(\frac{\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\left[\frac{\mu
h_{ij}}{2}-\frac{\sigma_{0}^{2}\sigma_{s}^{2}h_{ij}+\sigma_{0}^{4}}{2\sigma_{0}^{2}\mu
h_{ij}}\ln\left(\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)\right]}{N_{a}^{2}\sigma_{0}^{2}}\right).$
(28)

####  IV-A2 Ideal Spatial Tracking

In this scenario, we assume the case that the Rx only knows the instantaneous
channel coefficients $h_{1}$ and the vector of transmitted signal
$\underline{s}$ and needs to estimate position of the beam centre on the
detector array. To do so, one can obtain the number of bits “1” in the
observation window of length $L$ as $m=\sum\limits_{k=1}^{{L}}{s[k]}$. Under
this ideal scenario, the AoA of the received beam can be obtained based on the
ML criterion as

$\displaystyle\hat{\theta}_{x},\hat{\theta}_{y}=\underset{\theta_{x},\theta_{y}}{\operatorname{\text{\leavevmode\nobreak\
arg\leavevmode\nobreak\ max\leavevmode\nobreak\ }}}\leavevmode\nobreak\
\prod_{i=1}^{N_{a}}\prod_{j=1}^{N_{a}}p(r^{\prime\prime}_{ij}|h_{ij},m;\theta_{x},\theta_{y})$
(29)

where

$\displaystyle r^{\prime\prime}_{ij}$
$\displaystyle=\frac{1}{m}\sum_{k=1}^{L}r_{ij}[k]s[k]=\mu
h_{ij}+n_{ij}^{\prime\prime}$ (30)

and $n^{\prime\prime}$ is a zero-mean Gaussian noise with variance

$\displaystyle\sigma_{n_{ij}^{\prime\prime}}^{2}=\frac{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}{m}.$
(31)

Accordingly, from (III-A2), (30), and (31), the ideal spatial tracking method
of (29) can be simplified as

$\displaystyle\hat{\theta}_{x},\hat{\theta}_{y}$
$\displaystyle=\underset{\theta_{x},\theta_{y}}{\operatorname{\text{arg\leavevmode\nobreak\
min\leavevmode\nobreak\
}}}\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\Bigg{[}\frac{1}{m}\ln\left(\frac{\sigma_{s}^{2}h_{1}h_{2,ij}(\theta_{x},\theta_{y})+\sigma_{0}^{2}}{m}\right)$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
+\frac{|r_{ij}^{\prime\prime}-\mu
h_{1}h_{2,ij}(\theta_{x},\theta_{y})|^{2}}{\sigma_{s}^{2}h_{1}h_{2,ij}(\theta_{x},\theta_{y})+\sigma_{0}^{2}}\Bigg{]}.$
(32)

###  IV-B Optimal Data Detection Method Without Knowing CSI

Here, we consider a general condition where the Rx lacks any information
regarding the instantaneous channel coefficients. Consequently, the ML
detector utilizes the outputs of all APDs in the array to estimate the
transmitted signal vector as

$\displaystyle\hat{\underline{s}}=\underset{\underline{s}}{\operatorname{\text{\leavevmode\nobreak\
arg\leavevmode\nobreak\ max\leavevmode\nobreak\ }}}\leavevmode\nobreak\
\prod_{i=1}^{N_{a}}\prod_{j=1}^{N_{a}}p(\underline{r}_{ij}|\underline{s})$
(33)

where

$\displaystyle
p(\underline{r}_{ij}|\underline{s})=\prod_{k=1}^{L}p(r_{ij}[k]|s[k])$ (34)

and

$\displaystyle p(r_{ij}[k]|s[k])$
$\displaystyle=\int_{0}^{\infty}\int_{0}^{\infty}\int_{0}^{\infty}p(r_{ij}[k]|s[k],h_{1},h_{2,ij}(\theta_{x},\theta_{y}))$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\times
f_{h_{1}}(h_{1})f_{\theta_{x}}(\theta_{x})f_{\theta_{y}}(\theta_{y})\textrm{d}h_{1}\textrm{d}\theta_{x}\textrm{d}\theta_{y}.$
(35)

As we realize from (33), across an observation window of length $L$, the ML-
based method can detect the transmitted signal vector $\underline{s}$.
Specifically, data detection entails searching over $2^{L}$ potential states
of the metric $\eqref{f9}$, necessitating the computation of a three-level
integral. Evidently,this approach incurs relatively high computational
complexity, which contradicts the limited energy consumption constraints of
CubeSats. In the sequel, we seek to find sub-optimal detection methods having
lower computational complexity.

###  IV-C Equal Gain Combining Method

To reduce the computational complexity of the ML-based detection method, one
can employ the EGC technique at the Rx. More precisely, this approach obviates
the necessity to decide on a particular APD by combining the output of
$N_{a}\times N_{a}$ APDs with the same factor. In theory, the output photo-
current of the EGC method for an $N_{a}\times N_{a}$ APD array is similar to
that of a single APD with a large active area $A_{s}=N_{a}\times
N_{a}\times{w_{a}}^{2}$. However, unlike their terrestrial counterparts that
utilize very small APDs, for the ground-to-CubeSat FSO systems, the area of
the single APD must be large enough to overcome the effect of AoA
fluctuations. Nevertheless, as the size of an APD increases, its electrical
bandwidth and its quantum efficiency decrease. Given the importance of
maximizing detector sensitivity to minimize the required laser power for a
specified link margin, rather than employing a single large APD, we opt for an
array of APDs and implement the EGC technique to aggregate the total photo-
current generated by the array.

Accordingly, we denote by $r[k]$ the total photo-current generated by all APDs
in the detector array, which can be expressed as

$\displaystyle r[k]=\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}{r_{ij}[k]}=\mu
s[k]h_{t}+n[k]$ (36)

where

$\displaystyle h_{t}=\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}h_{ij}$ (37)

and $n[k]=\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}{n_{ij}[k]}$ is an additive
Gaussian noise with zero-mean and variance
$\sigma_{k}^{2}=\sigma_{s}^{2}h_{t}s[k]+N_{a}^{2}\sigma_{0}^{2}$. When $h_{t}$
is known at the Rx side, transmitted data can be detected symbol-by-symbol as

$\displaystyle r[k]\begin{subarray}{c}\hat{s}[k]=1\\\ >\\\ <\\\
\hat{s}[k]=0\end{subarray}h_{\textrm{th,EGC}}$ (38)

where $h_{\textrm{th,EGC}}$ is the detection threshold for EGC method. From
(36), the BER of the EGC method conditioned on $h_{t}$ can be obtained as

$\displaystyle\mathbb{P}^{\textrm{EGC}}_{e|h_{t}}=\frac{1}{2}Q\left(\frac{\mu
h_{t}-h_{\textrm{th,EGC}}}{\sqrt{N_{a}^{2}\sigma_{0}^{2}+\sigma_{s}^{2}h_{t}}}\right)+\frac{1}{2}Q\left(\frac{h_{\textrm{th,EGC}}}{N_{a}\sigma_{0}}\right).$
(39)

The optimal value for $h_{\textrm{th,EGC}}$ is obtained by differentiating
(39) with respect to $h_{\textrm{th,EGC}}$ and setting the results equal to
zero. By doing so, the optimal value for $h_{\textrm{th,EGC}}$ is obtained as

$\displaystyle h_{\textrm{th,EGC}}=$ $\displaystyle\frac{-\mu
h_{t}N_{a}^{2}\sigma_{0}^{2}}{\sigma_{s}^{2}h_{t}}+\frac{\sqrt{(N_{a}^{2}\sigma_{0}^{2}+\sigma_{s}^{2}h_{t})N_{a}^{2}\sigma_{0}^{2}}}{\sigma_{s}^{2}h_{t}}$
(40)
$\displaystyle\times\sqrt{\mu^{2}h_{t}^{2}+2\sigma_{s}^{2}h_{t}\ln{\displaystyle{\left(\frac{N_{a}^{2}\sigma_{0}^{2}+\sigma_{s}^{2}h_{t}}{N_{a}^{2}\sigma_{0}^{2}}\right)}}}.$

However, as we will observe in the numerical results section, despite its
simplicity, the EGC method fails to achieve ideal performance, even with
perfect channel estimation. Specifically, since the receiver in this method
incorporates all captured background noise with equal gain in the detection
process, it becomes more susceptible to variations in background noise
levels222In this setup, each detector in the array adds background noise to
the total noise power. Additionally, the combining algorithm does not provide
an optimal detection threshold to effectively mitigate the noise introduced by
each detector.. Moreover, by summing the output of the APDs in the array,
performing the fine beam tracking method is not directly possible. This
process requires additional hardware, which may not always be available due to
SWaP limitations.

In the sequel, we propose a sub-optimal method capable of data detection and
beam tracking, achieving performance close to the ideal receiver without
requiring knowledge of the channel. Simultaneously, this method aims to reduce
the computational load associated with the optimal ML detection method.

###  IV-D Sub-Optimal Detection Method and Spatial Beam Tracking

In this subsection, we propose a sub-optimal detection method to reduce the
computational complexity of the optimal ML detection method proposed in (33).
We also aim to achieve performance close to the ideal case. To this end, we
collect the received vector $r_{ij}$ of all APDs over an observation window of
length $L$. We then use the GLRT principle to estimate jointly the
instantaneous channel coefficient of each APD $h_{ij}$, and detect the
transmit data sequence. The proposed GLRT-based detection method uses the past
detected data to improve the accuracy of the channel estimation. For the
considered system model, the transmitted vector $\underline{s}$ and the
instantaneous channel coefficient $h_{ij}$ can be jointly estimated based on
the GLRT criterion as [41]

$\displaystyle\hat{\underline{s}},\hat{h}_{ij}$
$\displaystyle=\underset{\underline{s},h_{ij}\textrm{s}}{\operatorname{\text{arg\leavevmode\nobreak\
max\leavevmode\nobreak\
}}}\prod_{i=1}^{N_{a}}\prod_{j=1}^{N_{a}}p(\underline{r}_{ij}|\underline{s},h_{ij})$
(41)
$\displaystyle=\underset{\underline{s},h_{ij}}{\operatorname{\text{arg\leavevmode\nobreak\
max\leavevmode\nobreak\
}}}\prod_{i=1}^{N_{a}}\prod_{j=1}^{N_{a}}\prod_{k=1}^{L}\frac{1}{\sqrt{2\pi(\sigma_{s}^{2}h_{ij}s[k]+\sigma_{0}^{2})}}$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\times\exp\left(-\frac{|r_{ij}[k]-\mu
h_{ij}s[k]|^{2}}{2(\sigma_{s}^{2}h_{ij}s[k]+\sigma_{0}^{2})}\right).$

From (41), when $\underline{s}$ is known at the Rx side, the RV $h_{ij}$ for
the $(i,j)$th APD can be estimated as

$\displaystyle\hat{h}_{ij}$
$\displaystyle=\underset{h_{ij}}{\operatorname{\text{arg\leavevmode\nobreak\
min\leavevmode\nobreak\
}}}\underbrace{\sum_{k^{\prime\prime}\in\mathcal{K}^{\prime\prime}}\left(\ln(\sigma_{0}^{2})+\frac{|r_{ij}[k^{\prime\prime}]|^{2}}{\sigma_{0}^{2}}\right)}_{\text{A1}}$
(42)
$\displaystyle+\underbrace{\sum_{k^{\prime}\in\mathcal{K}^{\prime}}\left(\ln(\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2})+\frac{|r_{ij}[k^{\prime}]-\mu
h_{ij}|^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right)}_{\textrm{A2}}$

where $\mathcal{K}^{\prime}$ and $\mathcal{K}^{\prime\prime}$ are the subsets
of the set $\mathcal{K}=\\{k\\}$ where $k\in\\{1,\ldots,L\\}$ for which
$s[k]=1$, and $s[k]=0$, respectively. As evident from (42), the term $A1$ is
independent of the RV $h_{ij}$. Hence, (42) can be simplified as

$\displaystyle\hat{h}_{ij}$
$\displaystyle=\underset{h_{ij}}{\operatorname{\text{arg\leavevmode\nobreak\
min\leavevmode\nobreak\
}}}\sum_{k^{\prime}\in\mathcal{K}^{\prime}}\left(\ln(\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2})+\frac{|r_{ij}[k^{\prime}]-\mu
h_{ij}|^{2}}{\sigma_{s}^{2}h_{ij}+\sigma_{0}^{2}}\right).$ (43)

For any given $\underline{s}$, the optimal value of $h_{ij}$ is obtained by
differentiating (43) with respect to $h_{ij}$ and setting the resulting
expression to zero. By doing so, we have

$\displaystyle
h_{ij}^{\textrm{opt}}=\frac{-\sigma_{s}^{2}}{2\mu^{2}}-\frac{\sigma_{0}^{2}}{\sigma_{s}^{2}}+\sqrt{\frac{\sigma_{s}^{4}}{4\mu^{4}}+\frac{\sigma_{0}^{4}}{\sigma_{s}^{4}}+\mathbb{R}_{ij}}.$
(44)

where
$\mathbb{R}_{ij}=\frac{\sum_{k=1}^{L}(2\mu\sigma_{0}^{2}r_{ij}[k]s[k]+\sigma_{s}^{2}r_{ij}[k]^{2}s[k])}{m\mu^{2}\sigma_{s}^{2}}$.
Now, by substituting (44) in (41) and after some manipulations, the GLRT-based
method for data detection (41) is simplified as (45).

$\displaystyle\hat{\underline{s}}$
$\displaystyle=\underset{\underline{s}}{\operatorname{\text{arg\leavevmode\nobreak\
min\leavevmode\nobreak\
}}}N_{a}^{2}(L-m)\ln(\sigma_{0}^{2})+N_{a}^{2}\sum_{k^{\prime\prime}\in\mathcal{K}^{\prime\prime}}\frac{|r_{ij}[k^{\prime\prime}]|^{2}}{\sigma_{0}^{2}}+m\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\ln\left(\sqrt{\frac{\sigma_{s}^{8}}{4\mu^{4}}+\sigma_{0}^{4}+\sigma_{s}^{4}\mathbb{R}_{ij}}-\frac{\sigma_{s}^{4}}{2\mu^{2}}\right)$
(45) $\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
+\sum_{i=1}^{N_{a}}\sum_{j=1}^{N_{a}}\frac{1}{\sqrt{\frac{\sigma_{s}^{8}}{4\mu^{4}}+\sigma_{0}^{4}+\sigma_{s}^{4}\mathbb{R}_{ij}}-\frac{\sigma_{s}^{4}}{2\mu^{2}}}\sum_{k^{\prime}\in\mathcal{K}^{\prime}}\left(r_{ij}[k^{\prime}]+\frac{\sigma_{s}^{2}}{2\mu}+\frac{\mu\sigma_{0}^{2}}{\sigma_{s}^{2}}-\sqrt{\frac{\sigma_{s}^{4}}{4\mu^{2}}+\frac{\mu\sigma_{0}^{4}}{\sigma_{s}^{4}}+\mu\mathbb{R}_{ij}}\right)^{2}.$

It is noteworthy that unlike the metric in (33), which necessitates the
calculation of a three-level integral, the detection metric in (45) involves
straightforward additions and multiplications, and does not rely on channel
knowledge. The computational complexity of the metric in (45) is
$\mathcal{O}(LN_{a}^{2})$. Moreover, to detect the transmitted vector
$\underline{s}$ using the proposed detection method in (45), one needs to
search among $2^{L}$ possible received sequences and select the sequence that
minimizes the metric. Hence, the total computational complexity of the
proposed detection method is $\mathcal{O}(LN_{a}^{2}\times 2^{L})$.

We note that for an all-zero transmitted sequence, i.e., when $m=0$ &
$\mathbb{R}_{ij}=\frac{0}{0}$, our proposed method is unable to estimate
properly the channel, and subsequently, perform data detection. For an
observation window of length $L$, the occurrence probability of an all-zero
sequence is equal to $\frac{1}{2^{L}}$. Hence, the value of $L$ must be large
enough to ensure that the occurrence probability of the all-zero sequence is
lower than the desired BER. As we will show in the numerical result section,
the proposed GLRT-based method achieves an acceptable performance compared
with the ideal receiver when $L>20$. For such values of $L$, detection is done
by searching over $2^{20}=1048576$ possible states of metric (45). To avoid
this exhaustive search and to reduce further the Rx complexity, we propose a
sub-optimal method for implementing $\eqref{d5}$, which reduces the search
space to $L$. Particularly, by considering the nature of the OOK modulation,
one can reasonably assume that the maximum received signal levels are due to
the transmitted bit “1”, and the minimum received signal levels are due to the
transmitted bits “0”. Hence, the $m$ maximum received signal levels during the
observation interval $L$ can be selected as bit “1”, and the remainder can be
selected as bit “0”. Since the value of RV $m$ varies from $1$ to $L$, the Rx
searches among $L$ possible received sequences (i.e., sequences with which the
value of $m$ is equal to $1,2,\ldots,L$) to find the sequence that minimizes
$\eqref{d5}$. This way, the search space for finding $\hat{\underline{s}}$ is
reduced from $2^{L}$ to $L$, and thus, the computational complexity is
obtained as $\mathcal{O}(L^{2}N^{2}_{a})$. In the numerical result section, we
show that searching over this sub-space achieves a performance close to the
case that the detection is done over all possible states.

Finally, to perform fine beam tracking by determining the position of the
incident optical baeam on the APD array, we resort to the following
proposition.

Proposition: Suppose the center of received optical beam is located in the
($i^{\prime},j^{\prime}$)th APD in the array, i.e., we have
$\Big{\\{}(i^{\prime}-1-\frac{N_{a}}{2})w_{a}+\frac{w_{f}}{2}<\theta_{x}f_{c}<(i^{\prime}-\frac{N_{a}}{2})w_{a}-\frac{w_{f}}{2}\Big{\\}}$,
and
$\Big{\\{}(j^{\prime}-1-\frac{N_{a}}{2})w_{a}+\frac{w_{f}}{2}<\theta_{y}f_{c}<(j^{\prime}-\frac{N_{a}}{2})w_{a}-\frac{w_{f}}{2}\Big{\\}}$.
In this case, the channel coefficient corresponding to that APD,
$h_{2,i^{\prime}j^{\prime}}$, is greater than the other channel coefficients
corresponding to the other APDs in the array, i.e.,
$h_{2,i^{\prime}j^{\prime}}\geq h_{2,ij}$, where $(i,j)\in\\{1,...,N_{a}\\}$.

###### Proof:

Let rewrite (III-A2) as $h_{2,ij}=h_{2,ij}^{1}h_{2,ij}^{2}$, where

$\displaystyle
h_{2,ij}^{1}=\Bigg{[}Q\left(\frac{(i-1-N_{a}/2)w_{a}+w_{f}/2-f_{c}\theta_{x}}{\sigma_{I}}\right)$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
-Q\left(\frac{(i-N_{a}/2)w_{a}-w_{f}/2-f_{c}\theta_{x}}{\sigma_{I}}\right)\Bigg{]}$
(46) $\displaystyle
h_{2,ij}^{2}=\Bigg{[}Q\left(\frac{(j-1-N_{a}/2)w_{a}+w_{f}/2-f_{c}\theta_{y}}{\sigma_{I}}\right)$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\
-Q\left(\frac{(j-N_{a}/2)w_{a}-w_{f}/2-f_{c}\theta_{x}}{\sigma_{I}}\right)\Bigg{]}.$
(47)

From (IV-D) and (IV-D), it can be easily realized that
$\Big{\\{}(i^{\prime}-1-\frac{N_{a}}{2})w_{a}+\frac{w_{f}}{2}<\theta_{x}f_{c}<(i^{\prime}-\frac{N_{a}}{2})w_{a}-\frac{w_{f}}{2}\Big{\\}}$,
implies that $h^{1}_{2,i^{\prime}j^{\prime}}\geq h^{1}_{2,ij}$, and,
similarly,
$\Big{\\{}(j^{\prime}-1-\frac{N_{a}}{2})w_{a}+\frac{w_{f}}{2}<\theta_{y}f_{c}<(j^{\prime}-\frac{N_{a}}{2})w_{a}-\frac{w_{f}}{2}\Big{\\}}$,
implies that $h^{2}_{2,i^{\prime}j^{\prime}}\geq h^{2}_{2,ij}$. Thereby,
$h^{1}_{2,i^{\prime}j^{\prime}}\geq h^{1}_{2,ij}$ and
$h^{2}_{2,i^{\prime}j^{\prime}}\geq h^{2}_{2,ij}$, together imply that
$h_{2,i^{\prime}j^{\prime}}\geq h_{2,ij}$. ∎

Now, let define $r^{s}_{ij}=\frac{1}{L}\sum_{k=1}^{L}r_{ij}[k]$. Using (19),
we can rewrite $r^{s}_{ij}$ as

$\displaystyle r^{s}_{ij}=\frac{m\mu h_{1}h_{2,ij}}{L}+n^{s}_{ij}$ (48)

where $n^{s}_{ij}$ is a zero-mean Gaussian noise with variance

$\displaystyle\sigma^{2}_{s,ij}=\frac{mh_{1}h_{2,ij}\sigma_{s}^{2}}{L^{2}}+\frac{\sigma^{2}_{0}}{L}.$
(49)

As we can observe from (49), by increasing $L$, the variance of noise
$n^{s}_{ij}$ approaches zero, and thus, for large values of $L$, one can
reasonably neglect the effect of $n^{s}_{ij}$ and approximate $r^{s}_{ij}$ as
$r^{s}_{ij}\simeq\frac{m\mu h_{1}h_{2,ij}}{L}$. Given these assumptions and
also according to the proposition, the vertical and horizontal intervals on
the APD array at which the center of the received optical beam is placed, can
be estimated as

$\displaystyle(\hat{i}^{\prime},\hat{j}^{\prime})=\underset{i,j\in\\{1,...,N_{a}\\}}{\operatorname{\text{arg\leavevmode\nobreak\
max\leavevmode\nobreak\ }}}r^{s}_{ij}.$ (50)

By feeding this information back to the CubeSat’s mechanical subsystem, such
as a fast steering mirror, through a control message, the receiver can correct
any orientation errors and realign itself with the direction of the incoming
beam.

We note that (50) gives an estimate of the interval that the RVs $\theta_{x}$
and $\theta_{y}$ are located in, i.e.,
$\Big{\\{}\frac{(i^{\prime}-1-\frac{N_{a}}{2})w_{a}}{f_{c}}+\frac{w_{f}}{2f_{c}}<\theta_{x}<\frac{(i^{\prime}-\frac{N_{a}}{2})w_{a}}{f_{c}}-\frac{w_{f}}{2f_{c}}\Big{\\}}$,
and
$\Big{\\{}\frac{(j^{\prime}-1-\frac{N_{a}}{2})w_{a}}{f_{c}}+\frac{w_{f}}{2f_{c}}<\theta_{y}<\frac{(j^{\prime}-\frac{N_{a}}{2})w_{a}}{f_{c}}-\frac{w_{f}}{2f_{c}}\Big{\\}}$.
Nevertheless, the accurate estimation of $\theta_{x}$ and $\theta_{y}$
requires the perfect knowledge of the channel coefficients $h_{1}$ and
$h_{2,ij}$, which is beyond the scope of this work and can be an interesting
future direction.

##  V  Numerical Results And Discussion

In this section, numerical results are presented in terms of BER to assess the
performance of the proposed methods for both data detection and fine beam
tracking. Additionally, Monte-Carlo simulations are conducted to validate the
accuracy of the derived analytical expressions across various parameter values
related to CubeSat’s angular instabilities $\sigma_{x}^{2}$ and
$\sigma_{y}^{2}$, transmit power $P_{t}$, the number of APDs in the array
$N=N_{a}\times N_{a}$, and the observation window length $L$. Simulations are
performed based on the practical values of the parameters outlined in Table I
[40].

Table I: System Parameters Used Throughout Simulations Name | Parameter | Value
---|---|---
Altitude difference | $H_{s}-H_{g}$ |  $400$ km
Zenith Angle | $\xi$ | $30^{\circ}$
APD Gain | $G$ | $100$
Quantum Efficiency | $\eta$ | $0.7$
Wavelength | $\lambda$ |  $1550$ nm
Receiver Load | $R_{l}$ | $1\leavevmode\nobreak\ k\Omega$
Receiver Temperature | $T_{r}$ | $300\degree\leavevmode\nobreak\ K$
Optical Filter Bandwidth | $B_{o}$ |  $1$ nm
Spectral Radiance | $N_{b}(\lambda)$ |  $10^{-4}$ Watts/${\rm cm}^{2}$-$\micro$m-srad
Ground Turbulence Level | $C_{n}^{2}(0)$ |  $10^{-13}$ $\textrm{m}^{-2/3}$
Bit Time | $T_{b}$ | $10^{-9}$
Focal length | $f_{c}$ |  $3$ cm
Aperture Radius | $r$ |  $3$ cm
Beam Width | $w_{z}$ |  $30$ m
Background Power | $P_{b}$ |  $1$ nW
Detector width | $w_{a}$ | 250 $\micro$m
Dead-space width | $w_{f}$ | 5 $\micro$m

(a)

(b)

Figure 4: BER of the proposed GLRT-based detection method versus $P_{t}$ for
two different angular instabilities. The case of employing the EGC method and
the ideal Rx are considered as benchmarks.

We first investigate the performance of the proposed detection method in terms
of BER. For evaluation, we also consider the BER results when the Rx perfectly
knows the channel information (i.e., the ideal Rx), and when it employs the
EGC and MRC techniques as benchmarks. Fig. 4 demonstrates BER versus $P_{t}$
for two different angular instabilities, i.e., when
$\sigma_{x}^{2}=\sigma_{y}^{2}=2$ mrad (Fig. 4a), and when
$\sigma_{x}^{2}=\sigma_{y}^{2}=10$ mrad (Fig. 4b). The results of Fig. 4 are
obtained for $L=16$ and $N_{a}=6$. By comparing the results of Figs. 4a and
4b, one can readily observe that by increasing the CubeSat’s angular
instabilities, the link performance degrades significantly. Particularly, an
error floor can be noticed in case of large angular instabilities and small Rx
FoV due to the dominant effect of AoA fluctuations on the system performance.
Moreover, the results clearly prove the superiority of the proposed GLRT-based
detection method compared with the EGC and MRC methods. Indeed, our proposed
GLRT-based method outperforms the aforementioned combination techniques since
it can exploit the continuity of the channel, and thus efficiently and
accurately estimate the instantaneous channel coefficient to perform data
detection. In addition, it is obvious that the MRC technique outperforms the
EGC method at the expense of high computational complexity. Following the
higher computational complexity and relatively poor performance, it can be
concluded that the MRC technique is not a suitable method for data detection
of the considered scenario. Meanwhile, we notice a perfect match between the
analytical and simulation-based results which validates the accuracy of our
derived analytical expressions.

Figure 5: BER of the proposed GLRT-based detection method versus $L$ for two
different values of $P_{t}$. The case of employing the EGC method and the
ideal Rx are considered as benchmarks. Figure 6: BER of the proposed GLRT-
based detection method versus $L$ for two different values of $N_{a}$. The
case of employing the EGC method and the ideal Rx are considered as
benchmarks.

To have a deeper understanding about the effect of the length of observation
window, $L$, on the link performance, we have shown the BER curves of the
considered detection methods for $P_{t}$ of 5 dB and 10 dB versus $L$ in Fig.
5. The results of this figure are obtained for $N_{a}=10$ and
$\sigma_{x}^{2}=\sigma_{y}^{2}=5$ mrad. It is seen from Fig. 5 that the
performance of the both proposed GLRT-based and EGC-based detection methods
improve by increasing $L$. We note that this performance improvement comes at
the cost of increasing computational load as well as the detection delay.
Hence, choosing an optimal value for $L$ involves balancing a tradeoff between
tolerable complexity/delay and desirable BER. For instance, in the considered
setup in Fig. 5 and when transmit power is equal to 5 dB, the lowest value of
$L$ with which the system can achieve performance close to that of the ideal
receiver is equal to $18$. Meanwhile, from the results of Fig. 5 we can
observe that the optimal value for $L$ changes by varying $P_{t}$. For
instance, by increasing $P_{t}$ from 5 to 10 dB, the optimal value for $L$
changes from $18$ to $20$ to attain a lower BER. This can be justified by the
fact that increasing $P_{t}$ leads to higher receiver SNR, subsequently
reducing the BER of an ideal receiver. Therefore, for the proposed detection
method, a larger value for $L$ is required to attain lower values of BER by
mitigating detection and channel estimation errors. Besides $P_{t}$, another
important parameter that can affect the optimal value of $L$ is the number of
APDs in the detector array. To obtain more insight on this effect, we have
presented the BER plots versus $L$ for two different values of $N_{a}$ in Fig.
6. It can be seen form Fig. 6 that by increasing $N_{a}$ from $4$ to $8$, the
optimal value of $L$ for the GLRT-based method increases from $13$ to $20$,
and consequently, the BER decreases.

Figure 7: BER of the proposed GLRT-based detection method versus $N_{a}$ for
two different values of angular instability. The case of employing the EGC
method and the ideal Rx are considered as benchmarks.

The performance degradation due to the AoA fluctuations can be improved by
increasing the Rx FoV via increasing the number of APDs (or equally $N_{a}$)
in the detector array. Indeed, increasing $N_{a}$ makes a compromise between
increased the Rx background noise level on one hand, and the decrease in the
link interruption probability due to AoA fluctuations on the other hand.
Moreover, since the computational complexity of the proposed detection method
increases exponentially by increasing $N_{a}$, finding the optimal values for
$N_{a}$ is essential to design such ground-to-CubeSat FSO links. Similar to
$L$, the optimal value for $N_{a}$ is the minimum value that achieves
performance close to the ideal receiver. To this end, we have plotted the BER
curves of the different detection methods versus $N_{a}$ for two different
values of angular instabilities in Fig. 7. First, it can be seen from the
figure that the performance of the system degrades when the degree of
instabilities is varied form 5 to 12  mrad, since, as expected, higher degree
of instabilities results in the larger value of AoA fluctuations (or equally,
the incident optical beam is more likely to lie outside the Rx FoV). Second,
by increasing the degree of instabilities form 5 to 12  mrad, the optimal
value for $N_{a}$ increases from $7$ to $18$ to compensate the effect of AoA
fluctuations by employing a wider Rx FoV. Moreover, beyond an optimal point,
further enlarging the detector size does not necessarily improve link
reliability. Increasing the receiver FoV by enlarging the detector size
entails capturing more desired transmit power along with undesired background
noise. Consequently, beyond the optimal detector size (i.e., an optimal Rx
field-of-view), background noise becomes dominant over the signal level,
resulting in an increased BER.

Finally, as reducing computational complexity is often crucial to meet the
SWaP requirements of CubeSats, we assess the significance of the proposed
algorithm in terms of its processing time and computational load. It is
demonstrated that the optimal value for $L$ is approximately 20 to achieve
performance close to that of the ideal receiver at the target BER of
$10^{-4}$. Notably, achieving a lower BER (i.e., lower than $10^{-4}$) will
inevitably necessitate an increase in the optimal value for $L$. For the
target BER of $10^{-4}$, the required time for searching through the entire
search space in optimal detection is $2^{20}\times t_{r}=73400320t_{r}$, where
$t_{r}$ represents the processing time unit for each realization and is
directly related to the processing power of the processor unit in the
receiver. However, our proposed method, which can achieve performance close to
that of the ideal receiver, has a search complexity on the order of $L$,
requiring only $20t_{r}$ of processing time. This is significantly lower than
the processing time of the exhaustive search method, resulting in reduced
energy consumption and the ability to utilize simpler processing units in the
payload of the CubeSat.

##  VI Summary and Conclusion

In our study, we tackled the problem of data detection and fine beam tracking
for ground-to-CubeSat FSO links using an array of APDs at the receiver. We
developed a channel model tailored to this link and explored practical
scenarios where both channel information and instantaneous beam position are
unknown at the receiver. We proposed an efficient and practical data-aided
channel estimation method based on the GLRT criterion and evaluated its
performance under various conditions. Furthermore, we determined the center
position of the beam spot on the APD array by comparing the output signals
from the APDs. Our simulation results demonstrated that the proposed GLRT-
based method achieves performance close to the ideal receiver while
maintaining significantly lower complexity. This makes our method particularly
suitable for CubeSats, which operate under stringent SWaP constraints. The
reduced computational load and energy consumption allow for the use of simpler
processing units, enhancing the feasibility of implementing advanced optical
communication systems on CubeSats.

## References

  * [1] ITU-R, “M.2160 : Framework and overall objectives of the future development of IMT for 2030 and beyond,” Nov. 2023, accessed: May 1, 2024. [Online]. Available: https://www.itu.int/rec/R-REC-M.2160/en
  * [2] H.-B. Jeon _et al._ , “Free-Space Optical Communications for 6G Wireless Networks: Challenges, Opportunities, and Prototype Validation,” _IEEE Communications Magazine_ , vol. 61, no. 4, pp. 116–121, Apr. 2023.
  * [3] O. Kodheli _et al._ , “Satellite communications in the new space era: A survey and future challenges,” _IEEE Communications Surveys & Tutorials_, vol. 23, no. 1, pp. 70–109, Oct. 2021.
  * [4] H. Kaushal and G. Kaddoum, “Optical communication in space: Challenges and mitigation techniques,” _IEEE Communications Surveys & Tutorials_, vol. 19, no. 1, pp. 57–96, Aug. 2017.
  * [5] A. U. Chaudhry and H. Yanikomeroglu, “Free space optics for next-generation satellite networks,” _IEEE Consumer Electronics Magazine_ , vol. 10, no. 6, pp. 21–31, Oct. 2021.
  * [6] N. U. Hassan, C. Huang, C. Yuen, A. Ahmad, and Y. Zhang, “Dense small satellite networks for modern terrestrial communication systems: Benefits, infrastructure, and technologies,” _IEEE Wireless Communications_ , vol. 27, no. 5, pp. 96–103, Sep. 2020.
  * [7] N. Saeed, A. Elzanaty, H. Almorad, H. Dahrouj, T. Y. Al-Naffouri, and M.-S. Alouini, “Cubesat communications: Recent advances and future challenges,” _IEEE Communications Surveys & Tutorials_, vol. 22, no. 3, pp. 1839–1862, Apr. 2020.
  * [8] Skyquest, “Global CubeSat Market Insights,” Mar. 2024, accessed: May 1, 2024\. [Online]. Available: https://www.skyquestt.com/report/cubesat-market#:~:text=Global%20CubeSat%20Market%20Insights,period%20(2024%2D2031).
  * [9] J. Vanreusel, “Launching a cubesat: Rules, laws, and best practice,” in _Cubesat Handbook_. Elsevier, pp. 391–414, Jan. 2021.
  * [10] L. Li, Z. Xuejiao, Z. Jianhua, X. Changzhi, and J. Yi, “Advanced space laser communication technology on cubesats,” _ZTE communications_ , vol. 18, no. 4, pp. 45–54, Jan. 2021.
  * [11] R. Li, B. Lin, Y. Liu, M. Dong, and S. Zhao, “A survey on laser space network: terminals, links, and architectures,” _IEEE Access_ , vol. 10, pp. 34 815–34 834, Mar. 2022.
  * [12] B. Mathason _et al._ , “Cubesat lasercom optical terminals for near-earth to deep space communications,” in _Free-Space Laser Communications XXXI_ , vol. 10910. SPIE, pp. 24–29, Mar. 2019.
  * [13] D. R. Kolev _et al._ , “Latest developments in the field of optical communications for small satellites and beyond,” _Journal of Lightwave Technology_ , vol. 41, no. 12, pp. 3750–3757, June 2023.
  * [14] M. Dresscher _et al._ , “Key challenges and results in the design of cubesat laser terminals, optical heads and coarse pointing assemblies,” in _2019 IEEE International Conference on Space Optical Systems and Applications (ICSOS)_ , pp. 1–6, Oct. 2019.
  * [15] L. C. Andrews and R. L. Phillips, _Laser Beam Propagation Through Random Media_. SPIE press Bellingham, WA, 2005\.
  * [16] S. Janson _et al._ , “The NASA optical communications and sensor demonstration program: initial flight results.” in  Proc. Annu. AIAA/USU Conf. Small Satell., 2016.
  * [17] K. Cahoy _et al._ , “The CubeSat laser infrared crosslinK mission (CLICK),” in _International Conference on Space Optics—ICSO_ , vol. 11180\. SPIE, pp. 358–369, July 2019.
  * [18] P. Serra _et al._ , “Optical communications crosslink payload prototype development for the cubesat laser infrared crosslink (click) mission,” 2019.
  * [19] P. Grenfell, A. Aguilar, K. Cahoy, and M. Long, “Pointing, acquisition, and tracking for small satellite laser communications.” in  Proc. Annu. AIAA/USU Conf. Small Satell., 2018\.
  * [20] K. M. Riesing _et al._ , “Pointing, acquisition, and tracking for the TBIRD CubeSat mission: system design and pre-flight results,” in _Free-Space Laser Communications XXXIV_ , vol. 11993. SPIE, Mar. 2022, pp. 207–216.
  * [21] R. Samy, H.-C. Yang, T. Rakia, and M.-S. Alouini, “Ergodic capacity analysis of satellite communication systems with SAG-FSO/SH-FSO/RF transmission,” _IEEE Photonics Journal_ , vol. 14, no. 5, pp. 1–9, Oct. 2022.
  * [22] ——, “Hybrid SAG-FSO/SH-FSO/RF transmission for next-generation satellite communication systems,” _IEEE Transactions on Vehicular Technology_ , vol. 72, no. 11, pp. 14 255–14 267, Nov. 2023.
  * [23] J. Liang _et al._ , “Latency Versus Transmission Power Trade-Off in Free-Space Optical (FSO) Satellite Networks With Multiple Inter-Continental Connections,” _IEEE Open Journal of the Communications Society_ , vol. 4, pp. 3014–3029, Oct. 2023.
  * [24] T. V. Nguyen, T. V. Phamz, A. T. Phamx, and D. T. Ngoc, “Secrecy performance analysis of space-to-ground optical satellite communications,” in _International Conference on Advanced Technologies for Communications (ATC)_. IEEE, pp. 284–289, 2023.
  * [25] H. D. Le, P. V. Trinh, T. V. Pham, D. R. Kolev, A. Carrasco-Casado, T. Kubo-Oka, M. Toyoshima, and A. T. Pham, “Throughput analysis for tcp over the fso-based satellite-assisted internet of vehicles,” _IEEE Transactions on Vehicular Technology_ , vol. 71, no. 2, pp. 1875–1890, Dec. 2021\.
  * [26] P. G. Madoery, J. A. Fraire, J. M. Finochietto, H. Yanikomeroglu, and G. K. Kurt, “A Novel Non-Terrestrial Networks Architecture: All Optical LEO Constellations with High-Altitude Ground Stations,” _arXiv preprint arXiv:2403.15659_ , 2024.
  * [27] M. T. Dabiri, M. Hasna, S. Althunibat, and K. Qaraqe, “Modulating retroreflector-based satellite-to-ground optical communications: Sensing and positioning,” _IEEE Trans. Commun._ , 2024 (early access).
  * [28] H. Safi, A. Dargahi, and J. Cheng, “Beam tracking for UAV-assisted FSO links with a four-quadrant detector,” _IEEE Communications Letters_ , vol. 25, no. 12, pp. 3908–3912, Sep. 2021.
  * [29] M.-C. Tsai, M. S. Bashir, and M.-S. Alouini, “Improved angle-of-arrival estimation of narrow gaussian beams for mobile fso platforms,” _IEEE Transactions on Aerospace and Electronic Systems_ , (early access), Feb. 2024.
  * [30] B. Epple, “Using a GPS-aided inertial system for coarse-pointing of free-space optical communication terminals,” in _Free-Space Laser Communications VI_ , vol. 6304, pp. 306–315, Sep. 2006.
  * [31] B. L. Wilkerson, D. Giggenbach, and B. Epple, “Concepts for fast acquisition in optical communications systems,” in _Free-Space Laser Communications VI_ , vol. 6304, pp. 67–78, Sep. 2006.
  * [32] Y. Kaymak, R. Rojas-Cessa, J. Feng, N. Ansari, M. Zhou, and T. Zhang, “A survey on acquisition, tracking, and pointing mechanisms for mobile free-space optical communications,” _IEEE Communications Surveys & Tutorials_, vol. 20, no. 2, pp. 1104–1123, Feb. 2018.
  * [33] H. Safi, A. Dargahi, J. Cheng, and M. Safari, “Analytical channel model and link design optimization for ground-to-HAP free-space optical communications,” _IEEE/OSA J. Lightw. Technol._ , vol. 38, no. 18, pp. 5036–5047, May 2020.
  * [34] J. C. Ricklin and F. M. Davidson, “Atmospheric turbulence effects on a partially coherent gaussian beam: implications for free-space laser communication,” _JOSA A_ , vol. 19, no. 9, pp. 1794–1802, 2002.
  * [35] M. T. Dabiri, S. M. S. Sadough, and M. A. Khalighi, “Channel modeling and parameter optimization for hovering UAV-based free-space optical links,” _IEEE J. Sel. Areas Commun._ , vol. 36, no. 9, pp. 2104–2113, Sep. 2018.
  * [36] A. Khadjavi, “Calculation of solid angle subtended by rectangular apertures,” _JOSA_ , vol. 58, no. 10, pp. 1417–1418, Oct. 1968.
  * [37] R. Gagliardi and M. Sheikh, “Pointing error statistics in optical beam tracking,” _IEEE Trans. Aerosp. Electron. Syst._ , no. 5, pp. 674–682, Sep. 1980.
  * [38] M. S. Bashir and M. R. Bell, “Optical beam position estimation in free-space optical communication,” _IEEE Trans. Aerosp. Electron. Syst._ , vol. 52, no. 6, pp. 2896–2905, July 2016.
  * [39] F. M. Davidson and X. Sun, “Gaussian approximation versus nearly exact performance analysis of optical communication systems with PPM signaling and APD receivers,” _IEEE Trans. Commun._ , vol. 36, no. 11, pp. 1185–1192, Nov. 1988.
  * [40] Z. Ghassemlooy, W. Popoola, and S. Rajbhandari, _Optical Wireless Communications: System and Channel Modelling with MATLAB®_. CRC Press, 2012\.
  * [41] T. Song and P.-Y. Kam, “A robust GLRT receiver with implicit channel estimation and automatic threshold adjustment for the free space optical channel with IM/DD,” _IEEE/OSA J. Lightw. Technol._ , vol. 32, no. 3, pp. 369–383, Feb. 2014.