11institutetext:  LESIA, Observatoire de Paris, 5 place Jules Janssen, 92195
Meudon, France. 11email<EMAIL_ADDRESS>22institutetext:
Princeton Gravity Initiative, Princeton University, Princeton NJ 08544, USA
33institutetext: Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69,
D-53121 Bonn, Germany

# Photon ring test of the Kerr hypothesis: Variation in the ring shape

H. Paugnat 11 A. Lupsasca 22 F. H. Vincent 11 M. Wielgus 33

###### Abstract

Context. The Event Horizon Telescope (EHT) collaboration recently released
horizon-scale images of the supermassive black hole M87*. These images are
consistently described by an optically thin, lensed accretion flow in the Kerr
spacetime. General relativity (GR) predicts that higher-resolution images of
such a flow would present thin, ring-shaped features produced by photons on
extremely bent orbits. Recent theoretical work suggests that these “photon
rings” produce clear interferometric signatures that depend very little on the
astrophysical configuration and whose observation could therefore provide a
stringent consistency test of the Kerr hypothesis.

Aims. We wish to understand how the photon rings of a Kerr black hole vary
with its surrounding emission. Gralla, Lupsasca, and Marrone (GLM) found that
the shape of high-order photon rings follows a specific functional form that
is insensitive to the details of the astrophysical source, and proposed a
method for measuring this GR-predicted shape via space-based interferometry.
We wish to assess the robustness of this prediction by checking it for a
variety of astrophysical profiles, black hole spins, and observer
inclinations.

Methods. We use the ray tracing code Gyoto to simulate images of thin
equatorial disks accreting onto a Kerr black hole. We extract the shape of the
resulting photon rings from their interferometric signatures using a
refinement of the method developed by GLM. We repeat this analysis for
hundreds of models with different emission profiles, black hole spins, and
observer inclinations.

Results. We identify the width of the photon ring and its angular variation as
a main obstacle to the method’s success. We qualitatively describe how this
width varies with the emission profile, black hole spin, and observer
inclination. At low inclinations, our improved method is robust enough to
confirm the shape prediction for a variety of emission profiles; however, the
choice of baseline is critical to the method’s success. At high inclinations,
we encounter qualitatively new effects that are caused by the ring’s non-
uniform width and require further refinements to the method. We also explore
how the photon ring shape could constrain black hole spin and inclination.

###### Key Words.:

Physical data and processes: Gravitation – Accretion, accretion discs – Black
hole physics – Relativistic processes – Galaxies: individual: M87

##  1 Introduction

The existence of black holes (BH) is a key prediction of general relativity
(GR) in the strong-field regime. More specifically, the theory posits that the
spacetime geometry around these compact objects is described by the Kerr
metric. This “Kerr hypothesis” underlies a considerable amount of astrophysics
and also plays a driving role in theoretical physics, but it has yet to be
directly tested. A number of experiments are now poised to change this. The
recent successes obtained in gravitational-wave and radio interferometry by
LIGO-Virgo, GRAVITY, and EHT (LIGO Scientific Collaboration & Virgo
Collaboration, 2016; Gravity Collaboration, 2018; Event Horizon Telescope
Collaboration, 2019a) herald the advent of a new era of black hole astronomy,
and their future extensions promise to deliver observations that will enable
tests of the Kerr hypothesis with unprecedented precision.

In April 2019, the EHT collaboration published the very first images resolving
the nucleus of the galaxy M87 (Event Horizon Telescope Collaboration, 2019a,
b, c, d, e, f). Reconstructed from 1.3 mm very-long-baseline interferometry
(VLBI) observations, these images achieved an angular resolution comparable to
the expected size of M87*—the supermassive compact object at the center of
M87—and revealed a thick, non-uniformly bright ring with a diameter of
$\sim\\!40\,\mu$as encircling a central brightness deficit. These observed
features are roughly compatible with those appearing in simulated black hole
images (e.g., Luminet, 1979; James et al., 2015; Event Horizon Telescope
Collaboration, 2019e), and indeed M87* may be consistently modeled as a Kerr
BH surrounded by an accretion flow whose emission originates within a few
Schwarzschild radii of the event horizon. Direct (weakly lensed) images of
such a flow present an asymmetric bright ring that is consistent with the 2017
EHT observations (e.g., Gralla et al., 2019; Johnson et al., 2020; Chael et
al., 2021).

Due to the limited resolution of these first EHT images of M87*, this
consistency test of the Kerr hypothesis remains weak; instead, the Event
Horizon Telescope Collaboration (2019f) has thus far assumed the Kerr nature
of M87*, so as to constrain its mass from observations. However, as technology
continues to improve, increasingly better images will be obtained, and these
should enable correspondingly sharper probes of the spacetime geometry around
M87*. Nonetheless, it is not immediately clear how to best extract
relativistic signatures from such probes, which raises the question: how could
sharper images of M87* in principle lead to more stringent tests of the Kerr
hypothesis?

Addressing this question is made particularly difficult by the inherent
asymmetry of the imaging problem. Thanks to curved-spacetime ray tracing codes
such as Gyoto (Vincent et al., 2011) and others (Gold et al., 2020), solving
the forward problem has become a fairly routine task, and determining the
observational appearance of a compact object with specified properties
(spacetime metric, emission model) is now relatively straightforward. The
inverse problem, on the other hand, remains nigh intractable: inferring the
nature of an astrophysical source from its visual appearance requires one to
disentangle effects of the spacetime geometry from properties of the emission,
a problem subject to enormous degeneracy. A veritable zoo of theoretically
proposed objects (besides the Kerr BH) have been shown to also produce
simulated images with features compatible with the recent EHT
observations,111A non-exhaustive list includes: exotic compact objects (e.g.,
Vincent et al., 2016; Herdeiro et al., 2021), wormholes (e.g., Wielgus et al.,
2020), non-singular black holes (Lamy et al., 2018), etc.  and many of these
alternative models cannot yet be ruled out by observations with current EHT
resolution (e.g., Vincent et al., 2021).

This present blurriness is not the only limitation: even with higher
resolution, images that resolve only the direct emission will likely never
separate astrophysical and geometrical effects. (This partly explains the
present lack of robust constraints on the spin of M87*.) As such, these images
may always be compatible with models of a Kerr BH surrounded by some (possibly
exotic) emission, and hence they could never lead one to definitively conclude
a violation of the Kerr hypothesis, as this would require achieving greater
confidence in the astrophysics than in GR—this could of course break the
degeneracy in principle, but seems very unlikely in practice (e.g., Gralla,
2020; Bauer et al., 2022).

Recent work suggests a promising path forward. The strong gravity of a Kerr BH
creates around it a “photon shell”—a region outside its event horizon in which
its gravity is so strong that light rays can become trapped on unstably bound
orbits (Bardeen, 1973; Teo, 2003). As a result, sources in the vicinity of the
black hole can produce multiple relativistic images in the observer sky:
provided that the emission is optically thin, the primary image (consisting of
weakly lensed photons that travel directly to the observer after being
emitted) is generically superimposed with a series of mirror images arising
from highly bent photons that circumnavigated the photon shell multiple times
on their way from source to observer (Johnson et al., 2020; Gralla & Lupsasca,
2020a). If the emission is not spherical (so that it does not fully surround
the hole), then these successive images form a discrete sequence of “photon
rings” that are labeled by photon half-orbit number $n$ and converge to a
“critical curve” in the observer sky corresponding to the light rays that
asymptote to unstably bound orbits in the photon shell. Crucially, GR predicts
that these rings must display a self-similar structure governed by Kerr
critical exponents $\gamma$, $\delta$, and $\tau$ that respectively control
the successive demagnification, rotation, and time delay of successive images
(Johnson et al., 2020; Gralla & Lupsasca, 2020a). This suggests that observing
this characteristic lensing pattern may provide a way to extract information
about the black hole geometry, and to discriminate between different
spacetimes (e.g., Wielgus, 2021).

The photon ring has not yet been conclusively observed, but its GR-predicted
substructure has already been numerically confirmed with state-of-the-art
(general relativistic magneto-hydrodynamics) GRMHD-simulated images of M87*,
in which at least $n=3$ subrings can be distinguished (Johnson et al., 2020).
Although these rings (which are observable) are distinct from the theoretical
critical curve (which is not observable), they closely track its shape, which
is insensitive to the astrophysical details of the surrounding emission and
entirely fixed by the Kerr metric (Bardeen, 1973; Gralla & Lupsasca, 2020c).
Critical curves have been extensively discussed in the literature and
spacetimes other than Kerr are known to produce different shapes (e.g.,
Johannsen & Psaltis, 2010; Johannsen, 2013; Cunha & Herdeiro, 2018; Okounkova
et al., 2019; Medeiros et al., 2020). Thus, measuring the shape of a high-$n$
photon ring of M87* may be a good proxy for probing the shape of its critical
curve, and hence its geometry.

This idea was recently made precise by Gralla et al. (2020), who established
that, in a certain class of models, the shape of the $n=2$ subring is
extremely close to that of the critical curve. Moreover, they showed that on
very long baselines $u\sim 300$ G$\lambda$, which could be accessible with
space-based222The Earth’s atmosphere limits the frequency, and hence the
baseline length, of observations made with Earth-bound elements of the EHT. As
such, a photon ring shape measurement likely requires an interferometer with a
space leg. The ngEHT—the currently envisioned next-generation EHT (Doeleman et
al., 2019)—is slated to augment its coverage of Earth baselines with
additional ground-based stations; this will improve the image fidelity but not
resolve narrow features like the photon ring. Other proposals to deploy a
radio telescope to space include Haworth et al. (2019); Pesce et al. (2019);
Gurvits et al. (2022).  extensions of the EHT, this subring presents a clean
interferometric signature, as expected from Fourier theory (Johnson et al.,
2020; Gralla, 2020). They argued that a space mission targeting M87* could
resolve the shape of its $n=2$ photon ring with sufficient precision to enable
a sharp comparison with the GR prediction: simulations of their experiment
achieved sub-percent precision, forecasting a stringent consistency test of
the Kerr hypothesis. We will refer to their proposed procedure as the GLM
method.

Of course, this first study is subject to multiple limitations, as it relies
on a restricted class of emission models consisting of stationary,
axisymmetric disks composed of circular-equatorial orbiters. More precisely:
1) it only investigated a choice of 3 toy emission profiles, viewed from
observer inclinations of 10${}^{\circ}-30^{\circ}$ (the range of expected
relevance to M87*); 2) it considered only geometrically thin disks; 3) it
ignored astrophysical fluctuations, which amounts to studying (coherently)
time-averaged images.

It is crucial to assess the robustness of the GLM method in greater
generality. The aim of this paper is to extend the GLM analysis to many more
configurations; that is, to broaden its scope by performing a parameter survey
of emission profiles, as well as BH spins and inclinations, in order to remove
the first limitation outlined above. This work remains confined to the
framework of thin, equatorial accretion disks; a parallel study of the impact
of disk thickness is ongoing and will soon address the second limitation
(Vincent et al., 2022). In the meantime, it is worth noting that recent work
by Chael et al. (2021) suggests that geometrically thin models yield effective
approximations to numerically modeled GRMHD accretion flows.

Finally, another objective of this paper is to initiate the study of black-
hole-spin estimates based on astrophysics-independent shape measurements of
the $n=2$ photon ring. A first step is taken in this direction by establishing
that photon rings arising from astrophysically plausible configurations
constitute only a small subset of the full range of theoretically allowed
shapes.

The paper is organized as follows. First, Sec. 2 reviews the theoretical
foundations of the GLM method: the critical curve is defined and contrasted
with the concept of “black hole shadow,” the notion of “lensing bands” is
introduced, and the expected interferometric signature of the photon ring is
described. Next, Sec. 3 presents a numerical implementation of the GLM method
that improves upon the original technique. Then, Sec. 4 describes how the
choice of baseline may affect a shape measurement of the photon ring,
especially if its thickness varies substantially with angle. After that, Sec.
5 presents the results of a survey over emission profiles while Sec. 6
examines the influence of BH spin and inclination. We find that a measurement
is always possible at low inclinations, while new complications that may arise
at high inclination necessitate further refinements to the method. We provide
further discussion in Sec. 7, before ending with our conclusions and
perspectives in Sec. 8.

##  2 Theoretical overview

###  2.1 Shadow, critical curve, and photon rings

In the introduction, we asserted that the photon ring generically decomposes
into a stack of discrete subrings. For completeness, here we list the (mild)
astrophysical assumptions required for the presence of this substructure; this
also provides us with an opportunity to distinguish the photon ring from the
“black hole shadow” and critical curve, terms deserving of disambiguation.

First, the emission region must be optically thin, so that light may cross it
multiple times before being reabsorbed; otherwise, highly lensed ($n>0$)
images may be obscured (e.g., Beckwith & Done, 2005). The validity of this
assumption depends on both the source and the observing frequency. As
indicated in Fig. 3 of Johnson et al. (2020), GRMHD predicts it to likely hold
for 230 GHz photons collected from M87* by present-day EHT, and it is even
more likely true at the higher 345 GHz frequency that ngEHT will target; by
contrast, the core of Centaurus A is only expected to be optically thin past
1–5 THz (Janssen et al., 2021).

Second, the emission region must have a gap rather than fully surround the
hole; otherwise, the discrete structure visible in the intensity cross-
sections displayed in Johnson et al. (2020); Chael et al. (2021), which
resemble the layers of a wedding cake, will not be present. Instead, one will
see a continuous cross-section with a mild (logarithmic) divergence near the
the critical curve: the image-plane curve, first derived by Bardeen (1973)
(who used the term “apparent boundary”), corresponding to light rays that
asymptote to unstably bound spherical photon orbits in the photon shell. The
resulting image will then display a distinctive feature known as the “black
hole shadow”: a central brightness deficit precisely bounded by the critical
curve. This effect was first highlighted by Falcke et al. (2000) in the
context of a BH surrounded by a spherical, radially infalling accretion flow,
and later revisited by Johannsen & Psaltis (2010).

However, as recently pointed out by Gralla et al. (2019) and Narayan et al.
(2019), this spherical-infall scenario is highly fine-tuned, and in fact no
longer viewed as relevant for M87*: the latest EHT constraints derived from
the polarimetric image of M87* (Event Horizon Telescope Collaboration, 2021a,
b) strongly favor “magnetically arrested disk” (MAD) models (Narayan et al.,
2003; Igumenshchev et al., 2003) with emission concentrated near the midplane.
Such models do not display the traditional “shadow” feature (i.e., a central
brightness depression filling the critical curve), but rather a different
“inner shadow” feature associated with the inner edge of the emission, which
produces a smaller but even darker central brightness depression contained
well within the critical curve (Chael et al., 2021). For equatorial disks
extending to the horizon, the inner edge of the emission in principle
coincides with the lensed equatorial event horizon (though it may appear
larger due to redshift effects).

To summarize, the boundary of the black hole shadow is a mathematical curve
which is not in itself observable, and only happens to coincide with an image
feature—the central brightness depression—in special configurations currently
disfavored by EHT data on M87*. On the other hand, the photon ring is a
visible feature that always arises in images of a BH surrounded by an
optically thin emission region. Moreover, if this emission is gapped, then
this photon ring decomposes into subrings that converge (exponentially in $n$)
to the theoretical critical curve, which can be thought of as the $n\to\infty$
subring. The most easily accessible $n=1$ and $n=2$ subrings are close to, but
nonetheless still distinct from, the exact critical curve. Therefore their
shape is constrained to be similar, but not exactly equal, to that of the
shadow. This observation is key for the GLM method.

###  2.2 Bound geodesics form the photon shell

Before turning to the study of the observable photon ring shape, we first
elucidate the shape of the theoretical image-plane curve comprising the null
geodesics that asymptote to unstably bound orbits around a Kerr BH (Bardeen,
1973). The region of spacetime containing these orbits—the photon shell—is
extensively reviewed in Teo (2003); Johnson et al. (2020). Here we describe it
with Boyer-Lindquist coordinates $(t,r,\theta,\phi)$ on the spacetime of a
Kerr BH with mass $M$ and angular momentum $J=a_{*}M$:

$\displaystyle
ds^{2}=\frac{\Delta}{\Sigma}\left(\mathop{}\\!\mathrm{d}t-a_{*}\sin^{2}{\theta}\mathop{}\\!\mathrm{d}\phi\right)^{2}+\frac{\Sigma}{\Delta}\mathop{}\\!\mathrm{d}r^{2}+\Sigma\mathop{}\\!\mathrm{d}\theta^{2}$
$\displaystyle+\frac{\sin^{2}{\theta}}{\Sigma}\left[\left(r^{2}+a_{*}^{2}\right)\mathop{}\\!\mathrm{d}\phi-
a_{*}\mathop{}\\!\mathrm{d}t\right]^{2},$ (1)

with $\Delta(r)=r^{2}-2Mr+a_{*}^{2}$ and
$\Sigma(r,\theta)=r^{2}+a_{*}^{2}\cos^{2}{\theta}$. We will often use a
dimensionless spin parameter $a\equiv a_{*}/M\in[-1,1]$.

A Kerr photon with four-momentum $p^{\mu}$ has a conserved energy-rescaled
angular momentum and Carter constant

$\displaystyle\lambda=\frac{p_{\phi}}{-p_{t}},\qquad\eta=\frac{p_{\theta}^{2}}{p_{t}^{2}}-a_{*}^{2}\cos^{2}{\theta}+\lambda^{2}\cot^{2}{\theta}.$
(2)

Bound photon orbits in Kerr are “spherical”: they evolve at a constant Boyer-
Lindquist radius in the range $\tilde{r}_{-}\leq r\leq\tilde{r}_{+}$, where

$\displaystyle\tilde{r}_{\pm}=2M\left[1+\cos\left(\frac{2}{3}\arccos\left(\pm
a\right)\right)\right]$ (3)

is the radius of the retrograde (upper sign) or prograde (lower sign)
circular-equatorial orbit. As they rotate around the BH, these bound orbits
undergo polar librations between angles
$[\tilde{\theta}_{-},\tilde{\theta}_{+}]$, where

$\displaystyle\tilde{\theta}_{\pm}$
$\displaystyle=\arccos\left(\mp\sqrt{\tilde{u}_{+}}\right)\gtrless\frac{\pi}{2},$
(4) $\displaystyle\tilde{u}_{\pm}$
$\displaystyle=\frac{r}{a_{*}^{2}(r-M)^{2}}\left[-r^{3}+3M^{2}r-2a_{*}^{2}M\right.$
$\displaystyle\qquad\qquad\qquad\left.\pm
2\sqrt{M\Delta\left(2r^{3}-3Mr^{2}+a_{*}^{2}M\right)}\right],$ (5)

and their characteristic conserved quantities are given by

$\displaystyle\tilde{\lambda}=\frac{M\left(r^{2}-a_{*}^{2}\right)-r\Delta}{a_{*}(r-M)},\qquad\tilde{\eta}=-a_{*}^{2}\tilde{u}_{+}\tilde{u}_{-}.$
(6)

Bound orbits are either prograde or retrograde, depending on the sign of
$\tilde{\lambda}$. Orbits at generic radii $r\in[\tilde{r}_{-},\tilde{r}_{+}]$
are not closed, as they do not return to the same location after a full
rotation by $\Delta\phi=2\pi$, but rather explore the entire allowed region
within their shell of fixed $r$; resonant orbits which are closed form an
exceptional set of measure zero in the photon shell (Wong, 2021). The photon
shell occupies a 3D region of space given at any time $t$ by
$\tilde{r}_{-}\leq r\leq\tilde{r}_{-}$, $\theta_{-}\leq\theta\leq\theta_{+}$,
and $0\leq\phi<2\pi$. In the Schwarzschild case ($a=0$), it reduces to the
sphere at $r=3M$.

###  2.3 Asymptotically bound geodesics form the critical curve

A photon that is initially outside the photon shell but has exactly the same
conserved quantities $(\lambda,\eta)=(\tilde{\lambda},\tilde{\eta})$ as a
photon bound at orbital radius $r$ asymptotes to that orbit in its far future
or past. The critical curve is defined as the theoretical curve in the sky of
an observer corresponding to these asymptotically bound photons. It delineates
the region of photon capture by the BH (inside the curve) from that of photon
escape (outside the curve).

For the case of an observer at large distance $D\gg M$ and inclination $i\neq
0$ from the spin axis of a BH, Bardeen (1973) showed that a photon with
conserved quantities $(\lambda,\eta)$ appears in the sky at a position given
by333Only photons with vanishing angular momentum $\lambda=0$ can pass over
the poles. Thus, for the special case $i=0$ of an on-axis observer,
$\rho=D^{-1}\sqrt{\eta+a_{*}^{2}}$ and $\varphi_{\rho}=\phi_{o}$ is the
photon’s angle of arrival, while the critical curve is the circle
$\tilde{\rho}=D^{-1}\sqrt{a_{*}^{2}\left(1-\tilde{u}_{+}\tilde{u}_{-}\right)}$
corresponding to the bound orbit at the unique radius $\tilde{r}_{0}$ such
that $\tilde{\lambda}=0$, which is given in (Gralla & Lupsasca, 2020a), where
the “BH line of sight” is also defined.

$\displaystyle\rho$
$\displaystyle=\frac{1}{D}\sqrt{\eta+\lambda^{2}+a_{*}^{2}\cos^{2}{i}},\qquad\cos{\varphi_{\rho}}=-\frac{\lambda}{D\rho\sin{i}},$
(7)

where, following Johnson et al. (2020), we replaced Bardeen’s original
Cartesian coordinates $(\alpha,\beta)$ with (dimensionless) polar coordinates
$(\rho,\varphi_{\rho})$ that are centered about the “line of sight to the
black hole” and such that $\varphi_{\rho}=90^{\circ}$ corresponds to the BH
spin axis projected onto the plane perpendicular to this line of sight. Hence,
for such an observer, the critical curve is the parametric curve
$C_{\gamma}=(\tilde{\rho},\tilde{\varphi}_{\rho})$ obtained by tracing

$\displaystyle\tilde{\rho}$
$\displaystyle=\frac{1}{D}\sqrt{\tilde{\lambda}^{2}+a_{*}^{2}\left(\cos^{2}{i}-\tilde{u}_{+}\tilde{u}_{-}\right)},$
(8a) $\displaystyle\cos{\tilde{\varphi}_{\rho}}$
$\displaystyle=-\frac{\tilde{\lambda}}{D\tilde{\rho}\sin{i}},$ (8b)

over the radial extent $r\in[\tilde{r}_{-},\tilde{r}_{+}]$ of the photon
shell. Counter-intuitively, the angle $\tilde{\varphi}_{\rho}$ is
parameterized not by the angle $\phi$ around the BH, but rather by the radius
$r$ within the photon shell. This effect highlights the warped nature of the
Kerr spacetime.

Since Eq. (8b) admits two solutions for every choice of $r$, each orbit in the
photon shell corresponds to two angles $\tilde{\varphi}_{\rho}$ around the
critical curve: this corresponds to the fact that an asymptotically bound
photon can reach the observer on either an upwards or downwards libration.
Only the equatorial observer with $i=90^{\circ}$ can see the full set of
orbits within the photon shell; other observers see only the subset of orbits
that have sufficient inclination to reach them (that is, such that
$\tilde{\rho}^{2}\geq 0$). Examples of critical curves for various BH spins
and inclinations are plotted in many papers (e.g., Fig. 1 in Farah et al.,
2020).

Finally, we emphasize that the notion of a critical curve is not limited to
observers far from the BH. For instance, the critical curve for observers on
timelike circular-equatorial orbits around a Kerr BH is analyzed in Gates et
al. (2021).

###  2.4 Nearly bound photons form the photon ring

A photon with conserved quantities
$(\lambda,\eta)\approx(\tilde{\lambda},\tilde{\eta})$ very close to those of a
photon bound at orbital radius $\tilde{r}$ can describe multiple oscillations
in $\theta\in[\tilde{\theta}_{-},\tilde{\theta}_{+}]$ near the photon shell
before leaving it. Since Kerr bound orbits are unstable, their nearby, nearly
bound geodesics display an exponential rate of radial deviation with respect
to their (fractional) number $N$ of polar half-orbits,

$\displaystyle\left|r-\tilde{r}\right|\propto e^{\gamma N},$ (9)

where one full orbit is defined as a complete oscillation from a turning point
$\theta_{\pm}$ back to itself, while fractional $N$ is precisely defined in
Eq. (36) of Gralla & Lupsasca (2020a). The Lyapunov exponents governing the
orbital instability depend on both the BH spin and photon shell radius
(Johnson et al., 2020):

$\displaystyle\gamma=\frac{4}{a_{*}}\sqrt{r^{2}-\frac{Mr\Delta}{(r-M)^{2}}}\int_{0}^{1}\frac{\mathop{}\\!\mathrm{d}t}{\sqrt{\left(1-t^{2}\right)\left(\tilde{u}_{+}t^{2}-\tilde{u}_{-}\right)}}.$
(10)

By definition, when a nearly bound photon reaches a distant observer, it
appears at a position close to the critical curve (8):
$(\rho,\varphi_{\rho})\approx(\tilde{\rho}(\tilde{r}),\tilde{\varphi}_{\rho}(\tilde{r}))$
where $\tilde{r}$ is the radius of the nearby bound orbit. More precisely, its
arrival position can be labeled by the coordinates $(\tilde{r},d)$, where $d$
denotes the perpendicular distance from the nearest point
$(\tilde{\rho}(\tilde{r}),\tilde{\varphi}_{\rho}(\tilde{r}))$ on the critical
curve. Johnson et al. (2020) heuristically showed, and then Gralla & Lupsasca
(2020a) rigorously proved, that a light ray shot backwards into the geometry
from position $(\tilde{r},d)$ in the sky executes

$\displaystyle
N(\tilde{r},d)\approx-\frac{1}{\gamma(\tilde{r})}\log\left|d\right|$ (11)

half-orbits around the BH before either returning to asymptotic infinity, if
it was outside the critical curve ($d>0$), or crossing the event horizon, if
it was inside ($d<0$). Gralla & Lupsasca (2020a) also compute subleading
$\mathcal{O}(d^{0})$ corrections to $N(\tilde{r},d)$.

While Eq. (11) is a purely geometric statement about the lensing behavior of a
Kerr BH, it also constrains the structure of BH images. To describe these
constraints, we will from now on restrict our attention to emission regions
that are optically thin, stationary and axisymmetric. The assumption of
stationarity and axisymmetry amounts to considering only time-averaged images
(averaged over sufficiently long timescales), and the assumption of infinite
optical depth allows one to neglect absorption effects; as discussed in Sec.
2.1, this approximation becomes exact in the limit of large observation
frequency, and could already provide a good approximation at the frequencies
of current or near-future observations of M87*.

Since the emission region around the BH is optically thin, a light ray that
passes through it $2N$ times (each time at essentially the same inclination)
collects roughly twice as many photons as a neighboring light ray that passes
through it $N$ times. That is,

$\displaystyle I(\tilde{r},d)\propto
N(\tilde{r},d)\sim-\frac{1}{\gamma(\tilde{r})}\log\left|d\right|.$ (12)

As such, the lensing formula (11) implies a mild (logarithmic) divergence in
the observed intensity $I(\rho,\varphi_{\rho})$ near the critical curve. The
photon ring is defined as the intensity bump caused by this divergence; in
practice, it is cut off by absorption effects (finite optical depth) as $N$
grows large. In summary: GR predicts that embedded within every image of a BH
surrounded by an optically thin emission region, there lies a narrow photon
ring.

Due to Eq. (11), any axisymmetric ring emitting isotropically produces an
infinite sequence of lensed images within the photon ring (Gralla & Lupsasca,
2020a). Successive images within this sequence are produced by photons
describing increasingly many half-orbits around the black hole: if the first
(direct) image arises from photons with fractional half-orbit number $0\leq
N_{0}<1$, then the $n^{\text{th}}$ lensed image will arise from photons with
fractional half-orbit number $N\approx N_{0}+n$. Since the entirety of the
emission surrounding the black hole can be decomposed into such source rings,
it follows that the photon ring consists of a series of lensed images of the
main emission labeled by half-orbit number $n\in\mathbb{N}$.

Moreover, Eq. (11) implies that these successive images are increasingly
lensed towards the critical curve by an exponential demagnification factor: a
source ring whose $n^{\text{th}}$ image appears at perpendicular distance
$d_{n}(\tilde{r})$ from the critical curve will produce its
$(n+1)^{\text{th}}$ image at an exponentially smaller distance

$\displaystyle d_{n+1}(\tilde{r})=e^{-\gamma(\tilde{r})}d_{n}(\tilde{r}).$
(13)

Depending on whether the emission region is gapped or not, this lensing
pattern can result in one of two very distinctive image substructures within
the photon ring.

If the BH is entirely immersed in its surrounding emission (as in the case of
a spherically symmetric accretion flow onto the hole, for instance), then
successive images of the emission region continuously blend together and the
resulting intensity profile logarithmically rises to a smooth peak centered
about the critical curve. Such an intensity profile is displayed in Fig. 1 of
Narayan et al. (2019), with the corresponding image shown in Fig. 2 therein: a
bright photon ring encircling a central brightness deficit (the interior of
the critical curve). Moreover, as shown in that paper, if the emitting matter
is infalling, then beaming effects strongly attenuate the intensity inside the
critical curve, creating the distinctive shadow effect reviewed in Sec. 2.1
above. The reason is that light rays appearing inside the critical curve
cannot encounter a radial turning point, so any photons collected along such
rays had to have been emitted opposite to the infalling matter’s direction of
motion, thereby incurring a strong redshift.

On the other hand, if the emission is gapped (as in the case of emission
localized around the equator with gaps near the poles, for instance), then
each successive image of the emission region displays a corresponding gap, and
the resulting intensity profile no longer consists of a smooth logarithmic
peak (12), but rather a sequence of distinct peaks accumulating near the
critical curve. Such an intensity profile is displayed in Fig. 1 of Gralla et
al. (2019), with the corresponding image shown in Fig. 5 therein. In such
images, the photon ring decomposes into a stack of distinct photon subrings,
each of which is labeled by $n$ and corresponds to a single lensed image of
the main emission, demagnified by $e^{-\gamma n}$. As reviewed in Sec. 2.1,
the latest EHT observations favor M87* models with near-equatorial emission,
and thus we expect this subring structure to be present in actual images of
M87*.

High-resolution, time-averaged, 230 GHz images of M87* produced using state-
of-the-art GRMHD simulations display at least $n=3$ subrings and numerically
confirm this GR-predicted substructure (Johnson et al., 2020), which is also
present in all the images that we obtained from our toy models of M87*.

###  2.5 Photon subrings lie in lensing bands

Figure 1:  $n=1$ and $n=2$ lensing bands for a distant observer at inclination
$i=45^{\circ}$ from a BH with spin $a=0.99$. The critical curve is given as a
reference.

To study the subring structure of the photon ring, it is helpful to introduce
the (purely geometric) notion of a lensing band (Gralla et al., 2020; Chael et
al., 2021): a Kerr observer’s $n^{\text{th}}$ lensing band is defined as the
region in the sky corresponding to geodesics that cross the equatorial plane
at least $n+1$ times before escaping to asymptotic infinity (if $d>0$) or
crossing the event horizon $r=r_{+}$ of the BH (if $d<0$), where the inner and
outer horizon radii are

$\displaystyle r_{\pm}=M\pm\sqrt{M^{2}-a^{2}}.$ (14)

We present a general method for computing Kerr lensing bands in App. A. Figure
1 illustrates the $n=1$ and $n=2$ lensing bands for an observer at inclination
$i=45^{\circ}$ from a BH of spin $a=0.99$.

Each $n>0$ lensing band consists of a region bounded by two concentric, closed
curves, and therefore takes an annular shape. The inner and outer edges of the
$n^{\text{th}}$ lensing band correspond to the $n^{\text{th}}$ lensed images
of the equatorial circles of radius $r=r_{+}$ and $r\to\infty$, respectively.
The inner edge is always contained within the critical curve, while the outer
edge always encircles it on the outside. According to Eq. (13), these edges
converge to the critical curve exponentially fast in $n$, with widths scaling
as

$\displaystyle w_{n+1}(\tilde{r})\approx
e^{-\gamma(\tilde{r})}w_{n}(\tilde{r})\approx
e^{-n\gamma(\tilde{r})}w_{1}(\tilde{r}),$ (15)

where $w_{1}(\tilde{r})$ is the angle-dependent width of the first ($n=1$)
lensing band (here parameterized by photon shell radius) and the first
relation becomes exact in the limit of large $n\to\infty$. In short, the
lensing bands form a stacked sequence of annular shapes, each straddling the
critical curve, with the $(n+1)^{\text{th}}$ lensing band strictly contained
within the $n^{\text{th}}$ one.

From now on, we further restrict our attention to a thin disk of isotropic
emitters filling the equatorial plane. In this setting, the lensing bands are
particularly useful because they coincide with the photon subrings: the
$n^{\text{th}}$ lensed image of the disk, which produces the $n^{\text{th}}$
photon subring, exactly fills the $n^{\text{th}}$ lensing band. It follows
that the photon subrings also display a stacked structure with widths scaling
as in Eq. (15). Hence, measuring the angle-dependent ratio
$w_{n+1}(\tilde{r})/w_{n}(\tilde{r})\approx e^{-\gamma(\tilde{r})}$ of
successive subrings could (in this setting) allow for the Lyapunov exponents
$\gamma(\tilde{r})$ to be experimentally determined (Johnson et al., 2020).

###  2.6 Interferometric signature of a narrow curve in the sky

Having described the subring structure of BH images, we next turn to a
description of the corresponding subring structure in the visibility (Fourier)
domain. We begin by reviewing some facts about the interferometric signature
of a narrow feature in the sky.

We consider an infinitely thin, closed light curve in the sky. Assuming this
plane curve is convex, we may parameterize it by normal angle
$\varphi\in[0,2\pi)$ as $\@vec{r}(\varphi)=(\alpha(\varphi),\beta(\varphi))$,
where we temporarily revert to Cartesian coordinates $(\alpha,\beta)$. (Non-
convex curves decompose into multiple such segments $\@vec{r}_{i}(\varphi)$,
but we will not need this generalization here.) The inward-pointing unit
normal is $\hat{\@vec{n}}(\varphi)=-(\cos{\varphi},\sin{\varphi})$. Following
Gralla & Lupsasca (2020c), we define the “projected position” $f(\varphi)$ by

$\displaystyle
f(\varphi)\equiv-\@vec{r}(\varphi)\cdot\hat{\@vec{n}}(\varphi)=\alpha(\varphi)\cos{\varphi}+\beta(\varphi)\sin{\varphi}.$
(16)

This function encodes the whole curve
$\@vec{r}(\varphi)=(\alpha(\varphi),\beta(\varphi))$, which may be
reconstructed as

$\displaystyle\alpha(\varphi)$
$\displaystyle=f(\varphi)\cos{\varphi}-f^{\prime}(\varphi)\sin{\varphi},$
(17a) $\displaystyle\beta(\varphi)$
$\displaystyle=f(\varphi)\sin{\varphi}+f^{\prime}(\varphi)\cos{\varphi}.$
(17b)

The projected position $f(\varphi)$ thus encodes all other properties of the
curve, such as its radius of curvature
$\mathcal{R}(\varphi)=f(\varphi)+f^{\prime\prime}(\varphi)$. The real
motivation for introducing the quantity $f(\varphi)$ is that it naturally
connects to the interferometric signature of the curve. To see this, we first
note that the projected position can be uniquely decomposed into parity-even
and parity-odd parts

$\displaystyle f(\varphi)=\frac{d_{\varphi}}{2}+C_{\varphi},$ (18)

where

$\displaystyle d_{\varphi}$ $\displaystyle=f(\varphi)+f(\varphi+\pi),$ (19a)
$\displaystyle C_{\varphi}$
$\displaystyle=\frac{1}{2}\left[f(\varphi)-f(\varphi+\pi)\right].$ (19b)

The even part $d_{\varphi}=d_{\varphi+\pi}$ is the “projected diameter” of the
curve at angle $\varphi$, while the odd part $C_{\varphi}=-C_{\varphi+\pi}$ is
its centroid motion. Gralla & Lupsasca (2020c) illustrate several examples of
this shape decomposition.

We now arrive at the key point: if the curve is not infinitely thin, but still
very narrow relative to its diameter ($w\ll d_{\varphi}$), then its Fourier
transform $V(u,\varphi)$ (where $(u,\varphi)$ are polar coordinates in the
Fourier plane) exhibits a simple behavior in the regime

$\displaystyle\frac{1}{d}\ll u\ll\frac{1}{w},$ (20)

where it adopts the “universal” form, valid to leading order in an expansion
in $0<w/d\ll 1$ (Gralla, 2020),

$\displaystyle V(u,\varphi)\approx\frac{e^{-2\pi
iC_{\varphi}u}}{\sqrt{u}}\left[\alpha_{L}(\varphi)e^{-\frac{i\pi}{4}}e^{i\pi
d_{\varphi}u}+\alpha_{R}(\varphi)e^{\frac{i\pi}{4}}e^{-i\pi
d_{\varphi}u}\right],$ (21)

with $\alpha_{L,R}(\varphi)=\alpha_{R,L}(\varphi+\pi)>0$ encoding the
intensity profile of the curve and ensuring that
$V(u,\varphi+\pi)=V^{*}(u,\varphi)$, as required by definition [Eq. (29)
below]. In the context of interferometry, the complex function $V(u,\varphi)$
is known as the radio visibility, with the visibility plane $(u,\varphi)$
usually referred to as the baseline plane. Intuitively, short baselines $u\ll
1/d$ are not sufficient to resolve the shape of a bright curve, which then
just appears as a blob, while on very long baselines $u\gg 1/w$, the smooth
profile of the curve has been resolved and its signal therefore dies off
exponentially; on the other hand, in the regime (20), the curve appears
infinitely thin and as a result, its signal follows a very weak $u^{-1/2}$
power-law fall-off (21).

In the context of realistic high-frequency VLBI observations, it is in
practice significantly easier to measure the amplitude $|V|$ of the complex
visibility—the “visamp”—than its absolute phase, as the latter is more
susceptible to contamination by a variety of atmospheric and instrumental
effects. While the full complex visibility (21) encodes both $d_{\varphi}$ and
$C_{\varphi}$, the visibility amplitude depends on $d_{\varphi}$ only, and
obeys $|V(u,\varphi+\pi)|=|V(u,\varphi)|$:

$\displaystyle\left|V(u,\varphi)\right|$
$\displaystyle\approx\sqrt{\frac{\alpha_{L}^{2}(\varphi)+\alpha_{R}^{2}(\varphi)+2\alpha_{L}(\varphi)\alpha_{R}(\varphi)\sin\left(2\pi
d_{\varphi}u\right)}{u}}.$ (22)

To summarize, the clearest interferometric signature of a bright curve in the
sky is its visibility amplitude (22), which is strong in the universal regime
(20) and only encodes partial information about the shape, namely its
projected diameter $d_{\varphi}$ (but not the full shape information, which
also requires the centroid motion $C_{\varphi}$ that is only available in the
harder to measure visibility phase).

Figure 2:  Example of a visibility amplitude profile at fixed angle
$\varphi=0^{\circ}$ in the baseline plane, obtained from a BH image simulated
with Gyoto and exhibiting the clean ringing expected from the narrow photon
ring. Three regimes (23) are visible: in the range $[0,20]$ G$\lambda$, the
direct ($n=0$) image dominates; for $u\in[20,120]$ G$\lambda$, the first
($n=1$) subring prevails; beyond 120 G$\lambda$, the $n=2$ ring signature
takes over.

###  2.7 Interferometric signature of the photon ring

Combining the results of the previous sections, we are now in a position to
describe the interferometric signature produced by the photon ring. We expect
the photon subrings to be the only narrow features present in time-averaged BH
images, since other fine features (such as emission ropes and other flares)
should be transient and wash out after sufficient time-averaging. Thus, we
expect that the $n^{\text{th}}$ subring will dominate the visibility amplitude
in the regime (Johnson et al., 2020)

$\displaystyle\frac{1}{w_{n-1}}\ll u\ll\frac{1}{w_{n}},$ (23)

where, according to Eq. (15),
$w_{n}(\tilde{r})=e^{-\gamma(\tilde{r})}w_{n-1}(\tilde{r})$. Beyond, its
signal starts to die off exponentially and the $(n+1)^{\text{th}}$ subring
begins to dominate. This results in a distinctive cascade structure
illustrated in Fig. 2. In each range (23), the $n^{\text{th}}$ subring
produces a damped oscillation (22) governed by its projected diameter
$d_{\varphi}^{(n)}$. All that remains for us to do is to describe the
diameters $d_{\varphi}^{(n)}$, beginning with the projected diameter
$\tilde{d}_{\varphi}\equiv d_{\varphi}^{(\infty)}$ of the critical curve,
which they exponentially converge to in the limit $n\to\infty$.

###  2.8 Theoretical shape of the critical curve

A “phoval” (for “photon ring oval”) is a geometric shape with

$\displaystyle\frac{d_{\varphi}}{2}$
$\displaystyle=R_{0}+\sqrt{R_{1}^{2}\sin^{2}{\varphi}+R_{2}^{2}\cos^{2}{\varphi}},$
(24a) $\displaystyle C_{\varphi}$
$\displaystyle=\left(X-\chi\right)\cos{\varphi}+\arcsin\left(\chi\cos{\varphi}\right).$
(24b)

More precisely, this is a family of shapes parameterized by five parameters
$R_{0}$, $R_{1}$, $R_{2}$, $X$, and $\chi$, with the interesting property that
the Kerr critical curve is extremely well-approximated by some member of the
family for any choice of observer inclination and BH spin (Gralla & Lupsasca,
2020c).

We briefly summarize the geometric interpretation of these parameters. First,
the projected position $f_{\rm circ}(\varphi)=R_{0}$ describes a circle
centered at the origin with radius $R_{0}$. Next, the projected position
$f_{\rm
ellipse}^{2}(\varphi)=R_{1}^{2}\sin^{2}{\varphi}+R_{2}^{2}\cos^{2}{\varphi}$
describes an ellipse centered at the origin with radii $R_{1}$ and $R_{2}$.
Adding these two functions results in the projected position $f_{\rm
circlipse}(\varphi)=d_{\varphi}/2$ given in Eq. (24a), which defines a shape
known as a “circlipse.” It follows from the last paragraph that the projected
diameter of the critical curve is always exquisitely close to that of a
circlipse,

$\displaystyle\tilde{d}_{\varphi}\equiv d_{\varphi}^{(\infty)}\approx
2R_{0}+2\sqrt{R_{1}^{2}\sin^{2}{\varphi}+R_{2}^{2}\cos^{2}{\varphi}},$ (25)

with the best-fit parameters $R_{0}$, $R_{1}$, and $R_{2}$ depending in a
rather intricate way on the mass-to-distance ratio $M/D$, BH spin $a$, and
observer inclination $i$. Finally, adding the centroid motion (24b) to the
projected function produces the full phoval shape, which is translated along
the $\alpha$ axis (orthogonal to the spin axis projected onto the plane
perpendicular to the line of sight) by a horizontal offset $X\in\mathbb{R}$,
and carries a left-right asymmetry induced by a parameter $\chi\in[-1,1]$ that
matches the warping of the critical curve at high inclinations and high spins.
With these additional parameters, the best-fit phoval reproduces the Kerr
critical curve to a part in $10^{5}$ over the vast majority of parameter
space, and to a part in $10^{3}$ in the extremal limit $a\to 1$ (Gralla &
Lupsasca, 2020c). Away from extremality, the (parity-odd) centroid motion
gives a minor correction to the shape of the critical curve, which is almost
entirely encoded in its (parity-even) circlipse diameter.

A critical curve with best-fit circlipse parameters $R_{0}$, $R_{1}$, and
$R_{2}$ has diameters $d_{\perp}\equiv d_{0}$ and $d_{\parallel}\equiv
d_{\pi/2}$ along the directions perpendicular ($\alpha$ axis) and parallel
($\beta$ axis) to the projected spin

$\displaystyle d_{\perp}=2\left(R_{0}+R_{2}\right),\qquad
d_{\parallel}=2\left(R_{0}+R_{1}\right),$ (26)

which define a fractional asymmetry (Johnson et al., 2020)

$\displaystyle f_{A}\equiv
1-\frac{d_{\perp}}{d_{\parallel}}=\frac{R_{1}-R_{2}}{R_{0}+R_{1}}\in[0,1].$
(27)

The range $0\leq f_{A}\leq 1$ follows from the observation that $d_{\perp}\leq
d_{\parallel}$ over all parameter space, which also implies that $R_{2}\leq
R_{1}$. Moreover, we can see from Eq. (24a) that $d_{\perp}$ and
$d_{\parallel}$ are the minimal and maximal diameters of the circlipse,
respectively.

Finally, we point out that the two quantities $d_{\parallel}$ and $f_{A}$
encode the BH spin $a$ and observer inclination $i$, as illustrated in Fig. 7
of Johnson et al. (2020). The map $(a,i)\to(d_{\parallel},f_{A})$ is not
exactly bijective, however, but rather two-to-one, as $(a,i)$ and $(-a,i+\pi)$
both map to the same $(d_{\parallel},f_{A})$. Putting aside this discrete
twofold degeneracy (which can be broken by other means), this suggests that
one could in principle infer a BH’s spin and inclination from the projected
diameter $\tilde{d}_{\varphi}$ of its critical curve [Eq. (25)], or even just
its minimum and maximum $d_{\perp}$ and $d_{\parallel}$ (without any need for
the parameters $X$ and $\chi$ encoding its centroid).

However, this suggestion is misleading because the critical curve is purely
theoretical and not in itself observable. Nonetheless, the photon
subrings—which are observable—converge to it exponentially fast and moreover,
while the first subring is still noticeably different to the naked eye, the
second subring tracks it quite closely. It is therefore tempting to extract
the BH spin and inclination from the projected diameter $d_{\varphi}^{(2)}$.
Unfortunately, we have found this approach to be only moderately successful,
as the residual astrophysics-dependence of the $n=2$ ring shape introduces
significant uncertainties (see Sec. 5.3 below).

###  2.9 Testing GR with the photon ring shape

As we have explained, the projected diameter of the Kerr critical curve
depends only on the spacetime geometry: $\tilde{d}_{\varphi}$ encodes the BH
spin and inclination, and is fully determined by them. However, this
theoretical curve, and hence its diameter $\tilde{d}_{\varphi}$, are not
directly observable. Instead, what we may (in principle) observe is a sequence
of photon subrings labeled by half-orbit number $n$ and approaching the
critical curve as $n\to\infty$, with projected diameters $d_{\varphi}^{(n)}\to
d_{\varphi}^{(\infty)}=\tilde{d}_{\varphi}$. Thus, in the large-$n$ limit, the
rings enter a “universal” regime in which their dependence on astrophysical
conditions drops out and they converge to the GR-predicted, astrophysics-
independent shape $\tilde{d}_{\varphi}$. This suggests that a shape
measurement of a “large-$n$” subring could provide a test of GR and a means to
infer the BH parameters.

The first ($n=1$) subring is certainly not yet in the universal regime, as the
photons comprising it are not sufficiently bent and its shape still carries a
substantial imprint of the astrophysical details of the surrounding source. On
the other hand, the shape of the higher $n\geq 3$ subrings is much less
sensitive on astrophysics, but these rings are also much harder to detect
experimentally, for at least two reasons: first, their signature only comes to
dominate on exponentially long baselines [Eq. (23)], which (given present
observation frequencies) could only be accessed using a space element at
unrealistic separation from the Earth (e.g., Fig. 5 of Johnson et al., 2020);
second, this problem is compounded by the need for exceptional sensitivity, as
the flux in each subring is also exponentially suppressed in $n$ [Eq. (15)].
As such, the second subring seems to be the most promising target for a shape
test of GR, as it may sit in the “sweet spot” where it remains accessible with
current or near-future technology (Gralla et al., 2020), while $n=2$ is
“sufficiently large” that $d_{\varphi}^{(2)}\approx\tilde{d}_{\varphi}$.

Gralla et al. (2020) initiated a quantitative investigation of this idea and
found that in a wide class of toy models of M87* (subject to the limitations
listed in the introduction), the $n=2$ subring indeed produced the universal
visibility amplitude (22). Moreover, they showed that it is possible to
recover $d_{\varphi}^{(2)}$ from this signature, even in the presence of
simulated instrument noise. They found that, while the $n=2$ ring is still
noticeably different from its $n\to\infty$ limit (the critical curve), it
nonetheless adopts the same shape: more precisely, even though the $n=2$ ring
and the critical curve have different projected diameters
$d_{\varphi}^{(2)}\neq\tilde{d}_{\varphi}$, they both take the shape of a
circlipse. That is, $d_{\varphi}^{(2)}$ must also follow the functional form
(25).

We emphasize that, while the best-fit radii $R_{0}$, $R_{1}$, and $R_{2}$ for
$\tilde{d}_{\varphi}$ are determined only by the BH parameters $M/D$ and $a$,
as well as the observer inclination $i$, the best-fit radii $R_{0}$, $R_{1}$,
and $R_{2}$ for $d_{\varphi}^{(2)}$ additionally depend on astrophysical
conditions as well. This is illustrated in Fig. 7 of Gralla et al. (2020),
which shows how these best-fit parameters vary with the astrophysical model,
keeping the BH parameters fixed. Thus, the $n=2$ ring carries some
astrophysics-dependence (which makes it harder to extract the BH parameters
from its diameter), but it is weak enough that $d_{\varphi}^{(2)}$ does not
break out of the functional form

$\displaystyle\frac{d_{\varphi}^{(2)}}{2}\approx
R_{0}+\sqrt{R_{1}^{2}\sin^{2}\left(\varphi-\varphi_{0}\right)+R_{2}^{2}\cos^{2}\left(\varphi-\varphi_{0}\right)},$
(28)

which includes an additional parameter $\varphi_{0}$ to account for the fact
that the orientation of the image is generically not aligned with the Bardeen
coordinate system. This observation led Gralla et al. (2020) to propose a
novel GR test based on the $n=2$ ring shape: a space-based VLBI mission
targeting M87* could extract $d_{\varphi}^{(2)}$ and check to what extent it
follows the form (28).

If the deviation of $d_{\varphi}^{(2)}$ from its best-fit circlipse is large,
then the Kerr hypothesis fails the test; if it is small, then one can report a
test of the Kerr hypothesis and GR with a precision given by the root-mean-
square deviation of $d_{\varphi}^{(2)}$ from its best fit. With their
simulated experimental forecast, Gralla et al. (2020) attained a stringent,
sub-percent (0.04%) level of precision.

Finally, we wish to clarify in what sense we regard the GLM method as a test
of the Kerr hypothesis. Such a test can either be:

  * •

a consistency test aiming to confirm that some observable is consistent with
the prediction of the Kerr spacetime, in which case there is no notion of
model comparison; or,

  * •

a model-comparison test aiming to compare predictions of the Kerr hypothesis
with alternative predictions derived from other spacetime geometries, so as to
determine which one best explains the data (in the Bayesian sense).

In this article, we discuss the GLM method as a consistency test of the Kerr
hypothesis and thus do not consider any other spacetime but that of Kerr.
Performing an unambiguous model-comparison test of the Kerr hypothesis is a
difficult challenge, as one must be careful to properly take into account all
complex astrophysical degeneracies (e.g., Bauer et al., 2022).

Additional work is needed to establish the viability of this GLM method for a
wider range of configurations. In this paper, we explore whether the method
remains viable when some of the limitations of the original GLM analysis are
removed; this will be the object of our parameter survey in Secs. 5 and 6
below.

##  3 Implementation

To check the robustness and precision of the Kerr hypothesis test reviewed in
Sec. 2, we simulated high-resolution images of Kerr BHs using various types of
emission and confirmed that they all display a photon ring with $n=1$ and
$n=2$ subrings. Then, we inferred the $n=2$ ring diameter $d_{\varphi}^{(2)}$
from its characteristic damped oscillation (22) in the visibility amplitude on
long baselines, and verified that it follows the predicted functional form
(28), thereby establishing the viability of the GLM method as a test of the
Kerr hypothesis. We also studied the discrepancy between the $n=2$ subring and
the critical curve, and checked that a phoval with the minimal and maximal
diameters measured from the $n=2$ ring could fit within the $n=2$ lensing
band. Our simulations were performed with the relativistic ray tracing code
Gyoto (General relativitY Orbit Tracer of the Observatoire de Paris), whose
details are presented in Vincent et al. (2011).

###  3.1 Image simulation

Gyoto simulates an image, that is, a map $I(\rho,\varphi_{\rho})$ of the
specific intensity in the observer sky. To reproduce the complex visibility
$V(u,\varphi)$ that is directly sampled via VLBI observations, we had to
additionally compute its 2D Fourier transform,

$\displaystyle
V(u,\varphi)=\int_{0}^{\infty}\int_{0}^{2\pi}I(\rho,\varphi_{\rho}+\varphi)e^{-2\pi
iu\rho\cos{\varphi_{\rho}}}\rho\mathop{}\\!\mathrm{d}\rho\mathop{}\\!\mathrm{d}\varphi_{\rho},$
(29)

which by definition satisfies $V(u,\varphi+\pi)=V^{*}(u,\varphi)$. Instead of
computing this 2D FT, we followed Gralla et al. (2020) and made use of the
projection-slice theorem to directly compute $V(u,\varphi)$ along slices of
fixed angle $\varphi$ in the Fourier plane: for each $\varphi$, we computed
the Radon projection (i.e., the integrals of $I(\rho,\varphi)$ along lines
perpendicular to the slice of constant $\varphi$ across the image), and then
applied a 1D fast Fourier transform (FFT) to obtain the visibility
$V(u,\varphi)$ evaluated at that angle $\varphi$.

We made the following assumptions in our computations:

  * •

BH mass: $M=6.2\times 10^{9}M_{\odot}$, in the range of values for M87*
favored by stellar dynamical measurements (e.g., Gebhardt et al., 2011; Event
Horizon Telescope Collaboration, 2019f);

  * •

distance between the BH and observer: $D=16.9$ Mpc, in the range of values for
M87* (e.g., Event Horizon Telescope Collaboration, 2019f);

  * •

observation wavelength: $\lambda=1.3$ mm, matching the 230 GHz frequency of
current observations (Event Horizon Telescope Collaboration, 2019b);

  * •

we neglected absorption effects, a suitable approximation for optically thin
accretion flows, as discussed in Sec. 2 above;

  * •

we used Bardeen’s coordinates (7) in the sky (such that the projected spin
axis points in the direction $\varphi_{\rho}=90^{\circ}$);

  * •

we terminated the ray tracing of geodesics after their third equatorial
crossing, so as to avoid parasitic pixels produced by an under-resolved,
partially imaged $n=3$ ring.

From these values of $M$ and $D$, one obtains a conversion ratio from
microarcseconds to units of $M$ consistent with EHT priors (Event Horizon
Telescope Collaboration, 2019f):

$\displaystyle(M/D)_{\rm M87^{*}}=3.6212\,\mu{\rm as}.$ (30)

We took our field of view to be of 180 $\mu{\rm as}\approx 50M$ and ray traced
most images with a resolution of $10000\times 10000$ pixels, though to resolve
the $n=2$ ring, we found it necessary in some cases (when the ring is highly
peaked) to double the linear resolution.

###  3.2 Image treatment and analysis

####  3.2.1 Apodization

We find that the visibility amplitude profile $|V(u,\varphi)|$ obtained from
the FFT of a Radon slice indeed presents (on sufficiently long baselines) the
expected periodicity linked to the $n=2$ ring. Without further treatment,
however, this ringing signature does not quite take the “clean” form (22)
because it is “polluted” by some other oscillation with higher period and
lesser amplitude.

This effect arises from an artificial discontinuity in the Radon projection,
which is set to zero outside our field of view but does not exactly vanish at
its edges, where it instead decays to a tiny (but still nonzero) intensity
that is $\lesssim 10^{-4}$ times the maximum intensity. Though small, this
step function nonetheless produces its own ringing on long baselines, with
periodicity set by the size of the field of view, and this interferes with our
signal.

To eliminate this spurious artifact, one solution would be to increase the
field of view until the intensity at the edges is so small that this effect
vanishes. However, in order to still resolve the $n=2$ ring, these larger
images would have to be ray traced at a correspondingly higher resolution.
This would be too costly, especially since the bulk of the computation would
be dedicated to exterior pixels of very little interest to us. An ingenious
way to overcome this difficulty would be to use adaptative ray tracing, either
by employing algorithms designed for this purpose (e.g., Gelles et al., 2021),
or else by decomposing the total image into a sum of multiple layers, each
with its own resolution. The latter approach was adopted by Gralla et al.
(2020), who divided the image into layers labeled by half-orbit number $n$:
since the $n^{\text{th}}$ layer only has nonzero pixels in the $n^{\text{th}}$
lensing band, which is exponentially small, they were able to exponentially
increase the resolution in each layer while keeping fixed the number of pixels
computed in each layer.

Instead of computing images with an unnecessarily large field of view or using
adaptive ray tracing, we dealt with this effect by simply multiplying the
Radon projection by a window function. This method is known as apodization.

More specifically, we used the window function (also known as an apodization
or tapering function)

$\displaystyle W(x)=p(r(1-x))\,p(r(1+x)),$ (31)

where $x\in[-1,1]$ ranges over the Radon projection, $\epsilon=1/r$ is the
range of the cutoff, and $p$ is a “smooth plateau” function,

$\displaystyle p(x)=\frac{f(x)}{f(x)+f(1-x)},\qquad
f(x)=\left\\{\begin{array}[]{ll}e^{-1/x^{2}}&\mbox{if }x>0,\\\ 0&\mbox{if
}x\leq 0.\end{array}\right.$ (34)

This function is infinitely smooth ($C^{\infty}$), exactly equal to unity for
$x\in[-1+\epsilon,1-\epsilon]$, and exactly vanishes at the edges $x=\pm 1$.
As a result, multiplying the Radon transform with it smoothens the FFT while
retaining the periodicity that we are interested in.

####  3.2.2 Visibility amplitude profile fit

Once the visibility amplitude of a given image has been obtained, the next
step in the GLM test of the Kerr hypothesis is to extract the diameter
$d_{\varphi}^{(2)}$ of the $n=2$ ring at every angle $\varphi_{\rho}$ around
the image. This diameter can be inferred from the ringing signature displayed
by $|V(u,\varphi)|$ at the corresponding angle $\varphi=\varphi_{\rho}$ in the
Fourier plane, which is predicted to take the universal form (22) in the
baseline regime (23) appropriate to the $n=2$ ring. Naïvely, one would simply
determine the best-fit parameters $\alpha_{L}$, $\alpha_{R}$, and
$d_{\varphi}$, and take the latter to provide the desired measurement of
$d_{\varphi}^{(2)}$.

In practice, however, a simple fitting method—such as the FindFit routine
implemented in Mathematica—works only when the signal is extremely “clean,” in
the sense that $|V(u,\varphi)|$ displays no artifacts and follows the precise
functional form (22). This happens in ideal cases where the image perfectly
resolves the photon ring and the baselines chosen for the fit correspond
exactly to the domain (23) over which the contributions from the $n=2$ subring
are dominant. Unfortunately, these cases are quite rare and the signal is
typically less clean. As a result, the naïve fitting methods fail, even if—as
in Fig. 3—the periodicity is still manifestly present. A more robust approach
is clearly called for.

Markov chain Monte-Carlo (MCMC) methods offer such an approach and prove
effective even when the signal is buried in noise. For instance, the data
simulated by Gralla et al. (2020) for their experimental forecast contained so
much noise so as to render the periodicity invisible to the eye (see panels
(a) & (c) in their Fig. 8), and yet an MCMC analysis enabled the true ring
diameter to be extracted with great precision (see panels (b) & (d) in their
Fig. 8). While extremely robust, MCMC methods are also computationally
intensive and thus a “middle-ground” approach is more desirable, especially
for models without noise.

We devised such a procedure, which proved highly effective. We now describe it
using Fig. 3 as an illustrative example:

  1. 1.

For each $\varphi$, choose a baseline window in which to perform the fit (in
Fig. 3, we chose the range $u\in[2500,2600]$ G$\lambda$).

  2. 2.

Identify the local extrema (both maxima and minima) of the visibility
amplitude profile within this window.

  3. 3.

Interpolate between maxima to obtain the superior envelope $e_{\rm max}(u)$ of
the amplitude (green dashed curve in Fig. 3).

  4. 4.

Interpolate between minima to obtain the inferior envelope $e_{\rm min}(u)$ of
the amplitude (blue dashed curve in Fig. 3).

  5. 5.

Using a simple fitting routine, obtain the ring diameter $d_{\varphi}$ as the
parameter $d$ for which the refined model

$\displaystyle V_{\rm
fit}(u;d)=\sqrt{\alpha_{L}^{2}(u)+\alpha_{R}^{2}(u)+2\alpha_{L}(u)\alpha_{R}(u)\sin\left(2\pi
du\right)}$ (35)

best fits the data, where the functions $\alpha_{L/R}(u)$ are defined as

$\displaystyle\alpha_{L/R}(u)=\frac{e_{\rm max}(u)\pm e_{min}(u)}{2}.$ (36)

The precise interpolation method used to obtain the envelopes $e_{\rm
max/min}(u)$ is not important; we resorted to cubic interpolation with the
interp1d routine provided with SciPy. Likewise, the particular fitting method
used in the last step is not important either, as the model is by construction
very close to the data; we made use of the curve_fit routine, also provided
with SciPy.

Some words of explanation are in order. First, the use of the refined model
(35) allows us to consider signals that display the same periodic behavior as
that predicted by the formula (22), but with variable extrema $e_{\rm min}$
and $e_{\rm max}$. This relaxation significantly improves the method’s
robustness and enables it to determine the periodicity of even moderately
resolved rings, which can display non-monotonic modulations in their
envelopes. We thus expect the method to remain effective even at lower
resolution, thereby enabling a noticeable speed-up of parameter surveys. A
trade-off of the method is that, in order to attain sufficient precision on
the determination of the periodicity, we are required to examine the signal
over a relatively larger interval of baselines; in practice, we used a window
of 100 G$\lambda$, corresponding to $\sim 15$ periods.

This method also derives significant power from the flexible fall-off rate
built into the visamp model (35), which generalizes the $u^{-1/2}$ power-law
fall-off exhibited by the analytic formula (22) to arbitrary damping rates.
Such a step was already taken by Johnson et al. (2020), who noticed that their
signal presented damped oscillations with envelope $e_{\rm
max}(u)=u^{-1/2}e^{-(uw)^{\zeta}}$, with two additional parameters $w$ and
$\zeta$ needed to obtain a good fit (see their Fig. 4). The present approach
is more general still. The use of an envelope tailored to the oscillation
damping provides more robustness and allows one to ignore other ringing
patterns that could arise from pixel effects, or even beating produced by
interactions with other signals of larger periodicity: as shown in Fig. 3,
this method works even when such beats produce a local increase of the visamp.
Its drawback is that astrophysical noise (ignored here) may result in ill-
defined envelope functions (36).

We also wish to emphasize that Eq. (21) is only the leading term in an
expansion in powers of $0<w/d\ll 1$. Interestingly, preliminary analysis
suggests that subleading corrections set by the angle-dependent, nonzero
thickness of the ring $w(\varphi)$ could also predict a beating pattern;
analytic study of these corrections will be the subject of future work.

Finally, we note that while the model (35) is symmetric under
$L\leftrightarrow R$ interchange, this symmetry is broken by the definitions
(36), according to which $\alpha_{L}$ measures the oscillations’ mean and
$\alpha_{R}$ their amplitude, with $\alpha_{R}<\alpha_{L}$ by convention.

Figure 3:  Example of a visibility amplitude profile with its envelope and
best-fit model curve. Here, the naïve fitting method would fail because the
envelope has complex, non-monotonic variations; yet, the periodicity is
clearly visible and amenable to extraction via our refined fitting method.
Figure 4: Example of a model requiring a multi-peak fitting technique. (Left)
Fit comparison: naïvely fitting the model $V_{\rm fit}(u;d)$ [Eq. (35)] to the
visamp profile $|V(u)|$ (solid red curve) yields an excellent but suboptimal
fit (dashed yellow curve), since a better one exists (dashed black curve).
(Right) Goodness of fit $g(d_{\varphi})$ as a function of fitting diameter
[Eq. (38)]. The absolute best-fit corresponds to the global maximum of
$g(d_{\varphi})$, while other local maxima correspond to good but subopitmal
fits. These peaks are approximately periodic with a separation $\Delta
d_{\varphi}\approx 0.2\,\mu$as that exactly matches $1/u_{w}$, since our
baseline window is at $u_{w}\sim 1000\,$G$\lambda$.

####  3.2.3 Multi-peaked distribution for the diameter

The above fitting procedure returns a “best-fit” diameter $d_{\varphi}$ that
minimizes the normalized root-mean-square deviation (RMSD)

$\displaystyle\operatorname{RMSD}_{u}(d)=\frac{\sqrt{\left\langle\left[V_{\rm
fit}(u;d)-\left|V(u)\right|\right]^{2}\right\rangle_{u}}}{\left\langle V_{\rm
fit}(u;d)\right\rangle_{u}},$ (37)

where $\left\langle\cdot\right\rangle_{u}$ denotes an average over the chosen
baseline window.

Often, this $d_{\varphi}$ is only a local (rather than global) minimum of
$\operatorname{RMSD}_{u}(d)$, in which case it is not the absolute best-fit
diameter; even so, it still gives a good fit so long as
$0<\operatorname{RMSD}_{u}(d_{\varphi})\ll 1$. This is illustrated with the
example in Fig. 4, for which the above procedure yields a diameter
$d_{\varphi}=37.312\,\mu$as that is a local (but not global) minimum of Eq.
(37). Nevertheless, the resulting model $V_{\rm fit}(u;d_{\varphi})$ (dashed
yellow curve) closely tracks the visamp $|V(u)|$ (solid red curve), as
indicated by their small deviation (shown in the left panel inset) and as
measured by the correspondingly small numerical value of
$\operatorname{RMSD}_{u}(d_{\varphi})$.

Still, the global minimum of $\operatorname{RMSD}_{u}(d)$ is
$d_{\varphi}=37.115\,\mu$as, and with this diameter the model $V_{\rm
fit}(u;d_{\varphi})$ (dashed black curve) provides an even better fit to
$|V(u)|$: their deviation (left inset) is even narrower, as measured by the
smaller value of $\operatorname{RMSD}_{u}(d_{\varphi})$.

By definition, the global minimum of $\operatorname{RMSD}_{u}(d)$ is always
the absolute best-fit diameter. Equivalently, it is also determined as the
global maximum of the “goodness-of-fit” measure

$\displaystyle g(d_{\varphi})\equiv
e^{-\operatorname{RMSD}_{u}(d_{\varphi})}.$ (38)

We show $g(d_{\varphi})$ for the above example in the right panel of Fig. 4,
where we recognize the tallest peak to be the global maximum
$d_{\varphi}=37.115\,\mu$as, with the next local maximum corresponding to the
naïve best-fit parameter $d_{\varphi}=37.312\,\mu$as. Surprisingly, we also
observe an entire periodic sequence of local maxima (with several providing a
good fit) separated by a gap of $\Delta d_{\varphi}\approx 0.2\,\mu$as.

This multi-peak structure was already observed by Gralla et al. (2020) in
their experimental forecast. Its physical origin can be intuitively understood
as follows.

Over a sufficiently narrow baseline window, one may view $u\approx u_{w}\gg
1/d$ as fixed. Within such a window, the universal form (22) of the visibility
amplitude is periodic in the diameter, as $|V(u)|$ is approximately invariant
under shifts $d_{\varphi}\to d_{\varphi}+k/u_{w}$ for integer
$k\in\mathbb{Z}$. As a result, the probability distribution for $d_{\varphi}$
is multi-peaked, with peaks separated roughly by $\Delta d_{\varphi}\approx
1/u_{w}$. The example in Fig. 4 confirms this argument: the observed gap of
$\Delta d_{\varphi}\approx 0.2\,\mu$as between peaks corresponds to a baseline
length

$\displaystyle u_{w}\approx\frac{1}{2\times
10^{-7}\,\text{arcsec}}\approx\frac{1.03\times 10^{12}}{\text{rad}}\approx
1030\,\text{G}\lambda,$ (39)

which, as expected, lies well within the chosen baseline window of
$u\in[1000,1045]\,$G$\lambda$.

Another intuitive way to understand this degeneracy comes from the observation
that the number of damped oscillations (or “hops”) from the origin $u=0$ to
the baseline window $u\approx u_{w}$ is

$\displaystyle N_{w}\approx u_{w}d_{\varphi}.$ (40)

Hence, shifts $d_{\varphi}\to d_{\varphi}+k/u_{w}$ are equivalent to $N_{w}\to
N_{w}+k$, and this degeneracy in $d_{\varphi}$ implies that the number $N_{w}$
of hops from the origin cannot be measured with perfect precision using only
the visibility amplitude on very long baselines $u\approx u_{w}$: adding or
substracting a few units to this period number $N_{w}$ would barely shift the
extrema within our fitting window, and as a result the fit $V_{\rm
fit}(u,d_{\varphi}+k/u_{w})$ remains acceptable for several values of
$k\in\mathbb{Z}$.

In practice, we accounted for this degeneracy by adjuting our fitting
procedure to keep track of several of the most prominent peaks of the multi-
peak distribution, as follows:

  1. 1.

For each $\varphi$, obtain a first estimate of $d_{\varphi}$ by computing the
mean distance between local maxima of the visamp $|V(u)|$.

  2. 2.

Pick an interval of some width $\Delta d_{\varphi}\gtrsim 10/u_{w}$ centered
around this estimate and plot the goodness of fit $g(d_{\varphi})$ defined by
Eqs. (37)–(38), as in the right panel of Fig. 4.

  3. 3.

Find the peaks of $g(d_{\varphi})$ in this interval and keep the values of
$d_{\varphi}$ corresponding to local maxima above a chosen peak height,
together with the associated goodness of fit $g(d_{\varphi})$.

The resulting list should contain several diameters $d_{\varphi}$ that locally
maximize $g(d_{\varphi})$, including the global maximum. It is important to
note, however, that the latter is not always the “true” diameter $d_{\varphi}$
of the ring as measured from its image, which can sometimes correspond to one
of the smaller peaks in $g(d_{\varphi})$. There are also cases where the
multi-peaked distribution $g(d_{\varphi})$ is relatively flat, so we do not
obtain a good measurement at that angle $\varphi$. Despite such failures, one
can usually still infer the true diameter $d_{\varphi}$ by carrying out this
analysis at multiple angles, as we next explain.

####  3.2.4 Circlipse fit and test of the Kerr hypothesis

Figure 5:  Example of a successful circlipse fit and Kerr hypothesis test.
(Left) Every $5^{\circ}$, we determine the possible diameters $d_{\varphi}$ of
the $n=2$ ring by maximizing the goodness-of-fit measured $g(d_{\varphi})$
[Eq. (38)]. For each angle $\varphi$, we obtain a periodic set of diameters
separated by $\Delta d_{\varphi}\approx 1/u_{w}$, with likelihood proportional
to $g(d_{\varphi})$ and indicated by the coloring of the data point. (Right)
Circlipse fit for each of the possible rings. The darkest circlipse is most
likely to be the true shape of the image ring, since the corresponding “joint
goodness of fit” far exceeds that of the other solutions.

Having determined the absolute best-fit ring diameter $d_{\varphi}$ at every
angle $\varphi$ around the image, the final step of the GLM test is to fit it
to the GR-predicted functional form (28) of the ring,

$\displaystyle d_{\rm
fit}(\varphi;\@vec{R})=2R_{0}+2\sqrt{R_{1}^{2}\sin^{2}\left(\varphi-\varphi_{0}\right)+R_{2}^{2}\cos^{2}\left(\varphi-\varphi_{0}\right)},$
(41)

with parameters $\@vec{R}=\left\\{R_{0},R_{1},R_{2},\varphi_{0}\right\\}$. In
practice, we carried out this fitting by examining 36 angles in the range
$[0^{\circ},180^{\circ})$ spaced at regular intervals of $5^{\circ}$. That is,
we used $\varphi\in\left\\{0^{\circ},5^{\circ},\ldots,175^{\circ}\right\\}$.
The best-fit parameters $\@vec{R}_{\rm fit}$ are then obtained by minimizing

$\displaystyle\operatorname{RMSD}_{\varphi}(\@vec{R})=\frac{\sqrt{\left\langle\big{[}d_{\rm
fit}(\varphi;\@vec{R})-d_{\varphi}\big{]}^{2}\right\rangle_{\varphi}}}{\left\langle
d_{\rm fit}(\varphi;\@vec{R})\right\rangle_{\varphi}},$ (42)

where $\left\langle\cdot\right\rangle_{\varphi}$ denotes an average over the
chosen baseline angles.

If $0<\operatorname{RMSD}_{\varphi}(\@vec{R}_{\rm fit})\ll 1$, then the
circlipse fit is good and the Kerr hypothesis passes the GLM test at a level
of precision given by the normalized root-mean-square deviation
$\operatorname{RMSD}_{\varphi}(\@vec{R}_{\rm fit})$.444 Gralla et al. (2020)
use the width of the peak in $g(d_{\varphi})$ as a measure of the uncertainty
$\sigma_{\varphi}$ in the inferred diameter $d_{\varphi}$, enabling them to
use a standard chi-squared metric [their Eq. (15)] to assess goodness of fit.
Otherwise, the fit is poor and the Kerr hypothesis fails the test.

As explained in the last section, this naïve approach is not robust because of
the multi-peaked distribution for $d_{\varphi}$: keeping only the absolute
best-fit diameter $d_{\varphi}$ may yield a suboptimal fit, since it may force
us to use a wrong diameter at some angles $\varphi$.

However, we found a simple way to remedy this by using the data collected at
all angles simultaneously, as follows:

  1. 1.

For each angle $\varphi$, compute $g(d_{\varphi})$, identify several peak
values of $d_{\varphi}$, and display them all together as in Fig. 5.

  2. 2.

Separate these peaks into multiple circlipse-shaped subsets $\mathcal{C}_{i}$,
each with “joint goodness-of-fit” measure

$\displaystyle g(\mathcal{C}_{i})\equiv\prod_{\varphi}g(d_{\varphi}).$ (43)

  3. 3.

Finally, fit each $\mathcal{C}_{i}$ to a circlipse shape (41) with normalized
$\operatorname{RMSD}(\mathcal{C}_{i})\equiv\operatorname{RMSD}_{\varphi}(\@vec{R}_{\rm
fit})$ given by Eq. (42).

This “multi-fit method” results in several circlipses $\mathcal{C}_{i}$ that
are not equally likely, as measured by their different values of
$g(\mathcal{C}_{i})$. In favorable examples such as the one in Fig. 5, the
likeliest $\mathcal{C}_{i}$ (i.e., the one with maximal $g(\mathcal{C}_{i})$,
shown in dark blue in the right panel) does turn out to be the “true”
circlipse shape of the actual $n=2$ ring in the image.555However, we note that
even in such examples, the true circlipse does not minimize
$\operatorname{RMSD}(\mathcal{C}_{i})$, which tends to monotonically decrease
with the scale of $d_{\varphi}$ due to its normalization. That is, larger ring
diameters provide better fits to the circlipse shape according to the measure
(42).

That said, there is no reason for this to be the case in general, and in some
cases the “true” circlipse in the underlying image may well differ from the
likeliest circlipse $\mathcal{C}^{\rm max}$ as determined by this analysis. Of
course, there is no other way to infer the true circlipse in a real
experiment, so in principle we must take it to be the likeliest one
$\mathcal{C}^{\rm max}$, even if it does not provide the best fit to the
circlipse shape, that is, even if $\operatorname{RMSD}(\mathcal{C}^{\rm
max})\neq\min_{i}\operatorname{RMSD}(\mathcal{C}_{i})$.

In practice, the least favorable examples present a handful of circlipses
$\mathcal{C}_{i}^{\rm max}$ sharing similar near-maximal likelihoods
$g(\mathcal{C}_{i})$. In such cases, we may only infer the true ring diameter
up to a degeneracy of a few periods $\Delta d_{\varphi}\approx 1/u_{w}$, but
we may still report a test of the Kerr hypothesis at a level of precision
given by the maximal $\operatorname{RMSD}(\mathcal{C}_{i}^{\rm max})$ within
the set of best-fit circlipses $\mathcal{C}_{i}^{\rm max}$.

##  4 Importance of the choice of baseline window

Before reporting the results of our parameter surveys, we first describe some
new and important features of the GLM method that we noticed in the course of
our investigation. In order to determine the diameter $d_{\varphi}^{(2)}$ of
the $n=2$ ring, it is crucial to sample the visamp $|V(u,\varphi)|$ in the
appropriate regime (23), that is, on baselines long enough to resolve the
width of the $n=1$ ring (but not that of the $n=2$ ring) at image angle
$\varphi_{\rho}=\varphi$:

$\displaystyle\frac{1}{w_{1}(\varphi_{\rho})}\ll
u\ll\frac{1}{w_{2}(\varphi_{\rho})}.$ (44)

It is only in this regime that the signature of the $n=2$ ring can dominate
the visamp and produce damped periodic oscillations (22) that encode the ring
diameter $d_{\varphi}^{(2)}$. In the example of Fig. 2, this requires us to
choose a baseline window $u\gtrsim 100\,$G$\lambda$.

The angle-dependence of condition (44) is important. At low inclinations, this
dependence is relatively weak, and the baseline threshold
$b_{2}(\varphi)\approx 1/w_{1}(\varphi)$ past which the $n=2$ ring dominates
$|V(u,\varphi)|$ is approximately constant in $\varphi$. As a result, there
exist baseline windows for which a ring diameter measurement works for all
angles $\varphi$. On the other hand, at high inclinations, $b_{2}(\varphi)$
may vary so much that no choice of baseline can satisfy the condition (44) for
all angles simultaneously; worse, the visamp may even cease to encode the
$n=2$ ring diameter altogether.

These difficulties can be quantitatively fleshed out as follows. The GLM
method can only succeed for a fixed choice of baseline window $u\approx u_{w}$
that satisfies Eq. (44) for all angles, so that

$\displaystyle\max_{\varphi}b_{2}(\varphi)\ll
u_{w}\ll\min_{\varphi}b_{3}(\varphi).$ (45)

Here, $b_{n}(\varphi)\approx 1/w_{n-1}(\varphi)$ denotes the threshold past
which the visamp contribution from the $n^{\text{th}}$ ring overtakes its
predecessor’s. These thresholds are well-defined when there is a clear
transition between the baseline regimes in which each subring’s signature
dominates; however, these are not always cleanly delineated. To understand
why, we recall that $|V(u,\varphi+\pi)|=|V(u,\varphi)|$, which implies that
the baseline thresholds must obey (when they exist)

$\displaystyle b_{n}(\varphi+\pi)=b_{n}(\varphi).$ (46)

On the other hand, the subring widths $w_{n}(\varphi_{\rho})$ are only subject
to the periodicity constraint
$w_{n}(\varphi_{\rho}+2\pi)=w_{n}(\varphi_{\rho})$ and may display significant
variation with image angle $\varphi_{\rho}\in[0,2\pi)$. In particular,
$w_{1}(\varphi_{\rho})$ and $w_{1}(\varphi_{\rho}+\pi)$ may be very different
at high inclinations, in which case $b_{2}(\varphi_{\rho})$ fails to be
sharply defined. Worse, while it must always be the case that
$w_{2}(\varphi_{\rho})<w_{1}(\varphi_{\rho})$, it may happen that
$w_{2}(\varphi_{\rho})\sim w_{1}(\varphi_{\rho}+\pi)$ for some angles
$\varphi_{\rho}$, in which case the periodic ringing of the visamp
$|V(u,\varphi)|$ at $\varphi=\varphi_{\rho}$ is affected by both the $n=1$ and
$n=2$ rings but encodes neither’s diameter. We now describe these inclination-
dependent effects in detail.

###  4.1 Low to moderate inclinations

At low-to-moderate inclinations, each subring produces a distinct signature in
its own baseline regime, and so the threshold $b_{2}(\varphi)$ exists. This is
the case illustrated in Figs. 2 and 6. According to Eq. (15), the next
thresholds are related to it by

$\displaystyle b_{n+1}(\varphi)\approx e^{\gamma(\varphi)}b_{n}(\varphi),$
(47)

an approximate Kerr-lens equation that becomes exact as $n\to\infty$ and is
already an excellent approximation for $n\gtrsim 2$.

On the other hand, one cannot determine $b_{2}(\varphi)$ from first
principles, as the first subring’s width $w_{1}(\varphi_{\rho})$ must be
computed for each model separately. We may thus rewrite Eq. (45) as

$\displaystyle\max_{\varphi}b_{2}(\varphi)\ll
u_{w}\ll\min_{\varphi}e^{\gamma(\varphi)}b_{2}(\varphi),$ (48)

where $b_{2}(\varphi)$ is model-determined (astrophysics-dependent) and the
demagnification factor $e^{\gamma(\varphi)}$ is universal (GR-predicted).

Letting $\gamma_{0}=\min_{\varphi}\gamma(\varphi_{\rho})$ denote the minimal
demagnification factor, we see that a sufficient (though not necessary)
condition for the existence of a baseline window $u_{w}$ satisfying Eq. (48)
is

$\displaystyle\max_{\varphi}b_{2}(\varphi)\ll
e^{\gamma_{0}}\min_{\varphi}b_{2}(\varphi).$ (49)

At low inclinations, where $b_{2}(\varphi)\approx b_{0}$ is relatively flat
and $e^{\gamma_{0}}$ takes values of approximately $10\sim 20$, this
sufficient condition is manifestly satisfied, guaranteeing the existence of a
suitable choice of baseline window over which to carry the measurement. At
higher inclinations, however, $e^{\gamma_{0}}\sim 1$ (it is in fact exactly
unity for equatorial inclination and maximal spin $a=1$). Moreover,
$b_{2}(\varphi)$ varies significantly, so this condition is clearly violated.
This does not necessarily imply that Eq. (45) cannot be satisfied, but it does
mean that some choices of baseline window $u_{w}$ that are suitable at a given
angle $\varphi$ may not be suitable at other angles.

This effect is illustrated in Fig. 6, which displays a model from our
parameter survey with BH spin $a=0.5$ and moderate inclination $i=45^{\circ}$.
Panels B–D illustrate the angle-dependence of the threshold $b_{2}(\varphi)$,
which varies from $b_{2}(60^{\circ})\approx 200\,$G$\lambda$ to
$b_{2}(140^{\circ})\approx 700\,$G$\lambda$. If one were to extract a diameter
$d_{\varphi}$ from the periodicity of the visamp in the baseline window
$u_{w}\approx 600\,$G$\lambda$, then one would be measuring the $n=2$ ring
diameter $d_{\varphi}^{(2)}$ at $\varphi=60^{\circ}$ (panel B) and the $n=1$
ring diameter $d_{\varphi}^{(1)}$ at $\varphi=140^{\circ}$ (panel D). At an
intermediate angle of $\varphi=102^{\circ}$, one would land in a transition
region between the $n=1$ and $n=2$ regimes in which the visamp periodicity
encodes neither ring’s diameter, but rather some “superposition” of the two.

Of course, at moderate inclinations, this issue can always be avoided by
selecting a larger baseline window. For the example of Fig. 6, any choice
$u_{w}\gtrsim 1000\,$G$\lambda$ is sufficiently large: over such baseline
windows, the visamp periodicity always encodes the $n=2$ ring diameter (until
the $n=3$ signal takes over on even larger baselines $u\gtrsim
2000\,$G$\lambda$). Meanwhile, measurements of $d_{\varphi}$ on insufficiently
long baselines lead to tell-tale transitions between $d_{\varphi}^{(1)}$ and
$d_{\varphi}^{(2)}$ that inevitably result in the characteristic, non-periodic
pattern of data points shown in panel A of Fig. 6.

The reason is the following. For the angles in the blue zone,
$d_{\varphi}=d_{\varphi}^{(2)}$ is the “true” diameter of the $n=2$ ring’s
image, while for the angles in the gray zone, $d_{\varphi}=d_{\varphi}^{(1)}$
is the “true” diameter of the $n=1$ ring’s image. These regions are separated
by green zones in which neither ring dominates the signal, so that
$d_{\varphi}$ is equal to neither $d_{\varphi}^{(2)}$ nor $d_{\varphi}^{(1)}$:
in such regions, $d_{\varphi}$ does not track the diameter of any actual image
feature, and in particular it need no longer be a continuous curve. As a
result, the absolute best-fit diameter $d_{\varphi}$—defined to maximize Eq.
(38)—is allowed to “jump” across circlipse bands (and indeed, one can see two
such jumps in the second green zone $150^{\circ}<\varphi<175^{\circ}$ of panel
A). In fact, such jumps are necessary to ensure the periodicity of the
physical ring diameter $d_{\varphi}^{(2)}$ in the blue zone, which must
necessarily remain invariant under shifts $\varphi\to\varphi+\pi$.

On the other hand, the continuous curves $\mathcal{C}_{i}$ that the human eye
naturally joins points into need not be periodic as $\varphi\to\varphi+\pi$;
it is only the totality $\left\\{C_{i}\right\\}$ of such curves that needs to
respect that periodicity. It is therefore perfectly acceptable for these
curves to continuously wrap into their successor as one looks up across panel
A of Fig. 6. While these curves $C_{i}$ are still separated by a gap $\Delta
d_{\varphi}\approx 1/u_{w}$, they are manifestly not circlipses, and their
lack of periodicity is an indicator that the baseline window is inadequate for
a full measurement of $d_{\varphi}^{(2)}$.

To summarize, at low inclinations $i\lesssim 20^{\circ}$, the baseline regimes
corresponding to each subring are sharply delineated; in particular, the
baseline threshold $b_{2}(\varphi)$ exists and is relatively flat, so that if
a baseline window $u_{w}$ is suitable for a measurement of the $n=2$ ring
diameter $d_{\varphi}^{(2)}$ at some angle $\varphi$, then it is almost surely
suitable at every other angle around the baseline plane. At higher but still
moderate inclinations $i\lesssim 50^{\circ}$, the baseline threshold
$b_{2}(\varphi)$ still exists but now exhibits substantial variation. Because
of this, some baseline windows $u_{w}$ may be suitable for a measurement of
the $n=2$ ring diameter $d_{\varphi}^{(2)}$ at certain angles $\varphi$, but
not others, where they instead pick up the $n=1$ ring diameter
$d_{\varphi}^{(1)}$ or some unphysical “superposition” of the two.
Nonetheless, it should be possible to tell if this is happening from the shape
of $d_{\varphi}$, which would display the general pattern shown in panel A of
Fig. 6. Moreover, at these inclinations, this difficulty can always be avoided
by selecting a sufficiently large baseline window.

###  4.2 High inclinations

At higher inclinations $i\gtrsim 60^{\circ}$, subrings do not always produce
clearly delineated signatures in individual baseline regimes, and as a result
the baseline threshold $b_{2}(\varphi)$ does not exist at all angles
$\varphi$. Physically, this can occur because the baseline threshold must be
$\pi$-periodic [Eq. (46)], while the subring images need only be
$2\pi$-periodic. In particular, while
$w_{2}(\varphi_{\rho})<w_{1}(\varphi_{\rho})$ for every angle $\varphi_{\rho}$
around the image, at high inclinations it is usually the case that
$w_{2}(\varphi_{\rho})\sim w_{1}(\varphi_{\rho}+\pi)$ for angles
$\varphi_{\rho}$ opposite the spin axis. In that case, the visamp
$V(u,\varphi)$ at baseline angle $\varphi=\varphi_{\rho}$ (or equivalently,
$\varphi=\varphi_{\rho}+\pi$) is sensitive to the diameters of both the $n=1$
and $n=2$ subrings.

This phenomenon is clearly visible in the example shown in Fig. 7, which
corresponds to a BH of spin $a=0.5$ observed from a near-equatorial
inclination $i=89^{\circ}$. Since $w_{2}(\varphi_{\rho})\approx
w_{1}(\varphi_{\rho}+\pi)$ for $\varphi_{\rho}=90^{\circ}$, it follows that
the visamp profile $|V(u,90^{\circ})|$, which is plotted in red in Fig. 8,
cannot display a sharply delineation between regimes dominated by the $n=1$
and $n=2$ subrings. Nonetheless, this profile still displays oscillations.
According to Eq. (21), their periodicity ought to be set by the diameter of a
ring of fixed width in the image, which in this case must be the composite
ring consisting of the upper half of the $n=1$ ring joined together with the
lower half of the $n=2$ ring. In other words, we would expect the visamp for
this model to still adopt the universal form (22) on long baselines, but with
a projected diameter $d_{\varphi}$ corresponding to a hybrid $n=1\&2$ ring of
uniform width. In practice, this expectation is not quite realized in this
model because of the direct ($n=0$) image of the near-side of the equatorial
disk, which cuts across the photon ring. Since this horizontal line is itself
thin and bright like the hybrid $n=1\&2$ ring, the vertical intensity cross-
section consists of three spikes; as a result, the visamp profile at angles
$\varphi\approx 90^{\circ}$ is sensitive to multiple projected diameters,
including that between the upper half of the $n=1$ ring and the $n=0$ emission
line, and that between the $n=0$ emission line and the lower half of the $n=2$
ring. These diameters are roughly half that of the hybrid ring, and indeed the
visamp consists of a beating pattern with roughly twice the period expected
for the signature of the photon ring.

Figure 6:  Measurements of the photon ring shape carried out over
insufficiently long baseline windows display characteristic transition
effects. In this example, the BH spin is $a=0.5$, the observer inclination is
$i=45^{\circ}$, and the emission profile (52) has parameters $\mu=3r_{+}/2$,
$\gamma=-1$, and $\sigma=0.5M$. (A) Projected diameters $d_{\varphi}$ inferred
from the periodicity of the visamp in the range $[550,650]\,$G$\lambda$ via
the multi-fit method (Sec. 3). With this baseline choice, $d_{\varphi}$ does
not follow the circlipse shape (28) because it only tracks $d_{\varphi}^{(2)}$
within the blue zone, which is separated by green transition zones from the
gray zone where it corresponds to $d_{\varphi}^{(1)}$. The visamp profile
$|V(u,\varphi)|$ over this baseline window (shown as an orange band) is
dominated by the $n=2$ ring at $\varphi=60^{\circ}$ in the blue zone (B), by
the $n=1$ ring at $\varphi=140^{\circ}$ in the gray zone (D), and by neither
at $\varphi=100^{\circ}$ (C).  Figure 7:  (Left) Gyoto log-scale image of a BH
with spin $a=0.5$ and inclination $i=89^{\circ}$, surrounded by a thin disk
with emission profile (52) with parameters $\mu=r_{-}$, $\gamma=-2$, and
$\sigma=1.5M$. The bright horizontal line is the direct $n=0$ image of the
disk viewed almost edge-on. (Right) Zoom-in around a region with $\beta\approx
0$. The width of the $n=1$ subring for $\beta>0$ is of the same order of
magnitude as the width of the $n=2$ ring for $\beta<0$.  Figure 8:  Visibility
amplitude profiles $|V(u,\varphi)|$ of the image in Fig. 7. At
$\varphi=90^{\circ}$, the ringing signature is sensitive to the projected
diameter $d_{\varphi}$ of a hybrid ring composed of the upper half of the
$n=1$ ring and lower half of the $n=2$ ring, as well as to the distance
between this hybrid ring and the horizontal line of direct $n=0$ emission.
These visamps are not exact as they neglect emission from outside the field of
view of Fig. 7.

##  5 Parameter survey for the emission profile

###  5.1 Description of the models

In this paper, we restrict our attention to the “thin-disk” emission models
presented in Sec. III.A of Gralla et al. (2020) and Sec. 3.2 of Chael et al.
(2021), which consist of emission profiles that are stationary, axisymmetric,
and confined to the equatorial plane. The emitters forming the disk are
assumed to describe stable, circular-equatorial orbits all the way down to the
radius of the innermost stable circular orbit (ISCO), past which they follow
the prescription of Cunningham (1975) for infall to the horizon. This
specifies the redshift factor $g$ of light rays emitted from the disk, whose
observed (bolometric) intensity is then computed as

$\displaystyle
I(\alpha,\beta)=\sum_{n=0}^{N(\alpha,\beta)}f_{n}J\left(r_{s}^{(n)}\right)g^{4}\left(r_{s}^{(n)},\alpha\right),$
(50)

where $J(r)$ is an arbitrary radial emission profile,
$r_{s}^{(n)}=r_{s}^{(n)}(\alpha,\beta)$ denotes the radius at which a light
ray crosses the equatorial plane for the $n^{\text{th}}$ time after being shot
back into the geometry from position $(\alpha,\beta)$ in the observer sky, and
$N+1$ is the total number of times the light ray intersects the equatorial
plane.666Analytic formulas for $r_{s}^{(n)}$ and $N$ are derived by Gralla &
Lupsasca (2020a, b) and are succintly summarized in App. A.  Here, we also
introduced a “fudge factor”

$\displaystyle f_{n}=\begin{cases}1&\quad n=0,\\\ 1.5&\quad
n\in\left\\{1,2\right\\},\\\ 0&\quad n\geq 3,\end{cases}$ (51)

whose effect is to remove $n\geq 3$ subrings from the image (which would be
underresolved), while enhancing the intensity of the $n=1$ and $n=2$ subrings
relative to the direct $n=0$ image. This enhancement can account for the
effects of geometric thickness and results in images whose intensity cross-
sections match those observed in time-averaged GRMHD-simulated images (Johnson
et al., 2020; Chael et al., 2021). In this sense, the equatorial model (50)
provides a good approximation for realistic accretion flows, and its
simplicity and relatively low computational cost make it ideal for the survey
of emission profiles undertaken here.

Figure 9:  Impact of the emission profile on the observed ring thickness and
baseline domains. We show two extreme cases: (Left) A broad emission profile,
leading to large rings and an early transition between $n=1$ to $n=2$ regimes
in the visibility amplitude profile with $b_{2}(145^{\circ})\approx
200\,$G$\lambda$. (Right) A highly peaked emission profile, leading to thin
rings (the $n=2$ subring is not even visible in the image) and a late
transition between the $n=1$ and $n=2$ regimes with $b_{2}(145^{\circ})\gtrsim
1800\,$G$\lambda$.

While the equatorial source profile $J(r)$ can in principle take any form, in
this survey we considered only profiles of the form

$\displaystyle
J(r;\mu,\gamma,\sigma)=\frac{e^{-\frac{1}{2}\left[\gamma+\operatorname{arcsinh}\left(\frac{r-\mu}{\sigma}\right)\right]^{2}}}{\sqrt{\left(r-\mu\right)^{2}+\sigma^{2}}}.$
(52)

This is a three-parameter subset of Johnson’s $S_{U}$-distribution (the fourth
parameter $\delta$ is set to unity), already adopted by Gralla et al. (2020)
for its versatility: this functional form provides a simple way to generate
smooth profiles with emission typically concentrated near the event horizon of
the BH, and displaying substantial variation with changes in $\mu$, $\gamma$,
and $\sigma$. In short, these three parameters respectively determine a
location (the peak of the distribution), a shape (its asymmetry), and a scale
(its width). We now describe in detail how their variation affects ring width.

First, $\mu$ influences the location of the emission profile’s peak. We
considered values in the range $[r_{-},2r_{+}]$, where $r_{\pm}$ denote the
horizon radii given in Eq. (14). At the lower end of this range (left column
of Fig. 9), the peak is very close to or behind the event horizon, but it
moves farther out near the upper end of the range (right column of Fig. 9).
The observed thickness of the photon subrings is not significantly affected by
this variation.

Second, $\gamma$ controls the profile’s asymmetry. We considered values in the
range $[-2,2]$. When $\gamma<0$ (left column of Fig. 9), the steep end of the
distribution is to the left and (mostly) cut off by the horizon, leaving only
a gentle-sloping tail in the emission profile. When $\gamma>0$ (right column
of Fig. 9), the reverse is true: the profile’s steep end lies to the right and
remains visible outside the horizon, while its gentle-sloping side is cut off.
Hence, the rings tend to become thinner as $\gamma$ increases.

Third, $\sigma$ encodes the profile’s width. We considered values in the range
$[0.25M,1.5M]$. At the upper end of this range (left column of Fig. 9), the
profile is broad enough to produce thick rings, which makes the GLM test easy
to carry out on a baseline window $u\in\left[450,550\right]\,$G$\lambda$ for
all our considered values of $\mu$ and $\gamma$. At the lower end (right
column of Fig. 9), the profile is a steep spike (provided that the
distribution’s peak lies outside the horizon) that produces narrower rings: as
a result, the $n=1$ ring dominates on longer baselines, and a measurement of
the $n=2$ ring diameter may require (in the most extreme cases) a baseline
window $u_{w}\gtrsim 2000\,$G$\lambda$.

Images with thinner rings have a larger baseline threshold $b_{2}(\varphi)$
for carrying out a measurement of the $n=2$ ring diameter $d_{\varphi}^{(2)}$.
This trend is illustrated in Fig. 9 with two extreme examples that display
maximally thick rings (left column) and maximally thin rings (right column)
within our considered parameter range.

###  5.2 Choice of parameters for the survey

In order to assess the robustness of the GLM method described in Secs. 2 and
3, we conducted a survey over the emission-profile parameter space
$\left\\{\mu,\gamma,\sigma\right\\}$. We used “intermediate” values for the BH
spin and inclination of $a=0.5$ and $i=45^{\circ}$, and examined a total of
100 emission models corresponding to all the possible combinations of
parameters with values presented in Table 1.

Table 1: Parameter values considered in the emission profile survey.

Parameter | Considered values
---|---
$\mu$ | $\left\\{r_{-},r_{+}/2,r_{+},3r_{+}/2,2r_{+}\right\\}$
$\gamma$ | $\left\\{-2,-1,0,1,2\right\\}$
$\sigma/M$ | $\left\\{0.25,0.5,1.,1.5\right\\}$
777$r_{\pm}$ denote the radii (14) of the outer and inner horizons.

As a cross-check, we repeated the survey using the currently favored values of
BH spin and inclination for M87*. The power in the relativistic jet produced
by M87* suggests that this BH is highly spinning, since most low-spin GRMHD
models do not give rise to a sufficiently powerful jet (Event Horizon
Telescope Collaboration, 2019a). We therefore took $a=0.94$, the highest value
considered in EHT simulations (Event Horizon Telescope Collaboration, 2019e).
We also assumed that the spin axis of the BH and of its accretion disk are
aligned with its forward jet, which has a measured inclination of $i\approx
17^{\circ}$ (Walker et al., 2018). Since $r_{-}\approx r_{+}/2$ for this value
of BH spin, we examined the same combinations of parameters except for those
with $\mu=r_{+}/2$; this amounted to a total of 80 emission profiles in our
survey.

###  5.3 Results

We investigated these models to address the following questions:

  1. Q1.

Is the GLM test possible in principle? In particular, does $d_{\varphi}^{(2)}$
follow the GR-predicted functional form (41) of a circlipse?

  2. Q2.

If so, how do the best-fitting circlipse parameters $\@vec{R}_{\rm fit}$ vary
with the astrophysical profile? In particular, what can we learn about the BH
parameters from the photon ring shape?

Our survey led us to the following answers:

  1. A1.

Yes. In all the models considered here, we identified a regime in which the
$n=2$ ring dominates the visamp and extracted the projected ring diameter
$d_{\varphi}^{(2)}$ from the signal’s periodicity. For reasons discussed in
Sec. 4, this required using different baseline lengths, sometimes as low as
$u\gtrsim 150\,$G$\lambda$ (Fig. 2) or as high as $u\gtrsim 2000\,$G$\lambda$
(bottom right panel of Fig. 9). We always found an excellent fit to the
circlipse shape (41) as measured by the normalized RMSD (42): in all cases,

$\displaystyle\operatorname{RMSD}_{\varphi}(\@vec{R}_{\rm fit})<10^{-4}.$ (53)

  2. A2.

By analogy with Eq. (26), it is helpful to repackage the best-fit parameters
$\@vec{R}_{\rm fit}=\left\\{R_{0},R_{1},R_{2},\varphi_{0}\right\\}$ describing
the projected diameter of the $n=2$ ring into maximal/minimal diameters888 We
use the notation $(d_{+},d_{-})$ instead of $(d_{\parallel},d_{\perp})$ to
emphasize that these diameters describe the $n=2$ ring rather than the
critical curve. Another difference is that $d_{+}\neq d_{\pi/2}$ and
$d_{-}\neq d_{0}$ because of the rotation parameter $\varphi_{0}$ in Eq. (41),
though this angle is typically extremely small.

$\displaystyle
d_{+}=d_{\varphi_{0}+\pi/2}^{(2)}=2\left(R_{0}+R_{1}\right),\quad
d_{-}=d_{\varphi_{0}}^{(2)}=2\left(R_{0}+R_{2}\right)$ (54)

of the associated circlipse shape ($d_{-}\leq d_{+}$ since $R_{1}\geq R_{2}$
by definition). We plot the set of $(d_{+},d_{-})$ obtained from all our
emission models in Figs. 10 and 11, and provide mean values and standard
deviations for $(d_{+},d_{-})$ in Tables 2 and 3. These results give a good
overview of the variabilty in the circlipse shape obtained from different
astrophysical source profiles. We return to the question of spin inference in
Table 4 below.

Table 2:  Mean values and standard deviations for the inferred $n=2$ ring
diameters $d_{\pm}$ in the “intermediate” case $a=0.5$ and $i=45^{\circ}$.

| Mean value | Standard deviation
---|---|---
$d_{+}/M$ | 10.338 | 0.014
$d_{-}/M$ | 10.250 | 0.009
999The critical curve has $d_{\parallel}=10.319M$ and $d_{\perp}=10.236M$.

Table 3:  Mean value and standard deviations for the inferred $n=2$ ring
diameters $d_{\pm}$ for the “best guess for M87*” case $a=0.94$ and
$i=17^{\circ}$.

| Mean value | Standard deviation
---|---|---
$d_{+}/M$ | 9.849 | 0.022
$d_{-}/M$ | 9.751 | 0.022
101010The critical curve has $d_{\parallel}=9.834M$ and $d_{\perp}=9.733M$.

As discussed in Sec. 2.9, one expects the $n=2$ photon ring to closely track
the critical curve. This is verified in practice, though we find that the
choice of astrophysical source profile still has a significant impact on the
shape of the $n=2$ ring. In most of our models, the $n=2$ ring is slightly
larger in diameter than the critical curve, so that $d_{+}>d_{\parallel}$ and
$d_{-}>d_{\perp}$. In particular, we find that broader emission profiles
produce thicker $n=2$ rings and therefore a larger difference in these
diameters. Conversely, narrow profiles that produce thinner rings lead to
diameters that are closer to (and sometimes, even smaller than) those of the
critical curve. Either way, conflating the shape of the $n=2$ ring with that
of the critical curve can lead to erroneous conclusions.

Nevertheless, it is possible to approximate the critical curve diameters
$(d_{\parallel},d_{\perp})$ by $(d_{+},d_{-})$, and to then use these
parameters to estimate the BH spin and inclination (assuming a mass prior) as
described in Sec. 2.8. However, the precision of this spin and inclination
determination is significantly limited by the spread of results for the $n=2$
ring diameters at a given $(a,i)$ value. To illustrate this, Figs. 10 and 11
also display diameters of critical curves that are compatible with the
measured $n=2$ ring shape, but which correspond to $(a,i)$ values that are
different from the “true” values of the BH spin and inclination used in the
model (green and pink stars). The constraints on these BH parameters derived
from the inferred $n=2$ ring diameters have substantial error bars, for which
we provide lower bounds in Table 4.

Table 4:  Lower bounds on the uncertainties that we can expect from a measurement of BH spin and inclination using the $n=2$ ring diameters.  | Spin | Inclination
---|---|---
Intermediate case ($a=0.5$, $i=45^{\circ}$) | $\pm 0.07$ | $\pm 22^{\circ}$
M87* “best guess” ($a=0.94$, $i=17^{\circ}$) | $\pm 0.04$ | $\pm 1.5^{\circ}$
Figure 10:  Results of the survey of emission profiles in the “intermediate”
case with spin $a=0.5$ and inclination $i=45^{\circ}$.  Figure 11:  Results of
the survey of emission profiles in the best-guess case for M87* with spin
$a=0.94$ and inclination $i=17^{\circ}$.

###  5.4 Range of allowed phovals within the $n=2$ lensing band

These results highlight that the critical curve is not necessarily the most
relevant concept for the analysis of the $n=2$ photon ring shape. For our
thin-disk models, a more fruitful approach is to consider the $n=2$ lensing
band described in Sec. 2.5, since it contains all the possible geometric
shapes of the $n=2$ photon ring for a given BH spin and inclination. The
critical curve lies within this band as well, but so do the critical curves of
BHs with neighboring spins and inclinations, as Figs. 10 and 11 indicate.

Since Sec. 2.9, we have used “photon ring shape” to refer to the simplest
observational target for a measurement of the $n=2$ ring: its projected
diameter $d_{\varphi}^{(2)}$, which matches that of a circlipse [Eq. (24a)].
In this section, we study the full geometric shape of the $n=2$ ring, which
typically includes centroid motion effects [Eq. (24b)] and is therefore not a
circlipse in general. Instead, we expect this geometric shape to be close to a
phoval [Eq. (24)]; that is, the $n=2$ ring’s projected position should closely
match

$\displaystyle f_{\rm ph}(\varphi;\@vec{R})=R_{0}$
$\displaystyle+\sqrt{R_{1}^{2}\sin^{2}\left(\varphi-\varphi_{0}\right)+R_{2}^{2}\cos^{2}\left(\varphi-\varphi_{0}\right)}$
(55)
$\displaystyle+\left(X-\chi\right)\cos\left(\varphi-\varphi_{0}\right)+\arcsin\left[\chi\cos\left(\varphi-\varphi_{0}\right)\right],$

for some choice of parameters
$\@vec{R}=\left\\{R_{0},R_{1},R_{2},X,\chi,\varphi_{0}\right\\}$. (Here, as in
Eq. (28), we have included a rotation angle $\varphi_{0}$ to account for the
arbitrary orientation of the BH spin axis in the image.) This expectation is
based on the observation (reviewed in Sec. 2.8) that the Kerr critical curve
is always a phoval to high accuracy, together with the fact that the photon
rings converge very rapidly (exponentially fast in $n$) to the critical curve
[Eq. (13)].

Since we expect the geometric shape of the $n=2$ photon ring to be that of a
phoval, it is interesting to consider the range of all phoval shapes that can
possibly fit inside the $n=2$ lensing band. Since phovals are described by six
parameters, it is convenient to project this 6D region of allowed phoval
shapes down onto the 2D plane of diameters $(d_{+},d_{-})$ via the projection
(54). The resulting region in the $(d_{+},d_{-})$ plane gives a good sense of
the allowed range of possible $n=2$ ring shapes.

For the two choices of BH spin and inclination in our survey of emission
profiles, this projection results in the blue regions displayed in Figs. 10
and 11. Their numerical computation is nontrivial and we relegate the details
to App. B. We note that since $d_{-}\leq d_{+}$, these regions must lie below
the diagonal.111111We also note that these regions are in general different
from the set of $(d_{+},d_{-})$ corresponding to the circlipse shapes that fit
within the $n=2$ lensing band. These two notions only coincide at low
inclination, where the centroid motion is negligible and the allowed phovals
essentially reduce to circlipses—this explains the square shape of the blue
region in the $i=17^{\circ}$ case of Fig. 11. At higher inclinations, the
lensing band is no longer circlipse-shaped, so the set of circlipses it can
contain becomes an increasingly smaller subregion of the blue region of
allowed phovals.

As expected, all the emission profiles in our survey produce $n=2$ photon
rings that lie within the $n=2$ lensing band, with diameters $(d_{+},d_{-})$
lying within the blue region of allowed phoval shapes that fit inside the
band. The band’s outer and inner edges can themselves be fitted to the phoval
shape (55) with excellent agreement, and their diameters $(d_{+},d_{-})$—shown
in Figs. 10 and 11 as blue and orange dots, respectively—lie on opposite
corners of this blue region. Remarkably, all of our astrophysical source
profiles produce phovals whose diameters $(d_{+},d_{-})$ lie within a
surprisingly narrow subregion of the full allowed region.

More precisely, these diameters (the gold dots in Figs. 10 and 11) form a
relatively straight line connecting the corners of the blue region that
correspond to the outer and inner edges of the lensing band. In addition, they
are not uniformly scattered along this line; rather, they are tightly grouped
near its center. Both of these facts may be intuitively understood as follows.
We recall from Sec. 2.5 that the inner and outer edges of the $n^{\text{th}}$
lensing band are the $n^{\text{th}}$ lensed images of the equatorial circles
of radius $r=r_{+}$ (the event horizon) and $r\to\infty$, respectively. Images
of intermediate radii $r_{+}<r<\infty$ form curves that fill the lensing band
and interpolate between its inner and outer edges. From Fig. 6 of Gralla &
Lupsasca (2020a), these curves roughly appear to be dilated versions of each
other (with $\gamma(\varphi)$ controlling the angular-dependent dilation) and
indeed, we find that they map to a series of points in the $(d_{+},d_{-})$
plane tracing a diagonal line that connects the blue and orange dots (shown in
gray in Figs. 10 and 11): the same line that the gold dots are concentrated
around (see App. B for details of its computation).

Since an equatorial disk is a superposition of equatorial rings, the $n=2$
photon ring is also a superposition of the $n=2$ images of said rings. That
is, its position in the $(d_{+},d_{-})$ plane should be a weighted average of
points along the diagonal line, and must therefore also lie on this line.
Moreover, this position may be associated to an “effective” radius of
equatorial emission lying somewhere in the line’s middle, away from its
extremities. We intend to further explore this connection in future work.

##  6 Parameter survey for BH spin and inclination

###  6.1 Choice of parameters for the survey

To better assess the robustness of the GLM method described in Secs. 2 and 3,
we also conducted a parameter survey over BH spins and inclinations $(a,i)$.
Simulating 100 different emission models for each $(a,i)$, as we did in Sec.
5.2, would have been prohibitively expensive here; instead, we focused on a
handful of “representative” emission profiles of the form (52) that spanned
the full range of outcomes observed in that survey:

  * •

a “typical 1” profile for which the diameters $(d_{+},d_{-})$ of the $n=2$
ring tended to be near their average value (that is, near their mean value
averaged over all emission models);

  * •

an “overestimated” profile for which the diameters $(d_{+},d_{-})$ of the
$n=2$ ring were among the largest;

  * •

an “underestimated” profile for which the diameters $(d_{+},d_{-})$ of the
$n=2$ ring were among the smallest;

  * •

a “typical 2” profile for which they were slightly subaverage;

  * •

a “narrow” profile with a sharp peak (small $\sigma$) for which the circlipse
fit was more difficult (required longer baselines) and led to the smallest
values for $(d_{+},d_{-})$ when successful.

Table 5 lists the parameters $\left\\{\mu,\gamma,\sigma\right\\}$ of these
emission profiles.

Table 5:  Parameters of the five representative emission profiles selected for the survey over BH spins and inclinations.  Profile | $\mu$ | $\gamma$ | $\sigma/M$
---|---|---|---
Typical 1 | $r_{+}$ | $0$ | $1$
Overestimated | $r_{-}$ | $-2$ | $1.5$
Underestimated | $r_{+}$ | $2$ | $0.5$
Typical 2 | $r_{-}$ | $2$ | $1$
Narrow | $\frac{3}{2}r_{+}$ | $2$ | $0.25$

In this survey, we considered all the combinations $(a,i)$ with
$a\in\\{0,0.25,0.5,0.75,0.99\\}$ and
$i\in\\{1^{\circ},22^{\circ},45^{\circ},67^{\circ},89^{\circ}\\}$ (avoiding
the fine-tuned cases $0^{\circ}$ and $90^{\circ}$). In the Schwarzschild case
$a=0$, the inclination is defined relative to the equatorial disk, whose spin
is always aligned with the BH spin.

###  6.2 Results

We investigated these models to address the two questions posed in Sec. 5.3.
Our survey led us to the following answers:

  1. A1.

Yes, up to $i\lesssim 45^{\circ}$. In all the models considered here, we
identified a regime in which the $n=2$ ring dominates the visamp and extracted
the projected ring diameter $d_{\varphi}^{(2)}$ from the signal’s periodicity.
For reasons discussed in Sec. 4, the diameters inferred from images at
inclinations of $i=67^{\circ}$ and $i=89^{\circ}$ sometimes presented features
that did not allow for a full circlipse fit as described in Secs. 2 and 3.
However, for low-to-moderate inclinations $i\lesssim 45^{\circ}$, we always
found an excellent fit to the circlipse shape (41) as measured by the
normalized RMSD (42): in all such cases,

$\displaystyle\operatorname{RMSD}_{\varphi}(\@vec{R}_{\rm fit})<5\times
10^{-4}.$ (56)

For most configurations, these best fits were obtained in a window
$u\in\left[1000,1100\right]\,$G$\lambda$. However, the “narrow” profile
produced $n=2$ rings that were so thin that they only became dominant on much
longer baselines, requiring the fit to be carried out in a baseline window
$u_{w}\approx 2000\,$G$\lambda$. This also became necessary at $i=45^{\circ}$
for the “underestimated” profile.

  2. A2.

As in Sec. 5.3, we extracted a pair of diameters $(d_{+},d_{-})$ for each
emission model. This enabled us to produce, for each BH spin and inclination,
a plot like the one shown in Fig. 12 for the case of spin $a=0.5$ and
inclination $i=45^{\circ}$. Since our models were handpicked to form a
“representative” sample spanning the range of possible diameters
$(d_{+},d_{-})$ that one can obtain for each spin and inclination, it is
informative to draw a bounding box containing them (like the green box in Fig.
12) and then measure its size. The dimensions of these bounding boxes are
listed for each $(a,i)$ in Table 6 below. As we explain in Sec. 6.3, this
analysis bears relevance for the “inverse problem” of how to infer BH spin and
inclination from the $n=2$ photon ring shape.

Figure 12:  Parametrization of the region in the $(d_{+},d_{-})$ plane
containing all the astrophysical photon ring shapes (here, for $a=0.5$ and
$i=45^{\circ}$).  Table 6:  Triplets $(s_{\rm min},s_{\rm max},|\delta|_{\rm
max}/M)$ describing the region of the $(d_{+},d_{-})$ plane containing all the
models computed with Gyoto.

Spin | $i=1^{\circ}$ | $i=22^{\circ}$ | $i=45^{\circ}$
---|---|---|---
$a=0$ | $(0.23,0.59,3\times 10^{-3})$ | $(0.19,0.57,9\times 10^{-4})$ | $(0.18,0.56,5\times 10^{-3})$
$a=0.25$ | $(0.21,0.58,2\times 10^{-3})$ | $(0.19,0.58,5\times 10^{-4})$ | $(0.19,0.57,4\times 10^{-3})$
$a=0.5$ | $(0.23,0.60,5\times 10^{-3})$ | $(0.20,0.59,9\times 10^{-4})$ | $(0.19,0.58,6\times 10^{-3})$
$a=0.75$ | $(0.19,0.61,1\times 10^{-4})$ | $(0.20,0.60,2\times 10^{-3})$ | $(0.19,0.43,3\times 10^{-3})$
$a=0.99$ | $(0.23,0.73,4\times 10^{-4})$ | $(0.23,0.75,8\times 10^{-3})$ | $(0.06,0.94,2\times 10^{-3})$
121212$s_{min}$ and $s_{max}$ are affine parameters along the line connecting
the edges of the lensing band, while $\delta$ is a distance perpendicular to
this line.

We first note that this survey over BH spins and inclinations confirms one of
our crucial observations from Sec. 5.4: for each choice of $(a,i)$, all our
emission models produce $n=2$ rings with diameters $(d_{+},d_{-})$ that are
very tightly clustered near the line connecting the inner and outer edges of
the $n=2$ lensing band. In other words, the five observed ring shapes
corresponding to the models in Table 5 occupy a very small subregion (green
box in Fig. 12) of the full region of allowed phoval shapes (blue box). As
explained in Sec. 5.4, diameters along this line correspond to $n=2$ images of
equatorial circles of different radii $r_{s}\in[r_{+},\infty)$, so we may
equivalently say that the $n=2$ rings produced by our emission models have the
same diameters $(d_{+},d_{-})$ as images of equatorial rings of emission
concentrated near a fixed radius $r_{s}$.

Rather than pursue this association with “effective” emission radii, we
employed a more agnostic parameterization of this line in the $(d_{+},d_{-})$
plane and used an affine parameter $s$ such that $s=0$ corresponds to the
inner edge of the lensing band and $s=1$ to its outer edge. We also introduced
$\delta$, a perpendicular distance from the line, thereby defining a map
$(d_{+},d_{-})\leftrightarrow(s,\delta)$. Then, the diameters produced by our
models are all contained in a bounding box extending from some $s_{\rm min}$
to some $s_{\rm max}$ in the direction parallel to the line, with a
perpendicular width of $|\delta|_{\rm max}$ on each side on the line. For each
$(a,i)$, the resulting triplet $(s_{\rm min},s_{\rm max},|\delta|_{\rm max})$
characterizes the spread of the ring diameters obtained from our models. These
triplets are listed in Table 6.

We note that, aside from the extreme case with $a=0.99$, $i=45^{\circ}$, we
typically have $s_{\rm min}\gtrsim 0.18$ and $s_{\rm max}\lesssim 0.75$.
Hence, our models produce rings that are concentrated away from the edges of
the lensing band, as already observed in Sec. 5.4. This is consistent with an
“effective” emission radius $2M\lesssim r_{s}\lesssim 7M$, in line with the
choice of emission peak locations in Table 5.

###  6.3 Inferring BH spin and inclination from photon ring shape

In this survey, we have solved the “forward problem” of how to compute the
diameters $(d_{+},d_{-})$ of the $n=2$ photon ring given a choice of emission
profile, and a BH spin and inclination $(a,i)$. Purely geometric
considerations constrain the observed photon ring to lie within the $n=2$
lensing band, which corresponds to the blue regions in Figs. 10, 11, and 12.
Perhaps suprisingly, we have found that the ring shapes produced by our
emission models occupy a much smaller subregion (the green bounding box in
Fig. 12) of this full region of geometrically allowed shapes.

This observation is highly relevant to the “inverse problem” of how to infer
BH spin and inclination from a measurement of the $n=2$ ring diameters
$(d_{+},d_{-})$. Given such a data point, one could try to determine which
blue regions ($n=2$ lensing bands) it is contained within. The answer would
provide the set of $(a,i)$ for which lensing would allow the observed photon
ring shape.

Evidently, given the considerable size of these blue regions, these purely
geometric constraints on $(a,i)$ (which are derived from lensing alone) are
relatively weak. However, if one were to make additional assumptions about the
astrophysics of the source (such as the class of emission models being
considered), then one would obtain narrower subregions for each $(a,i)$ (the
smaller green boxes), and a given data point would therefore be compatible
with a smaller range of BH spins and inclinations.

Since our five models were picked to span the full range of $s$ observed in
our survey over emission profiles in Sec. 5, we expect $s_{\rm min}$ and
$s_{\rm max}$ to be close to the values they would take in more general
situations, at least so long as the emission peaks near an equatorial source
ring at small radius outside the horizon. On the other hand, our values of
$|\delta|_{\rm max}$ are probably less reliable, and would likely increase if
we were to consider geometrically thick models in which the emission is not
well-approximated by an equatorial source ring.

Thus, we do not claim here to have definitively computed the bounding boxes
for “realistic” choices of astrophysical profiles. Rather, we hope to have
provided useful organizing principles for tackling the problem of spin
inference: the choice between blue and green boxes clarifies the interplay
between geometric and astrophysical assumptions, with the additional
constraining power provided by the latter made manifest in the much smaller
size of the green regions relative to the blue ones.

An additional benefit of this framework is that it enables us to quantify the
uncertainty in a given spin estimate by computing the range of $(a,i)$ that
produce green or blue boxes containing a given measurement $(d_{+},d_{-})$ of
the $n=2$ photon ring shape. We leave a more detailed examination of this
problem to the future.

##  7 Discussion

###  7.1 Limitations of this work

We have revisited the proposal made by Gralla et al. (2020) for a test of the
Kerr hypothesis based on the photon ring shape. We have established the
robustness of their proposed test in a large class of models, provided that
the BH inclination is not too high. This confirms that M87* holds great
promise as an observational target for this test. However, as pointed out in
the introduction, this conclusion remains subject to some caveats.

First, we have only considered simple equatorial emission models. It is
important to test whether the GLM method remains viable in the presence of
more realistic accretion flows (such as fluctuating, geometrically thick
structures), for which effects such as absorption can become very significant
(for instance, the interferometric signature of the photon ring is suppressed
as the emission region grows thicker and more opaque). Encouraging preliminary
results in this direction will be reported elsewhere (Vincent et al., 2022).

Second, we have not presented a complete, quantitative study of the range of
baseline lengths for which the method works. In particular, we have not
identified an “ideal” baseline window in which the GLM test succeeds for the
largest set of configurations (astrophysical profiles, BH spins, and
inclinations), and we have largely ignored practical concerns in our
theoretical study of the baseline domain. With the current wavelength of
observations (1.3 mm), a space antenna would have to be almost 1.5 million km
away from Earth—close to $L_{2}$!—to sample the visibility at $u_{w}\sim
1000\,$G$\lambda$ (Johnson et al., 2020).

Third, we have not included any noise in our simulations. Generating mock data
with both instrumental and astrophysical noise is crucial, not only for
demonstrating the feasibility of this proposed test, but also for developing a
more robust, Bayesian approach to the multi-fit method described in Sec. 3. In
their original analysis, Gralla et al. (2020) included instrument noise but
not did not model astrophysical fluctuations, which can have a considerable
impact on the photon ring signatures present in instantaneous snapshots. A new
code tailored to include such fluctations—AART, for “Analytic Adaptive Ray
Tracing”—has recently been developed (Cárdenas-Avendaño et al., 2022) to
tackle this specific problem and will inform the requirements for (incoherent)
time-averaging.

While amplitude gain errors can be problematic for VLBI data analyses (e.g.,
Event Horizon Telescope Collaboration, 2019d, f), it is worth noting that they
may not pose a significant problem for a measurement of the photon ring
diameter, which would only require measuring the periodicity of oscillations
in the visibility amplitude (rather than their absolute scale).

###  7.2 Perspectives: Inferring black hole spin and inclination

If a future space mission targeting M87* successfully measured the shape of
its $n=2$ photon ring and confirmed that it matched the Kerr prediction (41),
then the natural next step would be to assume the Kerr nature of the source
and try to infer the BH spin and inclination $(a,i)$ from the measured photon
ring shape.

Our survey over emission profiles in Sec. 5 shows that this shape displays a
residual dependence on the astrophysical source profile (even for the $n=2$
ring). Though small, this dependence nonetheless introduces significant
uncertainty in estimates of the spin based on the $n=2$ ring shape alone, as
evidenced by the lower bounds in Table 4. (This uncertainty is undoubtedly
much greater for spin estimates based on the $n=1$ ring shape.)

Spin estimates can be tightened with additional assumptions about the nature
of the emission. As discussed in Secs. 5.4 and 6.2, the photon ring shapes
produced by our emission models inhabit only a small subset of the full range
of geometric shapes allowed by BH lensing. Thus, assuming that an observed
photon ring shape is produced by a model within this class considerably
narrows down the set of BH spins and inclinations compatible with the
measurement, excluding many more values of $(a,i)$ than could be ruled out
based on the Kerr geometry alone.

It would be extremely interesting to quantify the reduction in uncertainty
enabled by different astrophysical assumptions. This requires one to tackle
the inverse problem of how to determine the values of $(a,i)$ that are
compatible with a given photon ring shape; we hope that our framing of this
problem in Sec. 6.3 will prove useful for future discussions of spin
inference. We note that other approaches based on the identification of the
photon ring with the critical curve either fail to quantify the significant
systematic error in their spin estimates (see Figs. 10 and 11), or else they
effectively incorporate some version of our analysis.

Finally, spin estimates can also be substantially improved by combining the
constraints derived from multiple image features. Early work suggests that
measuring both the photon ring shape and that of the central brightness
depression could significantly constrain BH spin and inclination (Chael et
al., 2021). Of course, simultaneously resolving more than one ring (such as
the $n=1$ and $n=2$ subrings, for instance) should also enable sharper spin
estimates, though this question has so far only been explored in the simplest
cases with circular subrings (Broderick et al., 2022).

###  7.3 Perspectives: Other photon ring tests

The prospect of measuring multiple subrings naturally opens the door to even
more stringent tests of the Kerr hypothesis (or more precise techniques for
determining the BH parameters).

For instance, comparing the widths of successive subrings at different angles
around the image could in principle provide a measurement of the angle-
dependent demagnification factors $e^{-\gamma(\varphi)}$, and therefore of the
Lyapunov exponents $\gamma(\tilde{r})$ that offer a direct probe of the Kerr
geometry (Johnson et al., 2020; Gralla & Lupsasca, 2020a). If the
demagnification factors were not $e^{-\gamma(\varphi)}$, then this might
suggest that the spacetime is not Kerr (Wielgus, 2021). Conversely, if
$\gamma(\varphi)$ did match the form predicted from Kerr, then these exponents
would provide a function’s worth of constraints on the BH spin and inclination
$(a,i)$, which could allow these parameters to be inferred with great
precision.

Likewise, polarimetric imaging of successive subrings could also be used to
test the Kerr hypothesis, or conversely, to yield an estimate of the BH
parameters $(a,i)$ (Himwich et al., 2020), while also providing a probe of the
plasma and magnetic field around the BH (Event Horizon Telescope
Collaboration, 2021b). Lastly, the photon ring substructure predicts a
characteristic pattern of autocorrelations whose detection could also provide
a window into the Lyapunov exponents $\gamma$ (and other Kerr critical
exponents $\tau$ and $\delta$), even if the rings themselves remain unresolved
(Hadar et al., 2021). A related promising idea involves timing signatures such
as quasi-periodic oscillations (Johannsen & Psaltis, 2010).

Finally, a wide variety of other techniques for testing the Kerr hypothesis
and strong-field GR via VLBI observations have been proposed (e.g., Psaltis,
2019), many of which do not rely upon the photon ring. Similarly, some
constraints on BH inclination can be derived from other image features, such
as its brightness asymmetry (Medeiros et al., 2022).

###  7.4 Open questions and future directions

We have thoroughly investigated the shape of the $n=2$ photon ring and its
corresponding interferometric signature in a certain class of equatorial
models. We can now confirm the prediction that the ring should follow a
circlipse shape at low-to-moderate inclinations, but have yet to obtain a
precise characterization of its shape at high inclinations, where its image