# Hotspots and Photon Rings in Schwarzschild Black Hole Spacetimes

Prashant Kocherlakota,1,2,3 Luciano Rezzolla,3,4,5 Rittick Roy,6 and Maciek
Wielgus7
1Black Hole Initiative at Harvard University, 20 Garden St., Cambridge, MA
02138, USA
2Center for Astrophysics, Harvard & Smithsonian, 60 Garden St., Cambridge, MA
02138, USA
3Institut für Theoretische Physik, Goethe-Universität, Max-von-Laue-Str. 1,
60438 Frankfurt, Germany
4School of Mathematics, Trinity College, Dublin 2, Ireland
5Frankfurt Institute for Advanced Studies, Ruth-Moufang-Str. 1, 60438
Frankfurt, Germany
6Anton Pannekoek Institute for Astronomy, University of Amsterdam, Science
Park 904, 1098 XH, Amsterdam, The Netherlands
7Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, D-53121 Bonn,
Germany

(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Future black hole (BH) imaging observations with better sensitivity and
angular resolution are expected to resolve finer features corresponding to
higher-order images of both hotspots that are produced in the accretion flow,
as well as of the entire emitting region. In spherically symmetric spacetimes,
the image order is determined by the maximum number of half-loops executed
around the BH by the photons that form it. Due to the additional half-loop,
consecutive-order images arrive after a delay time of approximately $\pi$
times the BH shadow radius. Furthermore, the deviation of the diameters from
that of the shadow, widths, and flux-densities of consecutive-order images are
exponentially demagnified by the lensing Lyapunov exponent, a characteristic
of the spacetime. We compare the exact time delay between the appearance of
the zeroth and first-order images of a hotspot to our best analytic estimate
and find an error $\lesssim 50\%$ for hotspot locations within $\approx 10M$
from a Schwarzschild BH of mass $M$. We also explore the possibilities for
inferring kinetic properties of hotspots, and also our inclination, from
future BH movies. Furthermore, since the targets of such observations host
geometrically-thick accretion flows (also jets), we obtain simple theoretical
estimates of the variation in the diameters and widths of their first-order
images, for varying disk scale-heights. We find realistically that the
deviation of the former from that of the shadow is $\lesssim 30\%$ and that
the latter remain $\lesssim 1.3M$. Finally, we estimate the error in
recovering the lensing exponent, when using the first and second-order images,
to be $\lesssim 20\%$. This provides further evidence that future observations
could yield new and independent estimates of the shadow size and, in
principle, of the lensing exponent, allowing us to robustly learn about the
spacetimes of astrophysical BHs.

###### keywords:

black-hole physics –- accretion, accretion disks – relativistic processes –
methods: analytical

††pubyear: 2021††pagerange: Hotspots and Photon Rings in Schwarzschild Black
Hole Spacetimes–D

##  1 Introduction

The recent images of M87⋆ (EHT Collaboration et al., 2019a) and Sgr A⋆ (EHT
Collaboration et al., 2022a) obtained by the Event Horizon Telescope (EHT)
Collaboration show that these supermassive objects can be accurately described
as Kerr black holes (BHs) (EHT Collaboration et al., 2019c, 2022d), the
vacuum, spinning BHs of general relativity (GR).

The Kerr BH spacetime contains a photon shell $\mathrm{S}$, which is the union
of a set of spherical surfaces, each of which admits bound photon orbits (Teo,
2003). Of these, there are also two planar (circular) orbits at the equator
(Bardeen, 1973) that are fundamentally tied to the symmetries of the spacetime
(cf. Cunha & Herdeiro 2018). The (exterior; cf. Grenzebach 2015) photon shell
of a Kerr BH encloses its horizon and demarcates the fate of photons
approaching it into those that (a) fall into its interior and eventually end
up inside the BH, (b) remain on one of the spheres in $\mathrm{S}$, or (c)
escape to faraway observers, depending on their four-momenta or, equivalently,
their impact parameters.

Photon orbits of type (b) are the bound orbits, and of particular interest
since they determine the BH shadow boundary curve or the $n\\!=\\!\infty$
critical curve $\mathscr{C}_{\infty}$ on the image plane (Synge, 1966;
Bardeen, 1973). This is because $\mathscr{C}_{\infty}$ is the gravitationally-
lensed projection of $\mathrm{S}$ on the observer’s image plane. Since photons
orbits are described by null geodesics in GR, the Kerr shadow boundary curve
$\mathscr{C}_{\infty}$ is determined by the Kerr metric as well as the
observer’s viewing angle.

When a BH is lit up by a source of emission such as hot inflowing gas, the
central intensity depression that appears in the image can typically be
related to the BH shadow (region interior to $\mathscr{C}_{\infty}$): Photon
orbits that cross the horizon [type (a)] have shorter paths through the
spacetime than those that do not [type (c)] and thus pick up lesser emission
from the hot plasma present outside the BH, leading to much smaller
intensities on the image plane in the pixels that they arrive in (Jaroszynski
& Kurpiewski, 1997; Narayan et al., 2019; Bronzwaer et al., 2021; Bronzwaer &
Falcke, 2021; Bauer et al., 2022; Özel et al., 2022; Younsi et al., 2023;
Kocherlakota & Rezzolla, 2022; EHT Collaboration et al., 2022d). This does not
mean that the intensity maximum in the image coincides exactly with the shadow
boundary but that there is a strong association between the two.

Indeed the diameter of the emission ring in the EHT images depends not just on
the spacetime geometry of M87⋆ or of Sgr A⋆ but also on the arrangement,
absorption, and emissivity of the accreting material in their vicinity (the
“gastrophysics”). Nevertheless, the set of sizes of the observed bright
emission ring and that of the shadow boundary can be empirically related
through synthetic images obtained from simulations, thereby allowing for an
inference of the shadow diameter of M87⋆ (EHT Collaboration et al., 2019c;
Psaltis et al., 2020; Kocherlakota et al., 2021) and of Sgr A⋆ (EHT
Collaboration et al., 2022d) from the measured emission ring diameter. An
excellent analysis of this point is reported in Sec. 3. of EHT Collaboration
et al. (2022d).

This discussion illustrates how the photon shell and the shadow boundary curve
play a fundamental role in determining the properties of BH images.
Furthermore, it also elucidates how several confounding effects must be
carefully accounted for before extracting information regarding the shadow
boundary curve, such as its median diameter.

One of the exciting prospects for future black hole imaging observations
performed at higher angular resolutions and flux-sensitivities, with next-
generation microarcsecond resolution instruments achieved through space-based
very long baseline interferometry (Gralla et al., 2020; Gurvits et al., 2022;
Kurczynski et al., 2022), is the direct detection of the “photon ring” of a BH
(Gralla et al., 2019; Johnson et al., 2020), which forms the focus of our
study here.

Photons that appear in a region close to the shadow boundary
$\delta\mathscr{C}_{\infty}$ can have orbits that access the region close to
the photon shell $\delta\mathrm{S}$. Since photon orbits close on themselves
in $\mathrm{S}$, i.e., they have divergent angular deflections, it is natural
to expect that $\delta\mathrm{S}$ is also a region of strong gravitational
lensing (Luminet, 1979; Ohanian, 1987; Virbhadra & Ellis, 2000; Bozza &
Scarpetta, 2007). More specifically (Sec. 2.5), photons that appear in
$\delta\mathscr{C}_{\infty}$ at an image plane radius of $\eta$, and which
access $\delta\mathrm{S}$, have orbits whose deflection angles diverge
logarithmically in the limit of approach to the shadow boundary radius, i.e.,
$\lim_{\eta\rightarrow\eta_{\mathrm{PS}}}\not{\Delta}\vartheta(\eta)\propto\ln{|\eta-\eta_{\mathrm{PS}}|}$
(Luminet, 1979; Ohanian, 1987; Bozza & Scarpetta, 2007; Gralla et al., 2019;
Gralla & Lupsasca, 2020; Johnson et al., 2020). The slope of this divergence
is given by the lensing Lyapunov exponent $\gamma_{\mathrm{PS}}$, which is
determined purely by the spacetime metric.

The region $\delta\mathscr{C}_{\infty}$ on the image plane is referred to as
the photon ring111Not to be confused with the light ring, which is a bound,
planar null geodesic in the bulk of space (Cunha & Herdeiro, 2018). and
exhibits a rich substructure (barring when the emitting region is perfectly
spherical; Vincent et al. 2022). For an extended source of emission in the
bulk, such as an accretion disk, the photon ring is comprised of a series of
discrete – and often overlapping – “subrings” that are indexed by the number
of half-loops executed by the photons that arrive in them (see Sec. 2.3 for a
precise definition). Furthermore, these subrings are organized self-similarly
on the image plane, with the critical exponent $\gamma_{\mathrm{PS}}$
governing this self-similarity. More specifically, the median (over the image
plane polar angle) radii of consecutive order images,
$\langle\eta_{n+1}\rangle$ and $\langle\eta_{n}\rangle$, exhibit a scaling
relation in their deviation from the shadow radius as
$\langle\eta_{n+1}-\eta_{\mathrm{PS}}\rangle\approx\mathrm{e}^{-\gamma_{\mathrm{PS}}}\langle\eta_{n}-\eta_{\mathrm{PS}}\rangle$.
Similarly, their widths, $w_{n+1}$ and $w_{n}$, are related as
$w_{n+1}\approx\mathrm{e}^{-\gamma_{\mathrm{PS}}}w_{n}$, thus leading to a
diminishing subring flux density with increasing subring-order (Gralla et al.,
2019; Gralla & Lupsasca, 2020; Johnson et al., 2020). These approximate
scaling relations apply to the case of equatorial emitters viewed by an
observer on the spin-axis.

Therefore, the lensing Lyapunov exponent is, in principle, a new observable
that can be used to gain additional information about the spacetime geometry
(see, e.g., Wielgus 2021; Broderick et al. 2023; Staelens et al. 2023; Salehi
et al. 2023). The implications of a measurement of this exponent for the black
hole “no hair” conjecture have been explored in Kocherlakota et al. (2023), a
companion to this paper. This exponent is tied to a fundamental “phase space”
Lyapunov exponent (Appendix A), which governs the expansion of the meridional
null congruence at the photon sphere.

For the EHT targets, M87⋆ and Sgr A⋆, their photon subrings are the higher-
order images of the horizon-scale accretion flow. Naturally, their structures
are determined both by the properties of the hot plasma as well as by their
respective spacetime geometries. Thus, in future improved experiments, the
close vicinity of the photon shell in the bulk as well as of the critical
curve on the boundary of space will soon play a key role in determining BH
images, and even BH movies. With these experiments on the horizon, it is
imperative to gain a deeper understanding of how the characteristic features
of subrings – such as their diameters, widths, flux densities, asymmetries –
vary with the properties of the emitting region.

Significant progress has already been made in understanding the appearance of
photon subrings cast by geometrically-thin disks in the Kerr metric (see,
e.g., Gralla et al. 2020; Palumbo et al. 2022; Paugnat et al. 2022; Tsupko
2022; Bisnovatyi-Kogan & Tsupko 2022) as well as in non-Kerr spacetimes
(Wielgus, 2021; Ayzenberg, 2022; Cárdenas-Avendaño & Held, 2023), affording us
several valuable insights.

Astrophysical BHs that are presently relevant such as M87⋆ or Sgr A⋆, however,
do not host geometrically-thin accretion disks (EHT Collaboration et al.,
2019b, 2022c). Instead, the accretion flow extends all the way down to the
event horizon, and is geometrically-thick, with typical thicknesses of
$|h/r|\lesssim 0.4$ (Narayan et al., 2022). Here $h$ denotes the scale-height
of the flow and $r$ is the radial coordinate. Equivalently, this implies that
the “faces” of the volume harboring the accretion flow lie at colatitudes
$\vartheta=\pi/2\pm 0.2$.

Furthermore, the inflowing plasma is magnetized and drags magnetic field lines
to the BH and supports them there (Narayan et al., 2003). This can lead to
relativistic outflows, or jets, of highly-magnetised plasma along the spin
axis via the Blandford & Znajek (1977) mechanism, which describes an
electromagnetic Penrose process (Lasota et al., 2014). These jets have
generalized parabolic profiles (Narayan et al., 2022) and their boundaries
(“jet sheaths”) can produce a significant amount of emission (Sironi et al.,
2021) as well. This is also confirmed by observations (cf., e.g., Kim et al.
2018; Cruz-Osorio et al. 2022; Janssen et al. 2021; Lu et al. 2023). This
constitutes yet another nontrivial source of nonequatorial emission. In
particular, while a nonspinning Schwarzschild BH does not produce powerful
jets (due to the absence of an ergoregion), there is still a significant
nonequatorial outflow component (“disk wind”).

Motivated by the above considerations, one of our primary goals here is to
present a comprehensive survey of the impact of a varying source morphology on
the observed photon ring structure (cf. also Vincent et al. 2022). Another is
to estimate the error in inferring the critical lensing exponent
$\gamma_{\mathrm{PS}}$ from potential future joint measurements of diameters
or widths of a pair of subrings for varied source morphologies.

We achieve this by working with a simple yet flexible three-parameter model
for the shape of the emitting region, which can be used to interpolate
smoothly between a geometrically-thin-disk and a sphere. More specifically, we
employ conical surfaces, pictured in Fig. 1, to generate an axially-symmetric
wedge-shaped region – a “conical torus” – with the three parameters being used
to set the locations of its inner $r_{\mathrm{in}}$ and outer
$r_{\mathrm{out}}$ surfaces as well as its geometrical-thickness or half-
opening angle $\vartheta_{1/2}$ (this defines the latitudes of the bounding
surfaces). We note that this setup also allows one to easily model a jet-like
feature (cf. Papoutsis et al. 2023; Chang et al., in prep.). We will also
carefully consider the impact of the observer viewing angle, which may be
useful in preparation for observations of the photon ring of Sgr A⋆. Finally,
to cleanly isolate the impact of source morphology, we fix the spacetime
geometry here. For simplicity, in this demonstrative study, we choose it to
correspond to that of a Schwarzshild BH, i.e., of a vacuum nonspinning BH in
GR.

The restrictions above ensure that our model complexity is minimal but
sufficient for our present purposes, and, crucially, enable intuitive
connections between the morphological variations of the emitting region in the
bulk and that of the photon subrings on the image plane. For direct practical
applications to observations, one must additionally account for realistic
“gastrophysical” effects associated with the accreting plasma. Excellent
investigations considering the effects of the plasma synchrotron emissivity
profile, its velocity profile, and the increased optical depth at
submillimeter frequencies experienced by higher-order photons traversing
greater path-lengths through the plasma are well underway (Palumbo et al.,
2022; Paugnat et al., 2022; Vincent et al., 2022).

It is worth noting that our analysis corresponds essentially to using a
radially constant emissivity profile and ignoring all other effects. While it
is understood that it should fall off at least as $\sim r^{-4}$ (Narayan et
al., 2019; Kocherlakota & Rezzolla, 2022) from considerations of energy
conservation and that it likely falls off much faster (exponentially) in
realistic scenarios (Chael et al., 2021; Porth et al., 2019), we find this
simplification useful for the following reason. Given two emitting regions $A$
and $B$ of identical morphologies (and in particular of the same radial
extent), it is straightforward to see that the diameters and widths of the
photon subrings cast by $A$ are larger than those by $B$, if the emissivity
falls off more rapidly with radius in $B$ as compared to $A$ (see, e.g., Fig.
7 of Kocherlakota & Rezzolla 2022 for an analogous conclusion for their direct
images). Thus, the subring diameters and widths that we report here could be
seen as an upper bound for a particular morphology of the emitting region. The
gravitational and the Doppler redshifts will also introduce additional
asymmetry in the image, but it is hard to imagine how their median image sizes
or widths can beat this upper bound. It is in this way that we find the
simplifications adopted above to be reasonable for our present purposes.

Thus far, we have focussed on the higher-order images of the entire accretion
flow. However, compact flux eruption or flaring events associated with Sgr A⋆,
from within ten gravitational radii, have been observed across a multitude of
wavelengths (cf. Baganoff et al. 2001; Eckart et al. 2006; GRAVITY
Collaboration et al. 2018; Murchikova & Witzel 2021; Witzel et al. 2021).
Magnetic reconnection in an equatorial current sheet within the accreting
plasma can lead to local heating (Ripperda et al., 2022). Such hot transient
compact emission sources – or “hotspots” – can also move within the equatorial
plane (cf. Dexter et al. 2020). Propositions that the time delay between the
appearance of the primary ($n=0$) and secondary ($n=1$) images of a hotspot
can potentially be measured through observations have been forwarded (cf.
Tiede et al. 2020; Ball et al. 2021). This may be possible with future
upgrades and ground extensions to the EHT (Doeleman et al., 2019; Johnson et
al., 2023), and has been suggested as a means to infer the spin of Sgr A⋆
(Wong, 2021). In similar vein, Sahu et al. (2013) have proposed utilizing the
time delay between the appearance of Einstein rings to gather information
about the background spacetime geometry.

Below (see eq. 2.44), we find the characteristic delay time between the
appearance of higher-order ($n\geq 1$) images to be intimately linked to the
shadow size, in general spherically-symmetric spacetimes. Our last goal here
is to estimate the error in inferring the shadow size from a measurement of
the exact time delay between the appearance of the primary and secondary
images of a hotspot in a Schwarzschild BH spacetime. We will also discuss how
these images move on the image plane as the hotspot orbits a BH, and how we
can infer its kinetic properties (their locations, motion, etc.), which may
provide insights into their physical origin.

The broad outline of the paper is as follows. We introduce the mathematical
description of higher-order images in general spherically-symmetric spacetimes
in Sec. 2. We discuss the implications of an observation of a secondary
($n=1$) image of a hotspot orbiting a Schwarzschild BH (relevant for Sgr A⋆)
in Sec. 3. We discuss the qualitative as well as quantitative variations in
the morphological properties of photon subrings cast by geometrically-thin
disks for varying viewing angles, and by geometrically-thick disks for a face-
on observer (relevant for M87⋆) in a Schwarzschild spacetime in Sec. 4. We end
by presenting our conclusions in Sec. 5.

An expert reader already familiar with the theory of higher-order images may
still be interested in reading Sec. 2.4 and Sec. 2.6.1. To ensure readability,
we relegate overly technical comments and discussions to footnotes and
appendices respectively.

Finally, a note on convention: We will reserve $\nu$ to denote a frequency
throughout; It will never be used as a label for a tensor index. For example,
$I_{\nu}$ will be used to denote specific intensity
($\mathrm{W}~{}\mathrm{m}^{-2}\mathrm{sr}^{-1}\mathrm{Hz}^{-1}$) and not the
component of a one-form $I$. We employ geometrized units throughout, in which
$G=c=1$.

##  2 Higher-Order Images & the Universal Photon Ring Scaling Relations

The line element of an arbitrary static and spherically-symmetric spacetime
can be expressed in spherical-polar coordinates
$x^{\alpha}=(t,r,\vartheta,\varphi)$ as,

$\mathrm{d}s^{2}=\mathscr{g}_{\mu\nu}\mathrm{d}x^{\alpha}\mathrm{d}x^{\beta}=-f~{}\mathrm{d}t^{2}+\frac{g}{f}~{}\mathrm{d}r^{2}+R^{2}~{}\mathrm{d}\Omega_{2}^{2}\,,$
(2.1)

where the metric functions $f,g$, and $R$ are functions of $r$ alone, and
$\mathrm{d}\Omega_{2}^{2}=\mathrm{d}\vartheta^{2}+\sin^{2}{\vartheta}~{}\mathrm{d}\varphi^{2}$
is the standard line element on a unit 2-sphere. We will assume reasonably
that $g>0$ everywhere and that $R>0$ except at the center ($R=0$). Setting
$R(r)=r$ yields the metric in areal-polar coordinates. The metric above
describes a BH spacetime if $f(r)$ admits real, positive zeroes (with $R>0$),
the largest of which locates the event horizon, which we denote by
$r_{\mathrm{H}}$. For the Schwarzschild BH metric in particular, we have
$f(r)=1-2M/r$, $g(r)=1$, and $R(r)=r$. Clearly, its event horizon is located
at $r_{\mathrm{H}}=2M$.

The Lagrangian $\mathscr{L}$ that describes the motion along a null geodesic
$x^{\alpha}(\lambda)$ is defined as
$2\mathscr{L}:=g_{\alpha\beta}\dot{x}^{\alpha}\dot{x}^{\beta}$, where
$\dot{x}^{\alpha}=\mathrm{d}x^{\alpha}/\mathrm{d}\lambda=:k^{\alpha}$ is the
tangent and $\lambda$ an affine parameter along it. Of the associated momenta,
$p_{\alpha}=\partial_{\dot{x}^{\alpha}}\mathscr{L}=k_{\alpha}$, there are two
that are conserved (due to the Euler-Lagrange equations,
$\dot{p}_{\alpha}=\partial_{\alpha}\mathscr{L}$), namely the energy $E=-p_{t}$
and the azimuthal angular momentum $L=p_{\varphi}$ of the orbit. In addition
to the three conserved quantities $\\{2\mathscr{L}(=0),E,L\\}$, a fourth
constant of the motion exists, namely the (non-negative) Carter (1968)
constant $C=p_{\vartheta}^{2}+p_{\varphi}^{2}\csc^{2}{\vartheta}$. The Carter
constant defined in this way is simply the square of the total angular
momentum, and is associated with the existence of an irreducible symmetric
Killing tensor (cf., e.g., Kapec & Lupsasca 2020), which plays a role similar
to that of the familiar Laplace-Runge-Lenz vector in Newtonian gravity (Brill
& Goel, 1999; Gibbons & Warnick, 2007; Cariglia et al., 2014; Cariglia, 2014).
With four dynamical constants, we can separate the geodesic equation and write
the tangent to an arbitrary null geodesic as being (Kocherlakota & Rezzolla,
2022),222The tangents to arbitrary timelike ($2\mathscr{L}=-1$) and spacelike
($2\mathscr{L}=+1$) geodesics differ from eq. 2 only in the radial component,
which is given then instead as
$\dot{r}=\pm_{r}Eg^{-1/2}\sqrt{1+2\mathscr{L}E^{-2}f-\eta^{2}fR^{-2}}$.

$\displaystyle k^{\alpha}=$ $\displaystyle\
E\left(f^{-1},\pm_{r}g^{-1/2}\sqrt{1-\eta^{2}fR^{-2}},\right.$ $\displaystyle\
\left.\pm_{\vartheta}R^{-2}\sqrt{\eta^{2}-\xi^{2}\csc^{2}{\vartheta}},\pm_{\varphi}R^{-2}\xi\csc^{2}{\vartheta}\right)\,,$
(2.2)

where the indices $\pm_{r},\pm_{\vartheta}$, and $\pm_{\varphi}$ correspond to
the signs of the radial, polar, and azimuthal velocities respectively, and we
have introduced the useful (energy-rescaled, non-negative) quantities
$\xi=|L|/E$ and $\eta=\sqrt{C}/E$. Clearly, we must have $\eta\geq\xi$ along
any orbit so that $k^{\vartheta}$ remains real.

###  2.1 Planarity of Orbits and Meridional Photon Orbits

Due to spherical symmetry, all geodesics are spatially-planar. To see this,
let us first introduce

$\mathscr{i}=\pi/2-\arcsin{(\xi/\eta)}\,.$ (2.3)

The relevant Euler equation then is (for $0<\mathscr{i}<\pi/2$)

$\frac{\mathrm{d}\varphi}{\mathrm{d}\vartheta}=\frac{k^{\varphi}}{k^{\vartheta}}=\pm\frac{\cot{\mathscr{i}}\csc^{2}{\vartheta}}{\sqrt{1-\cot^{2}{\mathscr{i}}\cot^{2}{\vartheta}}}\,,$
(2.4)

the general solution to which is

$\varphi(\vartheta)=\mp\arcsin\left(\cot{\mathscr{i}}\cdot\cot{\vartheta}\right)+\varphi_{0}\,,$
(2.5)

where $\varphi_{0}$ is an integration constant. We can put this into a more
suggestive form,

$\displaystyle 0=$ $\displaystyle\
\mp(\sin{\varphi_{0}}\cdot\sin{\mathscr{i}})\sin{\vartheta}\cos{\varphi}$
(2.6) $\displaystyle\
\pm(\cos{\varphi_{0}}\cdot\sin{\mathscr{i}})\sin{\vartheta}\sin{\varphi}+(\cos{\mathscr{i}})\cos{\vartheta}\,,$

to immediately recognize that this is the equation of a plane passing through
the center ($r=0$), and that $\mathscr{i}$ is the inclination of the normal to
the orbital plane (i.e., the angle between the normal and the $z-$axis,
$\vartheta=0$). Clearly, orbits that differ only in the sign of the polar
($\pm_{\vartheta}$) or the azimuthal ($\pm_{\varphi}$) velocity lie in the
same plane.

Introducing the effective polar potential $\Theta$ as (compare with eq. 2.8
below),

$\Theta(\eta,\xi,\theta):=(k^{\vartheta}/E)^{2}=\left(\eta^{2}-\xi^{2}\csc^{2}{\vartheta}\right)/R^{4}\,,$
(2.7)

we see that the orbit must stay in
$\vartheta_{-}\leq\vartheta(\lambda)\leq\vartheta_{+}$, where
$\vartheta=\vartheta_{\pm}$ are the polar turning points ($\Theta=0$). For
$0\leq\mathscr{i}\leq\pi/2$, $\vartheta_{\pm}=\pi/2\pm\mathscr{i}$. Our
definition of the orbital inclination (2.3) restricts it to
$0\leq\mathscr{i}\leq\pi/2$ since orbits with larger inclinations, i.e., in
$(\pi/2,\pi]$, can simply be obtained by shifting $\varphi_{0}$ by $\pi$.

This helps to see that orbits with $\mathscr{i}=0$ (or $\eta=\xi\neq 0$) are
indeed equatorial orbits. Those with $\mathscr{i}=\pi/2$ (or $\eta\neq\xi=0$)
correspond to meridional orbits (i.e., constant-$\varphi$ orbits),
corresponding to motion along longitudes. Finally, orbits with $\eta=0$
necessarily require $\xi=0$ as well, i.e.,
$\dot{\vartheta}=\ddot{\vartheta}=\dot{\varphi}=\ddot{\varphi}=0$, which
correspond to radial geodesics. Null geodesics, in particular, of this kind
are members of the ingoing ($-_{r}$) and outgoing ($+_{r}$) principal null
congruences (PNCs) of these type-D static spacetimes.

Due to spherical symmetry, (1.) all geodesic orbits are simply copies of those
that lie, e.g., in the equatorial plane or a meridional plane, and (2.)
without loss of generality, we can set the observer to lie on the north pole
($\vartheta_{\mathrm{o}}=0$). The latter has the profoundly simplifying
consequence that all photons arriving on the screen move only on meridional
planes through the bulk of space. Thus, henceforth, we will restrict our
attention to meridional photon orbits. Our analytical methods here derive the
crux of their power from these simplifications, and are what attract us to a
complete consideration of our problem of interest (i.e., characterizing
higher-order images of hotspots and the entire accretion flow) first in
spherically-symmetric (nonspinning) spacetimes.

###  2.2 Circular Photon Orbits and Radial Turning Points

For any photon orbit, we can write an energy equation for the radial motion,
$\dot{r}^{2}-E^{2}\mathscr{R}=0$, with $\mathscr{R}$ playing the role of an
effective radial potential, and being given as (cf. 2),

$\mathscr{R}(\eta,r):=(k^{r}/E)^{2}=g^{-1}\left[1-\eta^{2}fR^{-2}\right]\,.$
(2.8)

Due to the strong gravity near ultracompact objects, it is generically
possible for photons to move on circular orbits (cf. Cunha & Herdeiro 2020;
Ghosh & Sarkar 2021), i.e., $\dot{r}=\ddot{r}=0$. These correspond to the
bound photon orbits of the spacetime, and photons on such orbits do not reach
faraway observers. These equations are respectively equivalent to
$\mathscr{R}=0$, i.e.,

$\eta_{\mathrm{PS}}:=R_{\mathrm{PS}}/\sqrt{f_{\mathrm{PS}}}\,,$ (2.9)

and $\partial_{r}\mathscr{R}=0$, which yields

$(\partial_{r}f)/f-(\partial_{r}R^{2})/R^{2}=0\,.$ (2.10)

A solution, $r=r_{\mathrm{PS}}$, to the above differential equation locates a
photon sphere in the bulk of space. In the above, and henceforth, the
subscript “$\mathrm{PS}$” indicates that the function is evaluated at
$r=r_{\mathrm{PS}}$. The radial stability of orbits of photons present nearby
is determined by the sign of
$\partial_{r}^{2}\mathscr{R}(\eta_{\mathrm{PS}},r_{\mathrm{PS}})$, i.e.,
unstable if positive and stable if negative.

The phase space of null geodesics is discussed below in Appendix A. The “phase
space Lyapunov exponent,” $\hat{\kappa}_{\mathrm{PS}}$, introduced there as

$\hat{\kappa}_{\mathrm{PS}}^{2}:=\frac{\partial_{r}^{2}\mathscr{R}(\eta_{\mathrm{PS}},r_{\mathrm{PS}})}{2}=-\frac{1}{2g_{\mathrm{PS}}}\left[\frac{\partial_{r}^{2}f_{\mathrm{PS}}}{f_{\mathrm{PS}}}-\frac{\partial_{r}^{2}R^{2}_{\mathrm{PS}}}{R^{2}_{\mathrm{PS}}}\right]\,,$
(2.11)

then, equivalently, determines the stability of the radial fixed point as
being unstable (stable) if $\hat{\kappa}_{\mathrm{PS}}$ is real (imaginary).

While there exist static spacetimes admitting multiple roots to eq. 2.10 (see,
e.g., Wielgus et al. 2020; Gan et al. 2021; Guo et al. 2023), we restrict our
focus here to those that admit a single unstable photon sphere outside the
horizon (cf. Cardoso et al. 2014; Keir 2016; Cunha et al. 2017).

The total angular momentum per unit energy of a meridional photon $\eta$ is
also its (apparent) impact parameter (Bardeen 1973; see also Ch. 4 of
Grenzebach 2015). That is, it is the radius on the image plane at which it
appears.

In particular, the radius of the shadow boundary curve, which is the
gravitationally-lensed projection of the (unstable) photon sphere on the image
plane, is simply $\eta=\eta_{\mathrm{PS}}$. It is worth noting that the
location of the photon sphere as well as the size of the shadow are
independent of the metric function $g$ (see Psaltis et al. 2020 for a
discussion on the implications for gravity tests). For the Schwarzschild BH
metric in particular, its photon sphere is located at $r_{\mathrm{PS}}=3M$,
its shadow size is $\eta_{\mathrm{PS}}=\sqrt{27}M$, and its phase space
Lyapunov exponent takes value $\hat{\kappa}_{\mathrm{PS}}=1/(\sqrt{3}M)$. For
a loose comparison, we note that the surface gravity $\kappa_{\mathrm{H}}$ of
the Schwarzschild horizon is $\kappa_{\mathrm{H}}=1/(4M)$.

From the effective potential (2.8), we can see that photons that initially
approach the BH ($\dot{r}<0$) will either reach the horizon or hit a “radial
turning point” at some $r=r_{\mathrm{tp}}(\eta)$, where their radial velocity
vanishes ($\dot{r}=0$), i.e.,

$\mathscr{R}(\eta,r_{\mathrm{tp}}(\eta))=0\,.$ (2.12)

The latter set of photons will then move outward and reach asymptotic
infinity. These are photons with large angular momenta $\eta$ (we set $E=1$
without loss of generality) that simply experience too much “centrifugal
force” at their closest point of approach. It is straightforward to see from
eqs. 2.8 and 2.9 that these are all photons with angular momenta
$\eta>\eta_{\mathrm{PS}}$. It is also clear that every point along a circular
photon orbit is a radial turning point,
$r_{\mathrm{tp}}(\eta_{\mathrm{PS}})=r_{\mathrm{PS}}$. For further discussion
on radial turning points, we direct the reader to see Appendix A of
Kocherlakota & Rezzolla (2022).

###  2.3 Angular Deflection for Higher-Order Images

The effective polar potential for meridional orbits ($\xi=0$) is

$\Theta(\eta,r):=(k^{\vartheta}/E)^{2}=\eta^{2}/R^{4}\,.$ (2.13)

It is clear then that meridional orbits do not admit any nontrivial polar
turning points (where $\dot{\vartheta}=0=\Theta$). However, the spherical-
polar coordinate system throws up two trivial turning points at $\vartheta=0$
and $\vartheta=\pi$, which are understood as follows. A photon approaching the
north pole has $\dot{\vartheta}<0$. After smoothly crossing it, it must have
$\dot{\vartheta}>0$, since the sense of rotation of an orbit about the center
($r=0$) remains invariant. We can ignore these trivial turning points but keep
track of the sign of the initial angular velocity of the photon, which
determines the sense of rotation or the “polarity” of the orbit.

Thus, the total angular deflection $\not{\Delta}\vartheta^{\pm}$ experienced
by a meridional photon with initial positive ($+$) or negative ($-$) polar
velocity ($\dot{\vartheta}$) is given as,

$\displaystyle\not{\Delta}\vartheta^{\pm}=\fint\mathrm{d}\vartheta=$
$\displaystyle\ \fint\dot{\vartheta}~{}\mathrm{d}\lambda=\pm
E\fint\sqrt{\Theta}~{}\mathrm{d}\lambda$ (2.14) $\displaystyle=$
$\displaystyle\
\fint\left(\dot{\vartheta}/\dot{r}\right)\mathrm{d}r=\pm\fint\sqrt{\Theta/\mathscr{R}}~{}\mathrm{d}r\,.$

In the above, the slash in $\not{\Delta}\vartheta^{\pm}$ is used to denote
that it is akin to an angular distance along the photon orbit, as opposed to
an angular displacement. The slash in $\fint$ denotes that the integrals are
path-dependent, and are evaluated along the worldline of the photon,
$x^{\mu}=x^{\mu}(\lambda)$ (see also Shaikh & Joshi 2019; Gralla & Lupsasca
2020). This convention makes it convenient to introduce the winding number
along individual photon orbits, as in eq. 2.18.

For photons traversing the photon sphere,
$\sqrt{\Theta}=\eta_{\mathrm{PS}}/R^{2}(r_{\mathrm{PS}})$, and it is clear
from eq. 2.14 that their angular deflection diverges. We can also anticipate
that it becomes arbitrarily large for some photons
($|\mathscr{R}(\eta,r)|\rightarrow 0$) that closely approach it.

We can “unslash” the integrals in eq. 2.14, and make them path-independent by
splitting the photon orbit up into pieces over which the map $\lambda\mapsto
r$ is bijective. The $r-$bijective pieces are naturally separated by the
radial turning points (2.12). Since we are only interested in photons that
make it to a faraway observer present at $r=\infty$, as discussed below, we
need only worry about a single radial turning point at most. Therefore, for
all observable photons, we can write,

$\displaystyle\not{\Delta}\vartheta^{\pm}(\eta,r_{\mathrm{e}})=\begin{cases}\pm\eta\left[+_{r}\int_{r_{\mathrm{e}}}^{\infty}1/(R^{2}\sqrt{\mathscr{R}})~{}\mathrm{d}r\right],&k^{r}_{\mathrm{e}}>0\,;\
r_{\mathrm{e}}>\\{r_{\mathrm{H}}~{}\mathrm{if}~{}\eta<\eta_{\mathrm{PS}},~{}r_{\mathrm{tp}}(\eta)~{}\mathrm{if}~{}\eta\geq\eta_{\mathrm{PS}}\\}\\\
\pm\eta\left[-_{r}\int_{r_{\mathrm{e}}}^{r_{\mathrm{tp}}(\eta)}1/(R^{2}\sqrt{\mathscr{R}})~{}\mathrm{d}r+_{r}\int_{r_{\mathrm{tp}}(\eta)}^{\infty}1/(R^{2}\sqrt{\mathscr{R}})~{}\mathrm{d}r\right],&k^{r}_{\mathrm{e}}\leq
0\,;\mathrm{if}~{}\eta>\eta_{\mathrm{PS}}~{}\mathrm{and}~{}\
r_{\mathrm{e}}\geq r_{\mathrm{tp}}(\eta)\end{cases}$ (2.15)

In the above, $r_{\mathrm{e}}$ and $k^{r}_{\mathrm{e}}$ denote the emission
location and radial velocity of the photon.

This equation can be understood simply as follows. The top row corresponds to
photons with positive radial velocity at emission ($k^{r}_{\mathrm{e}}>0$).
Since only photons with angular momenta $\eta\geq\eta_{\mathrm{ps}}$ possess
radial turning points, these must be emitted from outside their outermost
radial turning point to ensure they reach an asymptotic observer. Photons with
smaller angular momenta must simply be emitted from outside the horizon. The
bottom row corresponds to the case where photons are emitted with initial
nonpositive radial velocities. Photons emitted with zero radial velocity
($k^{r}_{\mathrm{e}}=0$) are, by definition, emitted from their orbital radial
turning points ($r=r_{\mathrm{tp}}(\eta)$). From the next instant onwards,
these photons must move away (“retreat”) from the BH to reach an asymptotic
observer. Photons with initially negative radial velocities
($k^{r}_{\mathrm{e}}<0$) can reach asymptotic observers only if they meet a
radial turning point. The two pieces of the integral, thus, simply correspond
to the angular deflection accumulated over the approaching ($-_{r}$) and
retreating ($+_{r}$) legs of the photon orbit. For an intuitive pictorial
representation of the above, we direct the reader to see the left panel of
Fig. A1 of Kocherlakota & Rezzolla (2022). Finally, the overall $\pm$ simply
makes explicit how photons that differ only in polarity experience identical
amounts of angular deflections. For equatorial orbits, one can simply replace
$\vartheta$ and $\eta$ by $\varphi$ and $\xi$ respectively in the discussion
above.

The polarity of the photon orbit determines its image plane polar angle
$\psi$. Two photons emitted from the same meridian (same $\varphi$) but with
opposite polarities, appear on opposite sides of the image plane origin.
Henceforth, we will use $\psi^{\mp}(\varphi)$ to denote keep track of this
feature (cf., e.g., Bardeen 1973; Luminet 1979; Grenzebach 2015; Tsupko 2022),

$\psi^{\mp}(\varphi)=\begin{cases}3\pi/2-\varphi\,,&-\ (\mathrm{i.e.,}\
k^{\vartheta}_{\mathrm{e}}<0)\\\ \pi/2-\varphi\,,&+\ (\mathrm{i.e.,}\
k^{\vartheta}_{\mathrm{e}}>0)\end{cases}\,.$ (2.16)

The image plane origin corresponds to the location on the image plane from
which the member of the ingoing PNC at the observer’s location enters their
past light cone.

Putting everything together, photons emitted from an initial location
$(r_{\mathrm{e}},0\leq\vartheta_{\mathrm{e}}\leq\pi,0\leq\varphi_{\mathrm{e}}<2\pi)$
that reach our observer, present on the $+z-$axis
($\vartheta_{\mathrm{o}}=0$), appear on their image plane at $(\eta,\psi)$.
The image plane polar angle $\psi$ is determined by eq. 2.16, which requires
knowledge also of the initial sign ($\pm$) of the photon polar velocity,
$k^{\vartheta}_{\mathrm{e}}=\dot{\vartheta}_{\mathrm{e}}$. Furthermore, the
image plane radius $\eta$ is determined via eq. 2.15 where

$\not{\Delta}\vartheta^{\pm}=\left[\vartheta_{\mathrm{o}}-\vartheta_{\mathrm{e}}\right]\mod{2\pi}=-\vartheta_{\mathrm{e}}\mod{2\pi}\,.$
(2.17)

Eq. 2.15 requires further knowledge of the initial sign of the photon radial
velocity $k^{r}_{\mathrm{e}}$. Quite wonderfully, as we shall see below in
Sec. 2.4, the source location
$(r_{\mathrm{e}},\vartheta_{\mathrm{e}},\varphi_{\mathrm{e}})$ as well as the
image order $n$ (introduced below) completely fix the photon orbit. Thus, we
can find the initial coordinate velocities (both $k^{\vartheta}_{\mathrm{e}}$
and $k^{r}_{\mathrm{e}}$) as well as the image location $(\eta,\psi)$ of that
photon uniquely.

The modulo piece in the equation above (2.17) is central to the notion of
higher-order images. There exist photons emitted from the same initial spatial
location and captured by an observer at the same final spatial location that
have orbits differing total angular deflection

$\not{\Delta}\vartheta^{\pm}_{m}=-\vartheta_{\mathrm{e}}\pm 2\pi m\,,\ \
0<\vartheta_{\mathrm{e}}<\pi\,,$ (2.18)

where $m$ is a positive integer, called the winding number (cf. also Luminet
1979; Ohanian 1987; Cunha & Herdeiro 2020; Wei 2020; Ye & Wei 2023), and each
distinct orbit is indexed by the pair $(m,\pm)$.

The total angular deflection for this class of orbits increases in the
sequence $\\{(0,-),(1,+),(1,-),(2,+),\cdots\\}$ taking values in
$\\{(0,\pi),(\pi,2\pi),(2\pi,3\pi),(3\pi,4\pi),\cdots\\}$ respectively. If we
index these sets by $n=\\{0,1,2,3,\cdots\\}$, then $n$ defines the order of
the image and we can rewrite eq. 2.18 as,

$\not{\Delta}\vartheta^{\pm}_{n}=\pi/2-\vartheta_{\mathrm{e}}+(-1)^{n+1}(2n+1)\pi/2\,,\
\ 0<\vartheta_{\mathrm{e}}<\pi\,.$ (2.19)

Thus, we find, naturally, that even$-n$ photons are associated with a negative
polarity whereas odd-$n$ photons are associated with a positive polarity.

As discussed above eq. 2.14, these polarities were introduced to track the
sign of $\dot{\vartheta}$ at emission. We can understand the reason for this
association as follows. Any “direct” or “primary” ($n=0$) photon that is
emitted from $0<\vartheta_{\mathrm{e}}<\pi$ and which reaches the observer
present at $\vartheta_{\mathrm{o}}=0$ must have had a negative polar angular
velocity, $\dot{\vartheta}<0$, along its entire orbit (thus,
$\dot{\vartheta}<0$ at emission). On the other hand, for the “indirect” or
“secondary” ($n=1$) photon, which is emitted from
$0<\vartheta_{\mathrm{e}}<\pi$ and which reaches the observer present at the
north pole, to experience an angular deflection in $(\pi,2\pi)$, it must first
head toward the south pole, cross it, and then move toward the north pole.
Thus, it must have had a positive polar velocity at emission, and so a
positive polarity. Similar reasoning applies to all even-order and odd-order
photons, and can be readily verified by Fig. 1.

Photons emitted from the same (static) source appear on an observer’s screen
on a single line. This is because photon orbits are planar in spherically-
symmetric spacetimes, and the source location, the observer location, and the
center of space define this plane. Furthermore, due to their alternating
polarities, consecutive-order images appear on opposite sides of the origin
(cf. eq. 2.16). This shows that in spherically-symmetric spacetimes the
rotation parameter $\delta_{0}$, introduced in eq. 72 of Gralla & Lupsasca
(2020), is trivial, i.e., $\delta_{0}=\pi$. By extension, we can define the
$n^{\mathrm{th}}-$order image of an extended source of emission as being
formed by the set of photons that have index $n$.

Since different-order photons experience different amounts of angular
deflections, they must appear at different radii on the image plane,
$\eta=\eta_{n}$. The latter are obtained simply by solving eq. 2.15 with eq.
2.19, i.e., the integral equation,

$\not{\Delta}\vartheta^{\pm}(\eta_{n},r_{\mathrm{e}})=\not{\Delta}\vartheta_{n}^{\pm}(\vartheta_{\mathrm{e}})\,.$
(2.20)

Furthermore, the total affine length $\not{\Delta}\lambda^{\pm}$ of a photon
orbit and the total time $\not{\Delta}t^{\pm}$ elapsed along it, which are
given respectively as

$\displaystyle\not{\Delta}\lambda^{\pm}=$ $\displaystyle\
\fint\mathrm{d}\lambda=\fint\left(1/\dot{r}\right)\mathrm{d}r=\fint(1/E)/\sqrt{\mathscr{R}}~{}\mathrm{d}r\,,$
(2.21) $\displaystyle\not{\Delta}t^{\pm}=$ $\displaystyle\
\fint\mathrm{d}t=\fint\left(\dot{t}/\dot{r}\right)\mathrm{d}r=\fint(1/f)/\sqrt{\mathscr{R}}~{}\mathrm{d}r\,,$
(2.22)

must also differ for different-order photons. We note that these (slashed)
path-integrals can be made (unslashed) path-independent integrals in precisely
the same way as in eq. 2.15. We can then find that photons forming
increasingly higher-order images have concomitantly longer paths and take
longer before getting to the observer, when emitted from the same event. Thus,
photon orbits of different orders connect the same initial and final spatial
locations but not the same two events in spacetime.

Finally, for the case when $\vartheta_{\mathrm{e}}=0$, pairs of photon orbits
of orders $2n+1$ and $2n+2$ have identical orbits under reflections across the
axis of the observer. In fact, due to the associated symmetries of this
configuration, a point source present at such a location (called a caustic)
does not form just two images but an entire ring, called an Einstein ring or a
critical curve (cf. Ref. Virbhadra & Ellis 2000). All of the photons that form
this ring on the image plane arrive at precisely the same time when emitted
from the same event due to identical orbits. Caustics are defined as locations
in the past light cone of the observer that have divergent magnifications and
critical curves are the maps of these points on the image plane (Virbhadra &
Ellis, 2000; Frittelli et al., 2000). There is another (half-)line of caustics
at $\vartheta_{\mathrm{e}}=\pi$ and pairs of photon orbits of orders $2n$ and
$2n+1$ are identical. Motivated by Fig. 2, we introduce the convention of
indexing the critical curves by $n$, where the total gravitational lensing
experienced by the photons that form it is given as
$\not{\Delta}\vartheta=n\pi$. Then odd-order critical curves are formed by the
caustics with $\vartheta_{\mathrm{e}}=\pi$ whereas the even-order critical
curves are formed by the other caustics at $\vartheta_{\mathrm{e}}=0$. With
this, it becomes clear that the $n=\infty$ critical curve (Gralla & Lupsasca,
2020) is simply the shadow boundary curve.

Figure 1: Photon orbits connecting point sources to an observer in a
Schwarzschild black hole spacetime. In the top-left panel, we show orbits of
photons emitted from just outside the horizon, $r_{\mathrm{e}}=2M^{+}$, which
reach an asymptotic observer, present on the $+z-$axis. The top edge of each
panel represents a radial slice of the observer’s screen. The different colors
indicate that they are emitted from different colatitudes
$\vartheta_{\mathrm{e}}=0^{\circ},15^{\circ},30^{\circ},\cdots,180^{\circ}$.
Across the top row, in each panel, lines of the same color represent photon
orbits that begin from the same initial spatial positions. Thus, a countably
infinite number of photon orbits connect any specific pair of emitter and
observer. The order ($n$) of photons with the same line color increases from
left to right across columns, as do their path lengths, travel times, and
total angular deflections. Across rows, the radial location of the emitters is
varied, $r_{\mathrm{e}}=2M^{+},2.9M,3.1M,6M$. The locations of the event
horizon, of the photon shell, and of the innermost stable circular orbit
(ISCO) are at $r_{\mathrm{H}}=2M$ (black circle), $r_{\mathrm{PS}}=3M$ (dashed
blue circle), and $r_{\mathrm{ISCO}}=6M$ (dotted green circle). Larger light
bending occurs for photons that appear closer to the shadow boundary on the
image plane $\eta_{\mathrm{PS}}=\sqrt{27}M$, shown as vertical red lines in
all panels. It is clear that these photons necessarily access the close
vicinity of the photon shell. Finally, we see that higher-order images (across
each row) occupy increasingly compact regions on the image plane. Figure 2:
Photon orbits and higher-order images of emitters in a Schwarzschild black
hole (BH) spacetime. The top-left panel presents the angular deflection for
photons emitted from different radii and with different angular momenta
($\eta$). The latter corresponds to their apparent impact parameters or screen
radii. Photons emitted towards and away from the BH are shown in dashed and
solid lines respectively. Photons emitted from the conical torus (see top-
right panel) occupy the darker-shaded regions. The dashed black line in the
$n=0$ region tracks the size of the “inner shadow” (cf. Chael et al. 2021). In
all panels, the vertical red lines show the size of the BH shadow. The top-
right panel shows our simple conical torus model. The bottom-left panel
displays a set of photon orbits reaching an observer on the $+z-$axis. The
event horizon and the photon sphere are shown as black and blue circles
respectively. The green and purple lines represent a cross-section of the
torus. The bottom-right figure shows the observer’s screen image of the torus,
with red and blue shading indicating regions collecting $n=0$ and $n=1$
photons, corresponding to the direct image and the first photon subring.

###  2.4 The Bulk $\leftrightarrow$ Boundary Lensing Maps

We will now consider carefully the properties of the lensing maps
$\eta_{n}(r_{\mathrm{e}},\vartheta_{\mathrm{e}})$, i.e., maps which relate
source locations in the bulk to order$-n$ image locations on the boundary.
Now, as noted above, these maps are obtained by solving the integral equation
(2.20) above. We will see below, when discussing the top left panel of Fig. 2,
how no other initial data (e.g., regarding the initial velocities) are
required. Note also that we can safely omit both $\varphi$ and $\psi$ from
this discussion since all photon orbits are planar in spherically symmetric
spacetimes, as a result of which we found the map $\varphi\mapsto\psi$ to be
trivial (2.16).

The discussion below revolves around two key competing observations:

  * I.

Photons that experience larger angular deflections possess larger angular
momenta, or, equivalently, larger impact parameters ($\eta$).

  * II.

The region on the image plane that collects photons that experience extremely
large angular deflections is a narrow ring near the shadow boundary
($|\eta-\eta_{\mathrm{ps}}|\ll 1$).

We will see below how the first observation applies mainly only for primary
image formation whereas the second plays the dominant role in shaping higher-
order image formation.

Figure 1 shows the first four order ($n=0-3$) images of emitters located at
different radii ($r_{\mathrm{e}}$) and colatitudes ($\vartheta_{\mathrm{e}}$)
in a Schwarzschild BH spacetime. It is immediately clear that images of
sources present inside the photon sphere in the bulk
($r_{\mathrm{e}}<r_{\mathrm{PS}}=3M$) always appear inside the shadow boundary
curve on the image plane ($\eta<\eta_{\mathrm{PS}}=\sqrt{27}M$). Moreover,
images of sources present outside the photon sphere can appear either inside
or outside the shadow boundary curve.333Most curiously, photons that are
emitted from $r_{\mathrm{e}}\gtrsim r_{\mathrm{PS}}$ can undergo large
deflections and their images (of arbitrary order) can still appear inside the
shadow boundary (cf. Sec. 4.1). Furthermore, it is also clear that
increasingly higher-order images (left to right panels) occupy increasingly
narrower regions on the image plane.

Photons that are emitted from the same emission radius and which are bent
through larger amounts typically appear closer to the shadow boundary on the
image plane (observation II). This is easily seen to be true for photons of
all orders when emitted from inside the photon sphere,
$r_{\mathrm{e}}<r_{\mathrm{PS}}$. However, photons that are emitted from the
same emission radius, $r_{\mathrm{e}}>r_{\mathrm{PS}}$, and which are bent
through larger amounts do not always appear closer to the shadow boundary on
the image plane.

To understand this, let us consider, e.g., the $n=0$ image of emitters located
at $r_{\mathrm{e}}=6M$ (bottom-left panel). First, although the source of
emission is always present outside the photon sphere in the bulk, its image
can appear inside the shadow boundary on the image plane, for small
$\vartheta_{\mathrm{e}}$. With increasing source-colatitude, photons
experience increasingly larger angular deflections and appear closer to the
shadow boundary, approaching it from the inside
($\eta\rightarrow\eta_{\mathrm{PS}}^{-}$). However, this trend reverses
briefly as photons begin appearing outside the shadow boundary, and their
images move away from the shadow boundary despite experiencing increasingly
greater angular deflections (observation I.). Beyond some large source-
colatitude, however, this trend reverses yet again, and the images appear
increasingly closer to the shadow boundary, now from the outside
($\eta\rightarrow\eta_{\mathrm{PS}}^{+}$; observation II.). Finally, notice
that these trends are universal, i.e., they apply for emission sourced from
any radius outside the photon sphere $r_{\mathrm{e}}>r_{\mathrm{PS}}$.

Thus, the primary lensing map
$\eta_{0}(r_{\mathrm{e}},\vartheta_{\mathrm{e}})$ is nonmonotonic
(nonbijective) w.r.t. variations in $\vartheta_{\mathrm{e}}$ but only for
photons that are emitted from outside the photon sphere $r>r_{\mathrm{PS}}$.
This is primarily because such radial locations can source photon orbits that
permit radial turning points.

Furthermore, we emphasize that this is a generic feature of all order ($n$)
lensing maps. However, the band of emission radii outside the photon sphere
for which this map is nonmonotonic shrinks exponentially with increasing image
order and therefore, as alluded to above, becomes less important.444 For an
image of order-$n$ this band is bounded from above
($r_{\mathrm{PS}}<r_{\mathrm{e}}\leq r_{n;\mathrm{C^{T}}}$) by the radius
$r_{n;\mathrm{C^{T}}}$ at which a photon emitted with zero radial velocity
($k^{r}_{\mathrm{e}}=0$) experiences an angular deflection of exactly $n\pi$,
i.e.,
$\not{\Delta}\vartheta(\eta_{n;\mathrm{C^{T}}},r_{n;\mathrm{C^{T}}})=n\pi$,
where $\eta_{n;\mathrm{C^{T}}}$ is a solution of the radial turning point
equation, $k^{r}(\eta,r_{n;\mathrm{C^{T}}}(\eta))=0$ (see the right panel of
Fig. 12).

This illustrates two important aspects of gravitational lensing. First,
photons emitted from the same radial location $r_{\mathrm{e}}$ in the bulk can
appear at different radii $\eta$ on the image plane, depending on their
angular locations $\vartheta_{\mathrm{e}}$. Second, and more importantly, two
photons that were emitted from the same radial location $r_{\mathrm{e}}$ in
the bulk and which appear at the same radii $\eta$ on the image plane can have
different trajectories, experiencing different amounts of angular deflections
as well as travel times (nonmonotonicity).

The top left panel of Fig. 2 (cf. also Gyulchev et al. 2021) neatly shows the
image radius on the observer’s sky ($\eta$; on the $x-$axis) of a source
located at a particular radius ($r_{\mathrm{e}}$; fixed along each line) and
at a particular colatitude ($\vartheta_{\mathrm{e}};y-$axis on the right
edge). The angular deflection experienced by these photons, of different
orders ($n$) is also shown
($\not{\Delta}\vartheta=|\not{\Delta}\vartheta^{\pm}|;y-$axis on the left
edge). Thus, horizontal lines in this plot correspond to emission coming from
conical surfaces, except for $\not{\Delta}\vartheta=(2n+1)\pi/2$ which
correspond to emission coming from the equatorial plane. The relation between
the source colatitude, $\vartheta_{\mathrm{e}}$, and the angular deflection
experienced by the photon, $\not{\Delta}\vartheta$, is straightforwardly given
by eq. 2.19.

Photons that were emitted in the radially-outward ($k^{r}_{\mathrm{e}}>0$) and
-inward ($k^{r}_{\mathrm{e}}<0$) directions are represented here in dotted and
dashed lines respectively. The meeting point of these two line types naturally
corresponds to the case when the photon is emitted with zero radial velocity
$k^{r}_{\mathrm{e}}=0$. That is, the emission radius in the bulk matches the
radial turning point of that photon orbit,
$r_{\mathrm{e}}=r_{\mathrm{tp}}(\eta)$. This photon appears at a radius
$\eta=\eta_{\mathrm{tp}}(r_{\mathrm{e}}):=R(r_{\mathrm{e}})/\sqrt{f(r_{\mathrm{e}})}$
on the image plane.

It is immediately clear from this panel that a combination of the source
location $(r_{\mathrm{e}},\vartheta_{\mathrm{e}})$ and the order ($n$) of the
image pick out a unique photon orbit, without needing any other initial
conditions.

It is also clear that photons emitted from any particular radius
$r=r_{\mathrm{e}}>r_{\mathrm{PS}}$ in the bulk all appear inside the radius
$\eta=\eta_{\mathrm{tp}}(r_{\mathrm{e}})$ on the boundary. By fixing
$\vartheta_{\mathrm{e}}$ (emission coming from a conical surface) and varying
$r_{\mathrm{e}}$, we can also see that photons emitted from larger radii in
the bulk $r_{\mathrm{e}}$ appear at larger radii on the image plane $\eta$.

Putting together findings from both figures, we find that while the lensing
map is nonmonotonic when fixing the source radius and varying the source
colatitude, it remains monotonic when fixing the source colatitude and varying
the source radius.

Since the peak of the emissivity profile in the bulk likely lies close to the
photon shell (cf. EHT Collaboration et al. 2019c, 2022b, 2022d), and since we
can only access a complicated superposition of all order images simultaneously
in the image domain, it is useful to keep in mind these subtle features of
gravitational lensing, to better infer properties of the accretion flow that
are made from observed patterns in future dynamical movies of M87⋆ and Sgr A⋆
(Conroy et al., 2023) and to set the potential scope for spacetime tomography
techniques (Tiede et al., 2020; Levis et al., 2022).

This also underscores the importance of considering nonequatorial emission
sources (i.e., sources located at $\vartheta_{\mathrm{e}}\neq\pi/2$) since the
total image is a combination of all-order images of the accretion disk as well
as of the forward and receding jets. Finally, it is also important to consider
the impact of non-face-on observer inclinations. This forms the motivation for
our analysis in Sec. 4.

To understand the impact of varying emitting region morphologies on inferences
of photon ring properties, we will consider below (in Sec. 4) the impact of
varying the morphological parameters of the emitting region (from that of a
geometrically-thin disk) on its image. Towards this end, and as noted above,
we model the emitting region as a conical torus. In the top-right panel of
Fig. 2 we show the bounding surfaces of our torus model.

In the bottom-left panel of Fig. 2, we show the spatial orbits of null
geodesics in a meridional plane in a Schwarzschild BH spacetime. If we take
this to be the $yz-$plane ($\varphi=\pi/2,3\pi/2$), then these photons appear
on the image plane Cartesian “$\alpha-$”axis (Bardeen 1973; See also eq.
2.16). Photons having large impact parameters ($\eta\gg\eta_{\mathrm{ps}}$)
naturally arrive at large radii on the image plane. Clearly, these do not
approach the BH closely and deviate only slightly from a straight line,
$\not{\Delta}\vartheta\approx\pi$. These are all $n=0$ photons. Similarly,
photons with small impact parameters ($\eta\approx 0$) also undergo nearly no
deflection. On the other hand, photons that appear close to the shadow
boundary $\delta\mathscr{C}_{\infty}~{}(\eta\approx\eta_{\mathrm{PS}})$ can be
strongly-lensed by the BH, $\not{\Delta}\vartheta\gg\pi$. This region
$\delta\mathscr{C}_{\infty}$ on the image plane ($|\bar{\eta}|\ll
1;\bar{\eta}:=\eta/\eta_{\mathrm{PS}}-1$) is called the “photon ring” (Gralla
& Lupsasca, 2020; Johnson et al., 2020). Photons that are strongly-lensed
necessarily access the close vicinity of the photon shell $\delta
S~{}(r_{\mathrm{e}}\approx r_{\mathrm{PS}})$ somewhere along their orbit (cf.
also Fig. 1).

In the bottom-right panel of Fig. 2, we show the image of the solid conical
torus as seen by an observer lying along the $+z-$axis, viewing the emitting
region face-on. We have also shown the region on the image plane occupied by
the $n=0$ or direct image as well as the $n=1$ or first-order image, which
collect photons that undergo deflections between $0$ and $\pi$ and between
$\pi$ and $2\pi$ respectively. It is clear to see how a higher-order image is
a demagnified (or thinner) version of a lower-order image.

###  2.5 Universal Photon Ring Scaling Relations

Now that we have qualitatively introduced the notion of a photon ring,
following Gralla & Lupsasca (2020) and Johnson et al. (2020), here we will see
how it can be identified quantitatively as a region on the image plane where
the total deflection angle increases logarithmically as
$\not{\Delta}\vartheta\propto\ln{|\eta-\eta_{\mathrm{PS}}|}\propto\ln{|\bar{\eta}|}$.
An excellent analysis for the analytic approximation of the angular deflection
of photons on orbits that approach the vicinity of the photon sphere in
arbitrary static and spherically-symmetric spacetimes was reported in Bozza &
Scarpetta (2007). We will now extend their analysis to include all of the
aforementioned quantities (elapsed affine time, elapsed coordinate time, total
deflection angle).

Towards this end, we first introduce the path-independent integral (with
$\dot{Q}(\eta,r)$ a regular function),

$\displaystyle\Delta
Q(\eta,r_{1},r_{2}):=\int_{r_{1}}^{r_{2}}\dot{Q}(\eta,r)/\sqrt{\mathscr{R}(\eta,r)}~{}\mathrm{d}r\,,$
(2.23) $\displaystyle\quad\quad\mathrm{where}\quad r_{2}\geq
r_{1}\geq\\{r_{\mathrm{H}}^{+}~{}\mathrm{if}~{}\eta<\eta_{\mathrm{PS}},~{}r_{\mathrm{tp}}(\eta)~{}\mathrm{if}~{}\eta\geq\eta_{\mathrm{PS}}\\}\,.$

Since by definition $\mathscr{R}(\eta,r_{\mathrm{tp}}(\eta))=0$, all integrals
of the type $\Delta Q(\eta,r_{\mathrm{tp}}(\eta),r_{2})$ have divergent
integrands with poles at precisely the same location. Depending on the radial-
falloff of $\dot{Q}$, these may or may not be finite but the dominant
contribution to these assuredly comes from the region close to
$r_{\mathrm{tp}}(\eta)$. We demonstrate in Appendix B how these dominant
pieces can be characterized by elliptic functions in general, for arbitrary
static and spherically-symmetric spacetimes. The path-independent integral
expression above (2.23) enables the following compact notation for path-
dependent integrals along photon orbits,

$\displaystyle\not{\Delta}Q(\eta,r_{\mathrm{e}}):=$ $\displaystyle\
\fint_{\lambda_{1}}^{\lambda_{2}}\dot{Q}(\eta,r(\lambda))~{}\mathrm{d}\lambda=\fint_{r_{\mathrm{e}}}^{\infty}\dot{Q}/\sqrt{\mathscr{R}}~{}\mathrm{d}r$
(2.24) $\displaystyle=$ $\displaystyle\ \begin{cases}\Delta
Q(\eta,r_{\mathrm{e}},\infty)\,,&k^{r}_{\mathrm{e}}>0\,;\
r_{\mathrm{e}}>\\{r_{\mathrm{H}}~{}\mathrm{if}~{}\eta<\eta_{\mathrm{PS}},~{}r_{\mathrm{tp}}(\eta)~{}\mathrm{if}~{}\eta\geq\eta_{\mathrm{PS}}\\}\\\
\Delta Q(\eta,r_{\mathrm{tp}}(\eta),r_{\mathrm{e}})+\Delta
Q(\eta,r_{\mathrm{tp}}(\eta),\infty)\,,&k^{r}_{\mathrm{e}}\leq
0\,;\mathrm{if}~{}\eta>\eta_{\mathrm{PS}}~{}\mathrm{and}~{}\
r_{\mathrm{e}}\geq r_{\mathrm{tp}}(\eta)\end{cases}\,.$ (2.25)

Equations 2.23 and 2.24 yield the total angular deflection
$\not{\Delta}\vartheta^{\pm}$ (2.15) for
$\dot{Q}=\dot{\vartheta}/E=\eta/R^{2}$, the total path length
$\not{\Delta}\lambda^{\pm}$ (2.21) for $\dot{Q}=1/E$, and the total elapsed
time $\not{\Delta}t^{\pm}$ (2.22) for $\dot{Q}=\dot{t}/E=1/f$.

In terms of a pair of conformal variables, namely the fractional deviations of
the radial coordinates in the bulk and on the image plane from their critical
values respectively,

$\bar{r}:=r/r_{\mathrm{PS}}-1\,,\ \ \bar{\eta}:=\eta/\eta_{\mathrm{PS}}-1\,,$
(2.26)

it is clear that the limit $\bar{r}\rightarrow 0$ takes us to the photon
sphere in the bulk and the limit $\bar{\eta}\rightarrow 0$ sends us to the
shadow boundary on the image plane.

Figure 3: Logarithmic divergence of the deflection angle with impact parameter
locates the photon ring on the image plane. The angular deflection
($\not{\Delta}\vartheta$) experienced by a photon depends on its emission
location ($r_{\mathrm{e}}$), the sign of its initial radial velocity
($k^{r}_{\mathrm{e}};\pm$), and its apparent impact parameter ($\eta$). The
latter is also the radius on the image plane at which it appears. A photon
that is emitted at the photon sphere in a Schwarzschild BH spacetime
($r_{\mathrm{e}}=r_{\mathrm{PS}}=3M$) with the critical impact parameter
($\eta=\eta_{\mathrm{PS}}=\sqrt{27}M$) undergoes infinite deflection due to
gravitational lensing. One that appears in the close vicinity of the shadow
boundary ($|\bar{\eta}|\ll 1$) experiences large but finite deflections and
necessarily accesses the close vicinity of the photon sphere ($|\bar{r}|\ll
1$) along its orbit. The exact angular deflection experienced by a photon
(2.15) and its approximation (2.5) are shown here in solid and dashed lines
respectively. The logarithmic divergence occurs for small $|\bar{\eta}|$
values. The left panel shows photons appearing inside the shadow boundary,
$\bar{\eta}<0$, i.e., in the “inner photon ring,” whereas the right panel
shows photons appearing in the “outer photon ring.” The lensing Lyapunov
exponent, $\gamma_{\mathrm{PS}}$, here takes value $\gamma_{\mathrm{PS}}=\pi$.

Photon orbits that terminate on the image plane in the close vicinity of the
critical curve ($|\bar{\eta}|\ll 1$) and which access the close vicinity of
the photon sphere somewhere along the orbit ($|\bar{r}(\lambda)|\ll 1$ for
some $\lambda$) experience extreme gravitational-lensing. More specifically,
these are photons that were either emitted from555We introduce the somewhat
technical definitions for these types of orbits in Appendix B. The reader can
verify from Fig. 1 and the top left panel of Fig. 2 that the properties of the
orbits listed here are the ones that yield the appropriate divergence
($\not{\Delta}Q\propto\ln{|\bar{\eta}|}$).

  * [A

] well inside the photon sphere ($r_{\mathrm{H}}<r_{\mathrm{e}}\lnsim
r_{\mathrm{PS}}$), and in the radially-outward direction
($k^{r}_{\mathrm{e}}>0;\bar{\eta}<0$),

  * [C

] close to the photon sphere ($r_{\mathrm{e}}\simeq r_{\mathrm{PS}}$) or a
radial turning point ($r_{\mathrm{e}}\simeq r_{\mathrm{tp}}$), with initially
nonnegative radial velocity ($k^{r}_{\mathrm{e}}\geq 0$), or

  * [E

] well outside the photon sphere ($r_{\mathrm{e}}\gnsim
r_{\mathrm{tp}}(\eta)$), and in the radially-inward direction
($k^{r}_{\mathrm{e}}<0;\bar{\eta}>0$).

For such extremely-lensed photon orbits, the leading-order behavior in
$\bar{\eta}$ of arbitrary path-integrals $\not{\Delta}Q$ for small
$|\bar{\eta}|$ is given as (see Appendix B),

$\not{\Delta}Q(\bar{\eta})\approx\not{\Delta}Q_{\mathrm{D}}(\bar{\eta})=\left\\{\begin{aligned}
&-\frac{\dot{Q}(0,0)}{\hat{\kappa}_{\mathrm{PS}}}\left[\ln{|\bar{\eta}|}+\tilde{c}\right]\,,\
&&\mathrm{[types\ A,\ E]}\\\
&-\frac{\dot{Q}(0,0)}{2\hat{\kappa}_{\mathrm{PS}}}\left[\ln{|\bar{\eta}|}+\tilde{c}_{\mathrm{o}}\right]\,,\
&&\mathrm{[type\ C]}\end{aligned}\right.$ (2.27)

where
$\dot{Q}(0,0):=\dot{Q}(\bar{\eta}=0,\bar{r}=0)=\dot{Q}(\eta_{\mathrm{ps}},r_{\mathrm{ps}})$
and we have absorbed the dependence on the locations of the emitter
($r_{\mathrm{e}}$) and the observer ($r_{\mathrm{o}}$) locations into some
constants $\tilde{c}$ (defined in Appendix B). These can be ignored for our
present purposes since we are primarily interested in the logarithmic scaling
and the scaling constant.

The scaling constant for some observable $Q$ is determined by a piece that is
specific to the observable $\dot{Q}(0,0)$ and a piece that is independent of
the observable. This latter universal scaling constant
$\hat{\kappa}_{\mathrm{PS}}$ is the null geodesic phase space Lyapunov
exponent that we introduced above in Sec. 2.1, and is defined purely by the
spacetime geometry. While the location of the photon sphere (2.10) as well as
the size of the shadow (2.9) are determined independently of the metric
function $g$, $\hat{\kappa}_{\mathrm{PS}}$ depends on all three metric
functions ($f,g,R$) but only on the derivatives of $f$ and $R$. Jacobson
(2007) discusses the physical implications underlying the metric function $g$.

Equation 2.27 presents a powerful closed-form expression for the logarithmic-
scaling behavior of arbitrary observables in the photon ring. This expression
was made general by being able to pull $\dot{Q}$ entirely outside the integral
in eq. 2.27. The universal behavior of various important quantities can now be
obtained trivially. For the total angular deflection
$\not{\Delta}\vartheta^{\pm}$, the total coordinate time
$\not{\Delta}t^{\pm}$, and the total affine time $\not{\Delta}\lambda^{\pm}$,
we find666In writing the equation below (2.5) and henceforth, we restrict to
photon orbits that are $\mathrm{type\ A}$ and $\mathrm{E}$ respectively for
brevity, since the analysis to cover $\mathrm{type\ C}$ orbits follows
straightforwardly (e.g., we have to account for the factor of $1/2$). All
qualitative statements however apply to all three cases.

$\displaystyle\not{\Delta}\vartheta^{\pm}\approx$ $\displaystyle\
\mp\frac{\eta_{\mathrm{PS}}}{R^{2}_{\mathrm{PS}}\hat{\kappa}_{\mathrm{PS}}}\left[\ln{|\bar{\eta}|}+\tilde{c}_{\vartheta}\right]\
=:\mp\frac{\pi}{\gamma_{\mathrm{PS}}}\left[\ln{|\bar{\eta}|}+\tilde{c}_{\vartheta}\right]\,;$
$\displaystyle\not{\Delta}t^{\pm}\approx$ $\displaystyle\
-\frac{1}{f_{\mathrm{PS}}\hat{\kappa}_{\mathrm{PS}}}\left[\ln{|\bar{\eta}|}+\tilde{c}_{t}\right]\
=:-t_{\ell;\mathrm{PS}}\left[\ln{|\bar{\eta}|}+\tilde{c}_{t}\right]\,;$
$\displaystyle\not{\Delta}\lambda^{\pm}\approx$ $\displaystyle\
-\frac{1}{E\hat{\kappa}_{\mathrm{PS}}}\left[\ln{|\bar{\eta}|}+\tilde{c}_{\lambda}\right]\,.$
(2.28)

Here $\tilde{c}$ are some constants that depend on the radial locations of the
observer and emitter. We remind the reader that, we had introduced the sign
convention in eq. 2.14 to track the polarity of the photon orbit, or,
equivalently, the sign of its polar velocity $\dot{\vartheta}$ at emission.

In the above, we have introduced the lensing Lyapunov exponent,
$\gamma_{\mathrm{PS}}$, as well as the Lyapunov time, $t_{\ell;\mathrm{PS}}$,
which are the critical exponents that control the logarithmic divergences of
the deflection angle and the elapsed time respectively, as

$\displaystyle\gamma_{\mathrm{PS}}$ $\displaystyle:=\ \frac{\pi
R^{2}_{\mathrm{PS}}}{\eta_{\mathrm{PS}}}\hat{\kappa}_{\mathrm{PS}}=\pi\left[\frac{f_{\mathrm{PS}}R^{2}_{\mathrm{PS}}}{2g_{\mathrm{PS}}}\left(\frac{\partial_{r}^{2}R^{2}_{\mathrm{PS}}}{R^{2}_{\mathrm{PS}}}-\frac{\partial_{r}^{2}f_{\mathrm{PS}}}{f_{\mathrm{PS}}}\right)\right]^{1/2}\,.$
(2.29) $\displaystyle t_{\ell;\mathrm{PS}}$ $\displaystyle:=\
\frac{1}{f_{\mathrm{PS}}\hat{\kappa}_{\mathrm{PS}}}=\left[\frac{f^{2}_{\mathrm{PS}}}{2g_{\mathrm{PS}}}\left(\frac{\partial_{r}^{2}R^{2}_{\mathrm{PS}}}{R^{2}_{\mathrm{PS}}}-\frac{\partial_{r}^{2}f_{\mathrm{PS}}}{f_{\mathrm{PS}}}\right)\right]^{-1/2}\,.$
(2.30)

Again, these are determined purely by the spacetime geometry. This shows how
the circular null geodesic, $(\eta,r)=(\eta_{\mathrm{PS}},r_{\mathrm{PS}})$
plays a key role in determining the behavior of a number of quantities
associated with photon orbits that arrive in the photon ring
$\delta\mathscr{C}_{\infty}$ and which access the region close to the photon
sphere $\delta\mathrm{S}$. We will soon see how these determine observable
features such as the sizes and widths of photon subrings, as well as the time
delay between the appearance of different-order images.

We now report in Fig. 3 the deflection angles experienced by photons in a
Schwarzschild BH spacetime with different impact parameters that access the
close vicinity of the photon sphere in the bulk, when computed exactly via eq.
2.15 as well as when computed approximately via eq. 2.5. This figure locates
the onset of the logarithmic scaling of $\not{\Delta}\vartheta$ with
$|\bar{\eta}|$, and thus the photon ring on the image plane as being the
region $|\bar{\eta}|\lesssim 10^{-1}$ for $\bar{\eta}>0$ and
$|\bar{\eta}|\lesssim 10^{-2}$ for $\bar{\eta}<0$. We also note that the
constant offset between the exact and approximate curves is accounted for by
the regular piece $\Delta Q_{\mathrm{R}}~{}(\dot{Q}=\eta/R^{2})$, as discussed
in Appendix B.

###  2.6 The Lensing Lyapunov Exponent, $\gamma_{\mathrm{PS}}$

The angular deflections $\not{\Delta}\vartheta^{\pm}_{n}$ experienced by two
photons that form the $n^{\mathrm{th}}$ and $(n+2)^{\mathrm{th}}$ order images
of any point source on the image plane of any observer differ by $2\pi$, i.e.,
(cf. eq. 2.19)

$\not{\Delta}\vartheta^{\pm}_{n+2}-\not{\Delta}\vartheta^{\pm}_{n}=\pm
2\pi\,.$ (2.31)

We remind the reader that, in our convention, even-order photons appear with a
negative sign whereas odd-order photons appear with a positive sign (see the
discussion below eq. 2.19).

With the equation above, and the universal angular deflection scaling relation
(2.5), we can immediately write,

$\not{\Delta}\vartheta^{\pm}_{n+2}-\not{\Delta}\vartheta^{\pm}_{n}=\pm
2\pi\approx\mp\frac{\pi}{\gamma_{\mathrm{PS}}}\ln{\left[\frac{\bar{\eta}_{n+2}}{\bar{\eta}_{n}}\right]}\,.$
(2.32)

This is equivalent to the “image radius scaling relation,”

$\frac{\bar{\eta}_{n+2}}{\bar{\eta}_{n}}=\frac{\eta_{n+2}-\eta_{\mathrm{PS}}}{\eta_{n}-\eta_{\mathrm{PS}}}\approx\mathrm{e}^{-2\gamma_{\mathrm{PS}}}\,.$
(2.33)

Notice how these equations are entirely independent of the spatial locations
of the emitter and observer.

For an extended stationary source of emission, if we denote the outer and
inner boundaries of the order$-n$ image by $\eta_{n;\mathrm{out}}(\psi)$ and
$\eta_{n;\mathrm{in}}(\psi)$ respectively, then its width is given as
$w_{n}:=\eta_{n;\mathrm{out}}-\eta_{n;\mathrm{in}}$. We can then apply eq.
2.33 to the outer and inner boundaries of the $(n+2)^{\mathrm{th}}$ and
$n^{\mathrm{th}}$ order images,777 Except when the outer boundary of the
emission zone is extremely close to the photon sphere such that the lensing
map becomes nonmonotonic (see Sec. 2.4). In this case, there is no guarantee
that the outer boundary of all order images is sourced by emission from the
same set of spatial points. to obtain the following “image width scaling
relation,”

$\frac{w_{n+2}}{w_{n}}=\frac{\eta_{n+2;\mathrm{out}}-\eta_{n+2;\mathrm{in}}}{\eta_{n;\mathrm{out}}-\eta_{n;\mathrm{in}}}\approx\mathrm{e}^{-2\gamma_{\mathrm{PS}}}\,.$
(2.34)

In this way, we understand why $\gamma_{\mathrm{PS}}$ is called the lensing
Lyapunov exponent. For the Schwarzschild BH spacetime,
$\gamma_{\mathrm{PS}}=\pi$. The image scaling relations, as well as the
lensing Lyapunov exponent, have been obtained also for the Kerr metric
(Johnson et al., 2020; Gralla & Lupsasca, 2020). This exponent has also been
discussed in various non-Kerr spacetimes (Wielgus 2021; Broderick et al. 2023;
Kocherlakota et al. 2023; Salehi et al. 2023).

The flux densities of the $(n+2)^{\mathrm{th}}$ and the $n^{\mathrm{th}}$
subrings also obey a scaling behavior as (see Appendix C)

$\frac{F_{\nu;{n+2}}}{F_{\nu;n}}\approx\frac{w_{n+2}}{w_{n}}\approx\mathrm{e}^{-2\gamma_{\mathrm{PS}}}\,.$
(2.35)

We can also explain in full generality how the relationship between the
locations of appearance of consecutive-order images (opposite polarities) is
dependent only on the lensing Lyapunov exponent and the relative inclination
of the source and emitter, but not on their radial locations. From (cf. eqs.
2.19 and 2.5)

$\not{\Delta}\vartheta^{\pm}_{n+1}+\not{\Delta}\vartheta^{\mp}_{n}=-2\vartheta_{\mathrm{e}}+\pi[1\pm
1]\approx\mp\frac{\pi}{\gamma_{\mathrm{PS}}}\ln{\left[\frac{\bar{\eta}_{n+1}}{\bar{\eta}_{n}}\right]}\,,$
(2.36)

we find the general image scaling relation,

$\frac{\bar{\eta}_{n+1}}{\bar{\eta}_{n}}\approx\mathrm{e}^{-\gamma_{\mathrm{PS}}}\cdot\mathrm{e}^{\pm\gamma_{\mathrm{PS}}(2\vartheta_{\mathrm{e}}/\pi-1)}\,.$
(2.37)

This is the central equation that can be used to obtain scaling relations for
any higher-order image observable, for arbitrary spatial locations of emitter
and observer in arbitrary spherically-symmetric and static spacetimes. In the
equation above, as well as in eqs. 2.38 and 2.43 below, we should choose the
upper sign when $n$ is odd (and $n\\!\\!+\\!\\!1$ is even) and the lower sign
when $n$ is even (and $n\\!\\!+\\!\\!1$ is odd).

We can recover eq. 2.33 from eq. 2.37 by seeing that, e.g.,
$(\not{\Delta}\vartheta^{+}_{n+2}+\not{\Delta}\vartheta^{-}_{n+1})-(\not{\Delta}\vartheta^{-}_{n+1}+\not{\Delta}\vartheta^{+}_{n})=(-2\vartheta_{\mathrm{e}}+2\pi)-(-2\vartheta_{\mathrm{e}})$.

One can also easily obtain the general width and flux density scaling
relations for extended sources at fixed colatitude ($\vartheta_{\mathrm{e}}$),
such as a cone or a ring, from eq. 2.37 above, as

$\displaystyle\frac{F_{\nu;n+1}}{F_{\nu;n}}$ $\displaystyle\approx\
\frac{w_{n+1}}{w_{n}}\approx\mathrm{e}^{-\gamma_{\mathrm{PS}}}\cdot\mathrm{e}^{\pm\gamma_{\mathrm{PS}}(2\vartheta_{\mathrm{e}}/\pi-1)}\,.$
(2.38)

For the special case of an equatorially-located emitter
($\vartheta_{\mathrm{e}}=\pi/2$) viewed by our observer on the $z-$axis
(“face-on”) in particular, eq. 2.37 simplifies to (cf. also Johnson et al.
2020; Bisnovatyi-Kogan & Tsupko 2022),

$\frac{\bar{\eta}_{n+1}}{\bar{\eta}_{n}}\approx\mathrm{e}^{-\gamma_{\mathrm{PS}}}\,.$
(2.39)

Concomitantly, for equatorially-located extended sources of emission, we can
write,

$\displaystyle\frac{F_{\nu;n+1}}{F_{\nu;n}}$ $\displaystyle\approx\
\frac{w_{n+1}}{w_{n}}\approx\mathrm{e}^{-\gamma_{\mathrm{PS}}}\,.$ (2.40)

####  2.6.1 Image Scaling Relations for Nonequatorial Sources

We now highlight the direct astrophysical implications of our general image
scaling relation (2.37). Emission from an astrophysical jet sheath can be
modeled approximately by a conical surface
(constant$-\vartheta_{\mathrm{e}}$). Clearly, for such a conical source of
emission, the scaling relations that relate image radii ($\eta_{n}$), widths
($w_{n}$), and flux densities ($F_{\nu;n}$) for consecutive order images
depend on its half-opening angle ($\vartheta_{\mathrm{e}}$). Let us now
consider four jet models specifically. $\mathtt{Model-I}$: A jet of opening
angle $\pi/10$ that is pointed towards our line of sight (forward or
approaching jet), i.e., $\vartheta_{\mathrm{e}}/\pi=0.1$. $\mathtt{Model-II}$:
A jet of opening angle $\pi/10$ that is pointed away our line of sight
(backward or retreating jet), i.e., $\vartheta_{\mathrm{e}}/\pi=0.9$.
$\mathtt{Model-III}$: A forward jet of opening angle $\pi/4$, i.e.,
$\vartheta_{\mathrm{e}}/\pi=0.25$. $\mathtt{Model-IV}$: A backward jet of
opening angle $\pi/4$, i.e., $\vartheta_{\mathrm{e}}/\pi=0.75$. The relations
between the $n=1$ and $n=2$ image radii, widths, and flux densities for each
of these jet models are given as,

$\displaystyle\frac{\eta_{2}-\eta_{\mathrm{PS}}}{\eta_{1}-\eta_{\mathrm{PS}}}\approx\frac{w_{2}}{w_{1}}\approx\frac{F_{\nu;2}}{F_{\nu;1}}\approx\mathrm{e}^{-1.8\gamma_{\mathrm{PS}}}\
$ $\displaystyle\ [\mathtt{Jet~{}Model-I}]$ (2.41)
$\displaystyle\frac{\eta_{2}-\eta_{\mathrm{PS}}}{\eta_{1}-\eta_{\mathrm{PS}}}\approx\frac{w_{2}}{w_{1}}\approx\frac{F_{\nu;2}}{F_{\nu;1}}\approx\mathrm{e}^{-0.2\gamma_{\mathrm{PS}}}\
$ $\displaystyle\ [\mathtt{Jet~{}Model-II}]$
$\displaystyle\frac{\eta_{2}-\eta_{\mathrm{PS}}}{\eta_{1}-\eta_{\mathrm{PS}}}\approx\frac{w_{2}}{w_{1}}\approx\frac{F_{\nu;2}}{F_{\nu;1}}\approx\mathrm{e}^{-1.5\gamma_{\mathrm{PS}}}\
$ $\displaystyle\ [\mathtt{Jet~{}Model-III}]$
$\displaystyle\frac{\eta_{2}-\eta_{\mathrm{PS}}}{\eta_{1}-\eta_{\mathrm{PS}}}\approx\frac{w_{2}}{w_{1}}\approx\frac{F_{\nu;2}}{F_{\nu;1}}\approx\mathrm{e}^{-0.5\gamma_{\mathrm{PS}}}\
$ $\displaystyle\ [\mathtt{Jet~{}Model-IV}]$

Since the image of an astrophysical BH could contain an imprint of a strong
jet component, our inference of the lensing Lyapunov exponent could be
corrupted, if not properly taken into account. This is a simple demonstration
of the importance of considering the properties of nonequatorial sources of
emission.

###  2.7 The Delay Time, $t_{d;\mathrm{PS}}$

From the universal scaling relations obtained above (2.5), we can also obtain
the time delay between the arrival of the $n^{\mathrm{th}}$ and
$(n+2)^{\mathrm{th}}$ order images on the image plane. Using eq. 2.33, we
find,

$\not{\Delta}t_{n+2}^{\pm}-\not{\Delta}t_{n}^{\pm}\approx-
t_{\ell;\mathrm{PS}}\ln{\left[\frac{\bar{\eta}_{n+2}}{\bar{\eta}_{n}}\right]}=2t_{\ell;\mathrm{PS}}\gamma_{\mathrm{PS}}=2\pi\eta_{\mathrm{PS}}\,.$
(2.42)

Similarly, we can also obtain the time delay between the arrival of the
consecutive order images on the image plane using eqs. 2.5 and 2.37 as being,

$\not{\Delta}t_{n+1}^{\pm}-\not{\Delta}t_{n}^{\mp}\approx\
\pi\eta_{\mathrm{PS}}\left[1\mp\left(\frac{2\vartheta_{\mathrm{e}}}{\pi}-1\right)\right]\,.$
(2.43)

Using the two equations above, we now introduce the characteristic delay time,
$t_{d;\mathrm{PS}}$, which yields an approximate measure of the time elapsed
between the appearance of consecutive order images on the image plane, as
being,

$t_{d;\mathrm{PS}}:=\pi\eta_{\mathrm{PS}}\,.$ (2.44)

For a Schwarzschild BH, $t_{d;\mathrm{PS}}=\pi\sqrt{27}M$ (cf. also the
parameter $\tau_{0}$ in eq. 72 of Gralla & Lupsasca 2020). We can understand
this delay time simply as being the half-orbital time of a photon moving on a
circular meridional orbit,

$\frac{t_{orb;\mathrm{PS}}}{2}=\frac{\pi}{\Omega_{\mathrm{PS}}}=\pi\eta_{\mathrm{PS}}\,,$
(2.45)

where $\Omega_{\mathrm{PS}}$ is its angular velocity (see eq. 2),

$\Omega_{\mathrm{PS}}=\frac{k^{\vartheta}}{k^{t}}=\left.\left[\frac{\mathscr{g}_{tt}}{\mathscr{g}_{\vartheta\vartheta}}\right]\right|_{r_{\mathrm{PS}}}\eta_{\mathrm{PS}}=\frac{1}{\eta_{\mathrm{PS}}}\,.$
(2.46)

This is a wonderful finding: A clean detection of the time delay between
higher-order images can yield an independent estimate of the shadow size
$\eta_{\mathrm{PS}}$ of any spherically-symmetric ultracompact object (see
Sec. 3 below). Current measurements of the shadow size are inferred through a
careful calibration procedure, as mentioned above, and an independent estimate
of the same – perhaps obtained “more directly” from the data – can yield a
vital sanity check.

Figure 4: Time delay between the primary and secondary images of a hotspot.
The delay time between the appearance of the $n^{\mathrm{th}}$ and
$(n+1)^{\mathrm{th}}$ order images on the image plane, of a static source, is
approximately given by $t_{d;\mathrm{PS}}=\pi\eta_{\mathrm{PS}}$, where
$\eta_{\mathrm{PS}}$ is the BH shadow size. We compute the exact time delay
$\Delta t_{0}$ between the appearance of the primary ($n=0$) and secondary
($n=1$) images of a hotspot (here, a point source) located at a radius
$r=r_{\mathrm{e}}$ and angular coordinate $\vartheta=\vartheta_{\mathrm{e}}$
in a Schwarzschild black hole spacetime ($\eta_{\mathrm{PS}}=\sqrt{27}M$). The
panel on the left shows the fractional error between the exact time delay and
the characteristic one ($t_{d;\mathrm{PS}}$). In the right panel, we use our
analytic prediction from eq. 3.1, which corrects for the viewing angle
($\vartheta_{\mathrm{e}}$), and find lower errors relative to the actual
value. We identify a hotspot located at $r_{\mathrm{e}}=11M$ and viewed at an
inclination of $22^{\circ}$ by the star (Wielgus et al., 2022a). By similarly
accounting for the (roughly linear) dependence of the time delay on the
hotspot distance from the BH, we expect a measurement of such a time delay to
yield an independent estimate of the BH shadow size.

###  2.8 The Lyapunov Time, $t_{\ell;\mathrm{PS}}$

We introduced above the characteristic delay time associated with the
appearance of higher-order images. Here we identify the existence of yet
another distinct fundamental timescale associated with the photon sphere.

The radial distance between two photons present close to the photon sphere
grows exponentially in affine time as given by eq. A.7. We can rewrite this
growth in terms of the coordinate time by starting with
$\dot{\bar{r}}/\dot{t}=-\pm_{r}\mathscr{g}_{tt}(r_{\mathrm{PS}})\hat{\kappa}_{\mathrm{PS}}\bar{r}$,
i.e.,

$\displaystyle\bar{r}(t)$ $\displaystyle=\
\bar{r}(0)\exp{\left[-\pm_{r}\mathscr{g}_{tt}(r_{\mathrm{PS}})\hat{\kappa}_{\mathrm{PS}}t\right]}$
$\displaystyle=\ \bar{r}(0)\exp{[\pm_{r}t/t_{\ell;\mathrm{PS}}]}\,.$ (2.47)

Thus, the Lyapunov time $t_{\ell;\mathrm{PS}}$, given in eq. 2.30, measures
the characteristic instability timescale for photons present close to the
photon sphere. In the language of dynamical systems, it is the time, as
measured by an asymptotic ($r\rightarrow\infty$) static observer
($u\propto\partial_{t}$), for the radial coordinate between nearby photon
orbits to increase by a factor of $\mathrm{e}\approx 2.72$. For the
Schwarzschild BH spacetime, $t_{\ell;\mathrm{PS}}=\sqrt{27}M$.

Recently, in Cardoso et al. (2021) (cf. also Ames & Thorne 1968), it was shown
that this time scale plays an important role in determining the late-time
characteristics of the observed luminosity evolution (light curve) of an
infalling star. Thus, the Lyapunov time can, in principle, be measured by
future black hole imaging measurements.

##  3 Hotspots in Schwarzschild Black Hole Spacetimes

Recent observations of flaring events associated with Sgr A⋆ have been modeled
in Wielgus et al. (2022a) (cf. also Vos et al. 2022). It is found that a
hotspot moving on a Keplerian orbit of radius $r_{\mathrm{K}}\approx 11M$ in a
Schwarzschild BH spacetime, and viewed at an inclination of
$\mathscr{i}\approx 22^{\circ}$, yields a good fit for the data (see Sec. 3.3
there). Theoretical aspects of hotspot production have been carefully
investigated in Ripperda et al. (2022).

In this section, we will explore the implications of a possible detection of
the secondary ($n=1$) image of such a hotspot. While hotspots are expected to
have sizes of up to $1M$ (cf. Sec. 3.4 of Ripperda et al. 2020), for our
present purposes, it suffices to consider the hotspot to be a point source.
This toy model greatly simplifies the problem and we can infer several useful
insights. It can be argued that our results may be directly useful for
observational purposes if the point source is interpreted as tracking the
location of the hotspot electron number density maximum (cf. Tiede et al.
2020).

###  3.1 Time Delay between Primary and Secondary Images

We begin with an investigation of whether the time delay between the
appearance of the primary ($n=0$) and secondary ($n=1$) images of a hotspot,
$\Delta t_{0}$, can also be used to infer the BH shadow size. As discussed
above in Sec. 2.7, the characteristic delay time $t_{d;\mathrm{PS}}$ is
related to the shadow size $\eta_{\mathrm{PS}}$ as
$t_{d;\mathrm{PS}}=\pi\eta_{\mathrm{PS}}$. For a Schwarzschild BH, we obtain
$t_{d;\mathrm{PS}}=\pi\sqrt{27}M$. For Sgr A⋆, this corresponds to
$t_{d;\mathrm{PS}}\approx\pi\sqrt{27}(GM_{\mathrm{SgrA}^{\star}}/c^{3})\approx
5~{}\mathrm{min}$, whereas for M87⋆, we find $t_{d;\mathrm{PS}}\approx
6~{}\mathrm{days}$.

In the left panel of Fig. 4 we show the fractional deviation of $\Delta t_{0}$
from the characteristic delay time $t_{d;\mathrm{PS}}=\pi\sqrt{27}M$ for
varying hotspot locations $(r_{\mathrm{e}},\vartheta_{\mathrm{e}})$ in a
Schwarzschild BH spacetime. As can be seen from the figure, this error can be
quite large depending on the hotspot location.

In eq. 2.43, we obtained an expression for the time delay between the
appearance of a pair of higher-order images (e.g., $n=1$ and $n=2$) that
corrected for the angular position of the hotspot. We now forward a proposal
that we can also use it to approximate $\Delta t_{0}$, i.e.,

$\Delta t_{0;\mathrm{approx}}\approx
2(\pi-\vartheta_{\mathrm{e}})\eta_{\mathrm{PS}}\,.$ (3.1)

The right panel of Fig. 4 shows that this approximation does a reasonably good
job of obtaining the exact time delay between the appearance of the primary
($n=0$) and secondary ($n=1$) images. As can be seen from the equation above,
the time delay $\Delta t_{0}$ has a linear dependence on the angular position
of the hotspot. We find that accounting for the (approximately linear)
dependence on the radial position of the hotspot yields an even better
approximation (not described here).

Measurements of such time delays at wavelengths other than the current EHT
observing one ($1.3\mathrm{mm}$) could lead to the possibility of measuring
the shadow size at multiple different frequencies. Since the shadow size of a
BH, in GR, is independent of the frequency at which these observations are
conducted, (nonsimultaneous) multifrequency observations of flaring events can
potentially be used to set up null tests of the achromaticity of the shadow
size.

###  3.2 Angle Offset between Primary and Secondary Images

In Appendix D of Wielgus et al. (2022a), the authors also reported the angular
offset on the image plane between the primary ($n=0$) and secondary ($n=1$)
images using full “slow-light” numerical simulations for their best-fit
hotspot model to be $\approx 140^{\circ}$. We will now attempt to offer simple
arguments to estimate this offset.

Let us consider a point source moving on a Keplerian orbit (i.e., a circular
timelike geodesic) of radius $r$ viewed by an observer face-on
($\mathscr{i}=0$). The $n=1$ photon appears on the screen an additional
$\Delta t_{0}\approx\pi\eta_{\mathrm{PS}}$ (3.1) amount of time after the
$n=0$ photon has appeared. Thus, at any particular observer time, the source
of the $n=0$ image has evolved by $\Delta t_{0}$ relative to the source of the
$n=1$ image. Since the Keplerian angular velocity, $\Omega_{\mathrm{K}}$, is
given by

$\Omega_{\mathrm{K}}(r)=[\partial_{r}f/\partial_{r}R^{2}]^{1/2}\,,$ (3.2)

the angle between the primary and secondary images $\Delta\psi_{0}(t)$ is
given as
$\Delta\psi_{0}(t):=\psi_{1}^{+}(t)-\psi_{0}^{-}(t)=\psi_{1}^{+}(\varphi_{\mathrm{e}}(t))-\psi_{0}^{-}(\varphi_{\mathrm{e}}(t+\Delta
t_{0}))=\psi_{1}^{+}(\varphi_{\mathrm{e}}(t))-\psi_{0}^{-}(\varphi_{\mathrm{e}}(t)+\Omega_{\mathrm{K}}(r)\Delta
t_{0})$, where $\varphi_{\mathrm{e}}(t)$ is the azimuthal coordinate of the
source. This reduces to

$\Delta\psi_{0}(r)=-\pi+\Omega_{\mathrm{K}}(r)\Delta
t_{0}=-\pi+\left[\frac{\partial_{r}f(r)}{\partial_{r}R^{2}(r)}\right]^{1/2}\pi\eta_{\mathrm{PS}}\,.$
(3.3)

We remind the reader that the meaning of the superscript decorations on $\psi$
is given in eq. 2.16. Notice also that for this special configuration of
source and observer, the offset angle becomes a function of the source orbit
radius $r$ alone. In a Schwarzschild BH spacetime, the Keplerian angular
velocity is given as $\Omega_{\mathrm{K}}=\sqrt{M/r^{3}}$.

For the best model parameters of Wielgus et al. (2022a) ($r=11M$), we find
$|\Delta\psi_{0}|\approx 154^{\circ}$, which is close to their value of
$140^{0}$. To compare, for a hotspot moving on the innermost stable circular
orbit ($r=r_{\mathrm{ISCO}}=6M$), we find $|\Delta\psi_{0}|\approx
116^{\circ}$. Thus, eq. 3.3 does a reasonable job, considering that we
neglected the observer inclination used there and have used an approximation
for the time delay.

Figure 5: Evolution of hotspots in a Schwarzschild black hole spacetime, and
their secondary images. Panel [a.] shows hotspot orbits of different radii
($6M,11M$) and inclinations ($22^{\circ},60^{\circ}$) relative to the
observer, located on the $+z-$axis. The normals to the orbital planes are
displayed in the upper corner. The relative orientation of the coordinate
($x,y$) and the image plane ($\alpha,\beta$) axes is also displayed. The time
evolution of the hotspot azimuthal and polar angles are shown in panels [b.]
and [c.] respectively, for one time period of each orbit. Panel [d.] shows the
orbits of their secondary ($n=1$) images on the observer’s screen. The
asymmetry is maximal along the projected normals since these correspond to
photons undergoing maximal ($\alpha<0$) and minimal ($\alpha>0$) angular
deflections. Panel [e.] shows the variation of the fractional radii of the
$n=1$ images, $\bar{\eta}=\eta/\eta_{\mathrm{PS}}-1$, with the image plane
polar angle $\psi$. Finally, panel [f.] shows the same as panel [e.] but for
the $n=0$ (primary) images. The $n=1$ images appear close to the Schwarzschild
shadow radius, $\eta_{\mathrm{PS}}=\sqrt{27}M$, due to strong lensing
($\bar{\eta}_{1}\lesssim 0.33$). The $y-$axis scales of panels [e.] and [f.]
demonstrate the level of image demagnification due to strong lensing. Figure
6: Evolution of the time delay and the angular offset between the primary and
secondary images of an orbiting hotspot. The panel on the left shows the time
for the secondary images to appear once the primaries have appeared, of
hotspots on Keplerian orbits around a Schwarzschild BH (see previous figure).
While the hotspot lights up at a time $t=0$, the primary image appears at a
time $t\approx 2000M$ since we have set here the observer to be located at
$r=2000M$. This travel time depends also on the radius of emission, as
indicated by the vertical lines. The delay time, $\Delta t_{0}$, itself is
periodic, depending strongly on the polar angle evolution of the source (see
panel [c.] of Fig. 5) as well as the source radius. The fluctuations in the
time delay are amplified for highly inclined orbits. The top right panel shows
the evolution of the image plane angle of both the primary and the secondary
images. The slopes of these curves are indicative of the angular velocity of
the hotspot and the amplitude of the fluctuations are indicative of the
inclination of the orbital plane. The bottom right panel shows the angular
offset, $\Delta\psi_{0}(t)=\psi_{1}(t)-\psi_{0}(t)$, between the primary and
secondary images, over one orbit, when both are visible on the screen. In
these panels on the right, the thicker vertical lines indicate the time of
appearance of the $n=1$ image.

###  3.3 Primary and Secondary Image Orbits

We will now explore properties of the trajectories, on the screen, of primary
and secondary images of hotspots on circular orbits around a Schwarzschild BH.

The time evolution of the azimuthal angle $\varphi_{\mathrm{e}}(t)$ of a
hotspot on a Keplerian orbit of inclination $\mathscr{i}(\neq\pi/2)$ can be
obtained from eq. 2 as (cf. also Sec. 2.1)

$\displaystyle\varphi_{\mathrm{e}}(t)=$ $\displaystyle\
\pm_{\varphi}\arctan{\left[\cos{\mathscr{i}}\cdot\tan{\left(\frac{\eta
f}{R^{2}}t\right)}\right]}+p\pi\,,$ $\displaystyle=$ $\displaystyle\
\pm_{\varphi}\arctan{\left[\cos{\mathscr{i}}\cdot\tan{\left(\Omega_{\mathrm{K}}t\right)}\right]}+p\pi\,,$
(3.4)

where $\Omega_{\mathrm{K}}$ depends on the radius of the orbit (3.2). Further,
since the polar angle of such an orbit is given in terms of the azimuthal
angle (2.5), we can write

$\vartheta_{\mathrm{e}}(t)=-\mathrm{arctan}\left[\cot{\mathscr{i}}\cdot\csc{\varphi_{\mathrm{e}}(t)}\right]+q\pi\,.$
(3.5)

Our choice of the integration constants above means that
$\varphi_{\mathrm{e}}(t=0)=0$ and $\vartheta_{\mathrm{e}}(t=0)=\pi/2$, and
that the normal to the orbit lies in the $yz-$plane. We will use the two
integers introduced above, $p$ and $q$, to ensure that the hotspot angular
coordinates lie in the principal sheet of the spherical-polar coordinate
system, i.e., $0\leq\varphi_{\mathrm{e}}(t)<2\pi$ and
$0<\vartheta_{\mathrm{e}}(t)<\pi$. The latter, in particular, is useful
because it enables us to directly employ eq. 2.19 to determine image formation
via eq. 2.20. We discuss this procedure in some more detail below, in Sec.
4.2.

Without loss of generality, we will consider below orbits at inclinations
$0\leq\mathscr{i}<\pi/2$888When the observer lies in the plane of the orbit,
$\mathscr{i}=\pi/2$, most of its image appears on a line through the image
plane origin. In addition, however, when the hotspot crosses the line of sight
(at the caustics), the image contains Einstein rings (or critical curves) of
all orders. For inclinations in $(\pi/2,\pi]$, the negative sign on the right-
hand side of eqs. 2.3 and 3.5 should be changed to positive signs. This is
equivalent to a passive rotation of the bulk and boundary coordinate systems
around the $+z-$axis by $\pi$. and having positive azimuthal angular
velocities ($+_{\varphi}$). Over a complete orbit of this kind, we can see
that the hotspot polar velocity ($\dot{\vartheta}_{\mathrm{e}}$) must be
positive for the first and last quarters (in time period) of the orbit (3.5).
The positive and negative polar velocity sections are separated by the polar
turning points, $\vartheta_{\pm}=\pi/2\pm\mathscr{i}$, at which
$\dot{\vartheta}_{\mathrm{e}}=0$. Along the entire orbit, we have
$\vartheta_{-}\leq\vartheta_{\mathrm{e}}(t)\leq\vartheta_{+}$.

In panel [a.] of Fig. 5, we show four hotspot orbits, of different radii
($r=6M,11M$) and orbital inclinations ($\mathscr{i}=22^{\circ},60^{\circ}$),
around a Schwarzschild BH (shown as a sphere). Panel [b.] shows the evolution
of the azimuthal angle (3.3) and panel [c.] similarly shows the evolution of
the polar angle (3.5) of each hotspot, over one respective orbital time
period, $T_{\mathrm{K}}=2\pi/\Omega_{\mathrm{K}}$ (cf. 3.2). The polar turning
points for each orbit are shown in horizontal lines of the same color in panel
[c.].

Panel [d.] shows the $n=1$ images of these orbits. For all of these orbits, an
$n=1$ photon emitted from the point closest to the observer (positive $z$,
negative $y$) is emitted in the negative $z$ direction, executes about one
full loop around the BH in the bulk, and appears on the image plane on the
negative-half of the $\alpha-$axis. One that is emitted from the point
furthest away from the observer (negative $z$, positive $y$) is also emitted
in the negative $z$ direction but appears on the positive half of the
$\alpha-$axis. The first of these two photons undergoes larger angular
deflection and thus appears closer to the shadow boundary curve, whose radii
are marked in horizontal and vertical red lines. In this way, we understand
that the asymmetry of the $n=1$ image orbits is maximal in the direction of
the normal projected onto the screen.

To see this better, panel [e.] shows the variation of the fractional radii of
the $n=1$ image, $\bar{\eta}_{1}=\eta_{1}/\eta_{\mathrm{PS}}-1$, with the
image plane polar angle, $\psi$. Since $\bar{\eta}_{1}=0$ implies that the
image lies on the shadow boundary curve, this quantity measures the distance
of the image from the latter. Indeed, since $n=1$ photons are strongly lensed,
they appear close to the shadow boundary, with $\bar{\eta}_{1}\lesssim 0.33$
across all orbits and times. As discussed above, we see that images appearing
on the negative $\alpha-$axis (vertical orange line here) are always the
closest to the shadow boundary curve. Furthermore, the $n=1$ photons that
appear on the $\beta-$axis, were emitted either from
$\varphi_{\mathrm{e}}=\pi$ or $0$ (2.16), and from the equatorial plane (3.5).
Thus, this image plane axis collects photons that undergo gravitational
deflections of exactly $\not{\Delta}\vartheta_{n}=n\pi+\pi/2$ (2.19), i.e.,
these are photons that execute precisely an integer number of half-loops
around the black hole. Finally, panel [f.] shows the equivalent of panel [e.]
but for the $n=0$ images. The $y-$axis scales of these two panels show the
effect of image demagnification due to strong lensing (2.37). We should
highlight also that the image asymmetry of the $n=0$ and $n=1$ orbits is
roughly aligned (see Fig. 9).

In the previous two sections, we discussed useful approximations of the time
delay between the appearance of the primary and secondary images as well as
the angular offset between them respectively. We now show, in Fig. 6, the
exact trends for these observables, for the hotspot orbits discussed above.

The left panel shows the time delay, i.e., the time for the secondary image to
appear once the primary has appeared. The hotspot lights up in the bulk at an
observer time $t=0$ (cf. panel [b.] of Fig. 5). However, due to finite travel
time, the primary image appears at a time $t\approx 2000M$ since we have set
the observer to be located at $r=2000M$. This $n=0$ photon travel time depends
weakly on the source radius but not on the orbital inclination here since, as
discussed above, we have chosen $\varphi_{\mathrm{e}}(t=0)=0$ and
$\vartheta_{\mathrm{e}}(t=0)=\pi/2$. This is indicated in the vertical lines.
The delay time itself depends strongly on the source colatitude (cf. panel
[c.] of Fig. 5), and is clearly periodic. We show also the characteristic
delay time, $t_{d;\mathrm{PS}}=\pi\sqrt{27}M$ (2.44) for reference.

In the top right panel, we show the time evolution of the image angular
coordinate for both the $n=0$ as well as the $n=1$ images. The slopes of these
curves are indicative of the orbital period and the amplitude of fluctuations
are a characteristic of the orbital inclination. The thin vertical lines show
the arrival time of the $n=0$ photons whereas the thicker ones show the same
for the $n=1$ photons. Notice how the $n=0$ photon emitted from $r=6M$ reaches
the screen after the one emitted from $r=11M$ and also how this trend is
reversed for the $n=1$ photon. We can understand this as follows. The first
$n=0$ photons are both emitted in the radially outward direction but the first
$n=1$ photons are both emitted toward the BH (see the top left panel of Fig.
2). The $n=0$ photon emitted from further away travels a smaller distance
whereas the $n=1$ photon emitted from further away travels a larger distance.
These trends hold reasonably well for emitter locations outside the photon
sphere. Finally, the bottom right panel shows the nontrivial evolution of the
angle offset between the primary and secondary images. For related “slow
light” descriptions of the entire accretion flow, see Bronzwaer et al. (2018).

In concluding our analysis of hotspots, we note that there is much that we can
learn about the source of flux eruption events (see, e.g., Ripperda et al.
2022), and potentially about spacetime geometry, from future images and
movies, resolved or otherwise.

Thus far, we have focussed our attention on point sources. We move now to a
discussion of the image formation of the entire accretion flow.

Figure 7: Higher order images of a thin disk and a solid sphere in a
Schwarzschild BH spacetime. The left panel shows the first four ($n=1-4$)
photon subrings (higher-order images) of an equatorial-thin disk when viewed
face-on. We choose the inner and outer boundaries of the disk to be located at
$6M$ and $\approx 2\times 10^{4}M$ respectively, to mimic an Novikov-Thorne
accretion disk, which successfully explains quasar spectra. In the panel on
the right, we show the subrings cast by a spherically-symmetric emitting
region that extends down to the horizon ($2M$) from the same outer radius, to
mimic a Bondi-Michel accretion process. While the universal self-similar
scaling of the subrings is clear to see, the structures of the two photon
rings are also significantly different. The subrings cast by the thin-disk are
well separated and appear outside the shadow boundary, which is shown as a red
line in all sectors. For the spherical emission zone, all order subrings
overlap and straddle the shadow boundary. In this case, due to emission from
the line of caustics $\vartheta_{\mathrm{e}}=0,\pi$, we also see the formation
of Einstein rings or critical curves, shown here as the bounding dashed and
solid lines. If sources of emission are present close to the photon sphere
($3M$), or closer to our line of sight, we should expect the photon ring to
look qualitatively similar to the panel on the right. Reducing the outer
boundary of the emission zone will cause the outer edges of the subrings to
move inwards.

##  4 Photon Rings in Schwarzschild Black Hole Spacetimes

The images of extended sources of emission such as accretion disks around
black holes can also be decomposed into the leading-order ($n=0$) and higher-
order ($n>0$) images. We will refer to the latter as (photon) subrings.

The properties of these subrings depend on the morphology of the emitting
region as well as the inclination $\mathscr{i}$ of the observer. In this
section, we will be interested in understanding, both qualitatively and
quantitatively, the variation in the diameters $d_{n}$, widths $w_{n}$, and
asymmetries $\mathscr{A}_{n}$ of the first few higher-order images (cf. also
Cárdenas-Avendaño & Lupsasca 2023). While the first-higher order ($n=1$) image
will likely be accessible in the near future, with space-based very long
baseline interferometry, our interest in the $n=2$ image is in exploring the
potential for a measurement of the lensing Lyapunov exponent (cf. Sec. 2.6).

To study the variations in these characteristic subring features, we will
employ a simple geometric model for the morphology of the emitting region.
This will be a conical torus concentric with the BH, as pictured in the top
right panel of Fig. 2. The inner and outer boundaries of the torus are located
at $r=r_{\mathrm{in}}$ and $r=r_{\mathrm{out}}$ respectively, and its scale
height is characterized by the half-opening angle
$0\leq\vartheta_{1/2}\leq\pi/2$. Finally, the angle between the axis of the
wedge and the $z-$axis determines the viewing inclination,
$0\leq\mathscr{i}<\pi/2$. Since this configuration is composed of a series of
conical surfaces with half-opening angles $\vartheta$ in
$\pi/2-\vartheta_{1/2}\leq\vartheta\leq\pi/2+\vartheta_{1/2}$, this simplistic
model is also suitable for analyzing the images of jets as well.

Figure 8: Higher-order image radii of equatorial emitters in a Schwarzschild
BH spacetime. The left panel shows the radial distance from the black hole
shadow boundary ($\eta_{\mathrm{PS}}=\sqrt{27}M$), $\bar{\eta}_{1}$, at which
the first-order image of an equatorial emitter appears. The $n=1$ image always
appears close to the shadow boundary $|\bar{\eta}|\lesssim 0.2$ and its image
radius exhibits an overall variation of $\approx 20\%$. In particular, for a
“maximal” thin-disk, extending from the horizon to infinity, the diameter and
width of its $n=1$ subring are $\approx 12.4M$ and $\approx 1.3M$
respectively. The scaling relation (2.39) obeyed by the image fractional radii
of consecutive-order images,
$\bar{\eta}_{n+1}=\bar{\eta}_{n}\mathrm{e}^{-\pi}$, is also clear to see. The
inset zooms the $y-$axis by $\sim\\!\mathrm{e}^{\pi}\approx 23$ to highlight
this further. In the panel on the right, we zoom in on the shaded regions in
the left panel and switch to the fractional radial coordinate,
$\bar{r}=r/r_{\mathrm{PS}}-1$, for the $x-$axis. This shows how while most
$n=1$ photons emitted from outside the photon sphere ($\bar{r}>0$) appear
outside the shadow boundary ($\bar{\eta}_{1}>0$), this is not true for all.
Those that are emitted from just outside the photon sphere can appear inside
the shadow boundary. The demarcating photon orbit satisfies the radial turning
point equation $\mathscr{R}(\bar{\eta}_{\mathrm{tp}},\bar{r}_{\mathrm{e}})=0$
and experiences an angular deflection $\not{\Delta}\vartheta$ of exactly
$3\pi/2$ (cf. Appendix B).

We begin by considering, qualitatively, the images of a geometrically-thin
disk ($\vartheta_{1/2}=0$), viewed face on (i.e., $\mathscr{i}=0$), and a
spherical emission zone ($\vartheta_{1/2}=\pi/2$), in a Schwarzschild BH
spacetime. We model the former following Novikov & Thorne (1973) and, thus,
choose $(r_{\mathrm{in}},r_{\mathrm{out}})=(6M,\infty)$. For the latter, we
set $(r_{\mathrm{in}},r_{\mathrm{out}})=(2M,\infty)$, to mimic Michel (1972)
accretion. We will take a closer look at the thin-disk in Sec. 4.1 below.

Fig. 7 shows the first four subrings ($n=1-4$) for the thin-disk and the
spherical zone in the left and right panels respectively. The outer (inner)
edges of the subring correspond to the order$-n$ images of the outer (inner)
boundaries of the emission zone. Since their images are circularly-symmetric
on the sky, we show only a quarter of each photon subring in each sector. As
we go clockwise from the top-left sector (which shows the $n=1$ subring) the
radial scale is successively zoomed by a factor of
$\approx\mathrm{e}^{\pi}\approx 23$. We show both the image plane radial
coordinate $\eta$ as well as the fractional or conformal radial coordinate
$\bar{\eta}=\eta/\eta_{\mathrm{PS}}-1$, which measures the radial distance
from the shadow boundary curve, $\eta=\eta_{\mathrm{PS}}$.

This immediately showcases the qualitative self-similar scaling exhibited by
the subrings as well as the quantitative scaling discussed above, in eqs. 2.37
and 2.38. The former applies to the conformal radii of the edges of the
subrings and the latter to the widths of the subrings. Since all of the
photons collected from the thin-disk execute an exact number of half-loops
before appearing on the image plane, the image self-similarity is explained by
even simpler scaling relations (2.39, 2.40).

We can also see clear qualitative differences in the higher-order images of
these two emission zone morphologies. There are gaps between the subrings cast
by the thin-disk whereas those cast by the spherical zone overlap. The shadow
boundary curve is shown as a red circle in all sectors. We can use this to see
that the photon ring lies entirely outside the shadow boundary curve
($\bar{\eta}>0$) for the thin-disk whereas, for the spherical model, it
straddles the shadow boundary curve. For the thin disk, the “inner photon
ring” ($-1\ll\bar{\eta}<0$) is empty because of the observer inclination as
well as absence of emission from close to the photon sphere (located at
$r=r_{\mathrm{PS}}=3M$).

We move now towards a quantitative estimation of the variation in subring
characteristics by considering, systematically, a sequence of configurations.
In Sec. 4.1 we will consider the variation in the subring characteristics for
geometrically-thin ($\vartheta_{1/2}=0$) emitting regions viewed face-on
($\mathscr{i}=0$) due to varying inner $r_{\mathrm{in}}$ and outer
$r_{\mathrm{out}}$ boundaries. In Sec. 4.2 we will consider the scenario of
thin-disks viewed at an inclination. Finally, in Sec. 4.3 we allow the
geometrical-thickness ($\vartheta_{1/2}\neq 0$) and the boundaries of the disk
to vary, while keeping the inclination fixed to zero, and analyze how it
affects the structure of the photon ring.

###  4.1 Geometrically-Thin Disk viewed Face-on

For this configuration of a thin-disk viewed from zero inclination, the image
is circularly symmetric, and the subrings are perfect rings (cf. also
Bisnovatyi-Kogan & Tsupko 2022). Furthermore, the image contains no Einstein
rings due to the absence of emission from the line of caustics (line of
sight).

The quantitative relation between the order$-n$ image radius
$\eta_{n}(r_{\mathrm{e}})$ and the source radius $r_{\mathrm{e}}$ for an
equatorial emitter is obtained straightforwardly by solving the integral
equation (2.20)

$|\not{\Delta}\vartheta^{\pm}(\eta_{n}(r_{\mathrm{e}}),r_{\mathrm{e}})|=(2n+1)\pi/2\,,$
(4.1)

In the above, we have used the understanding that photons of all orders ($n$)
execute an exact number of half-loops around the BH (2.19). Remember also that
$\not{\Delta}\vartheta^{\pm}$ is negative for even-order photons and positive
for odd-order photons (cf. eq. 2.19).

The radii of the outer and inner edges of the subrings,
$\eta_{n;\mathrm{out}}$ and $\eta_{n;\mathrm{in}}$, for a thin-disk can then
be obtained from the above by setting $r_{\mathrm{e}}=r_{\mathrm{out}}$ and
$r_{\mathrm{e}}=r_{\mathrm{in}}$ respectively. The former radius
($r_{\mathrm{out}}$) corresponds to the outer boundary of the disk and the
latter to its inner boundary. The width of a subring is naturally given as
$w_{n}=\eta_{n;\mathrm{out}}-\eta_{n;\mathrm{in}}$.

The left panel of Fig. 8 shows the image fractional radius,
$\bar{\eta}:=\eta/\eta_{\mathrm{PS}}-1$, of the first-order image ($n=1$) for
an equatorial source present at a radius $r_{\mathrm{e}}$ from a Schwarzschild
BH. We use the fractional radius here to conveniently measure the distance at
which the image appears from the BH shadow boundary,
$\eta_{\mathrm{PS}}=\sqrt{27}M$. Indeed, if we denote the $n=1$ subring
diameter by $d_{1}(=2\eta_{1})$ and the shadow diameter by
$d_{\mathrm{sh}}(=2\eta_{\mathrm{PS}})$, then the deviation of the former from
the latter is given as $d_{1}/d_{\mathrm{sh}}-1=\bar{\eta}_{1}$. Thus, we can
infer from the figure that, for this configuration, the $n=1$ subring diameter
appears very close to the shadow, is $-4\times
10^{-2}\lesssim\bar{\eta_{1}}\lesssim 19\times 10^{-2}$. The overall variation
in the $n=1$ subring diameter is about $20\%$.

The maximal diameter and the maximal width of the $n=1$ subring are given as
$d_{1}=2~{}\mathtt{max}[\eta_{1}]\approx 12.4M$ and
$w_{1}=\eta_{\mathrm{PS}}(\mathtt{max}[\bar{\eta}_{1}]-\mathtt{min}[\bar{\eta}_{1}])\approx
1.3M$ respectively. The maximal sets of the equatorial emission photon
subrings, corresponding to the emission region located anywhere between the
event horizon and infinity, are referred to as (equatorial) lensing bands
(see, e.g., Paugnat et al. 2022; Wielgus 2021).

Figure 9: Angular deflection experienced by photons emitted by a ring, and
image asymmetry. In the panel on the left, we show the angular deflection
experienced by photons emitted from a ring before reaching an observer at an
inclination of $\mathscr{i}$ w.r.t. the normal to the ring. Such rings have
been pictured in the top left panel of Fig. 5. The radius of the ring is
irrelevant for this panel. Order$-n$ photons undergo angular deflections in
the range $(n\pi,(n\\!+\\!1)\pi)$. Photons appearing at the same image plane
polar angle, $\psi$, but of different orders, $n_{1}$ and $n_{2}$, differ in
their angular deflections by precisely $(n_{1}-n_{2})\pi$. This is easy to
understand for a ring viewed face-on ($\mathscr{i}=0$). However, this remains
true for rings viewed from arbitrary inclinations because the $n$ and
$n\\!+\\!1$ photons appearing at the same $\psi$ were emitted from antipodal
points on the ring, i.e.,
$\varphi_{\mathrm{e};n+1}-\varphi_{\mathrm{e};n}=\pi$ (2.16). Thus, if we
denote their source colatitudes as $\vartheta_{\mathrm{e};n}$ and
$\vartheta_{\mathrm{e};n\\!+\\!1}$ respectively, then
$\vartheta_{\mathrm{e};n+1}$ \- $\vartheta_{\mathrm{e};n}=\pi$ (4.3). Thus,
their angular deflections are also simply offset by $\pi$ (2.19). The panel on
the right shows how the asymmetry, $\mathscr{A}_{1}$, (or ellipticity, axis-
ratio) of the $n=1$ image of a ring of radius $r_{\mathrm{e}}$ changes with
the observer inclination. For small ring radii, $r_{\mathrm{e}}\lesssim 10M$,
and moderate inclinations $0\leq\mathscr{i}\lesssim\pi/4$, the asymmetry
varies essentially with radius of the ring and remains quite small. Figure 10:
Higher-order images of a geometrically-thin disk in a Schwarzschild BH
spacetime, viewed at inclination. The panel on the left shows the fractional
median radii of the first three photon subrings whereas the one on the right
shows their median widths. The squiggly lines indicate the locations of the
event horizon ($2M$; black), the photon sphere ($3M$; blue), and the innermost
stable circular orbit ($6M$; green). This figure shows that the deviations
from circularity in the shapes of the higher-order images due to the
observer’s inclination do not impact their median diameters or widths. This is
reminiscent of the independence of the shadow boundary with the observer
inclination $\mathscr{i}$ in static and spherically-symmetric spacetimes. We
see that the $n=1$ subring diameter always remains close to the shadow
boundary remains $|\Braket{\bar{\eta}_{1}}_{\psi}|<0.2$. Thus, a measurement
of the subring diameter yields an accurate estimate of the shadow diameter.
The right panel estimates the widths of the $n=1$ subrings to be at most
$\approx 1M$. Finally, the color bar conveniently shows the error in
recovering the lensing Lyapunov exponent from a joint measurement of the first
two subrings, when using either their diameters (left) or their widths
(right), to be at most about $10\%$. Subring asymmetry due to observer
inclination is discussed in figures 5 and 9.

The fractional (or conformal) radius is also the perfect coordinate with which
to display the self-similar (or conformal) scaling symmetry of the photon ring
(2.39). To that end, we also show in this panel the scaled fractional radii,
i.e., $\mathrm{e}^{(n-1)\pi}\bar{\eta}_{n}$, of the next three higher-order
images ($n=2,3,4$).

The scaling between the $n=1$ and $n=2$ image fractional radii is captured
very well by the relations obtained above (2.39), i.e.,
$\bar{\eta}_{2}\approx\mathrm{e}^{-\pi}\bar{\eta}_{1}$, for small emission
radii $r_{\mathrm{e}}\lesssim 6M$. Nevertheless, we find a maximum error of
only $\approx 9\%$ when using the two lowest order subrings diameters to infer
the lensing Lyapunov exponent, over all radii. Thus, a measurement of the
image radii (or widths) of a pair of higher-order images for small thin-disks
will yield an accurate measurement of the lensing Lyapunov exponent.

We can also see that the scaling relations for images of order $n\geq 2$
continue to work very well for arbitrarily large source radii. To emphasize
this point, the inset shows the (unscaled) fractional radius of the $n=2$
image as well as the scaled radii of the next pair of higher-order images. We
can infer then that the variations in the subring diameters and widths are
suppressed by a factor of $\approx e^{\pi}\approx 23$ per increasing subring
order.

This figure also shows, as expected, how photons emitted from inside the
photon sphere ($r<r_{\mathrm{PS}}-3M$) in the bulk appear always inside the
shadow boundary ($\bar{\eta}<0$). Similarly, $n\geq 1$ photons emitted from
outside the photon sphere typically appear outside the shadow boundary curve.
However, and as discussed above in Sec. 2.4, this is not always true: Photons
emitted from just outside the photon shell can appear inside the shadow
boundary. This is clearly demonstrated for the $n=1$ image by the blue solid
line in the right panel here.

We emphasize that the converse of the observation above may be important in
the context of “spacetime tomography,” i.e., for attempts to map out the
spacetime geometry of astrophysical ultracompact objects by observing, e.g.,
multiple flaring events occurring near M87⋆ or Sgr A⋆ (see, e.g., Tiede et al.
2020). A holistic qualitative description of the different types of photon
orbits that participate in image formation is presented in Appendix B.

###  4.2 Geometrically-Thin Disk viewed from an Inclination

In this section, we will consider the impact of a nonzero observer inclination
on the properties of the photon subrings cast by a thin-disk (cf. also Tsupko
2022). The image and consequently the photon ring is no longer circularly-
symmetric. Since the thin-disk can be decomposed into a series of concentric
great-circles of different radii, an analysis of the formation of its higher-
order images is identical to the image formation of circular hotspot orbits,
which we discussed in Sec. 3.3.

The order$-n$ image, $\eta_{n}(r_{\mathrm{e}},\psi)$, of a ring of emitters of
radius $r_{\mathrm{e}}$ viewed from an inclination, $\mathscr{i}$, can be
determined by solving the integral equation (2.20)

$\displaystyle\not{\Delta}\vartheta^{\pm}\left(\eta_{n}(r_{\mathrm{e}},\psi),r_{\mathrm{e}}\right)=\not{\Delta}\vartheta_{n}^{\pm}(\vartheta_{\mathrm{e}}(\mathscr{i},\psi))\,,$
(4.2)

where, in the above we have (cf. eqs. 3.5, 2.16),

$\vartheta_{\mathrm{e}}(\mathscr{i},\psi)=\begin{cases}+\arctan{[\cot{\mathscr{i}}\cdot\sec{\psi}]}\,,&\mathrm{for\
even\ n}\\\ -\arctan{[\cot{\mathscr{i}}\cdot\sec{\psi}]}\,,&\mathrm{for\ odd\
n}\end{cases}$ (4.3)

for rings of inclinations $0\leq\mathscr{i}<\pi/2$. For rings with
inclinations in $(\pi/2,\pi]$, the signs on the right are flipped. For small
inclinations $\mathscr{i}\approx 0$, we can find that the emission is sourced
from colatitudes
$\vartheta_{\mathrm{e}}(\psi)\approx\pi/2-\mathscr{i}\cos{\psi}$ (cf. Sec. III
D of Gralla et al. 2019).

Fig. 9 shows the variation in the angular deviation
$\not{\Delta}\vartheta_{n}(\vartheta_{\mathrm{e}}(\mathscr{i},\psi))$
discussed above, with the image plane polar angle, $\psi$, for rings of
various inclinations. For photons emitted from a particular ring, the
additional angular deflection experienced by two that appear at the same angle
but which are of different orders differs by $\pi$, in spherically-symmetric
spacetimes. This may be slightly unintuitive but is true because (a) antipodal
points on the ring of emitters appear at antipodal points on the image plane
and (b) consecutive order images of the same point emitter appear at antipodal
points on the image plane.

This means that the orientation of the asymmetry of arbitrary higher-order
images is identical. If the normal to the plane of the ring $n_{\mathrm{d}}$
appears on the image plane as pointing along the negative $\alpha-$axis (see
Fig. 5), this is also the direction on the image plane along which a higher-
order image is also most stretched. In the right panel of Fig. 9, we show the
asymmetry,

$\mathscr{A}_{n}=1-d_{n}(\psi=\pi/2)/d_{n}(\psi=0)\,,$ (4.4)

of the secondary ($n=1$) image. This definition measures the deviation between
the maximum and minimum image diameters, in the directions parallel and
perpendicular to the projected normal $n_{\mathrm{d}}$ respectively. As
expected, viewing an axisymmetric emission source from higher inclinations
leads to higher image asymmetry. Interestingly, the asymmetry remains small
$\mathscr{A}_{1}\lesssim 0.10$ independently of the viewing angle for emission
coming from close to the BH $r\lesssim 10M$ (cf. also Medeiros et al. 2022).

We expect the shape asymmetry of the $n=1$ image to closely track the shadow
boundary curve (see Fig. 10 below), with the latter defined purely by the
spacetime geometry and the observer inclination. It should then be possible to
infer the spin of a black hole from such a measurement (cf., e.g., Johannsen &
Psaltis 2010, Johnson et al. 2020, Kocherlakota et al. 2021). This may be a
particularly good avenue to measure the spin of M87⋆ since we know its
inclination rather well (cf. Walker et al. 2018).

Figure 11: Higher-order images of a geometrically-thick disk in a
Schwarzschild BH spacetime viewed face-on. Same as Fig. 11 but for
geometrically-thick disks, which we model as a conical torus (see right panel
of Fig. 2). The scale-height of the disk is shown on the $x-$axis (i.e., the
surfaces of the disk are at latitudes $\pm\vartheta_{1/2}$). The left panel
shows the $n=1$ subring diameter to be sensitive to variations in the disk
scale-height, even for realistic values of the morphological parameters
$\vartheta_{1/2}\lesssim 0.2$ and $r_{\mathrm{out}}\lesssim 20M$ (shown in
bright green lines; cf. Narayan et al. 2022; Chatterjee & Narayan 2022). Over
these ranges, the $n=1$ subring diameter nevertheless closely tracks the
shadow diameter, $|\bar{\eta}_{1}|\lesssim 0.25$. In the panel on the right,
while we have chosen a large outer boundary radius to show the maximal
possible width variations, we can infer the impact of a varying outer boundary
from the left panel (a negligible change for $r_{\mathrm{out}}\gtrsim 100M$).
As with the left panel, the disk geometrical-thickness plays an important role
in determining subring widths and, for realistic morphologies, the maximal
$n=1$ subring width is approximately $\approx 1.3M$. Since our interferometers
measure angular sizes in practice, we note that a ring of width $w_{n}\approx
1M$ has an angular thickness of one angular gravitational radius
$\theta_{\mathrm{g}}=GM/c^{2}D=M/D$ on the sky. For M87⋆ and for Sgr A⋆, these
have been inferred by the EHT to be
$\theta_{\mathrm{g}}=3.8^{+0.4}_{-0.4}\mu\mathrm{as}$ (EHT Collaboration et
al., 2019a) and $\theta_{\mathrm{g}}=4.8^{+1.1}_{-0.7}\mu\mathrm{as}$ (EHT
Collaboration et al., 2022a) respectively. Finally, the color bars in each
panel show that the lensing Lyapunov exponent can be inferred, in principle,
with an error of $\lesssim 20\%$ by comparing the $n=1$ and $2$ subring
diameters (left) or widths (right).

The left panel of Figure 11 shows the change in the median fractional radii
(similar to Fig. 8) of the first three ($n=1\\!-\\!3$) photon subrings
respectively, with changing inclination of the observer $\mathscr{i}$ and
varying size $r_{\mathrm{out}}$ of a geometrically-thin emission disk. These
image radii correspond to the order$-n$ gravitationally-lensed sizes of the
outer boundary of the disk. Equivalently, these also describe the median image
radii of a single ring of radius $r=r_{\mathrm{out}}$.

The panel on the right shows the variation in the median widths of these
subrings for a varying disk inner boundary. For concreteness, we pick the disk
outer boundary to be located at $r_{\mathrm{out}}\approx 2\times 10^{4}M$.
Changing the outer boundary to smaller values naturally reduces the width of
the ring. We pick this large outer boundary radius to show the magnitude of
the maximal width possible in principle.

These panels demonstrate how the median diameters and widths of the subrings
depend primarily on the size of the thin disk $r_{\mathrm{out}}$ and are
independent of the viewing inclination $\mathscr{i}$. Thus, the maximal
diameters and widths reported in the previous section continue to indicate the
range of possible median diameters and widths, for thin-disks viewed from
arbitrary inclinations.

In addition to encoding the extent of the physical region that sources the
observed emission, the median subring diameters and widths capture the scaling
exponent, which is specific to spacetime geometry. This is seen from the error
($\lesssim 9\%$; shown in the color bars) in inferring the lensing Lyapunov
exponent from a joining measurement of the diameters or widths of the $n=1$
and $n=2$ subrings.

###  4.3 Geometrically-Thick Disk viewed Face-on

There is mounting evidence that the ultracompact objects M87⋆ and Sgr A⋆ host
geometrically-thick accretion flows (EHT Collaboration et al., 2019b, 2022c),
and thus have geometrically-thick emitting regions. The effective scale-
heights, $h/r$, for hot accretion flows around Kerr BHs were carefully studied
using general relativistic magnetohydrodynamics (GRMHD) simulations, and an
upper bound of $h/r\lesssim 0.4$ was obtained (see Fig. 7 of Narayan et al.
2022; See also Porth et al. 2019; Chatterjee & Narayan 2022). This can be
converted into the faces of the emission zone being located at
$(h/2)/r=\tan{[\pm(\pi/2-\vartheta_{\mathrm{e}})]}\approx\pm[\pi/2-\vartheta_{\mathrm{e}}]$,
or, equivalently, $\vartheta_{\mathrm{e}}\approx\pi/2\pm 0.2$.

With this in mind, we will now explore the impact of a varying disk scale-
height on potential inferences of the Schwarzschild photon subring
characteristics. Since we understand the impact of the observer inclination on
the subring characteristics from the previous section, we will set here, for
clarity, the inclination to vanish (face-on observer). The image morphology
for this configuration is also circularly symmetric, and higher-order images
are perfect annuli concentric with the shadow boundary curve.

As discussed above, we will employ a conical torus, which is formed out of the
intersection of cones and spheres, to model a thick-disk. It is parametrized
by three parameters, two that control its inner $r_{\mathrm{in}}$ and outer
$r_{\mathrm{out}}$ spherical surfaces, and a third $\vartheta_{1/2}$ that
modifies its conical faces, i.e., its scale-height. It is to be imagined that
photons are emitted from the region $r_{\mathrm{in}}\leq r_{\mathrm{e}}\leq
r_{\mathrm{out}}$ and between colatitudes
$\pi/2-\vartheta_{1/2}\leq\vartheta_{\mathrm{e}}\leq\pi/2+\vartheta_{1/2}$ (or
equivalently from latitudes between $\pm\vartheta_{1/2}$). The realistic upper
bound (from GRMHD simulations) on the scale height is equivalent to
$\vartheta_{1/2}\lesssim 0.2$.

In the top-right panel of Fig. 2, we show such a conical torus, with
parameters
$\\{r_{\mathrm{in}},r_{\mathrm{out}},\vartheta_{1/2}\\}=\\{2M,18M,\pi/10\\}$.
Its primary ($n=0$) and secondary ($n=1$) images, as seen by a face-on
observer, are shown in the lower right panel there. The regions with
relatively darker shading in the top left panel there capture the properties
of the photons that, as discussed below, form the photon ring ($n\geq 1$) for
this source morphology.

We can understand the higher-order image formation of our conical torus by
considering first the image formation of a conical surface. As can be seen
from the top left panel of Fig. 2, photons of any particular order that are
emitted from increasingly larger radii, $r_{\mathrm{e}}$, but from the same
conical surface (fixed $\vartheta_{\mathrm{e}}$; i.e., as we move along a
fixed horizontal line in this plot) always appear at increasingly larger
radii, $\eta$. On such surfaces, the lensing map, $r_{\mathrm{e}}\mapsto\eta$,
is neatly monotonic (see also related discussion in Sec. 2.4). Therefore, for
emission from any single cone, the lensed size of the emission outer (inner)
boundary radius determines the outer (inner) edge of its image. Understanding
the images of cones in black hole spacetimes finds an interesting application
in the analysis of jet shapes and jet opening angles (Papoutsis et al., 2023).

Since the image of a conical torus can be thought of as the sum of the images
of single cones, to determine the edges of the higher-order images of the
former, we must also understand how the images of different single cones
compete. As can be seen from the solid lines in the top left panel of Fig. 2,
as we move across a spherical surface (fixed $r_{\mathrm{e}}$) located inside
the photon sphere ($r_{\mathrm{e}}<r_{\mathrm{{PS}}}$), the lensing map
displays monotonic behaviour. The inner and outer edges of the order$-n(>0)$
image of smaller spheres correspond to the Einstein rings of its north and
south poles for odd-$n$, and of its south and north poles for even-$n$. Its
primary ($n=0$) image is the entire region inside the first Einstein ring or
the first critical curve.

As we move across a spherical surface (fixed $r_{\mathrm{e}}$) located outside