Vol.0 (20xx) No.0, 000–000

11institutetext: Yunnan Observatories, Chinese Academy of Sciences, 396
Yangfangwang, Guandu District, Kunming, 650216, P. R. China<EMAIL_ADDRESS>
22institutetext: Key Laboratory for the Structure and Evolution of Celestial
Objects, Chinese Academy of Sciences, 396 Yangfangwang, Guandu District,
Kunming, 650216, P. R. China<EMAIL_ADDRESS>
33institutetext: Center for Astronomical Mega-Science, Chinese Academy of
Sciences, 20A Datun Road, Chaoyang District, Beijing, 100012, P. R. China;
44institutetext: University of Chinese Academy of Sciences, Beijing, 100049,
P. R. China
Received 20xx month day; accepted 20xx month day

# New simulations of the X-ray spectra and polarizations of accretion-disc
corona systems with various geometrical configurations I. Model Description

Xiao-Lin Yang 11223344 Jian-Cheng Wang 11223344 Chu-Yuan Yang 112233

###### Abstract

Energetic X-ray radiations emitted from various accretion systems are widely
considered to be produced by Comptonization in the hot corona. The corona and
its interaction with the disc play an essential role in the evolution of the
system and are potentially responsible for many observed features. However
many intrinsic properties of the corona are still poorly understood,
especially for the geometrical configurations. The traditional spectral
fitting method is not powerful enough to distinguish various configurations.
In this paper we intent to investigate the possible configurations by modeling
the polarization properties of X-ray radiations. The geometries of the corona
include the slab, sphere and cylinder. The simulations are implemented through
the publicly available code, LEMON, which can deal with the polarized
radiative transfer and different electron distributions readily. The results
demonstrate clearly that the observed polarizations are dependent on the
geometry of the corona heavily. The slab-like corona produces the highest
polarization degrees, the following are the cylinder and sphere. One of the
interesting things is that the polarization degrees first increase gradually
and then decrease with the increase of photon energy. For slab geometry there
exists a zero point where the polarization vanishes and the polarization angle
rotates for $90^{\circ}$. These results may potentially be verified by the
upcoming missions for polarized X-ray observations, such as $IXPE$ and $eXTP$.

###### keywords:

polarization — radiative transfer — radiation mechanisms: non-thermal —
relativistic processes — scattering — X-rays: galaxies

##  1 Introduction

Active galactic nuclei (AGNs), $\gamma$-ray bursts (GRBs) and X-ray binaries
are the most powerful X-ray objects in the universe. Accreting ambient
materials and then releasing gravitational energy by the central compact
objects is one of the most energetic phenomena in astrophysics. The released
energies will eventually heat the accreting gases and result in the radiations
range from radio to gamma-rays. The X-rays are widely believed to be produced
by the inverse Comptonization of soft photons in the corona which is a hot
region close to the central objects (e.g. Haardt & Maraschi (1991, 1993);
Eardley et al. (1975); Thorne & Price (1975); Gilfanov (2010)). And the soft
photons are usually multi-temperature blackbody emissions and come from the
accretion disc (Shakura & Sunyaev (1973); Page & Thorne (1974); Abramowicz et
al. (1988); Yuan & Narayan (2014)).

The corona plays a very important role in the disc-corona system. While due to
the complicated physical processes involved in the accretion system, the
evolution, formation and heating of the corona and its geometrical
configurations are still under debate (e.g. see the discussions given by
Dreyer & Böttcher (2021); Ursini et al. (2022); You et al. (2021)). Various
physical processes can lead to the formation of a corona and they show similar
spectral profiles or spectral energy distribution (SED). For example the
corona with an extended slab-like geometry is usually sited above the disc and
may be a result of magnetic instabilities (Galeev et al. 1979; Di Matteo
1998). Materials accreted around a neutron star instead of a black hole will
accumulate and finally form into a transition layer which shows the
characteristics of a corona (Sunyaev & Revnivtsev 2000; Long et al. 2022). At
the vicinity of the black hole, the quite active magnetic processes could
release considerable energies by magnetic reconnection which can heat the
plasma into a very high temperature (e.g. Wilkins & Fabian (2012)). The corona
can even be formed by the evaporation of the inner part of an accretion disc,
or as the transfer region between the jet and the black hole, either as a
failed jet (Ghisellini et al. 2004) or as a standing shock wave (Miyamoto &
Kitamoto 1991; Fender et al. 1999; Done et al. 2007). The geometrical
configurations of these corona models are deeply connected with their physical
origins. Thus, to distinguish the geometry of corona from observable
quantities will provide significant constrains on the physics of the accretion
system (e.g. Ursini et al. (2022); Long et al. (2022); Dreyer & Böttcher
(2021)).

Former researches have put some constrains on the geometry of corona. For
example, the size and location of X-ray corona have been estimated to be
within a few gravitational radii by the microlensing observations (Kochanek
2004; Reis & Miller 2013; Chartas et al. 2016). Comparing the time lags
between the direct and reflected radiations (or radiations in different energy
bands) can provide further constrains on the geometrical parameters of the
disc-corona system (e.g. Ingram et al. (2019); Mastroserio et al. (2021)).
However, the polarization of X-ray radiations is an alternative and unique way
to give possible new constrains on the corona geometry (Dreyer & Böttcher
2021; Long et al. 2022; Ursini et al. 2022; You et al. 2021), since the
polarizations induced by the inverse Comptonization are intrinsically
dependent on geometry and electron distribution (Schnittman & Krolik 2010;
Laurent et al. 2011; Beheshtipour et al. 2017). The Compton scattering can be
simply divided as Thomson and Klein-Nishina regimes according to the energy of
incident photons and the cross sections are intrinsically polarization
dependent (Fano 1949; Chandrasekhar 1960). It can induce polarizations for
anisotropic and unpolarized photons that scatter off non-relativistic
electrons (Bonometto et al. 1970; Schnittman & Krolik 2009). For photons
scattering off energetic electrons the polarization will be suppressed due to
the beaming effect (e.g. Dreyer & Böttcher (2021)). Thus the polarized
radiations among the X-ray bands would be reasonably expected (Poutanen &
Vilhu 1993; Schnittman & Krolik 2010; Beheshtipour et al. 2017; Ursini et al.
2022) and they can be used as a useful probe to distinguish the geometrical
configurations of the corona-disc systems.

Following the previous studies, we are here motivated to provide constrains on
the corona geometry by modeling the observed X-ray spectra and polarizations.
Our paper is organized as follows. The model and method are introduced in
Section 2. In the Section 3, we present the results of our calculations. The
discussions and conclusions are finally provided in Section 4.

##  2 Model and Method

In this section, we give an introduction to the model and method used in this
paper, which are based on our public available code Lemon (Yang et al. 2021).
Here we mainly discuss how to generate photons effectively in these
geometrical configurations for the scattering, and show the estimation
procedures.

###  2.1 the Geometries of the Corona

Figure 1: The schematic diagrams for the geometries of the corona with a slab
(top panel) and cylinder (bottom panel) configurations. The red and blue
regions represent the accretion disc and the corona, respectively. The black
hole is represented by the central black disc.

Figure 2: The schematic descriptions for the geometrical relationships between
the spherical corona and the accretion disc for various heights $H$ of sphere
center. The red and blue regions represent the accretion disc and the corona,
respectively. Seed photons are emitted from the disc surface and only those
with their momentum direction $\mathbf{p}$ (the blue vector) falling into the
cone of $\beta$ can reach the corona atmosphere. Paths of photons with a red
momentum vector will not intersect with the corona and should be rejected.

In this paper, we will calculate the observed spectra and polarizations from
three kinds of coronas with slab, sphere and cylinder geometries,
respectively. Their geometrical configurations are shown in Figure 1 and 2.
The slab-like corona is a thin layer and sandwiches the accretion disc
completely. Its geomertry is determined by the inner and outer radii: $R_{\rm
in},\,\,R_{\rm out}$, and the height $H$. The cylinder-like corona is located
above the disc with the top and bottom surfaces putting at height $H$ and
$H+H_{c}$. The radius of the cylinder is $R_{c}$. As the top surface of the
cylinder is set to be sufficiently high, $H_{c}$ is also very large, thus the
effects of its minor changes on the results can be ignored. The corona with a
sphere geometry is shown in Fig. 2, which is described by the sphere radius
$R_{\rm sp}$ and the height $H$ of the sphere center with respect to the disc.

The physical parameters to describe these coronas are the electron temperature
$kT_{e}$, the number density $n_{e}$ and the optical depth $\tau$ of Thomson
scattering. For the sake of simplicity, both $kT_{e}$ and $n_{e}$ for all of
three kinds of coronas are set to be constants throughout the corona. The
values of the Thomson optical depth $\tau$ for three cases are fixed and given
by:
$\tau=\sigma_{T}n_{e}H,\tau=\sigma_{T}n_{e}R_{c},\tau=\sigma_{T}n_{e}R_{\rm
sp}$, respectively, where $\sigma_{T}$ is the Thomson cross section. Hence, as
$\tau$ and $H,R_{c},R_{\rm sp}$ are provided, one can calculate the electron
number density by

$n_{e}=\frac{\tau}{\sigma_{T}H},\,\,n_{e}=\frac{\tau}{\sigma_{T}R_{c}},\,\,n_{e}=\frac{\tau}{\sigma_{T}R_{\rm
sp}},$ (1)

and vice versa.

###  2.2 Photon Generation

For all the configurations, the low energy seed photons are emitted from the
standard geometrically thin and optically thick accretion disc (Shakura &
Sunyaev 1973). We assume that the accretion disc is in a multi-temperature
state and its surface temperature changes with the disc radius $R$. The
distribution function of the temperature is given by (Shakura & Sunyaev 1973;
Tamborra et al. 2018)

$\displaystyle\begin{split}T(R)=\left[\frac{3GM\dot{m}}{8\pi
R^{3}r_{g}^{3}\sigma_{\rm SB}}\left(1-\sqrt{\frac{R_{\rm
in}}{R}}\right)\right]^{1/4},\end{split}$ (2)

where $r_{g}$ is the gravitational radius, $M$ is the mass of the central
black hole, $G$ is the gravitational constant, $\sigma_{\rm SB}$ is the
Stefan-Boltzmann constant and $\dot{m}$ is the accretion rate. The radius $R$
is in unit of $r_{g}$. For convenience, we will use the energy release rate
$\eta$ and the Eddington luminosity $L_{\rm EDD}$ to replace $\dot{m}$ as
$\eta L_{\rm EDD}=\dot{m}c^{2}$, where $c$ is the light speed. The Eddington
luminosity is defined by $L_{\rm EDD}=4\pi GMm_{p}c/\sigma_{T}$, where $m_{p}$
is the proton mass (Yuan & Narayan 2014). Then the expression of $T(R)$ can be
rewritten as

$\displaystyle\begin{split}T(R)=\left[\frac{\eta
C}{R^{3}}\left(1-\sqrt{\frac{R_{\rm in}}{R}}\right)\right]^{1/4},\end{split}$
(3)

where $C=3L_{\rm EDD}/(8\pi r_{g}^{2}\sigma_{\rm SB})$. In all the
simulations, the values of $\eta$ and $M$ are given in advance and set to be
$0.1$ and 10 $M_{\sun}$, corresponding to a mass accretion rate of
$\dot{m}\approx 1.40\times 10^{18}$ g$/$s (see table 1).

In order to describe the Keplerian motions of the accretion disc, we define a
static reference frame, whose basis vectors are given by $\mathbf{e}_{x}$,
$\mathbf{e}_{y}$ and $\mathbf{e}_{z}$. With $\mathbf{e}_{i}$, we can further
define a local reference frame at radius $R$ and azimuth angle $\varphi$, the
basis vectors of which are constructed by

$\displaystyle\left\\{\begin{array}[]{rl}\mathbf{e}^{\prime}_{x}&=-(\cos\varphi\mathbf{e}_{x}+\sin\varphi\mathbf{e}_{y}),\\\
\mathbf{e}^{\prime}_{y}&=\sin\varphi\mathbf{e}_{x}-\cos\varphi\mathbf{e}_{y},\\\
\mathbf{e}^{\prime}_{z}&=\mathbf{e}_{z}.\end{array}\right.$ (7)

Then with respect to the static reference frame, the Keplerian velocity of the
disc at $R$ and $\varphi$ can be expressed as

$\displaystyle\begin{split}\mathbf{V}_{k}(R,\varphi)=-\frac{c}{\sqrt{R}}\mathbf{e}^{\prime}_{y}.\end{split}$
(8)

Within the comoving frame of the disc at $(R,\varphi)$, we assume that the
emissivity $\widetilde{j}$ along direction
$\widetilde{\mathbf{\Omega}}=(\widetilde{\theta},\widetilde{\phi})$ and at
frequency $\widetilde{\nu}$ is a modified black body radiation given by

$\displaystyle\begin{split}\widetilde{j}(R,\varphi,\widetilde{\mathbf{\Omega}},\widetilde{\nu})=\frac{2h\widetilde{\nu}^{3}}{c^{2}}\frac{1+a\widetilde{\mu}}{e^{h\widetilde{\nu}/k_{\rm
B}T(R)}-1},\end{split}$ (9)

where $a=2.06$ is a constant and $\widetilde{\mu}=\cos\widetilde{\theta}$.
Factor $1+a\widetilde{\mu}$ means that the emissivity obeys the limb-darkening
law (Tamborra et al. 2018). Then the corresponding emissivity $j$ in the
static reference can be obtained through a Lorentz transformation (Pomraning
1973)

$\displaystyle\begin{split}j(R,\varphi,\mathbf{\Omega},\nu)=(\gamma
D_{p})^{2}\widetilde{j}(R,\varphi,\widetilde{\mathbf{\Omega}},\widetilde{\nu}),\end{split}$
(10)

where $\gamma=1/\sqrt{1-V_{k}^{2}/c^{2}}$ is the Lorentz factor,
$D_{p}=1+\widetilde{\mathbf{\Omega}}\cdot\mathbf{V}_{k}/c$ and

$\displaystyle\begin{split}\widetilde{\nu}&=\nu\gamma\left(1-\frac{\mathbf{\Omega}\cdot\mathbf{V}_{k}}{c}\right),\\\
\widetilde{\mathbf{\Omega}}&=\frac{1}{\gamma
D_{m}}\left[\mathbf{\Omega}-\frac{\gamma(1+\gamma
D_{m})}{1+\gamma}\frac{\mathbf{V}_{k}}{c}\right],\end{split}$ (11)

where $D_{m}=1-\mathbf{\Omega}\cdot\mathbf{V}_{k}/c$. Notice that the
quantities related to the Lorentz transformation are defined with respect to
the frame $\mathbf{e}^{\prime}_{i}$ given by Equation (7), using Equation (8),
the above equations can be written explicitly as

$\displaystyle\left\\{\begin{array}[]{rl}\widetilde{\nu}&=\displaystyle\nu\gamma\left(1+\Omega_{y}\beta\right),\\\
\widetilde{\Omega}_{x}&=\displaystyle\frac{\Omega_{x}}{\gamma D_{m}},\\\
\widetilde{\Omega}_{y}&=\displaystyle\frac{1}{D_{m}}(\Omega_{y}+\beta),\\\
\widetilde{\Omega}_{z}&=\displaystyle\frac{\Omega_{z}}{\gamma
D_{m}},\end{array}\right.$ (16)

where $\beta=V_{k}/c=1/\sqrt{R}$, $D_{m}=1+\Omega_{y}\beta$ and
$D_{p}=1-\widetilde{\Omega}_{y}\beta$.

In our former paper (Yang et al. 2021) we have explained that the Monte Carlo
radiative transfer is actually equivalent to the evaluation of the Neumann
solution of the radiative transfer equation. Each term of Neumann solution is
a multiple integral which can be written as

$\begin{split}I_{m}=\int n(P_{0})K(P_{0}\rightarrow P_{1})\cdots
K(P_{m-1}\rightarrow P_{m})f(P_{m})dP_{0}dP_{1}\cdots dP_{m},\end{split}$ (17)

where $n(P_{0})$ is the emitted photon number density and related to the
emissivity $j(P_{0})$ through $n(P_{0})=j(P_{0})/h\nu$, where $h$ is the Plank
constant, $K(P\rightarrow P^{\prime})$ is the transfer kernel, $f(P_{m})$ is
the recording function and $P=(\mathbf{r},\mathbf{\Omega},\nu)$, which are the
position and momentum vectors and frequency of the photon, respectively. The
generation of photons are actually related to the calculation of the integral
in terms of $P_{0}$ by Monte Carlo method, which can be separately written as

$\begin{split}&\int
n(\mathbf{r},\mathbf{\Omega},\nu)\delta(z)d\mathbf{r}d\mathbf{\Omega}d\nu,\end{split}$
(18)

where a $\delta$-function $\delta(z)$ is inserted, since the seed photons are
emitted by the disc that is located on the equatorial plane. We assume that
$\nu_{1}\leq\nu\leq\nu_{2}$, since for sufficient small and large frequency
$\nu$, the contributions from the black body radiation can be ignored.
$\nu_{1}$ and $\nu_{2}$ are free parameters and set to be
$h\nu_{1}=10^{-6}m_{e}c^{2}$, $h\nu_{2}=10^{-2}m_{e}c^{2}$, respectively. For
convenience, we introduce a new variable $y$ to replace $\nu$ and
$\nu=10^{y}$, $y_{1}\leq y\leq y_{2}$, where $y_{1}=\log_{10}(\nu_{1})$,
$y_{2}=\log_{10}(\nu_{2})$. Then Equation (18) can be rewritten as

$\begin{split}&\int
n(\mathbf{r},\mathbf{\Omega},\nu)\delta(z)d\mathbf{r}d\mathbf{\Omega}d\nu,\\\
=&\ln 10\int n(R,\varphi;\mu,\phi;\nu)\nu RdRd\varphi d\mu d\phi
dy.\end{split}$ (19)

To utilize Monte Carlo method to evaluate the above integral, we split the
integrand into two parts, one is used as PDFs for $R,\varphi,\mu,\phi$ and
$y$, the other one is used as weight for the integral. For $R$, $\varphi$ and
$y$, we assign them with the PDFs given by

$\left\\{\begin{array}[]{lr}p(R)=\displaystyle\frac{R}{\mathcal{N}_{R}},&\,\,R_{\rm
in}\leq R\leq R_{\rm out},\\\
p(\varphi)=\displaystyle\frac{1}{2\pi},&\,\,0\leq\phi\leq 2\pi,\\\
p(y)=\displaystyle\frac{1}{y_{2}-y_{1}},&\,\,y_{1}\leq y\leq
y_{2},\end{array}\right.$ (20)

where $\mathcal{N}_{R}=(R^{2}_{\rm out}-R_{\rm in}^{2})/2$ is the
normalization factor. Then $R,\varphi$ and $y$ can be sampled directly by

$R_{e}=\sqrt{R_{\rm
in}^{2}+2\mathcal{N}_{R}\xi_{1}},\quad\varphi_{e}=2\pi\xi_{2},\quad
y_{e}=y_{1}+(y_{2}-y_{1})\xi_{3},$ (21)

where $\xi_{i}$ are random numbers, whose PDFs are $p(\xi)=1$ and
$0\leq\xi\leq 1$. From now on, we will use $\xi$ to represent random numbers,
unless otherwise stated. With $R_{e},\varphi_{e}$ and $y_{e}$, the position
vector of the emission site and the frequency of the photon can be obtained as

$\displaystyle\begin{split}\mathbf{r}_{e}&=R_{e}(\cos\varphi_{e}\mathbf{e}_{x}+\sin\varphi_{e}\mathbf{e}_{y}),\\\
\nu_{e}&=10^{y_{e}}.\end{split}$ (22)

The sampling of the initial direction $\mathbf{\Omega}_{e}=(\mu_{e},\phi_{e})$
is more complicated and will be discussed in the following sections for the
three geometrical configurations, respectively.

While in order to discuss the initial weight $w_{\rm ini}$, we suppose that
$(\mu_{e},\phi_{e})$ has already been obtained. Then $w_{\rm ini}$ equals the
remaining part of the integrand of Equation (19), i.e.,

$\displaystyle\begin{split}w_{\rm ini}&=\ln
10\frac{n(R_{e},\varphi_{e};\mu_{e},\phi_{e};\nu_{e})\nu_{e}R_{e}}{p(R_{e})p(\varphi_{e})p(y_{e})},\\\
&=\mathcal{C}j(R_{e},\varphi_{e};\mu_{e},\phi_{e};\nu_{e}),\end{split}$ (23)

where $\mathcal{C}=\mathcal{N}_{R}2\pi(y_{2}-y_{1})\ln 10/h$ and $j$ is given
by Equaton (10).

####  2.2.1 Slab Case

For the slab corona, the procedure is quite simple, since any photons emitted
by the disc will enter the corona automatically. We first construct a local
triad at $\mathbf{r}_{\rm e}$ by

$\displaystyle\begin{split}\mathbf{e}^{\prime}_{x}=-\hat{\mathbf{r}}_{\rm
e},\quad\mathbf{e}^{\prime}_{z}=\mathbf{e}_{z},\quad\mathbf{e}^{\prime}_{y}=\mathbf{e}^{\prime}_{z}\times\mathbf{e}^{\prime}_{x},\end{split}$
(24)

where $\hat{\mathbf{r}}_{e}=\mathbf{r}_{e}/R_{e}$ is the unit vector of
$\mathbf{r}_{\rm e}$. With this triad one can obtain
$\mathbf{\Omega}_{e}=(\mu_{e},\phi_{e})$ directly

$\displaystyle\begin{split}\mu_{e}=\xi_{1},\quad\phi_{e}=2\pi\xi_{2},\end{split}$
(25)

and

$\displaystyle\begin{split}\mathbf{\Omega}_{e}=\sqrt{1-\mu^{2}_{e}}\cos\phi_{e}\mathbf{e}^{\prime}_{x}+\sqrt{1-\mu^{2}_{e}}\sin\phi_{e}\mathbf{e}^{\prime}_{y}+\mu_{e}\mathbf{e}^{\prime}_{z}.\end{split}$
(26)

####  2.2.2 Cylinder Case

Figure 3: The region formed by the effective photon momentum directions in the
$\phi^{\prime}$-$\mu^{\prime}$ plane of the triad defined by Equation 24 for
the cylinder case. The top panel shows the case as $R_{e}>R_{c}$ and the
bottom panel as $R_{e}\leq R_{c}$. The functions of the boundary curves of
this region are given by Equation 30. A sample of momentum directions will be
accepted if it falls into the grey region, otherwise it will be rejected. The
sampling algorithms are given by Equations 28 and 31. The area $S(R_{e})$ of
this region is a function of $R_{e}$ and should be multiplied with $w_{\rm
ini}$ to regulate the weight.

Comparing to the slab, the sampling procedures of $\mathbf{\Omega}_{e}$ for
the cylinder and sphere cases are more complicated. This is because if we
sample the emission direction $\mathbf{\Omega}_{e}$ isotropically in
$\mathbf{e}^{\prime}_{xyz}$, the efficiency will be quite low, since many
photon samples will miss the cylinder (or sphere) directly. To increase the
efficiency, for the cylinder case, we need to get the region formed by the
effective directions of the photon in the
$\mu^{\prime}_{e}$-$\phi^{\prime}_{e}$ plane. Here, a direction denoted by
$\mathbf{\Omega}_{e}$ is effective, it means that a photon assigned with this
direction can reach the corona eventually. Using the geometrical definitions
given in Figure 1, we can obtain the region in the
$\mu^{\prime}_{e}$-$\phi^{\prime}_{e}$ plane directly, which is shown in
Figure 3. As $R_{e}>R_{c}$ (top panel of Figure 3) the blue and red curves are
the boundaries of this region and their function expressions are given by

$\displaystyle\begin{split}\mu^{\prime}_{\rm
low}(\phi^{\prime})&=\frac{H}{\sqrt{\left(R_{e}\cos\phi^{\prime}+\sqrt{R_{c}^{2}-R_{e}^{2}\sin^{2}\phi^{\prime}}\right)^{2}+H^{2}}},\,\,-\phi^{\prime}_{1}\leq\phi^{\prime}\leq\phi^{\prime}_{1},\\\
\mu^{\prime}_{\rm
up}(\phi^{\prime})&=\frac{H+H_{c}}{\sqrt{\left(R_{e}\cos\phi^{\prime}-\sqrt{R_{c}^{2}-R_{e}^{2}\sin^{2}\phi^{\prime}}\right)^{2}+(H+H_{c})^{2}}},\,\,-\phi_{1}\leq\phi^{\prime}\leq\phi_{1},\end{split}$
(27)

where $\phi^{\prime}_{1}=\arccos(\sqrt{R_{e}^{2}-R_{c}^{2}}/R_{e})$. Then the
direction is sampled by

$\displaystyle\begin{split}&{\rm do}\\\ &\quad\quad\mu^{\prime}=\mu_{\rm
low}+(\mu_{\rm up}-\mu_{\rm
low})\xi_{1},\quad\phi^{\prime}=-\phi^{\prime}_{1}+2\phi^{\prime}_{1}\xi_{2},\\\
&\quad\quad{\rm if(\mu^{\prime}_{\rm
low}(\phi^{\prime})\leq\mu^{\prime}\,\,and\,\,\mu^{\prime}\leq\mu^{\prime}_{\rm
up}(\phi^{\prime}))exit},\\\ &{\rm end\,\,do}\end{split}$ (28)

where

$\displaystyle\begin{split}\mu_{\rm low}=\mu^{\prime}_{\rm
low}(0)&=\frac{H}{\sqrt{\left(R_{e}+R_{c}\right)^{2}+H^{2}}},\\\ \mu_{\rm
up}=\mu^{\prime}_{\rm
up}(0)&=\frac{H+H_{c}}{\sqrt{\left(R_{e}-R_{c}\right)^{2}+(H+H_{c})^{2}}},\end{split}$
(29)

which are the minimum and maximum of two functions given by Equation 30. As
$R_{e}\leq R_{c}$ (bottom panel of Figure 3), the functions of boundary curves
become

$\displaystyle\begin{split}\mu^{\prime}_{\rm
low}(\phi^{\prime})&=\frac{H}{\sqrt{\left(R_{e}\cos\phi^{\prime}+\sqrt{R_{c}^{2}-R_{e}^{2}\sin^{2}\phi^{\prime}}\right)^{2}+H^{2}}},\,\,&-\pi\leq\phi^{\prime}\leq\pi,\\\
\mu^{\prime}_{\rm
up}(\phi^{\prime})&=1,\,\,&-\pi\leq\phi^{\prime}\leq\pi,\end{split}$ (30)

and the sampling procedure becomes

$\displaystyle\begin{split}&{\rm do}\\\ &\quad\quad\mu^{\prime}=\mu_{\rm
low}+(1-\mu_{\rm low})\xi_{1},\quad\phi^{\prime}=-\pi+2\pi\xi_{2},\\\
&\quad\quad{\rm if(\mu^{\prime}_{\rm
low}(\phi^{\prime})\leq\mu^{\prime})exit},\\\ &{\rm end\,\,do}\end{split}$
(31)

Once $\mu^{\prime}$ and $\phi^{\prime}$ are obtained, $\mathbf{\Omega}_{e}$
can be constructed from Equations (26).

There is a subtle thing that needs to be clarified, i.e., to make the final
results correct one must update the initial weight $w_{\rm ini}$ by
multiplying a factor $S(R_{e})$, namely, $w^{\prime}_{\rm ini}=w_{\rm
ini}S(R_{e})$. And $S(R_{e})$ is the area of the grey region in Figure 3 and
can be calculated by

$\displaystyle\begin{split}S(R_{e})=\left\\{\begin{array}[]{lll}&\displaystyle\int_{-\phi^{\prime}_{1}}^{\phi^{\prime}_{1}}\left[\mu^{\prime}_{\rm
up}(\phi^{\prime})-\mu^{\prime}_{\rm
low}(\phi^{\prime})\right]d\phi^{\prime},&R_{e}>R_{c},\\\ &\displaystyle
2\pi-\int_{-\pi}^{\pi}\mu^{\prime}_{\rm
low}(\phi^{\prime})d\phi^{\prime},&R_{e}\leq
R_{c}.\end{array}\right.\end{split}$ (32)

####  2.2.3 Sphere Case

Figure 4: The same as Figure 3, but for the corona with a sphere geometry. The
effective grey regions for panels from top to bottom correspond to the three
geometrical configurations shown in Figure 2. The blue boundary curve in the
top panel is given by Equation (42).

The procedure for the sphere case is similar to that of the cylinder case, but
the local triad at $\mathbf{r}_{\rm e}$ is constructed in a different way by
(see Figure 2)

$\displaystyle\begin{split}\left.\begin{array}[]{ll}\mathbf{e}^{\prime}_{x}=-\mathbf{e}_{z},\quad\mathbf{e}^{\prime}_{z}=-\hat{\mathbf{r}}_{\rm
e},\quad\mathbf{e}^{\prime}_{y}=\mathbf{e}^{\prime}_{z}\times\mathbf{e}^{\prime}_{x}.\end{array}\right.\end{split}$
(33)

In addition, we should construct another triad at $\mathbf{r}_{\rm e}$ defined
through

$\displaystyle\begin{split}\left\\{\begin{array}[]{ll}\mathbf{e}^{\prime\prime}_{x}&=\cos\alpha\mathbf{e}^{\prime}_{x}+\sin\alpha\mathbf{e}^{\prime}_{z},\\\
\mathbf{e}^{\prime\prime}_{z}&=-\sin\alpha\mathbf{e}^{\prime}_{x}+\cos\alpha\mathbf{e}^{\prime}_{z},\\\
\mathbf{e}^{\prime\prime}_{y}&=\mathbf{e}^{\prime}_{y},\end{array}\right.\end{split}$
(34)

where the definition of $\alpha$ is shown in Figure 2. In this triad frame,
all of the effective directions $(\mu^{\prime\prime},\phi^{\prime\prime})$
also form a grey region in the $\mu^{\prime\prime}$-$\phi^{\prime\prime}$
plane (see, Figure. 4). Then we can draw an effective direction
$(\mu^{\prime\prime},\phi^{\prime\prime})$ isotropically in the rectangle
region, $[-\pi,\pi]\times[\mu_{\beta},1]$, by the following algorithm (for
case (a))

$\displaystyle\begin{split}&{\rm do}\\\
&\quad\quad\mu^{\prime\prime}=\mu_{\beta}+(1-\mu_{\beta})\xi_{1},\quad\phi^{\prime\prime}=-\pi+2\pi\xi_{2},\\\
&\quad\quad{\rm if(|\phi^{\prime\prime}|>\phi_{1})exit}\\\ &\quad\quad{\rm
else\,\,\,if(\mu^{\prime\prime}>\mu^{\prime\prime}(\phi^{\prime\prime}))exit}\\\
&{\rm end\,\,do},\end{split}$ (35)

where $\mu_{\beta}=\cos\beta=\sqrt{1-R^{2}_{\rm sp}/(R^{2}_{e}+h^{2})}$, and
$h=OH$ is the height of the sphere center of the corona above the disc, and

$\displaystyle\begin{split}\phi_{1}=\arccos\left(\tan\alpha\cot\beta\right).\end{split}$
(36)

A direction $(\mu^{\prime\prime},\phi^{\prime\prime})$ will be accepted if it
falls into the grey region, otherwise it will be rejected. For case (b)

$\displaystyle\begin{split}\mu^{\prime\prime}&=\mu_{\beta}+(1-\mu_{\beta})\xi_{1},\quad\phi^{\prime\prime}=-\pi+2\pi\xi_{2}.\end{split}$
(37)

For case (c)

$\displaystyle\begin{split}&\mu^{\prime\prime}=\mu_{\beta}+(1-\mu_{\beta})\xi_{1},\quad\phi^{\prime\prime}=\pi/2+\pi\xi_{2},\\\
&{\rm
if(\phi^{\prime\prime}>\pi)\phi^{\prime\prime}=\phi^{\prime\prime}-2\pi}.\end{split}$
(38)

With $(\mu^{\prime\prime},\phi^{\prime\prime})$, the initial direction of a
photon can be similarly written as

$\displaystyle\begin{split}\hat{\mathbf{p}}_{e}&=p^{\prime\prime}_{x}\mathbf{e}^{\prime\prime}_{x}+p^{\prime\prime}_{y}\mathbf{e}^{\prime\prime}_{y}+p^{\prime\prime}_{z}\mathbf{e}^{\prime\prime}_{z}=p_{x}\mathbf{e}_{x}+p_{y}\mathbf{e}_{y}+p_{z}\mathbf{e}_{z},\end{split}$
(39)

where

$\displaystyle\begin{split}\left\\{\begin{array}[]{ll}p^{\prime\prime}_{x}&=\sqrt{1-\mu^{\prime\prime
2}}\cos\phi^{\prime\prime},\\\ p^{\prime\prime}_{y}&=\sqrt{1-\mu^{\prime\prime
2}}\sin\phi^{\prime\prime},\\\
p^{\prime\prime}_{z}&=\mu^{\prime\prime}.\end{array}\right.\end{split}$ (40)

From Equations (33), (34), (39), one can obtain that

$\begin{split}\left\\{\begin{array}[]{ll}p_{x}&=-(\sin\alpha\cos\phi_{\rm
e}p^{\prime\prime}_{x}-\sin\phi_{\rm
e}p^{\prime\prime}_{y}+\cos\alpha\cos\phi_{\rm e}p^{\prime\prime}_{z}),\\\
p_{y}&=-(\sin\alpha\sin\phi_{\rm e}p^{\prime\prime}_{x}+\cos\phi_{\rm
e}p^{\prime\prime}_{y}+\cos\alpha\sin\phi_{\rm e}p^{\prime\prime}_{z}),\\\
p_{z}&=-\cos\alpha p^{\prime\prime}_{x}+\sin\alpha
p^{\prime\prime}_{z}.\end{array}\right.\end{split}$ (41)

In panel $(a)$ of Figure 4, there is a boundary curve plotted with blue color.
The expression of this curve can be simply derived. From Figure 1, one can see
that a photon reaching the corona sphere must satisfy the condition that
$p_{z}\geq 0$ and $p_{z}=0$ yields the boundary curve. Using the expressions
of $p_{z}$ given by Equation (41) and $p^{\prime\prime}_{x}$,
$p^{\prime\prime}_{z}$ given by Equation (40), from $p_{z}=0$, we can obtain
the function of the curve

$\displaystyle\begin{split}\mu^{\prime\prime}(\phi^{\prime\prime})=\frac{\cos\alpha\cos\phi^{\prime\prime}}{\sqrt{1-\cos^{2}\alpha\sin^{2}\phi^{\prime\prime}}}.\end{split}$
(42)

Finally, we need to update the weight $w_{\rm ini}$ by using the areas
$S(R_{e})$ of the grey region in Figure 4 and the analytical expressions of
$S(R_{e})$ for three panels can be obtained as

$\displaystyle\begin{split}S(R_{e})=\left\\{\begin{array}[]{lll}&\displaystyle
2\pi(1-\mu_{\beta})-2\arcsin(\cos\alpha\sin\phi_{1})+2\mu_{\beta}\phi_{1},&(a),\\\
&\displaystyle 2\pi(1-\mu_{\beta}),&(b),\\\
&\displaystyle\pi(1-\mu_{\beta}),&(c).\end{array}\right.\end{split}$ (43)

###  2.3 Scattering distance sampling

As a photon is generated, we then discuss how to determine the distance
between any two scattering points randomly. As discussed in (Yang et al.
2021), the scattering distance is a random variable, and its probability
density distribution function (PDF) is given by

$\begin{split}p(s)ds=\frac{1}{N}\exp\left(-\int_{0}^{s}\sigma_{a}(s^{\prime})n_{e}(s^{\prime})ds^{\prime}\right)\sigma_{a}(s)n_{e}(s)ds,\end{split}$
(44)

where $\sigma_{a}(s)=\sigma_{a}[h\nu,T_{e}(s)]$ is the averaged scattering
cross section (Hua 1997), $T_{e}(s)$ and $n_{e}(s)$ are the temperature and
number density of hot electrons at $s$, respectively, $N$ is the normalization
factor and

$\begin{split}N=1-\exp\left(-\int_{0}^{s_{m}}\sigma_{a}(s^{\prime})n_{e}(s^{\prime})ds^{\prime}\right),\end{split}$
(45)

where $s_{m}$ is the maximum distance that the path can extend in the corona
region. Usually, $s_{m}$ can be obtained by solving an algebraic equation
derived from the geometrical conditions, i.e., with the provided vectors
$\hat{\mathbf{p}}_{\rm ini},\mathbf{r}_{\rm ini}$ of the initial position and
momentum direction, we have $\mathbf{r}_{\rm ini}+\hat{\mathbf{p}}_{\rm
ini}s_{m}=\mathbf{r}_{b}$, where $\mathbf{r}_{b}$ is a position vector that
falls on the boundary of the corona. Taking the square for both sides, we get
a quadratic equation of $s_{m}$ as

$\begin{split}s_{m}^{2}+2s_{m}\mathbf{r}_{\rm ini}\cdot\hat{\mathbf{p}}_{\rm
ini}+r^{2}_{\rm ini}=r^{2}_{b}.\end{split}$ (46)

In our model, we adopt the assumption that both $n_{e}$ and $kT_{e}(s)$ are
constants in the corona, then the above Equation becomes

$\begin{split}p(s)ds=\frac{1}{N}\exp\left(-\sigma_{a}n_{e}s\right)\sigma_{a}n_{e}ds,\end{split}$
(47)

where $N=1-\exp(-\sigma_{a}n_{e}s_{m})$, which can be sampled by the inverse
cumulated distribution function (CDF) method directly. That is

$\begin{split}s=-\frac{1}{n_{e}\sigma_{a}}\ln(1-N\xi).\end{split}$ (48)

When the scattering distance is determined, the weight $w_{\rm ini}$ should be
updated as well by multiplying the normalization factor $N$, i.e., $w_{\rm
ini}^{\prime}=w_{\rm ini}\cdot N$.

###  2.4 Scattering sampling

In our model, we will consider inverse Comptonization of the photon in the
sphere corona with the electrons assigned with various distribution functions.
As discussed in Yang et al. (2021), Lemon can incorporate any kind of
scattering readily. Here, we mainly consider three kinds of electron
distributions, i.e., the relativistic thermal, the $\kappa$ and the power law
distributions, which are respectively given by (Xiao 2006; Pandya et al. 2016)

$\displaystyle\frac{dn_{e}}{d\gamma
d\Omega}=\left\\{\begin{array}[]{ll}&\displaystyle\frac{n_{e}}{4\pi}\frac{\gamma\sqrt{\gamma^{2}-1}}{\Theta
K_{2}(\Theta)}\exp(-\gamma/\Theta),\\\
&\displaystyle\frac{n_{e}}{4\pi}\frac{\gamma\sqrt{\gamma^{2}-1}}{N_{\kappa}}\left(1+\frac{\gamma-1}{\kappa
w}\right)^{-(\kappa+1)},\\\
&\displaystyle\frac{n_{e}}{4\pi}\frac{1}{N_{p}}\gamma^{-\alpha_{e}},\end{array}\right.$
(52)

where $\gamma=1/\sqrt{1-v^{2}/c^{2}}$ is the Lorentz factor,
$\gamma_{1}\leq\gamma\leq\gamma_{2}$, $d\Omega=\sin\theta d\theta d\phi$ is
the solid angle, $\Theta=kT/(m_{e}c^{2})$ is the dimensionless temperature,
$\kappa,w$ are free parameters, $\gamma_{1},\,\,\gamma_{2}$ are boundary
values and $N_{\kappa},\,\,N_{p}$ are the normalization factors. In general,
we assume that $\gamma$ of these distributions is confined between
$\gamma_{1}$ and $\gamma_{2}$. $N_{p}$ can be written as
$N_{p}=(\gamma_{1}^{1-\alpha_{e}}-\gamma_{2}^{1-\alpha_{e}})/(\alpha_{e}-1)$,
while $N_{\kappa}$ will be evaluated numerically since its analytical
expression involves special functions and is quite complicated (Pandya et al.
2016). In Yang et al. (2021), we have discussed how to sample $dn_{e}/(d\gamma
d\Omega)$ by a method proposed by Hua (1997). Here we present a new method to
deal with these two PDFs in a unified way. This method is a combination of
inverse CDF and rejection method. For the thermal distribution, an auxiliary
function is introduced as

$\begin{split}f_{1}(\gamma)&=\frac{1}{A}\gamma^{2}\exp(-\gamma/\Theta),\end{split}$
(53)

where $A$ is the normalization factor. Then the procedure of the method is
given by (Schnittman & Krolik 2013)

$\displaystyle\begin{split}&①{\rm
get\,\,a\,\,\gamma_{f_{1}}\,\,by\,\,sampling\,\,f_{1}(\gamma)},\\\ &②{\rm
if\left(\gamma_{f_{1}}\xi\leq\sqrt{\gamma_{f_{1}}^{2}-1}\right)accept\,\,\gamma_{f_{1}}}\\\
&{\rm\quad else\,\,\,goto\,\,①}.\end{split}$ (54)

The algorithm of sampling $\gamma_{f_{1}}$ from $f_{1}(\gamma)$ is the inverse
CDF method, which is equivalent to solve the following algebraic equation
(Schnittman & Krolik 2013)

$\begin{split}g(u)=g(u_{1})-\xi[g(u_{1})-g(u_{2})],\end{split}$ (55)

where $u=\gamma/\Theta$, $g(u)=(u^{2}+2u+2)\exp(-u)$,
$u_{1}=\gamma_{1}/\Theta$ and $u_{2}=\gamma_{2}/\Theta$. It turns out that
this equation can be solved by iterative method numerically. To accomplish it,
we rewrite the above equation as

$\begin{split}u=\ln(u^{2}+2u+2)-\ln G_{0},\end{split}$ (56)

where $G_{0}=g(u_{1})-\xi[g(u_{1})-g(u_{2})]$. Then the root $u_{0}$ of
equation. (55) can be obtained by the following algorithm:

$\displaystyle\begin{split}&x_{1}=0,\\\ &{\rm do}\\\
&\quad\quad{x_{2}=\ln(x_{1}^{2}+2x_{1}+2)-\ln G_{0}}\\\ &\quad\quad{\rm
if}(|x_{2}-x_{1}|<\epsilon){\rm exit}\\\ &\quad\quad{x_{1}=x_{2}}\\\ &{\rm
end\,\,do}\end{split}$ (57)

and $u_{0}\approx x_{2}$, where $\epsilon$ is the tolerance value for the
accuracy of the root. Usually, we set it to be $\epsilon=10^{-10}$. Once the
root $u_{0}$ is obtained, the sample of the Lorentz factor is given by
$\gamma_{f_{1}}=u_{0}\Theta$.

For the $\kappa$ distribution, the algorithm is similar and the auxiliary
function is given by

$\begin{split}f_{\kappa}(\gamma)&=\frac{1}{A}\frac{\gamma^{2}}{(\gamma+a)^{\kappa+1}},\end{split}$
(58)

where $a=\kappa w-1$ and $A$ is the normalization factor. Also one needs to
solve the following equation to get a trial sample of $\gamma$,

$\begin{split}\gamma=\left[\frac{1}{d_{1}G_{0}}(a_{1}\gamma^{2}+b_{1}\gamma+c_{1})\right]^{1/\kappa}-a,\end{split}$
(59)

where $a_{1}=\kappa(\kappa-1)$, $b_{1}=2a\kappa$, $c_{1}=2a^{2}$,
$d_{1}=\kappa(\kappa-1)(\kappa-2)$ and
$G_{0}=g(\gamma_{1})-\xi[g(\gamma_{1})-g(\gamma_{2})]$, the function
$g(\gamma)$ is given by

$\displaystyle\begin{split}g(\gamma)=\frac{a_{1}\gamma^{2}+b_{1}\gamma+c_{1}}{d_{1}(\gamma+a)^{\kappa}}.\end{split}$
(60)

Equation (59) can be numerically solved by iterative method as well and the
algorithm is similar

$\displaystyle\begin{split}&x_{1}=1,\\\ &{\rm do}\\\
&\quad\quad{x_{2}=\left[\frac{1}{d_{1}G_{0}}(a_{1}x_{1}^{2}+b_{1}x_{1}+c_{1})\right]^{1/\kappa}-a}\\\
&\quad\quad{\rm if}(|x_{2}-x_{1}|<\epsilon){\rm exit}\\\
&\quad\quad{x_{1}=x_{2}}\\\ &{\rm end\,\,do}\end{split}$ (61)

and $\gamma_{f_{\kappa}}\approx x_{2}$. For the power law distribution, the
samples of $\gamma$ can be obtained directly by inverse CDF method as

$\displaystyle\begin{split}\gamma=\left[\gamma_{1}^{1-\alpha_{e}}(1-\xi)+\gamma_{2}^{1-\alpha_{e}})\right]^{\frac{1}{1-\alpha_{e}}}.\end{split}$
(62)

The momentum direction (MD) $(\theta,\phi)$ of the electron can be drawn
isotropically as follows

$\displaystyle\begin{split}\mu=-1+2\xi_{1},\quad\phi=2\pi\xi_{2},\end{split}$
(63)

where $\mu=\cos\theta$. Then we calculate the ratio of the total cross section
of Klein-Nishina and Thomson by (Hua 1997)

$\begin{split}\frac{\sigma_{\rm
KN}}{\sigma_{T}}=\frac{3}{4\epsilon}\left[\left(1-\frac{4}{\epsilon}-\frac{8}{\epsilon^{2}}\right)\ln(1+\epsilon)+\frac{1}{2}+\frac{8}{\epsilon}-\frac{1}{2(1+\epsilon)^{2}}\right],\end{split}$
(64)

where

$\begin{split}\epsilon=\frac{2h\nu}{m_{e}c^{2}}\gamma(1-v\mu),\end{split}$
(65)

and $h\nu$ is the energy of the incident photon, $v=\sqrt{1-1/\gamma^{2}}$ is
speed of the electron. A random number $\xi$ is generated to determine the
acceptence of $\gamma,\mu,\phi$ as the Lorentz factor and velocity direction
of the scattering electron. If $\xi\leq\sigma_{\rm KN}/\sigma_{T}$ they are
accepted, otherwise rejected.

###  2.5 Lorentz Transformation of Stokes parameters

The polarization states are described by the Stokes parameters (SPs):
$\mathbf{S}=(I,Q,U,V)^{T}$ and the polarization vector (PV) $\mathbf{f}$. The
Compton scattering will inevitably change the polarization states of the
photons. In this subsection we will discuss how to describe and trace these
changes in Lemon in a detailed way. These discussions however can be found in,
e.g., Krawczynski (2012). For the purpose of completeness, we shall give these
descriptions in a more consistent way as follows.

1\. As the MD $\mathbf{k}$ of the incident photon is given, we first construct
the photon triad as
$\mathbf{e}_{z}(p)=\mathbf{k},\,\mathbf{e}_{x}(p)=\mathbf{f},\,\mathbf{e}_{y}(p)=\mathbf{e}_{z}(p)\times\mathbf{e}_{x}(p)$.
For unpolarized radiations, the basis vector $\mathbf{e}_{x}(p)$ is set as
$\mathbf{e}_{x}(p)=\mathbf{k}\times\mathbf{z}/|\mathbf{k}\times\mathbf{z}|$,
where $\mathbf{z}=(0,0,1)$.

2\. Then in the photon triad $\mathbf{e}_{i}(p)$ we obtain the MD
$\mathbf{p}_{e}$ and Lorentz factor $\gamma_{e}$ of the scattering electron by
sampling the distribution function given by Equation (52) as discussed in the
above subsection. Equivalently, $\mathbf{p}_{e}$ can be expressed by
$(\mu_{e},\phi_{e})$, i.e.,
$\mathbf{p}_{e}=\sqrt{1-\mu_{e}^{2}}\cos\phi_{e}\mathbf{e}_{x}(p)+\sqrt{1-\mu_{e}^{2}}\sin\phi_{e}\mathbf{e}_{y}(p)+\mu_{e}\mathbf{e}_{z}(p)$.
With $\mathbf{p}_{e}$ we can construct the static electron triad with respect
to the photon triad $\mathbf{e}_{i}(p)$ as
$\mathbf{e}_{z}(e)=\mathbf{p}_{e},\,\mathbf{e}_{y}(e)=\mathbf{k}\times\mathbf{p}_{e}/|\mathbf{k}\times\mathbf{p}_{e}|,\,\mathbf{e}_{x}(e)=\mathbf{e}_{y}(e)\times\mathbf{e}_{z}(e)$.
We denote the MD of the incident photon as $(\mu_{e},\phi_{e})$ in this triad.

3\. To complete the Compton scattering, we need to transform the SPs into the
rest frame of the electron. To do so, we should first carry out a rotation and
get the SPs defined with respect to the electron triad $\mathbf{e}_{i}(e)$,
that is $\mathbf{S}_{e}=\mathbf{M}(\phi_{e})\mathbf{S}$, where

$\begin{split}\mathbf{M}(\phi_{e})=\left(\begin{array}[]{cccc}1&0&0&0\\\
0&\cos 2\phi_{e}&\sin 2\phi_{e}&0\\\ 0&-\sin 2\phi_{e}&\cos 2\phi_{e}&0\\\
0&0&0&1\end{array}\right)\end{split}$ (66)

is the rotation matrix (see Chandrasekhar (1960)). At the same time the PV
$\mathbf{f}$ has been rotated into the plane of $\mathbf{p}_{e}$ and
$\mathbf{k}$ as well. As demonstrated by Krawczynski (2012), $\mathbf{S}_{e}$
will keep invariant under the Lorentz transformation, which means that we can
get the SPs in the rest frame of the electron directly as:
$\widetilde{\mathbf{S}}_{e}=\mathbf{S}_{e}$ (from now on, the quantities in
the rest frame of the electron are denoted with a tilde $\tilde{}$). The MD
$(\mu_{e},\phi_{e})$ and frequency $\nu$ of the incident photon will be
transformed as:

$\begin{split}\widetilde{\mu}_{e}&=\frac{\mu_{e}-v_{e}}{1-\mu_{e}v_{e}},\,\,\widetilde{\phi}_{e}=\phi_{e},\\\
\widetilde{\nu}&=\nu\gamma(1-\mu_{e}v_{e}),\end{split}$ (67)

where $v_{e}$ is the velocity of the electron.

4\. With $(\widetilde{\mu}_{e},\widetilde{\phi}_{e})$ we could reconstruct the
MD $\widetilde{\mathbf{k}}$ of the photon in the rest frame of the electron.
Using $\widetilde{\mathbf{k}}$, we can construct another photon triad with
respect to the electron rest frame as:
$\widetilde{\mathbf{e}}_{z}(p)=\widetilde{\mathbf{k}},\,\widetilde{\mathbf{e}}_{y}(p)=\widetilde{\mathbf{k}}\times\widetilde{\mathbf{z}}/|\widetilde{\mathbf{k}}\times\widetilde{\mathbf{z}}|,\,\widetilde{\mathbf{e}}_{x}(p)=\widetilde{\mathbf{e}}_{y}(p)\times\widetilde{\mathbf{e}}_{z}(p)$,
where $\widetilde{\mathbf{z}}=(0,0,1)$ is the $z$-axis of the electron frame.
In frame $\widetilde{\mathbf{e}}_{i}(p)$, we redenote the MD of the electron
as $\widetilde{\mathbf{p}}_{e}$.

5\. Simulating the Compton scattering in the rest frame of the electron, we
sample the Klein-Nishina differential cross section $d\sigma_{\rm
KN}/(d\widetilde{\mu}^{\prime}d\widetilde{\phi}^{\prime})$ to get the
scattered frequency $\widetilde{\nu}^{\prime}$ and MD
$(\widetilde{\mu}^{\prime}_{e},\widetilde{\phi}^{\prime}_{e})$ (see, e.g.,
Pozdnyakov et al. (1983) and Hua (1997)), which is defined with respect to the
triad $\widetilde{\mathbf{e}}_{z}(p)$ and

$\begin{split}\widetilde{\nu}^{\prime}=\frac{\widetilde{\nu}}{\displaystyle
1+h\widetilde{\nu}/(m_{e}c^{2})(1-\widetilde{\mu}^{\prime}_{e})}.\end{split}$
(68)

From $(\widetilde{\mu}^{\prime}_{e},\widetilde{\phi}^{\prime}_{e})$ we can get
the scattered MD vector $\widetilde{\mathbf{k}}^{\prime}$ and use it to
construct the scattered photon triad as
$\widetilde{\mathbf{e}}^{\prime}_{z}(p)=\widetilde{\mathbf{k}}^{\prime},\,\widetilde{\mathbf{e}}^{\prime}_{y}(p)=\widetilde{\mathbf{k}}^{\prime}\times\widetilde{\mathbf{k}}/|\widetilde{\mathbf{k}}^{\prime}\times\widetilde{\mathbf{k}}|,\,\widetilde{\mathbf{e}}^{\prime}_{x}(p)=\widetilde{\mathbf{e}}^{\prime}_{y}(p)\times\widetilde{\mathbf{e}}^{\prime}_{z}(p)$.

Then we carry out a rotation given by
$\widetilde{\mathbf{S}}_{p}=\mathbf{M}(\widetilde{\phi}^{\prime}_{e})\widetilde{\mathbf{S}}_{e}$
to get the SPs defined with respect to the scattering plane determined by
$\widetilde{\mathbf{k}}^{\prime}$ and $\widetilde{\mathbf{k}}$. From
$\widetilde{\mathbf{S}}_{p}$, the scattered SPs
$\widetilde{\mathbf{S}}^{\prime}_{p}$ can be obtained by the Fano’s Matrix as
(Fano 1949, 1957; McMaster 1961):

$\begin{split}\widetilde{\mathbf{S}}^{\prime}_{p}=\mathbf{F}(\widetilde{\nu},\widetilde{\nu}^{\prime},\widetilde{\theta}^{\prime}_{e})\widetilde{\mathbf{S}}_{p},\end{split}$
(69)

where

$\begin{split}\mathbf{F}(\widetilde{\nu},\widetilde{\nu}^{\prime},\widetilde{\theta}^{\prime}_{e})=\left(\frac{\widetilde{\nu}^{\prime}}{\widetilde{\nu}}\right)^{2}\left(\begin{array}[]{cccc}F_{0}&F_{3}&0&0\\\
F_{3}&F_{33}&0&0\\\ 0&0&F_{11}&0\\\ 0&0&0&F_{22}\end{array}\right)\end{split}$
(70)

and

$\begin{split}F_{0}&=\displaystyle\frac{\widetilde{\nu}^{\prime}}{\widetilde{\nu}}+\frac{\widetilde{\nu}}{\widetilde{\nu}^{\prime}}-\sin^{2}\widetilde{\theta}^{\prime}_{e},\quad
F_{3}=\sin^{2}\widetilde{\theta}^{\prime}_{e},\quad
F_{11}=2\cos\widetilde{\theta}^{\prime}_{e},\\\
F_{22}&=\left(\frac{\widetilde{\nu}^{\prime}}{\widetilde{\nu}}+\frac{\widetilde{\nu}}{\widetilde{\nu}^{\prime}}\right)\cos\widetilde{\theta}^{\prime}_{e},\quad
F_{33}=1+\cos^{2}\widetilde{\theta}^{\prime}_{e}.\end{split}$ (71)

6\. Now we construct a triad by using $\widetilde{\mathbf{k}}^{\prime}$ and
$\widetilde{\mathbf{p}}_{e}$, i.e.,
$\widetilde{\mathbf{e}}^{\prime\prime}_{z}(p)=\widetilde{\mathbf{k}}^{\prime},\,\widetilde{\mathbf{e}}^{\prime\prime}_{y}(p)=\widetilde{\mathbf{k}}^{\prime}\times\widetilde{\mathbf{p}}_{e}/|\widetilde{\mathbf{k}}^{\prime}\times\widetilde{\mathbf{p}}_{e}|,\,\widetilde{\mathbf{e}}^{\prime\prime}_{x}(p)=\widetilde{\mathbf{e}}^{\prime\prime}_{y}(p)\times\widetilde{\mathbf{e}}^{\prime\prime}_{z}(p)$.
We then can obtain the angle between the plane
$\widetilde{\mathbf{k}}^{\prime}$-$\widetilde{\mathbf{p}}_{e}$ and plane
$\widetilde{\mathbf{k}}^{\prime}$-$\widetilde{\mathbf{k}}$ (or equivalently
the angle between basis vector $\widetilde{\mathbf{e}}^{\prime}_{x}(p)$ and
$\widetilde{\mathbf{e}}^{\prime\prime}_{x}(p)$) as (refer to Equation. (8) of
Krawczynski (2012)):

$\begin{split}\phi_{0}=-{\rm
sign}[\widetilde{\mathbf{e}}^{\prime}_{x}(p)\cdot\widetilde{\mathbf{e}}^{\prime\prime}_{y}(p)]\arccos[\widetilde{\mathbf{e}}^{\prime}_{y}(p)\cdot\widetilde{\mathbf{e}}^{\prime\prime}_{y}(p)].\end{split}$
(72)

With $\phi_{0}$ we can do a rotation to get the SPs defined with respect to
the triad $\widetilde{\mathbf{e}}^{\prime\prime}_{i}(p)$ as
$\widetilde{\mathbf{S}}^{\prime\prime}_{p}=\mathbf{M}(\phi_{0})\widetilde{\mathbf{S}}^{\prime}_{p}$,
which could be transformed back into the static frame directly.

7\. At this stage, the Compton scattering has been completed and we implement
another Lorentz transformation to bring all the quantities back into the
static reference. For the scattered MD $\widetilde{\mathbf{k}}^{\prime}$ we
first transform it from the $\widetilde{\mathbf{e}}_{y}(p)$ frame to the
electron rest frame through

$\begin{split}\widetilde{k}^{\prime m}(e)=\widetilde{k}^{\prime
i}(p)\widetilde{e}_{(i)}^{m}(p),\end{split}$ (73)

where $\widetilde{e}_{(i)}^{m}(p)$ are the components of the $i$-th basis
vector $\widetilde{\mathbf{e}}_{i}(p)$ with respect to the electron rest frame
and $\widetilde{k}^{\prime m}(p),\widetilde{k}^{\prime i}(e)$ are the
components of $\widetilde{\mathbf{k}}^{\prime}$ with respect to the photon and
electron triad respectively. Then the Lorentz transformation of
$\widetilde{\nu}^{\prime}$ and $\widetilde{k}^{\prime m}(e)$ can be written as

$\begin{split}&\nu^{\prime}=\widetilde{\nu}^{\prime}\gamma[1+\widetilde{k}^{\prime
z}(e)v_{e}],\,\,k^{\prime z}(e)=\frac{\widetilde{k}^{\prime
z}(e)+v_{e}}{1+\widetilde{k}^{\prime z}(e)v_{e}},\\\ &k^{\prime
x}(e)=\frac{\widetilde{k}^{\prime x}(e)}{\gamma[1+\widetilde{k}^{\prime
z}(e)v_{e}]},\,\,k^{\prime y}(e)=\frac{\widetilde{k}^{\prime
z}(e)}{\gamma[1+\widetilde{k}^{\prime z}(e)v_{e}]}.\end{split}$ (74)

The scattered direction in the static frame can be obtained from

$\begin{split}k^{\prime m}=k^{\prime
i}(e)e_{(i)}^{j}(e)e_{(j)}^{m}(p),\end{split}$ (75)

and the scattered vector $\mathbf{k}^{\prime}=k^{\prime m}\mathbf{e}_{m}$. The
scattered SPs in the static frame are given by
$\widetilde{\mathbf{S}}^{\prime\prime}_{p}$ and the scattered polarization
vector $\mathbf{f}^{\prime}$ is obviously in the
$\mathbf{p}_{e}$-$\mathbf{k}^{\prime}$ plane and can be expressed as
$\mathbf{f}^{\prime}=[\mathbf{p}_{e}-(\mathbf{p}_{e}\cdot\mathbf{k}^{\prime})\mathbf{k}^{\prime}]/N$,
where $N$ is the normalization factor. As $\mathbf{k}^{\prime}$,
$\widetilde{\mathbf{S}}^{\prime\prime}_{p}$ and $\mathbf{f}^{\prime}$ are
obtained, one can continue to trace the next transfer and scattering of the
photon. Simultaneously, one can record the contributions of the photon made to
the observed quantities at the scattering point, which will be discussed in
the next section.

###  2.6 Spectrum and polarization estimation

Lemon used a scheme that can improve the efficiency and accuracy of spectrum
and polarization evaluations, since the information at any scattering point
can perform contributions to the observed spectrum and polarization. Under the
Neumann expansion solution of differential-integral equation, the scheme is
equal to the introduction of a $\delta$-function and recording function (for
more detailed discussions, refer to Yang et al. (2021)). This scheme has
actually been applied widely and implemented in many codes dealing with
Ly$\alpha$ radiative transfer (e.g., see Seon et al. (2022), where the scheme
is named as “peeling-off technique”, also known as “next event estimation” or
“shadow rays” (Yusef-Zadeh et al. 1984; Laursen & Sommer-Larsen 2007; Yajima
et al. 2012), one may also refer to Whitney (2011) and Noebauer & Sim (2019)
for reviews). We will use this scheme to reduce the noise generated by the
Monte Carlo method and obtain the results with high signal-to-noise ratio.

Now we discuss the specific procedures of the estimation scheme using the
conventions and triads established in the last subsection. The observer is
assumed to be located at a direction $\mathbf{n}_{\rm obs}(\mu_{\rm
obs},\phi_{\rm obs})$. Due to the axial symmetry of the system, we can choose
$\phi_{\rm obs}$ randomly, i.e., $\phi_{\rm obs}=2\pi\xi$. Then using the
coefficients of the photon and electron triad, we can transform
$\mathbf{n}_{\rm obs}$ into the two triads directly by

$\begin{split}n_{\rm obs}^{(m)}(p)&=n_{\rm obs}^{i}e_{i}^{(m)}(p),\\\ n_{\rm
obs}^{(m)}(e)&=n_{\rm obs}^{i}e_{i}^{(j)}(p)e_{j}^{(m)}(e),\end{split}$ (76)

where $e_{i}^{(m)}$ is the inverse of matrix $e_{(m)}^{i}$. Since all of the
triads are orthonormal, we have $e_{i}^{(m)}e_{(m)}^{j}=\delta_{i}^{j}$,
$e_{i}^{(m)}e_{(n)}^{i}=\delta_{n}^{m}$ and $e_{i}^{(m)}={e^{T}}_{(m)}^{i}$,
where the superscript $T$ represents the matrix transpose. Through a Lorentz
transformation, we get the observer direction in the rest frame of the
electron as

$\begin{split}\widetilde{n}_{\rm obs}^{(x)}(e)&=\frac{n_{\rm
obs}^{(z)}(e)}{\gamma D},\,\,\widetilde{n}_{\rm obs}^{(y)}(e)=\frac{n_{\rm
obs}^{(z)}(e)}{\gamma D},\\\ \widetilde{n}_{\rm obs}^{(z)}(e)&=\frac{n_{\rm
obs}^{(z)}(e)-v_{e}}{D}.\end{split}$ (77)

where $D=1-n_{\rm obs}^{(z)}(e)v_{e}$. Using the triad matrix
$\widetilde{e}^{i}_{(m)}(p)$ of $\widetilde{\mathbf{e}}_{i}(p)$, we can
further transform the observer direction into the photon frame as

$\begin{split}\widetilde{n}_{\rm obs}^{(m)}(p)=\widetilde{n}_{\rm
obs}^{(n)}(e)\widetilde{e}_{n}^{(m)}(p),\end{split}$ (78)

or in the vector form $\widetilde{\mathbf{n}}_{\rm obs}(p)=\widetilde{n}_{\rm
obs}^{(m)}(p)\widetilde{\mathbf{e}}_{m}(p)$. Then the cosine of the scattering
angle in the rest frame of the electron is given by
$\widetilde{\mu}^{\prime}_{e}=\cos\widetilde{\theta}^{\prime}_{e}=\widetilde{n}_{\rm
obs}^{(z)}(p)$, from which we can get the scattered frequency in the electron
rest and static frame respectively as

$\begin{split}\widetilde{\nu}^{\prime}_{\rm
obs}&=\frac{\widetilde{\nu}}{\displaystyle
1+h\widetilde{\nu}/(m_{e}c^{2})(1-\widetilde{\mu}^{\prime}_{e})},\\\
\nu^{\prime}_{\rm obs}&=\widetilde{\nu}^{\prime}_{\rm
obs}\gamma[1+\widetilde{k}^{\prime z}(e)v_{e}],\end{split}$ (79)

where $\widetilde{\nu}$ is the frequency of the incident photon given by
Equation (67), $m_{e}c^{2}$ is the rest energy of the electron. To get SPs for
a given observer direction, we should rotate the SPs
$\widetilde{\mathbf{S}}_{e}$ into scattering plane by
$\widetilde{\mathbf{S}}_{p}^{\rm
obs}=\mathbf{M}(\widetilde{\phi}^{\prime}_{e})\widetilde{\mathbf{S}}_{e}$, and
obviously

$\begin{split}\cos\widetilde{\phi}^{\prime}_{e}=\frac{\widetilde{n}_{\rm
obs}^{(x)}(p)}{\sin\widetilde{\theta}^{\prime}_{e}},\quad\sin\widetilde{\phi}^{\prime}_{e}=\frac{\widetilde{n}_{\rm
obs}^{(y)}(p)}{\sin\widetilde{\theta}^{\prime}_{e}}.\end{split}$ (80)

Then the scattered SPs can be obtained by the Fano’s Matrix as

$\begin{split}\widetilde{\mathbf{S}}^{\prime\rm
obs}_{p}=\mathbf{F}(\widetilde{\nu},\widetilde{\nu}^{\prime}_{\rm
obs},\widetilde{\theta}^{\prime}_{e})\widetilde{\mathbf{S}}_{p}^{\rm
obs}.\end{split}$ (81)

With $\widetilde{\mathbf{n}}_{\rm obs}(p)$ we can construct two triads as
before:

$\begin{split}\widetilde{\mathbf{e}}^{\prime}_{z}(p)&=\widetilde{\mathbf{n}}_{\rm
obs}(p),\,\widetilde{\mathbf{e}}^{\prime}_{y}(p)=\widetilde{\mathbf{n}}_{\rm
obs}(p)\times\widetilde{\mathbf{k}}/|\widetilde{\mathbf{n}}_{\rm
obs}(p)\times\widetilde{\mathbf{k}}|,\\\
\widetilde{\mathbf{e}}^{\prime}_{x}(p)&=\widetilde{\mathbf{e}}^{\prime}_{y}(p)\times\widetilde{\mathbf{e}}^{\prime}_{z}(p),\end{split}$
(82)

and

$\begin{split}\widetilde{\mathbf{e}}^{\prime\prime}_{z}(p)&=\widetilde{\mathbf{n}}_{\rm
obs}(p),\,\widetilde{\mathbf{e}}^{\prime\prime}_{y}(p)=\widetilde{\mathbf{n}}_{\rm
obs}(p)\times\widetilde{\mathbf{p}}_{e}/|\widetilde{\mathbf{n}}_{\rm
obs}(p)\times\widetilde{\mathbf{p}}_{e}|,\\\
\widetilde{\mathbf{e}}^{\prime\prime}_{x}(p)&=\widetilde{\mathbf{e}}^{\prime\prime}_{y}(p)\times\widetilde{\mathbf{e}}^{\prime\prime}_{z}(p),\end{split}$
(83)

from which we can get rotation angle $\phi_{0}$ given by Eq. (72) and then the
scattered SPs defined with respect to the plane of
$\widetilde{\mathbf{n}}_{\rm obs}(p)$ and $\widetilde{\mathbf{p}}_{e}$ through
$\widetilde{\mathbf{S}}^{\prime\prime\rm
obs}_{p}=\mathbf{M}(\phi_{0})\widetilde{\mathbf{S}}^{\prime}_{p}$. Finally,
the polarization vector is given by $\mathbf{f}^{\prime}_{\rm
obs}=[\mathbf{p}_{e}-(\mathbf{p}_{e}\cdot\mathbf{n}_{\rm obs})\mathbf{n}_{\rm
obs}]/N$.

As $\nu^{\prime}_{\rm obs}$, $\widetilde{\mathbf{S}}^{\prime\prime\rm
obs}_{p}$ and $\mathbf{f}^{\prime}_{\rm obs}$ are obtained, we can eventually
record the contribution made by this scattering site to the observed spectrum
and polarization.. Before that, we should do a final rotation to get the SPs
defined with respect to the static triad of the observer, which is constructed
by

$\begin{split}\mathbf{e}_{z}(\rm obs)&=\mathbf{n}_{\rm
obs},\,\mathbf{e}_{x}({\rm obs})=(-n_{y}^{\rm obs},n_{x}^{\rm
obs},0),\mathbf{e}_{y}({\rm obs})=\mathbf{e}_{z}({\rm
obs})\times\mathbf{e}_{x}({\rm obs}).\end{split}$ (84)

By the definition of $\mathbf{e}_{x}({\rm obs})$, we have assumed that the
polarization angle is measured from the a direction parallel to the disc in
the sky plane. A polarization vector along the north-south direction
corresponds to a polarization angle of 90 degrees. This definition is
different from the convension adopted by Ursini et al. (2022) with a 90∘
rotation. The triad associated with $\mathbf{f}^{\prime}_{\rm obs}$ can be
obtained by

$\begin{split}\mathbf{e}_{z}(p)&=\mathbf{n}_{\rm
obs},\,\mathbf{e}_{x}(p)=\mathbf{f}^{\prime}_{\rm
obs},\mathbf{e}_{y}(p)=\mathbf{e}_{z}(p)\times\mathbf{e}_{x}(p).\end{split}$
(85)

Then the rotation angle can be calculated by

$\begin{split}\phi_{0}=-{\rm sign}[\mathbf{e}_{x}(p)\cdot\mathbf{e}_{y}({\rm
obs})]\arccos[\mathbf{e}_{y}(p)\cdot\mathbf{e}_{y}({\rm obs})].\end{split}$
(86)

With $\phi_{0}$, we get the SPs defined with respect to the observer triad:
$\mathbf{S}^{\prime\prime\rm
obs}_{p}=\mathbf{M}(\phi_{0})\widetilde{\mathbf{S}}^{\prime\prime\rm
obs}_{p}$. Then one can synthesize the spectrum $\nu L_{\nu}$, polarization
degree (PD) $\delta$ and polarization angle (PA) by

$\left\\{\begin{array}[]{rl}\nu L_{\nu}&=I_{\rm
obs}=\displaystyle\sum\frac{I^{\prime\prime}_{\rm obs}w\sigma_{\rm
KN}\exp(-\tau_{\rm obs})h\nu}{\Delta\ln\nu},\\\ Q_{\rm
obs}&=\displaystyle\sum\frac{Q^{\prime\prime}_{\rm obs}w\sigma_{\rm
KN}\exp(-\tau_{\rm obs})h\nu}{\Delta\ln\nu},\\\ U_{\rm
obs}&=\displaystyle\sum\frac{U^{\prime\prime}_{\rm obs}w\sigma_{\rm
KN}\exp(-\tau_{\rm obs})h\nu}{\Delta\ln\nu},\\\
\delta&=\displaystyle\frac{Q_{\rm obs}}{I_{\rm obs}},\\\
PA&=\displaystyle\frac{1}{2}\arctan\left(\frac{U_{\rm obs}}{Q_{\rm
obs}}\right),\end{array}\right.$ (87)

where $(I^{\prime\prime}_{\rm obs},Q^{\prime\prime}_{\rm
obs},U^{\prime\prime}_{\rm obs})^{T}=\mathbf{S}^{\prime\prime\rm obs}_{p}$,
$w$ is the weight, $\sigma_{\rm KN}$ is the Klein-Nishina differential cross
section and $\tau_{\rm obs}$ is the optical depth from the scattering site to
the boundary (e.g. Yang et al. (2021)).

##  3 Results

In this paper we mainly investigate the effects of geometries on the observed
spectrum and polarization of the disc-corona system. The corona is assumed to
be assigned with a slab, a sphere and a cylinder geometry composed of hot
electron gas. In the following, we will present the primary results of our
models.

###  3.1 Example demonstrations

Figure 5: The spectra of a corona with spherical configuration for various scattering number. The parameters are: $R_{\rm in}$ = 6.0 $r_{g}$, $R_{\rm out}$ = 100.0 $r_{g}$, $H=0.0$ $r_{g}$, $\tau=1.0$, $T_{e}=100$ keV. Figure 6: The same as Figure 5, but for different viewing inclination angles. The top and bottom panels are the spectra and PDs, respectively. Table 1: Parameters values for all simulations. parameters |  values
---|---
$R_{\rm out}$ |  $200$ $r_{g}$
$R_{\rm in}$ |  $6\,\,r_{g}$
$\eta$ |  0.1
$\dot{m}$ |  $1.40\times 10^{18}\rm g/s$
$M_{\rm BH}$ |  $10\,\,M_{\sun}$

We first demonstrate the global pictures of the spectra for various viewing
angles and scattering numbers. To accomplish this, the values of all other
parameters must be fixed. The results are shown in Figure 5 and 6. In Figure 5
we show the spectra with various scattering numbers (top panel) viewed at
$\mu_{\rm obs}=0.5$. From the figures, one can see the typical characteristics
for a Comptonized spectrum, i.e., at high energy bands, the spectrum is
composed by a power law followed with a steep exponential high-energy cut-off
due to the KN effect (e.g. Fabian et al. (2015); Mastroserio et al. (2021)).
Also as the higher the energy is, the louder the noise is. The green line
represents the spectrum formed by the multi-temperature black body radiations
escaping from the disc directly. As the scattering number increases, the
corresponding spectrum becomes harder. In Figure 6, we show the spectra and
polarizations with respect to the inclination angles, where the curves with
different colors represent different cosine of the viewing angles. From the
figure one can see that as the inclination angle increases the PD $\delta$
will increase, but the intensity $I$ will decrease. This is because that
photons will go through a larger optical depth and thus suffer from more
scatterings at higher inclination.

###  3.2 Parameter Settings

The parameters of the disc are set to be the same for the three kinds of
coronas (see Table 1). The inner and outer radii of the disc are set to be
$R_{\rm in}=6$ $r_{g}$ and $R_{\rm out}=200$ $r_{g}$, respectively. The mass
of the central black hole and the mass accretion rate of the disc are 10
$M_{\sun}$ and $1.4\times 10^{18}$ g/s.

The slab-like corona is composed of parallel planes sandwiching and covering
the disc completely (Poutanen & Vilhu 1993; Schnittman & Krolik 2010). Then,
the inner and outer radii of the slab are set to be the same as the disc. The
height and the temperature of the slab are given by $h$ and $kT_{e}$. Former
researches show that the dependence of results on the height $h$ is not
significant (e.g. Ursini et al. (2022)). Hence, in all the calculations, we
will set $h$ to be $1\,\,r_{g}$. The Thomson optical depth $\tau$ of the slab
is given as $n_{e}\sigma_{T}h$. We will simulate the cases with
$\tau=0.5,\,\,1.0$ and $2.0$, respectively.

The geometry of corona with a spherical configuration is determined by its
radius $R_{\rm sp}$ and height $H$. For appropriately setting values of these
two parameters, the configuration of the corona can either be extending that
can fully cover the disk, or be compact and located above the disc. Many
evidences have shown that the size of the hot corona is most likely very small
and close to the central compact (Fabian et al. 2015; Ursini et al. 2022)
(Fabian et al. 2015; Ursini et al. 2020), which is the well known lamppost
model (Zdziarski et al. 1996; Życki et al. 1999). In other models, the size of
the sphere can be as large as fully covering the whole disc (e.g. Tamborra et
al. (2018)). Thus in our simulations, we will set $R_{\rm sp}$ = 10 $r_{g}$
for a compact configuration and $R_{\rm sp}$ = 200.0 $r_{g}$ for a extending
configuration. The optical depth $\tau=n_{e}\sigma_{T}R_{\rm sp}$ of the
sphere is measured from the center to the boundary and set to be
$\tau=0.5,\,\,1.0$ and $2.0$ as well.

The corona with a cylindrical geometry is usually used to describe the
outflowing materials, or a jet (Ghisellini et al. 2004). Then the corona will
be assigned with a bulk velocity $\beta=v/c$. As long as the motion is not
extremely relativistic, its influence on the results will not be significant
(Ursini et al. 2022). Thus in our simulations we can set $\beta=0$ and focus
on the impact of other relevant geometrical parameters on the results. The
Thomson optical depth $\tau=n_{e}\sigma_{T}R_{c}$ is measured along the
horizontal direction and is also set to be $\tau=0.5,\,\,1.0,\,\,2.0$. The
electron temperature is assigned with the value of $kT_{e}=100.0$ keV for all
the situations.

###  3.3 Spectra and Polarizations

With the parameter settings and assumptions given above, we carry out the
simulations to study the dependencies of observed spectra and polarizations on
the geometries of the corona. The results will be presented as follows.

Figure 7: The synthesized spectra (top panel), PDs (middle panel) and PAs
(bottom panel) of a corona with slab geometry viewing from various inclination
angles. The parameters are: $T_{e}=100$ keV, $\tau=0.5,\,\,1.0,\,\,2.0$ for
panels from left to right, the height $h=1.0$ $r_{g}$, the inner and outer
radii of the disc $R_{\rm in}=6.0$ $r_{g}$ and $R_{\rm out}=200.0$ $r_{g}$.
Figure 8: The same as Figure 7, but for a corona with sphere geometry. The
parameters for the sphere are: the radius of the sphere $R_{\rm sp}=200.0$
$r_{g}$, the height of the sphere center above the disc $h=0.0$ $r_{g}$.
Figure 9: The same as Figure 7, but for a corona with cylinder geometry. The
parameters for the sphere are: the radius $R_{c}=10$ $r_{g}$, the height
$H=20$ $r_{g}$ and the intrinsic height $H_{c}=100$ $r_{g}$.

Figure 10: The synthesized spectra of coronas with a slab (top panel), sphere
(middle panel) and cylinder (bottom panel) geometries for various scattering
numbers. The cosine of viewing angle is $\mu_{\rm obs}=0.5$. Other parameters
are the same as that given in Figure 7, 8 and 9 for the three kinds of
coronas, respectively.

We first present the spectra, PDs and PAs of the radiations emerging from
three kinds of coronas. The results are plotted in Figures 7, 8 and 9, which
show that as the optical depth increases, the spectra become harder for all
the cases. The profiles of the spectra are similar for three cases and also
for different inclination angles. This supports the perspective of
distinguishing the geometries for the coronas with different configurations
through fitting the spectral energy distributions is not so effective (see
also Tamborra et al. (2018); Dreyer & Böttcher (2021)). The polarizations of
the three geometrical configurations, however, show significant differences
both in the profile and the magnitude. One can see that in the X-ray bands,
which we are mostly interested in, the slab corona has the biggest PD, whose
value can be up to around 10 $\%$ depending on the viewing inclination. And
followed are the magnitude of PD for the sphere and cylinder coronas, which
goes to 1-2 % and is less than 1 $\%$, respectively. These results are
consistent with those given by Ursini et al. (2022). However, in the high
energy bands, the results change into the opposite situation, where the
cylinder corona has the highest PD up to almost 15 $\%$. One can see that
there exists zero points in PDs both for the slab and sphere coronas, as shown
in Figures 7 and 8, where the $Q$ component of the SPs vanishes and changes
its sign simultaneously. From the parameters in Figures 7 and 8, we can see
that both of the two coronas have extending geometrical configurations. This
may be taken as a special feature for this kind of coronas.

From Figures 7 and 8, one can conclude that the trends and profiles of PD for
the slab and the sphere are similar to each other but different to that of the
cylinder case. This is because the former two coronas have similar geometrical
configurations due to the parameters we chosen, i.e., they both have the
extended configurations that fully cover the disc. On the contrary, the
cylinder corona has a compact configuration above the disc, which makes its
irradiation by disc to be less isotropical compared with the extending
configurations. For the extending corona, a higher optical depth yields a
higher PD, while for the compact corona, the conclusion seems opposite (see
Figure 9, where the maximum of PD decreases as optical depth increases).

According to our definition of the reference direction of the PA, for all the
cases, the PA oscillates around zero and the polarization direction is
horizontal in the low energy bands. This result is simply due to the fact that
all of the geometrical configurations considered here are the axial-symmetry,
which will yield an orientation of the polarization vectors that can either be
horizontal or vertical (Connors et al. 1980). However, in the high energy
bands, the PA oscillates turbulently to lead the results unreliable.

In Figure 10, we give the spectra of all configurations for various scattering
numbers and optical depths. For the spherical and cylindrical configurations,
the profiles of the spectra are quite similar for different optical depth. But
for the slab configuration, the spectra become harder as the optical depth
increases.

Figure 11: The spectra (top row panels) and PDs (bottom row panels) of a
corona with a sphere geometry for different heights $H$ and viewing angles
$\mu_{\rm obs}$. The parameters are: the radius $R_{sp}=10.0$ $r_{g}$, the
potical depth $\tau=1.0$ and temperature $kT_{e}=100$ keV. The height $H$ of
the spherical center varies from $10.0$ to $80.0$ $r_{g}$ and the
corresponding results are plotted with different color and line styles. The
cosine of the viewing angles $\mu_{\rm obs}$ is $0.2$, 0.5 and 0.9 for panels
in the left, middle and right columns, respectively. Figure 12: The same as
Figure 11, but for a corona with cylinder geometry. The parameters are: the
radius $R_{c}=10$ $r_{g}$ and the intrinsic height $H_{c}=100$ $r_{g}$.

For the corona with a spherical or a cylindrical geometry, its height $H$
measured from the disc is an important parameter. The changes of $H$ will not
only alter the geometry of the system, but also affect the efficiency that the
disc illuminates the corona. For a higher $H$, the illumination is less
isotropical. Hence we do simulations to investigate the impact of $H$ on the
spectra and PDs. The results are shown in Figures 11 and 12 for the spherical
and cylindrical cases, respectively. From these figures one can see that as
$H$ increases, the spectra become softer in the high energy bands, due to the
reduction of the luminosity and flux that could be captured by the corona. As
the viewing angles change, the spectra tend to keep the same. But for the PDs,
the magnitude of variations is considerable. Avaragely, the PDs of the
cylinder corona are bigger than that of the sphere corona. The profiles of the
PDs for the two configurations are analogous, due to the quite similar
geometrical parameters we choose. However, for the sphere case, with the
increase of $H$, the $Q$ component will change its sign from positive to
negative. But for the cylinder case, the sign of $Q$ always keeps negative.
Also one can see that with the increase of $H$, the PDs in the low energy
bands keep almost unchanged, just as the trend of the spectra. However, the
changes of the PDs in the high energy bands are significant. This difference
of PDs between the low and high energy bands may be used as a probe to
discriminate a corona with or without a compact geometry.

##  4 Discussion and Conclusions

The geometrical configurations of corona in the accreting systems are poorly
understood due to the degeneracy of spectroscopy differentiating the system.
The polarimetry provides a more effective and powerful option to overcome this
dilemma. With the advent of the missions of IXPE (Weisskopf et al. 2016), eXTP
(Zhang et al. 2016, 2019), the era of high-quality data of polarization
observations is arriving in the near future. Hence it is necessary and urgent
to make the theoretical studies in advance.. For this purpose, we carry out
the simulations of radiative transfer in corona with different geometries.
These simulations are based on the public available Monte Carlo code, Lemon
(Yang et al. 2021), which is based on the Neumann series expansion solution of
differential-integral equations. By using the code, one can increase the
signal to noise ratio dramatically and simplify the calculations when the
configuration of the system has geometric symmetries. The main contents and
results of this paper can concluded as follows:

We have discussed detailedly how to simulate the polarized radiative transfer
in the three geometry configurations, namely, the slab, sphere and cylinder.
We emphasized how to generate photons efficiently for the sphere and cylinder
cases. To accomplish this, we have derived the shapes of the regions formed by
the effective momentum directions in the $\phi$-$\mu$ plane of the triad. This
method can increase the efficiency of the simulation, especially for coronas
with a compact geometry, since the solid angle subtended by the corona with
respect to the emission site is small, which further reduces the probability
that a photon can be received by the corona. In our model, we considered the
effect due to the Keplerian motion of the disc on the spectra. This effect
works by affecting the photon generation. In our model, it can be taken into
account readily through a Lorentz transformation, which connects the
emissivities in the comoving reference frame of the disc and the static
reference frame. The Keplerian motion of the disc will make the frequencies of
the seed photons to have blue or red shifts, which will further affect the
final results. While due to the low speed for most part of the disc ($R\gg
1$), the Keplerian motion seems to have a minor impact on the spectra and
polarizations. However if we focus on the radiations emitted from the most
inner part of the disc, both the Keplerian motion and general relativity
effects should be taken into account necessarily. Then we discussed how to
obtain the scattering distance between any two scattering sites by the inverse
CDF method. Because the electron distribution has a very important impact on
the Comptonization spectra, we proposed a new scheme to deal with three often
used distribution functions, i.e., the thermal, $\kappa$ and power law, in a
uniformed way. Next, we demonstrated how to implement the polarized Compton
scattering in the Klein-Nishina regime consistently and detailedly, which
involves the complicated triads constructions, Lorentz boosts, SPs
transformation and rotations. Finally, we discussed the procedure for
evaluating the contributions made by any scattering site to the observed
quantities when the inclination angle of the observer is provided.

We use our model to simulate the radiative transfer in three kinds of coronas
with different geometrical configurations. The results demonstrate that the
polarizations of the observed radiations are significantly dependent on the
geometries of the corona. Different configurations will produce PDs with
different magnitudes and profiles in the X-ray bands. The corona with an
extending configuration, such as the slab, yields a higher PD while the
compact one yields a less polarized result. With the increase of the photon
energy, the PD will increase gradually as well untill a maximum is reached.
After that, the PD will decrease to zero due to the relativistic beaming
effect (Dreyer & Böttcher 2021). The maximum of PD for the extending
configurations increases with the increase of the optical depth, but for the
compacting configurations the conclusion is the opposite.

Our results are consistent with those of the former researches. However, our
model and code are flexible and can deal with different geometrical
configurations readily. One just needs to modify the photon generation and
tracing parts. However, the results presented here are quite theoretical and
not connected with the practical observations. Also, our model does not
include the effects of the general relativity, which inevitably plays an
important role in the radiative transfer around a black hole. These effects
will be included in the future work. Nonetheless, our model includes the
essetial ingredients of Comptonization in a hot electron corona. Thus,
hopefully, our model will provide some useful insights for the observations of
the upcoming X-ray missions, such as $XIPE$ and $eXTP$.

##  5 Acknowledgments

We thank the anonymous referee for helpful comments and suggestions that
improve the manuscript significantly. We thank the Yunnan Observatories
Supercomputing Platform, on which our code was partly tested. We acknowledge
the ﬁnancial support from the National Natural Science Foundation of China
U2031111, 11573060, 12073069, 11661161010 and 11673060.

##  6 code AVAILABILITY

The source code of this article is part of Lemon and can be downloaded from:
https://bitbucket.org/yangxiaolinsc/corona_geometry/src/main/.

## References

  * Abramowicz et al. (1988) Abramowicz, M. A., Czerny, B., Lasota, J. P., et al. 1988, ApJ, 332, 646
  * Beheshtipour et al. (2017) Beheshtipour, B., Krawczynski, H., & Malzac, J. 2017, ApJ, 850, 14
  * Bonometto et al. (1970) Bonometto, S., Cazzola, P., & Saggion, A. 1970, A&A, 7, 292
  * Chartas et al. (2016) Chartas, G., Rhea, C., Kochanek, C., et al. 2016, Astronomische Nachrichten, 337, 356
  * Chandrasekhar (1960) Chandrasekhar, S. 1960, New York Dover: 1960
  * Connors et al. (1980) Connors, P. A., Piran, T., & Stark, R. F. 1980, ApJ, 235, 224
  * Di Matteo (1998) Di Matteo, T. 1998, MNRAS, 299, L15
  * Done et al. (2007) Done, C., Gierliński, M., & Kubota, A. 2007, A&A Rev., 15, 1
  * Dreyer & Böttcher (2021) Dreyer, L. & Böttcher, M. 2021a, ApJ, 906, 18
  * Eardley et al. (1975) Eardley, D. M., Lightman, A. P., & Shapiro, S. L. 1975, ApJ, 199, L153
  * Fabian et al. (2015) Fabian, A. C., Lohfink, A., Kara, E., et al. 2015, MNRAS, 451, 4375
  * Fano (1949) Fano, U. 1949, Journal of the Optical Society of America (1917-1983), 39, 859
  * Fano (1957) Fano, U. 1957, Reviews of Modern Physics, 29, 74
  * Fender et al. (1999) Fender, R., Corbel, S., Tzioumis, T., et al. 1999, ApJ, 519, L165
  * Galeev et al. (1979) Galeev, A. A., Rosner, R., & Vaiana, G. S. 1979, ApJ, 229, 318
  * Ghisellini et al. (2004) Ghisellini, G., Haardt, F., & Matt, G. 2004, A&A, 413, 535
  * Gilfanov (2010) Gilfanov, M. 2010, Lecture Notes in Physics, Berlin Springer Verlag, 17
  * Hua (1997) Hua, X. M. 1997, Comput. Phys., 11, 660
  * Haardt & Maraschi (1991) Haardt, F. & Maraschi, L. 1991, ApJ, 380, L51
  * Haardt & Maraschi (1993) Haardt, F. & Maraschi, L. 1993, ApJ, 413, 507
  * Ingram et al. (2019) Ingram, A., Mastroserio, G., Dauser, T., et al. 2019, MNRAS, 488, 324
  * Kochanek (2004) Kochanek, C. S. 2004, ApJ, 605, 58
  * Krawczynski (2012) Krawczynski, H. 2012, ApJ, 744, 30
  * Laurent et al. (2011) Laurent, P., Rodriguez, J., Wilms, J., et al. 2011, Science, 332, 438
  * Laursen & Sommer-Larsen (2007) Laursen, P. & Sommer-Larsen, J. 2007, ApJ, 657, L69
  * Long et al. (2022) Long, X., Feng, H., Li, H., et al. 2022, ApJ, 924, L13
  * Mastroserio et al. (2021) Mastroserio, G., Ingram, A., Wang, J., et al. 2021, MNRAS, 507, 55
  * McMaster (1961) McMaster, W. H. 1961, Reviews of Modern Physics, 33, 8
  * Miyamoto & Kitamoto (1991) Miyamoto, S. & Kitamoto, S. 1991, ApJ, 374, 741
  * Noebauer & Sim (2019) Noebauer, U. M. & Sim, S. A. 2019, Living Reviews in Computational Astrophysics, 5, 1
  * Page & Thorne (1974) Page, D. N. & Thorne, K. S. 1974, ApJ, 191, 499
  * Pandya et al. (2016) Pandya, A., Zhang, Z., Chandra, M., et al. 2016, ApJ, 822, 34
  * Pomraning (1973) Pomraning, G. C. 1973, International Series of Monographs in Natural Philosophy, Oxford: Pergamon Press, 1973
  * Poutanen & Vilhu (1993) Poutanen, J. & Vilhu, O. 1993, A&A, 275, 337
  * Pozdnyakov et al. (1983) Pozdnyakov, L. A., Sobol, I. M., & Syunyaev, R. A. 1983, Astrophys. Space Phys. Res., 2, 189
  * Reis & Miller (2013) Reis, R. C. & Miller, J. M. 2013, ApJ, 769, L7
  * Schnittman & Krolik (2009) Schnittman, J. D. & Krolik, J. H. 2009, ApJ, 701, 1175
  * Schnittman & Krolik (2010) Schnittman, J. D. & Krolik, J. H. 2010, ApJ, 712, 908
  * Schnittman & Krolik (2013) Schnittman, J. D., & Krolik, J. H. 2013, ApJ, 777, 11
  * Seon et al. (2022) Seon, K.-. il ., Song, H., & Chang, S.-J. 2022, ApJS, 259, 3
  * Shakura & Sunyaev (1973) Shakura, N. I. & Sunyaev, R. A. 1973, A&A, 500, 33
  * Sunyaev & Revnivtsev (2000) Sunyaev, R. & Revnivtsev, M. 2000, A&A, 358, 617
  * Tamborra et al. (2018) Tamborra, F., Matt, G., Bianchi, S., et al. 2018, A&A, 619, A105
  * Thorne & Price (1975) Thorne, K. S. & Price, R. H. 1975, ApJ, 195, L101
  * Ursini et al. (2022) Ursini, F., Matt, G., Bianchi, S., et al. 2022, MNRAS, 510, 3674
  * Ursini et al. (2020) Ursini, F., Dovčiak, M., Zhang, W., et al. 2020, A&A, 644, A132
  * Weisskopf et al. (2016) Weisskopf, M. C., Ramsey, B., O’Dell, S., et al. 2016, Proc. SPIE, 9905, 990517
  * Whitney (2011) Whitney, B. A. 2011, Bulletin of the Astronomical Society of India, 39, 101
  * Wilkins & Fabian (2012) Wilkins, D. R. & Fabian, A. C. 2012, MNRAS, 424, 1284
  * Xiao (2006) Xiao, F. 2006, Plasma Physics and Controlled Fusion, 48, 203
  * Yajima et al. (2012) Yajima, H., Li, Y., Zhu, Q., et al. 2012, MNRAS, 424, 884
  * Yang et al. (2021) Yang, X.-L., Wang, J.-C., Yang, C.-Y., et al. 2021, ApJS, 254, 29
  * You et al. (2021) You, B., Tuo, Y., Li, C., et al. 2021, Nature Communications, 12, 1025
  * Yuan & Narayan (2014) Yuan, F. & Narayan, R. 2014, ARA&A, 52, 529
  * Yusef-Zadeh et al. (1984) Yusef-Zadeh, F., Morris, M., & White, R. L. 1984, ApJ, 278, 186
  * Zhang et al. (2016) Zhang, S. N., Feroci, M., Santangelo, A., et al. 2016, Proc. SPIE, 9905, 99051Q
  * Zhang et al. (2019) Zhang, S., Santangelo, A., Feroci, M., et al. 2019, Science China Physics, Mechanics, and Astronomy, 62, 29502
  * Zdziarski et al. (1996) Zdziarski, A. A., Johnson, W. N., & Magdziarz, P. 1996, MNRAS, 283, 193
  * Życki et al. (1999) Życki, P. T., Done, C., & Smith, D. A. 1999, MNRAS, 309, 561