# Impact of overlapping signals on parameterized post-Newtonian coefficients
in tests of gravity

Yixuan Dang Department of Physics, University of Oxford, Parks Road, Oxford
OX1 3PU, UK  Ziming Wang Department of Astronomy, School of Physics, Peking
University, Beijing 100871, China  Kavli Institute for Astronomy and
Astrophysics, Peking University, Beijing 100871, China  Dicong Liang Kavli
Institute for Astronomy and Astrophysics, Peking University, Beijing 100871,
China  Lijing Shao Kavli Institute for Astronomy and Astrophysics, Peking
University, Beijing 100871, China  National Astronomical Observatories,
Chinese Academy of Sciences, Beijing 100012, China

###### Abstract

Gravitational waves have been instrumental in providing deep insights into the
nature of gravity. Next-generation detectors, such as the Einstein Telescope,
are predicted to have a higher detection rate given the increased sensitivity
and lower cut-off frequency. However, this increased sensitivity raises
challenges concerning parameter estimation due to the foreseeable overlap of
signals from multiple sources. Overlapping signals (OSs), if not properly
identified, may introduce biases in estimating post-Newtonian (PN)
coefficients in parameterized tests of general relativity (GR). We investigate
how OSs affect $-1$PN to 2PN terms in parameterized GR tests, examining their
potential to falsely suggest GR deviations. We estimate the prevalence of such
misleading signals in next-generation detectors, and their collective
influence on GR tests. We compare the effects of OSs on coefficients at
different PN orders, concluding that overall the 1PN coefficient suffers the
most. Our findings also reveal that while a non-negligible portion of OSs
exhibit biases in PN coefficients that might individually prefer to conclude
deviations from GR, collectively, the direction to deviate is random and a
statistical combination will still be in favor of GR.

Gravitational Waves (678) — General Relativity (641) — Compact binary stars
(283)

††facilities: Einstein Telescope

##  1 Introduction

Gravitational waves (GWs), predicted by general relativity (GR), have provided
essential tests of gravity since their first direct detection by the
LIGO/Virgo Collaboration in 2015 (Abbott et al., 2016). Current detected GWs
are all generated from compact binary coalescences, such as the mergers of
binary black holes (BBHs) or neutron stars, and carry information about the
sources and the nature of gravity itself. To date, observatories including
LIGO, Virgo, and KAGRA have identified approximately 90 GW events, offering us
some of the most profound insights into GR and significantly enhancing our
grasp of gravitational dynamics in the strong-field regime (Abbott et al.,
2017a, b, 2019a, 2021a). Various tests of GR were also performed with current
detections with no conclusive sign of GR being falsified (Abbott et al.,
2019b, c, 2021b, 2021c; Ghosh, 2022).

Next generation (XG) detectors, for example, the Einstein Telescope (ET;
Punturo et al., 2010; Hild et al., 2011; Sathyaprakash et al., 2012) and the
Cosmic Explorer (CE; Reitze et al., 2019a, b), are expected to receive far
more signals than current detectors do due to the largely increased
sensitivity and a lower cut-off frequency (Sathyaprakash et al., 2019;
Kalogera et al., 2021). In particular, ET, whose sensitivity is expected to be
better by a factor of 10 than current-generation ground based detectors, has
been shown to detect overlapping signals (OSs) in the mock data challenges
(Regimbau et al., 2012; Meacher et al., 2016).

While the enhanced sensitivity of detectors is a significant achievement, it
introduces a new challenge in the parameter estimation (PE) process that
several studies have explored, notably in dealing with the potential overlap
of signals from multiple GW sources (Samajdar et al., 2021; Pizzati et al.,
2022; Relton & Raymond, 2021; Regimbau & Hughes, 2009). Such OSs can
complicate the PE process, introducing biases across parameters (Antonelli et
al., 2021; Janquart et al., 2022; Wang et al., 2023). In this work we
investigate if misidentifying OSs can distort PE on parameters in the post-
Newtonian (PN) waveform (Yunes & Pretorius, 2009), suggesting deviations from
GR. In the case where a weaker secondary signal overlaps with a primary
signal, erroneously classifying this as a single, pure signal can lead to
inaccuracies. It leads us to three crucial questions:

  * •

Is it possible for OSs, if inaccurately perceived as from a singular source,
to misguide us into concluding a deviation from GR?

  * •

Are such deceptive OSs frequent enough to pose a genuine concern for XG
detectors?

  * •

If a substantial number of such signals are identified, do they, as a group,
indicate a deviation from GR?

There exist abundant studies on tests of parameterized PN coefficients for GR
in general (Arun et al., 2006; Mishra et al., 2010; Cornish et al., 2011;
Abbott et al., 2021c; Wang et al., 2022). For OSs, Hu & Veitch (2023) have
shown that they may lead to incorrect measurement of deviation from the 0PN
term in the GR waveform.

However, a comprehensive assessment of the effects of OSs on PN coefficients
is lacking. Our work aims to explore the aforementioned questions related to
OSs regarding false-alarmed deviations in various PN coefficients, including
$-1$PN, 0PN, 1PN, 2PN coefficients, and offers a comparative analysis of how
each PN coefficient responds to misidentified OSs. We anticipate not only a
potential for misleading deviations across all PN terms but also that some PN
terms might be more adversely affected than others when leaving out the
secondary signals in the OSs.

Our results show that there exist OSs resulting in large bias in multiple PN
coefficients, and notably, these misleading signals represent a non-negligible
portion of all overlapping instances to be registered by XG detectors. We also
find that the relative size of the fractions of misleading OSs at different PN
orders depends on how we define OSs. If only two GWs with very small merger-
time difference can be identified as OSs, the fraction of misleading OSs at
higher PN orders is larger. But generally, the impact of overlapping on the
$1$PN term is most significant among $-1$PN to 2PN coefficients, that is, GR-
violation effects are most likely to be observed at 1PN when misidentifying an
OS as one single GW signal. Fortunately, these signals only prefer models with
deviation from GR individually, but not collectively, manifested in the rapid
decrease of Bayes factor with increased occurrences of misleading signals.
This aligns with the result in Hu & Veitch (2023), which considered the
deviation from the 0PN term when OSs exist.

This paper is organized as follows. We introduce PE methods used in this work
in Sec. 2, including how we calculate the Fisher matrix, biases, and Bayes
factor for OSs. In Sec. 3 we list our simulation setup for both the waveforms
and the population model for BBHs. Our results are given in Sec. 4, and Sec. 5
concludes the work.

##  2 Parameter estimation methods

To determine the parameters $\bm{\theta}$ of the GW signal from data $g(t)$
received by GW detectors, one can calculate the posterior distribution of the
parameters using Bayes’ theorem,

$P(\bm{\theta}|g(t))\propto P(\bm{\theta})P(g|\bm{\theta}),$ (1)

where $P(\bm{\theta})$ is the prior of the parameters $\bm{\theta}$ given
assumed background knowledge, and $P(g|\bm{\theta})$ is the conditional
probability of receiving data $g(t)$ given parameters $\bm{\theta}$. The data
$g(t)$ consist of the GW signal $h(t)$ and noise $n(t)$. Assuming the noise is
stationary and Gaussian, the likelihood reads (Finn, 1992)

$P(g|\bm{\theta})\propto e^{-\frac{1}{2}(g-h,g-h)},$ (2)

where the inner product $(u,v)$ is defined as

$(u,v):=2\Re\int_{-\infty}^{\infty}\frac{u^{*}(f)v(f)}{S_{n}(|f|)}{d}f,$ (3)

where $S_{n}$ is the power spectral density of the noise, $u(f)$ and $v(f)$
are the Fourier transforms of $u(t)$ and $v(t)$, respectively.

In principle, given the GW data, one can use some sampling techniques, such as
the Markov-Chain Monte Carlo (MCMC) method (Christensen & Meyer, 1998;
Christensen et al., 2004; Sharma, 2017) and Nested Sampling (Skilling, 2004,
2006), to obtain the posterior distributions of $\bm{\theta}$. However, the
full Bayesian inference takes considerable computational time and resources,
especially when the dimension of parameter space is high. In this work, we
adopt the Fisher matrix (FM) as a fast approximation of the full Bayesian
inference, which has been successfully applied to conducting PE for numerous
OSs in previous works (Antonelli et al., 2021; Hu & Veitch, 2023; Wang et al.,
2023). Mathematically, the FM is defined as (Cutler & Flanagan, 1994;
Vallisneri, 2008)

$F_{\alpha\beta}=(h_{,\alpha},h_{,\beta})\,,$ (4)

where $h_{,\alpha}\equiv\partial h/\partial\theta^{\alpha}$. Under the high
signal-to-noise ratio (SNR) limit, or equivalently the linear signal
approximation (LSA), the posterior can be well approximated by a Gaussian
distribution (Finn, 1992; Vallisneri, 2008). Denoting the inverse of the FM as
$\Sigma\equiv F^{-1}$ and assuming flat priors, the standard deviations and
correlation coefficients of parameters are given by (no summation over
repeated indices)

$\sigma_{\alpha}=\sqrt{\Sigma^{\alpha\alpha}}\,,\quad
c_{\alpha\beta}=\frac{\Sigma^{\alpha\beta}}{\sqrt{\Sigma^{\alpha\alpha}\Sigma^{\beta\beta}}}\,.$
(5)

Considering a superposition of two independent GW signals in the detector, the
data can be written as

$g(t)=h^{(1)}(t;\tilde{\bm{\theta}}^{(1)})+h^{(2)}(t;\tilde{\bm{\theta}}^{(2)})+n(t)\,,$
(6)

where $\tilde{\bm{\theta}}^{(i)}$ represents the true parameters of the $i$-th
signal. We consider the situation where an OS is mistakenly categorized as a
signal containing a single event, and therefore only one-signal model
$h(t;{\bm{\theta}})$ is used in the PE. Usually, the louder signal can be
identified, while the weaker one induces biases in the PE of the louder signal
(Pizzati et al., 2022). Without loss of generality, we assume that signal 1 is
louder than signal 2.

When conducting PE, it is more convenient to consider the ratio of the bias
relative to the statistical uncertainty of the parameters, that is, the
reduced bias (Wang et al., 2023)

$B_{\alpha}:=\frac{{\Delta\theta}^{\alpha}_{\rm bias}}{{\sigma_{\alpha}}}\,,$
(7)

where
${\Delta\theta}^{\alpha}_{\rm{bias}}={\Sigma^{\alpha\beta}\big{(}h_{,\beta}^{(1)},h^{(2)}\big{)}}$
refers to the absolute bias. Note that the index $\alpha$ is not summed. In
this formula, we have assumed LSA, meaning that the SNR of signal 1 is large
enough.

In PE, when $|B_{\alpha}|$ exceeds 1, the bias is significant. When
$\theta^{\alpha}$ is chosen as a PN-deviation parameter in the parameterized
PN waveform test, we call the OSs with $|B_{\alpha}|>1$ as “misleading OSs”,
which means that this kind of OSs may lead to a false conclusion of a
deviation from GR.

Bayes factor is another useful indicator for assessing whether the data show a
stronger preference for GR or for a deviation from GR. The Bayes factor
$K_{\mathrm{M_{1}},\mathrm{M_{0}}}$ is defined as the ratio of evidences,

$K_{\mathrm{M_{1}},\mathrm{M_{0}}}=\frac{P(D|\mathrm{M_{1}})}{P(D|\mathrm{M_{0}})},$
(8)

where $D$ stands for data, and $\mathrm{M_{1}}$, $\mathrm{M_{0}}$ are two
competing models. Bayes factor compares how well each model predicts the data.
Roughly speaking, when $K_{\mathrm{M_{1}},\mathrm{M_{0}}}$ tends to an extreme
large/small value, the data favor/disfavor $\mathrm{M_{1}}$ over
$\mathrm{M_{0}}$. Otherwise, the preference for $\mathrm{M_{1}}$ and
$\mathrm{M_{0}}$ is comparable.

In our calculation, $D$ is the GW signal $g(t)$, $\mathrm{M_{1}}$ represents a
model allowing the waveform deviating from GR at $-1$, 0, 1 and 2 PN orders,
respectively denoted as $\rm{M}_{-1\mathrm{PN}}$, $\rm{M}_{0\mathrm{PN}}$,
$\rm{M}_{1\mathrm{PN}}$ and $\rm{M}_{2\mathrm{PN}}$. The model
$\mathrm{M_{0}}$ represents GR. As an example, assuming flat priors and LSA,
the Bayes factor of $\rm{M}_{-1\mathrm{PN}}$ against ${\rm{M}_{0}}$ is given
by

$\displaystyle K_{\rm{M}_{-1\mathrm{PN}},\rm{M}_{0}}$
$\displaystyle=\frac{P\big{(}g(t)|\rm{M}_{-1\mathrm{PN}}\big{)}}{P\big{(}g(t)|\rm{M}_{0}\big{)}}$
(9)
$\displaystyle=p\big{(}{\delta\phi_{-1}}\big{)}\sqrt{2\pi}\frac{e^{-\frac{1}{2}\hat{\chi}_{\rm{M}_{-1\mathrm{PN}}}^{2}}}{e^{-\frac{1}{2}\hat{\chi}_{\rm{M}_{0}}^{2}}}\sqrt{\frac{|\Sigma_{\rm{M}_{-1\mathrm{PN}}}|}{|\Sigma_{\rm{M}_{0}}|}},$

where $p\big{(}{\delta\phi_{-1}}\big{)}$ is the prior probability for the
$-1$PN deviation parameter $\delta\phi_{-1}$. We incline to consider
$p\big{(}{\delta\phi_{-1}}\big{)}$ as a constant. The residual chi-square,
$\hat{\chi}^{2}_{{\rm M}_{i}}$, is defined as

$\hat{\chi}^{2}_{{\rm
M}_{i}}=\left(h^{(2)}-h_{,\alpha}^{(1)}\,\Delta{\theta}^{\alpha}_{{\rm
bias,M}_{i}},h^{(2)}-h_{,\alpha}^{(1)}\,\Delta{\theta}^{\alpha}_{{\rm
bias,M}_{i}}\right),$ (10)

where $\Delta{\theta}^{\alpha}_{{\rm bias,M}_{i}}$ is the absolute bias of
parameter $\alpha$ under hypothesis ${\rm M}_{i}$.

##  3 Simulation setup

We generate non-spinning, stellar-mass BBH signals using the IMRPhenomD
waveform model. Similarly to Wang et al. (2023), we average the waveform over
the sky position, inclination and polarization angles. This corresponds to a
conservative scenario, where two GW signals come from the same direction with
the same inclination angle and polarization angle. Using the angle-averaged
waveform will affect the fraction of misleading OSs. However, here we are more
interested in comparing the fractions of misleading OSs for different PN
deviations, which are hardly affected by the angle-averaging.

When generating OSs, we use the GR waveform for both signals 1 and 2, and the
free parameters for each signal are the chirp mass $\mathcal{M}$ in the
detector frame, the symmetry mass ratio $\eta$, the luminosity distance
$d_{L}$, and the time of coalescence $t_{c}$. When we recover with one single
GR waveform, the estimated parameters are simply $\\{{\cal
M},\eta,d_{L},t_{c}\\}$ of the primary signal.

For the GR-deviation model at the $i$-th PN order, we add a deviation
parameter $\delta\phi_{i}$, which leads to a phase correction in the Fourier-
domain waveform (Abbott et al., 2021c)

$\Delta\Phi_{i}(f)=\frac{3}{128}\mathfrak{u}^{2i-5}\delta\phi_{i}\,,$ (11)

where $\mathfrak{u}\equiv(\pi Mf)^{1/3}$, and $M$ is the redshifted total
mass. We only consider the modification in the inspiral stage, which cuts off
when the frequency is twice the innermost stable circular orbit frequency.
Therefore, for the $i$-th PN deviation model, the estimated parameters are
$\\{{\cal M},\eta,d_{L},t_{c},\delta\phi_{i}\\}$. Note that deviation
parameters in different PN orders are usually studied separately, assuming
only one PN deviation parameter is nonzero, which is more effective to pick up
possible deviations at different PN orders (Abbott et al., 2021c). The case
with multiple nonzero PN-deviation parameters was also investigated, but with
much less constraining power due to increased degeneracy of parameters.

As for the detector, we choose the ET as an example, which is expected to
detect thousands of OSs per year (Pizzati et al., 2022). We will use the
designed ET-D sensitivity curve (Punturo et al., 2010; Hild et al., 2011;
Maggiore et al., 2020).

Now we introduce the population model. When generating OSs, we inject signal 1
and signal 2 with GR waveform, which requires the binary masses, the
luminosity distance and the time of coalescence of the two signals. We also
request the SNR of signal 1 to be larger than 30 to ensure the validity of FM
approximation.

The BBH mass population is generated according to the POWER+PEAK
phenomenological model given in Abbott et al. (2019d, 2023). For the
distribution of the luminosity distance, we adopt the model given in Regimbau
et al. (2017), which gives the merger rate in the detector frame,

$R_{z}(z)=\frac{R_{m}(z)}{1+z}\frac{dV(z)}{dz}\,,$ (12)

where $V$ is the comoving volume. $R_{m}(z)$ is the merger rate per comoving
volume at redshift $z$,

$R_{m}(z)=\int_{t_{\min}}^{t_{\max}}R_{f}\left(\tilde{z}\left[\tilde{t}(z)-t_{d}\right]\right)P\left(t_{d}\right)dt_{d}\,.$
(13)

Here, $t_{d}$ is time delay between the formation of a binary and its merger,
$\tilde{t}(z)$ is the cosmic time when the merger happens, $\tilde{z}$ is the
redshift as a function of cosmic time, $P(t_{d})$ is the distribution of the
time delay, which is inversely proportional to $t_{d}$ and normalized between
$t_{\min}$ and $t_{\max}$. Following previous works, we choose
$t_{d}^{\mathrm{min}}=50\,\mathrm{Myr}$ (Meacher et al., 2015; Abbott et al.,
2018), and $t_{d}^{\mathrm{max}}$ equal to the Hubble time (Ando, 2004;
Belczynski et al., 2006; Nakar, 2007; O’Shaughnessy et al., 2008; Dominik et
al., 2012, 2013). The star formation rate $R_{f}\left(z\right)$ is given by

$R_{f}(z)=\nu\frac{ae^{b\left(z-z_{p}\right)}}{a-b+be^{a\left(z-z_{p}\right)}}\,,$
(14)

where $\nu=0.146\,\mathrm{M}_{\odot}\,\mathrm{yr}^{-1}\,\mathrm{Mpc}^{-3}$,
$a=2.80$, $b=2.46$ and $z_{p}=1.72$ (Vangioni et al., 2015; Behroozi & Silk,
2015).

The population model of $t_{c}$ needs more discussion. On the one hand, both
the biases and the Bayes factor only depend on the difference between
$t_{c}^{(1)}$ and $t_{c}^{(2)}$, denoted as $\Delta t_{c}$, so we only need to
consider the distribution of $\Delta t_{c}$. On the other hand, in this work
we are interested in the fraction of misleading OSs in all OSs, but there is
no unified criterion for OSs (Wang et al., 2023). Usually, this is done by
setting a threshold $\Delta t_{c,\rm{th}}$ for $\Delta t_{c}$ between two
signals, which ranges from 0.1 s to several seconds in the previous works
(Pizzati et al., 2022; Relton & Raymond, 2021; Himemoto et al., 2021). To
avoid the ambiguity of choosing $\Delta t_{c,\rm th}$, we vary $\Delta
t_{c,\rm{th}}$ in a wide range, from 0.1 s to ${\cal{O}}(10)$ s and generate
the corresponding $t_{c}$ samples for each threshold. Consider an observation
with a duration of $T$, and two GW events that are equally likely to occur at
any time in this observation, i.e., both $t_{c}^{(1)}$ and $t_{c}^{(2)}$ obey
a uniform distribution in the time interval $(0,T)$. We find that the
conditional distribution of $\Delta t_{c}$ given a threshold,
$P\Big{(}t_{c}^{(2)}-t_{c}^{(1)}\Big{|}|t_{c}^{(2)}-t_{c}^{(1)}|\leq\Delta
t_{c,\rm th}\Big{)}$, can be well approximated by a uniform distribution.

##  4 Results

Figure 1: Reduced bias of each PN deviation term as a function of $\Delta
t_{c}$, for an OS characterized by parameters in Eq. (15). Figure 2: Fraction
of misleading OSs for each PN-deviation model against $\Delta t_{c,\rm th}$.

As an illustration, we consider an OS with the following parameters,

$\displaystyle\mathcal{M}^{(1)}$
$\displaystyle=32\,\mathrm{M}_{\odot}\,,\quad$
$\displaystyle\mathcal{M}^{(2)}$ $\displaystyle=35\,\mathrm{M}_{\odot}\,,$
$\displaystyle\eta^{(1)}$ $\displaystyle=0.246\,,\quad$
$\displaystyle\eta^{(2)}$ $\displaystyle=0.248\,,$ $\displaystyle d_{L}^{(1)}$
$\displaystyle=10000\,\rm{Mpc}\,,\quad$ $\displaystyle d_{L}^{(2)}$
$\displaystyle=45000\,\rm{Mpc}\,,$ $\displaystyle\rm{SNR}^{(1)}$
$\displaystyle=41.74\,,\quad$ $\displaystyle\rm{SNR}^{(2)}$
$\displaystyle=9.96\,.$ (15)

In Figure 1 we plot the change of the reduced bias with respect to $\Delta
t_{c}$ for this OS. It answers our first question in Sec. 1, namely, it is
possible that, under the effects of a weaker signal 2, some PN terms of signal
1 possess a reduced bias larger than 1. It is also worth attention that large
biases of $-1$PN, 0PN, 1PN coefficients are obtained around $\Delta
t_{c}=1\,{\rm s}$, while the bias of $2$PN coefficient has peaked at a much
smaller $\Delta t_{c}$ around zero. This can be explained by the fact that the
higher order PN terms are mainly prominent in the late inspiral phase of the
signal, and we expect the deviation from a higher order PN term to be more
significant when $\Delta t_{c}$ between two signals is small.

Simply highlighting a single instance with parameters in Eq. (15) is not
sufficient to underscore the risks posed by all OSs, as such instances may be
too infrequent to a warrant general concern. To answer the second question in
Sec. 1, we generate populations for different $\Delta t_{c,\rm th}$ and
calculate the biases for every PN-deviation waveform under consideration.
Afterwards, we calculated the fraction of misleading OSs in all OSs in each
scenario. The findings are presented in Figure 2. The observation is that the
fractions of misleading OSs, defined as $|B_{\delta\phi}|>1$, are around ten
percent for populations with small $\Delta t_{c,\rm th}$. Such fractions are
not negligible. All fractions then decrease as $\Delta t_{c,\rm th}$
increases, because most signals with strong overlapping features happen at
small $\Delta t_{c}$. When $\Delta t_{c,\rm th}$ is large, e.g., larger than
10 s, all fractions are inversely proportional to $\Delta t_{c,\rm th}$. This
means that almost all misleading OSs have $\Delta t_{c}\leq 10\,{\rm
s}$.111The fraction of the misleading OSs is  $\frac{N({\rm
misleading\,\,OSs})}{N({\rm all\,\,OSs})}=\frac{N({\rm
misleading\,\,OSs})}{N({\rm OSs}|\Delta t_{c}\leq 10{\rm s})}\cdot\frac{N({\rm
OSs}|\Delta t_{c}\leq 10{\rm s})}{N({\rm all\,\,OSs})}\,.$ If all misleading
OSs have $\Delta t_{c}\leq 10\,{\rm s}$, the first factor on the right-hand
side of the above equation does not depend on $\Delta t_{c,\rm th}$, while the
second factor is inversely proportional to $\Delta t_{c,\rm th}$.

The behaviors of these fractions in Figure 2 at different PN orders are also
noteworthy. We note that the fraction of misleading OSs in the 2PN-deviation
case is larger than that in the 1PN-deviation case when $\Delta t_{c,\rm th}$
is within $0.5$ s, while for a larger $\Delta t_{c,\rm th}$ it is the
opposite. Since the 2PN term is more prominent in the late inspiral phase, it
is more likely to generate misleading OSs when $\Delta t_{c}$ is small (as
shown in Figure 1), corresponding to a relatively large fraction. As a trade-
off, almost all misleading OSs in the 2PN-deviation case have $\Delta
t_{c}\lesssim 0.5\,{\rm s}$, corresponding to the quick decrease after $\Delta
t_{c,{\rm th}}\simeq 0.5\,$s. In comparison, in the 0PN and 1PN deviation
cases, similar behavior of these fractions is observed at $\Delta t_{c,\rm
th}$ between 1 s and 10 s. The fraction of misleading OSs in the 1PN-deviation
case decreases faster than that at the 0PN order, meaning that the number of
misleading OSs in the 0PN deviation case with $1\,{\rm s}<\Delta
t_{c}<10\,{\rm s}$ is larger than that in the 1PN-deviation case. However,
since $\Delta t_{c}$ is relatively large when the 0PN term dominates the phase
evolution, it is difficult to generate large biases, and the total number of
misleading OSs in the 1PN-deviation case is still larger than that in the 0PN-
deviation case for a large $\Delta t_{c,\rm th}$. In general, for most $\Delta
t_{c,\rm th}$, the fraction of large bias in 1PN coefficient stays the highest
among the PN-deviation cases we investigate in this work.

Figure 3: Bayes factor divided by $p(\delta\phi)$ over accumulation of
misleading OSs, denoted as $N_{\mathrm{MOS}}$, with different $\Delta t_{c,\rm
th}$ for each panel.

From the above findings, we observe a non-negligible proportion of misleading
OSs to parameterized PN waveform tests. These OSs, when analyzed individually
in tests of GR, tend to favor their respective PN-deviation model more
strongly relative to the GR waveform. However, this is not the case when all
misleading OSs are considered collectively because of the disagreement on the
values of deviation parameters arising from each misleading OS. Figure 3 shows
the evolution of Bayes factors of each PN-deviation model against GR with the
increase of number of events. For all PN-deviation models, their corresponding
Bayes factors exceed 1 only for the first few misleading OSs, and then start
drastically decreasing as more misleading OSs enter the test of GR. It is
evident that biases in PN deviation terms caused by OSs are not constant as
they do not represent a “real” deviation. This is because these biases arise
solely from the misclassification of an OS as a single-event signal.

Figure 4: Distribution of large biases of each PN deviation parameter for
different $\Delta t_{c,\rm th}$.

However, in our PN-deviation models, the deviation parameters are assumed to
be constant. When combining two misleading OSs, where, say, one has the
reduced bias $B_{\alpha}>1$ and the other has $B_{\alpha}<-1$, the Bayes
factor will rapidly decrease due to the inconsistent biases. To support this,
we plot the distribution of reduced biases of the misleading OSs for each PN
deviation model in Figure 4, from which we can observe symmetric distributions
of reduced biases.

Another interesting aspect noted in Figure 3 is that the Bayes factor for the
2PN deviation model declines more rapidly than those for the $-1$PN, 0PN, and
1PN deviation models, especially at small $\Delta t_{c,\rm th}$. In Figure 4,
one can find that the distributions of biases of the $-1$PN, 0PN, and 1PN
deviation models cluster more densely near 1, while the distribution for the
2PN deviation model is broader, leading to a greater variance in individual
biased values. This is consistent with the argument above that it is more
likely to have larger biases in the deviation parameter at the 2PN order.
Therefore, the Bayes factor for the 2PN deviation model decreases more rapidly
with accumulation of misleading OSs.

##  5 Conclusions and Discussions

The XG ground-based detectors are anticipated to exhibit a significantly
enhanced detection rate of BBHs relative to their forerunners, indicating the
unavoidable appearance of OSs in future detection (Regimbau & Hughes, 2009;
Meacher et al., 2016; Kalogera et al., 2021). A primary concern regarding such
OSs is how they can potentially introduce biases if an OS is mistakenly
treated as a single-event signal, thereby misleading parameterized tests of
GR. Hu & Veitch (2023) studied the biases in PE when testing GR at the 0PN
order with these OSs, and left the comprehensive investigation at different PN
orders for future work.

In this work, we extend Hu & Veitch (2023) and investigate how OSs affect
$-1$PN to 2PN terms in parameterized tests of GR. We adopt the FM
approximation to calculate the biases introduced by OSs. FM is a fast and
effective method to study the biases of OSs (Finn, 1992; Antonelli et al.,
2021; Hu & Veitch, 2023; Wang et al., 2023). For each PN deviation parameter,
we calculate the proportion of signals exhibiting significant bias, and
further explore the Bayes factor over accumulation of such misleading OSs. We
find that a non-negligible fraction of OSs could produce a significant bias in
each PN deviation parameter, especially when the difference in coalescence
times, $\Delta t_{c}$, is small. It is worth noting that the fraction of
misleading OSs in the 1PN-deviation case is usually the largest, instead of
the 0PN-deviation case studied by Hu & Veitch (2023). This underscores the
importance of consideration for GR deviations at different PN orders when
discussing the impact of OSs on testing gravity.

However, different misleading OSs lead to inconsistent values for deviation
parameters, and if analyzed collectively, they do not show a preference for
non-GR deviations over GR. This is supported by the decreasing trend of the
Bayes factor with an increasing number of misleading OSs, shown in Figure 3.

In real detection, there exist other sources of PE biases when analyzing OSs,
such as the waveform modeling (Hu & Veitch, 2023), instrumental calibration
(Sun et al., 2020; Hall et al., 2019) and detector glitches (Pankow et al.,
2018). We leave a detailed assessment and comparison of all these effects on
testing gravity for future work. However, as the misleading OSs discussed, we
emphasize the importance of considering the effects at different PN orders to
obtain a global understanding.

## Acknowledgements

This work was supported by the National Natural Science Foundation of China
(11975027, 11991053), the China Postdoctoral Science Foundation (2021TQ0018),
the National SKA Program of China (2020SKA0120300), the Roy Skinner fund from
Balliol College at University of Oxford, the Max Planck Partner Group Program
funded by the Max Planck Society, and the High-Performance Computing Platform
of Peking University.

## References

  * Abbott et al. (2016) Abbott, B. P., et al. 2016, Phys. Rev. Lett., 116, 061102, doi: 10.1103/PhysRevLett.116.061102
  * Abbott et al. (2017a) —. 2017a, Phys. Rev. Lett., 119, 141101, doi: 10.1103/PhysRevLett.119.141101
  * Abbott et al. (2017b) —. 2017b, Phys. Rev. Lett., 118, 221101, doi: 10.1103/PhysRevLett.118.221101
  * Abbott et al. (2018) —. 2018, Phys. Rev. Lett., 120, 091101, doi: 10.1103/PhysRevLett.120.091101
  * Abbott et al. (2019a) —. 2019a, Phys. Rev. X, 9, 031040, doi: 10.1103/PhysRevX.9.031040
  * Abbott et al. (2019b) —. 2019b, Phys. Rev. Lett., 123, 011102, doi: 10.1103/PhysRevLett.123.011102
  * Abbott et al. (2019c) —. 2019c, Phys. Rev. D, 100, 104036, doi: 10.1103/PhysRevD.100.104036
  * Abbott et al. (2019d) —. 2019d, Astrophys. J. Lett., 882, L24, doi: 10.3847/2041-8213/ab3800
  * Abbott et al. (2021a) Abbott, R., et al. 2021a.  https://arxiv.org/abs/2111.03606
  * Abbott et al. (2021b) —. 2021b, Phys. Rev. D, 103, 122002, doi: 10.1103/PhysRevD.103.122002
  * Abbott et al. (2021c) —. 2021c.  https://arxiv.org/abs/2112.06861
  * Abbott et al. (2023) —. 2023, Phys. Rev. X, 13, 011048, doi: 10.1103/PhysRevX.13.011048
  * Ando (2004) Ando, S. 2004, JCAP, 06, 007, doi: 10.1088/1475-7516/2004/06/007
  * Antonelli et al. (2021) Antonelli, A., Burke, O., & Gair, J. R. 2021, Mon. Not. Roy. Astron. Soc., 507, 5069, doi: 10.1093/mnras/stab2358
  * Arun et al. (2006) Arun, K. G., Iyer, B. R., Qusailah, M. S. S., & Sathyaprakash, B. S. 2006, Phys. Rev. D, 74, 024006, doi: 10.1103/PhysRevD.74.024006
  * Behroozi & Silk (2015) Behroozi, P. S., & Silk, J. 2015, Astrophys. J., 799, 32, doi: 10.1088/0004-637X/799/1/32
  * Belczynski et al. (2006) Belczynski, K., Perna, R., Bulik, T., et al. 2006, Astrophys. J., 648, 1110, doi: 10.1086/505169
  * Christensen et al. (2004) Christensen, N., Dupuis, R. J., Woan, G., & Meyer, R. 2004, Phys. Rev. D, 70, 022001, doi: 10.1103/PhysRevD.70.022001
  * Christensen & Meyer (1998) Christensen, N., & Meyer, R. 1998, Phys. Rev. D, 58, 082001, doi: 10.1103/PhysRevD.58.082001
  * Cornish et al. (2011) Cornish, N., Sampson, L., Yunes, N., & Pretorius, F. 2011, Phys. Rev. D, 84, 062003, doi: 10.1103/PhysRevD.84.062003
  * Cutler & Flanagan (1994) Cutler, C., & Flanagan, E. E. 1994, Phys. Rev. D, 49, 2658, doi: 10.1103/PhysRevD.49.2658
  * Dominik et al. (2012) Dominik, M., Belczynski, K., Fryer, C., et al. 2012, Astrophys. J., 759, 52, doi: 10.1088/0004-637X/759/1/52
  * Dominik et al. (2013) —. 2013, Astrophys. J., 779, 72, doi: 10.1088/0004-637X/779/1/72
  * Finn (1992) Finn, L. S. 1992, Phys. Rev. D, 46, 5236, doi: 10.1103/PhysRevD.46.5236
  * Ghosh (2022) Ghosh, A. 2022, Summary of Tests of General Relativity with GWTC-3.  https://arxiv.org/abs/2204.00662
  * Hall et al. (2019) Hall, E. D., Cahillane, C., Izumi, K., Smith, R. J. E., & Adhikari, R. X. 2019, Class. Quant. Grav., 36, 205006, doi: 10.1088/1361-6382/ab368c
  * Hild et al. (2011) Hild, S., et al. 2011, Class. Quant. Grav., 28, 094013, doi: 10.1088/0264-9381/28/9/094013
  * Himemoto et al. (2021) Himemoto, Y., Nishizawa, A., & Taruya, A. 2021, Phys. Rev. D, 104, 044010, doi: 10.1103/PhysRevD.104.044010
  * Hu & Veitch (2023) Hu, Q., & Veitch, J. 2023, Astrophys. J., 945, 103, doi: 10.3847/1538-4357/acbc18
  * Janquart et al. (2022) Janquart, J., Baka, T., Samajdar, A., Dietrich, T., & Van Den Broeck, C. 2022.  https://arxiv.org/abs/2211.01304
  * Kalogera et al. (2021) Kalogera, V., et al. 2021.  https://arxiv.org/abs/2111.06990
  * Maggiore et al. (2020) Maggiore, M., et al. 2020, JCAP, 03, 050, doi: 10.1088/1475-7516/2020/03/050
  * Meacher et al. (2016) Meacher, D., Cannon, K., Hanna, C., Regimbau, T., & Sathyaprakash, B. S. 2016, Phys. Rev. D, 93, 024018, doi: 10.1103/PhysRevD.93.024018
  * Meacher et al. (2015) Meacher, D., Coughlin, M., Morris, S., et al. 2015, Phys. Rev. D, 92, 063002, doi: 10.1103/PhysRevD.92.063002
  * Mishra et al. (2010) Mishra, C. K., Arun, K. G., Iyer, B. R., & Sathyaprakash, B. S. 2010, Phys. Rev. D, 82, 064010, doi: 10.1103/PhysRevD.82.064010
  * Nakar (2007) Nakar, E. 2007, Phys. Rept., 442, 166, doi: 10.1016/j.physrep.2007.02.005
  * O’Shaughnessy et al. (2008) O’Shaughnessy, R. W., Kalogera, V., & Belczynski, K. 2008, Astrophys. J., 675, 566, doi: 10.1086/526334
  * Pankow et al. (2018) Pankow, C., et al. 2018, Phys. Rev. D, 98, 084016, doi: 10.1103/PhysRevD.98.084016
  * Pizzati et al. (2022) Pizzati, E., Sachdev, S., Gupta, A., & Sathyaprakash, B. 2022, Phys. Rev. D, 105, 104016, doi: 10.1103/PhysRevD.105.104016
  * Punturo et al. (2010) Punturo, M., et al. 2010, Class. Quant. Grav., 27, 194002, doi: 10.1088/0264-9381/27/19/194002
  * Regimbau et al. (2017) Regimbau, T., Evans, M., Christensen, N., et al. 2017, Phys. Rev. Lett., 118, 151105, doi: 10.1103/PhysRevLett.118.151105
  * Regimbau & Hughes (2009) Regimbau, T., & Hughes, S. A. 2009, Phys. Rev. D, 79, 062002, doi: 10.1103/PhysRevD.79.062002
  * Regimbau et al. (2012) Regimbau, T., et al. 2012, Phys. Rev. D, 86, 122001, doi: 10.1103/PhysRevD.86.122001
  * Reitze et al. (2019a) Reitze, D., et al. 2019a, Bull. Am. Astron. Soc., 51, 035.  https://arxiv.org/abs/1907.04833
  * Reitze et al. (2019b) —. 2019b, Bull. Am. Astron. Soc., 51, 141.  https://arxiv.org/abs/1903.04615
  * Relton & Raymond (2021) Relton, P., & Raymond, V. 2021, Phys. Rev. D, 104, 084039, doi: 10.1103/PhysRevD.104.084039
  * Samajdar et al. (2021) Samajdar, A., Janquart, J., Van Den Broeck, C., & Dietrich, T. 2021, Phys. Rev. D, 104, 044003, doi: 10.1103/PhysRevD.104.044003
  * Sathyaprakash et al. (2012) Sathyaprakash, B., et al. 2012, Class. Quant. Grav., 29, 124013, doi: 10.1088/0264-9381/29/12/124013
  * Sathyaprakash et al. (2019) Sathyaprakash, B., Buonanno, A., Lehner, L., et al. 2019, Bull. Am. Astron. Soc., 51, 251, doi: 10.48550/arXiv.1903.09221
  * Sharma (2017) Sharma, S. 2017, Ann. Rev. Astron. Astrophys., 55, 213, doi: 10.1146/annurev-astro-082214-122339
  * Skilling (2004) Skilling, J. 2004, in American Institute of Physics Conference Series, Vol. 735, Bayesian Inference and Maximum Entropy Methods in Science and Engineering, 395, doi: 10.1063/1.1835238
  * Skilling (2006) Skilling, J. 2006, Bayesian Analysis, 1, 833, doi: 10.1214/06-BA127
  * Sun et al. (2020) Sun, L., et al. 2020, Class. Quant. Grav., 37, 225008, doi: 10.1088/1361-6382/abb14e
  * Vallisneri (2008) Vallisneri, M. 2008, Phys. Rev. D, 77, 042001, doi: 10.1103/PhysRevD.77.042001
  * Vangioni et al. (2015) Vangioni, E., Olive, K. A., Prestegard, T., et al. 2015, Mon. Not. Roy. Astron. Soc., 447, 2575, doi: 10.1093/mnras/stu2600
  * Wang et al. (2023) Wang, Z., Liang, D., Zhao, J., Liu, C., & Shao, L. 2023.  https://arxiv.org/abs/2304.06734
  * Wang et al. (2022) Wang, Z., Zhao, J., An, Z., Shao, L., & Cao, Z. 2022, Phys. Lett. B, 834, 137416, doi: 10.1016/j.physletb.2022.137416
  * Yunes & Pretorius (2009) Yunes, N., & Pretorius, F. 2009, Phys. Rev. D, 80, 122003, doi: 10.1103/PhysRevD.80.122003


$0^{\prime}$ (such as in Figure 1). In this case one can find
$x\in\mathrm{supp}_{p}$, $y\in\mathrm{supp}_{p^{\prime}}$ and $\gamma:0\to 0$
passing at least once through every sites of $\\{1,\cdots,n\\}$ and with
$|\gamma|\leq K+1$. The contribution to the total weight of $\pi$ from such a
path is

$\Phi\big{(}\phi_{x}\nabla_{x_{\gamma_{1}}}^{e_{\gamma_{1}}}\phi\big{)}\Phi\big{(}\phi_{y}\nabla_{x_{\gamma_{M-1}}}^{e_{\gamma_{M-1}}}\phi\big{)}\prod_{k=1}^{M-2}\Phi\big{(}\nabla_{x_{\gamma_{k}}}^{e_{\gamma_{k}}}\phi\nabla_{x_{\gamma_{k+1}}}^{e_{\gamma_{k+1}}}\phi\big{)}.$

up to additional factors from the other pairs, bounded from above by $c^{K}$
for some universal $c<\infty$. Proceeding as before, we obtain the following
upper bound on the target quantity:

$c^{K}\sum_{\begin{subarray}{c}x\in\mathrm{supp}_{p}\\\
y\in\mathrm{supp}_{p^{\prime}}\end{subarray}}\sum_{\gamma:x\to
y}\sum_{\begin{subarray}{c}x_{1},\cdots,x_{n}\\\
e_{1},\cdots,e_{n}\end{subarray}}\Big{|}\Phi(\phi_{x}\nabla_{x_{\gamma_{1}}}^{e_{\gamma_{1}}}\phi)\Phi(\phi_{y}\nabla_{x_{\gamma_{M-1}}}^{e_{\gamma_{M-1}}}\phi)\prod_{k=1}^{M-2}\Phi(\nabla_{x_{\gamma_{k}}}^{e_{\gamma_{k}}}\phi\nabla_{x_{\gamma_{k+1}}}^{e_{\gamma_{k+1}}}\phi)\Big{|}.$

The contribution of a fixed path $\gamma$ is then bounded from above as before
by $C^{M}$ times

$\sum_{x_{1},\cdots,x_{M-1}}(|x-x_{1}|+1)^{1-d}(|y-x_{M-1}|+1)^{1-d}\prod_{k=1}^{M-2}(|x_{k}-x_{k+1}|+1)^{-d}.$

Once again, repeated use of Lemma C.4 gives the summability of the infinite
sum, implying the claim in the case that there is no pairing in $\pi$ linking
$0$ to $0^{\prime}$. This concludes the proof. ∎

## References

  * [1] T. Balaban.  A low temperature expansion for classicalN-vector models. I. A renormalization group flow.  Communications in Mathematical Physics, 167(1):103–154, January 1995.
  * [2] T. Balaban.  A low temperature expansion for classicalN-vector models. II. Renormalization group equations.  Communications in Mathematical Physics, 182(3):675–721, December 1996.
  * [3] T. Balaban.  A Low Temperature Expansion for Classical N-Vector Models III. A Complete Inductive Description, Fluctuation Integrals.  Communications in Mathematical Physics, 196(3):485–521, September 1998.
  * [4] T. Balaban and M. O’Carroll.  Low Temperature Properties for Correlation Functions in Classical N-Vector Spin Models.  Communications in Mathematical Physics, 199(3):493–520, January 1999.
  * [5] M. Biskup.  Reflection positivity and phase transitions in lattice spin models, volume 1970 of Methods of Contemporary Mathematical Statistical Physics, Lecture Notes in Mathematics.  Springer-Verlag Berlin Heidelberg, 2009.
  * [6] J. Bricmont, J.-R. Fontaine, J. L. Lebowitz, E. H. Lieb, and T. Spencer.  Lattice systems with a continuous symmetry. III. Low temperature asymptotic expansion for the plane rotator model.  Comm. Math. Phys., 78(4):545–566, 1980.
  * [7] J. Bricmont, J.-R. Fontaine, J. L. Lebowitz, and T. Spencer.  Lattice systems with a continuous symmetry. I. Perturbation theory for unbounded spins.  Comm. Math. Phys., 78(2):281–302, 1980.
  * [8] C.M. Fortuin and Kastelein P.W.  On the random-cluster model: I. introduction and relation to other models.  Physica, 57(4):536–564, 1972.
  * [9] J. Fröhlich and P.-F. Rodriguez.  Some applications of the Lee-Yang theorem.  J. Math. Phys., 53:1–15, 2012.
  * [10] J. Fröhlich and P.-F. Rodriguez.  On cluster properties of classical ferromagnets in an external magnetic field.  J. Stat. Phys., 166:828–840, 2017.
  * [11] J. Fröhlich, B. Simon, and Thomas Spencer.  Infrared bounds, phase transitions and continuous symmetry breaking.  Communications in Mathematical Physics, 50:79–95, 1976.
  * [12] J. Fröhlich, R. B. Israel, E. H. Lieb, and B. Simon.  Phase transitions and reflection positivity. i. general theory and long range lattice models.  Comm. Math. Phys., 62(1):1–34, 1978.
  * [13] H.-O. Georgii.  Gibbs Measures and Phase Transitions.  De Gruyter, 2011.
  * [14] G. F. Lawler and V. Limic.  Random walk: a modern introduction, volume 123 of Cambridge Studies in Advanced Mathematics.  Cambridge University Press, Cambridge, 2010.
  * [15] B. Lees and L. Taggi.  Site monotonicity and uniform positivity for interacting random walks and the spin O(N) model with arbitrary N.  Comm. Math. Phys., Dec 2019.
  * [16] E. Michta and G. Slade.  Asymptotic behaviour of the lattice green function.  arXiv:2101.04717, 2021.
  * [17] E. Seiler.  Gauge Theories as a Problem of Constructive Quantum Field Theory and Statistical Mechanics, volume 159 of Lecture Notes in Physics.  Springer-Verlag Berlin Heidelberg, 1982.
  * [18] R.P. Soni.  On an inequality for modified Bessel functions.  Journal of Mathematics and Physics, 44:406–407, 1965.


# Polaris: The Mathematics of Navigation and the Shape of the Earth

John P. Boyd
Department of Climate & Space Sciences and Engineering
University of Michigan, 2455 Hayward Avenue, Ann Arbor MI 48109
<EMAIL_ADDRESS>
###### Abstract

For millenia, sailors have used the empirical rule that the elevation angle of
Polaris, the North Star, as measured by sextant, quadrant or astrolabe, is
approximately equal to latitude. Here, we show using elementary trigonometry
that Empirical Law 1 can be converted from a heuristic to a theorem.

A second ancient empirical law is that the distance in kilometers from the
observer to the North Pole, the geodesic distance measured along the spherical
surface of the planet, is the number of degrees of colatitude multiplied by
111.1 kilometers. Can Empirical Law 2 be similarly rendered rigorous? No;
whereas as the _shape_ of the planet is controlled by trigonometry, the _size_
of our world is an accident of cosmological history. However, Empirical Law 2,
can be rigorously verified by _measurements_. The association of 111 km of
north-south distance to one degree of latitude trivially yields the
circumference of the globe as 40,000 km. We also extend these ideas and the
parallel ray approximation to three different ways of modeling a Flat Earth.
We show that photographs from orbit, taken by a very expensive satellite, are
unnecessary to render the Flat Earth untenable; simple mathematics proves
Earth a sphere just as well.

Table 1:  List of Symbols

$\Omega$ is the angular frequncy of the Earth’s rotation. When two angles add
up to $\pi/2$, the angles are said to be “complementary”; the cosine of one
member of the pair is equal to the sine of the other.

The “declination [angle]” is the angular distance of a point north or south of
the celestial equator.  $d$ length of a geodesic from the observer’s position
to the North Pole. $D$ distance from observer’s position to the Pole Star,
Polaris $D^{star}$ distance from observer’s position to target star $E^{disk}$
equatorial circumference $E^{globe}$ equatorial circumference $L$ thickness of
the layer of atmospheric refraction $W$ scaling for the nondimensional
distance-to-North Pole coordinate $W^{earth}$ width (radius) of spherical
earth $W^{star}$ width (radius) of star $\eta$ angle between the central ray
from star to earth and the least parallel ray $\epsilon$ elevation [angle],
angular distance between the horizon and the target star. More precisely,
$\epsilon$ is the angle between a ray from the center of the earth to the star
and a second ray of the same longitude but in the equatorial plane
$\epsilon=\pi/2-\zeta$ , synonym: “altitude angle” $\epsilon_{app}$ apparent
elevation $\zeta$ zenith [angle], angular distance between directly overhead
[observer’s zenith] and the target star $\zeta=\pi/2-\epsilon$ $\varphi$
latitude, $\varphi=\pi/2-\theta$ $\theta_{i}$ angle of incoming ray, relative
to the zenith $\theta_{o}$ angle of outgoing ray $\theta$ nondimensional
distance-to-North-Pole $\theta$ colatitude, $\theta=\pi/2-\varphi$, [sphere
only] $\psi(\theta_{0},\theta_{1})$ angular width of arc between polar
coordinates $\psi(0,\theta)=\theta$ $\Omega$ $(2\pi/86400)$ per second;
angular frequency of earth’s rotation. $\mathcal{M}$ North-Pole-to-equator
distance $\rho$ air density $r$ radial coordinate in spherical coordinates
$\nu$ Eratoshtenes parameter, $\,(W^{earth}+W^{star})/D^{star}$ $\digamma$
refraction correction to the elevation angle

##  1 Polaris and the True North Star

The key measurement for navigation is the following.

######  Definition 1 (elevation angle)

The elevation angle $\epsilon$ is the angular distance between the horizon and
the target star. More precisely, $\epsilon$ is the angle between a ray from
the center of the earth to the star and a second ray in the same longitude but
in the equatorial plane. Alternatively, the elevation can be defined as the
complementary angle to the zenith angle $\zeta$: $\epsilon\equiv\pi/2-\zeta$.
A Synonyms for “elevation angle include“altitude angle” and simply “altitude”.

######  Definition 2 (zenith angle)

The zenith angle $\zeta$ is the angular distance between directly overhead
[observer’s zenith] and the target star. A second equivalent definition
$\zeta\equiv\pi/2-\epsilon$.

The “Elevation angle” $\epsilon$ is ninety degrees when the Pole Star is at
zenith and zero degrees when the Pole Star is bisected by the horizon.

Both angles can be measured by a variety of “inclinometers” [angle-measuring
devices] including the mariner’s astrolabe, quadrant, sextant and the Arab
instrument known as the _Kamal_.

As explained below, angle measurements of the Pole Star are especially useful
for navigation. However, one difficulty is that Polaris is 3/4 of a degree
from the celestial North Pole today, and was more than three degrees away in
Elizabethan times. However, mariners soon realized that by looking at stars
near Polaris and using the set of empirical rules with the charming name of
the “Regiment of the Stars”, the distance-from-the-North-Pole errors could be
greatly reduced.

######  Definition 3 (Regiment of the North Star)

A set of rules that allow sightings of stars such as Kochab to correct for the
difference between the position of Polaris and the true north celestial pole.

We use “Pole Star” denote a star that is not Polaris, but rather Polaris-
after-the-Regiment-of-the-North-Star has been applied. Our Pole Star is thus
aligned with the earth’s rotation axis and is the true north celestial pole.

Practical star sightings incur additional errors because of the limitations of
the sextant, the mariner’s astrolable and other inclinometers. Again,
correction tables and instruments were developed many centuries ago [7, 8].
The ful twelve-minute video [11] is a careful discussion of all these
corrections. For the videographer’s example, a sextant sighting of 29 degrees,
53.5 minutes is altered to a best estimate of 29 degrees, 21.3 minutes, a
correction of about half a degree or about 30 nautical miles.

The importance of angles for navigation has driven the development of planar
and spherical trigonometry [16, 17, 18, 19].

##  2 The First Empirical Law: Latitude Is Elevation

> For at least two millennia, navigators have known how to determine their
> latitude. Knowing the latitude of the desired destination, the navigator
> could sail north or south to that latitude and then sail east or west to
> reach that destination. To do this, it was necessary to have a method of
> measuring angles. In early days Arabs used one or two fingers at arms length
> to measure the angle between the horizon and Polaris. Later they used a
> _Kamal_ , which is a piece of cord with knots tied in it. This could be used
> to measure the angle between the horizon and Polaris. A knot could be tied
> in the cord as a measure of the homeport latitude before leaving, so the
> desired latitude was premeasured. Arabs tied knots on the cord at intervals
> of one _issabah_ , Arabic for “finger”, whic denotes 1 degree and 36
> minutes.
>
> In the 10th century AD, Arabs introduced to Europe the astrolabe and the
> quadrant. The quadrant spans 90 degrees and is divided into whole degrees. A
> plumb bob establishes the vertical. The quadrant was popular with Portuguese
> explorers in the 15th century. In addition to Polaris, they used
> observations of the Sun for determining latitude, particularly in the
> southern hemisphere.

— P. Kenneth Seidelmann in Sec. 7 of [13]

Waters’ treatise on navigation in the sixteenth and seventeenth centuries [20]
discusses the Polaris-latitude connection on pgs. 41 to 50.

Royal Air Force navigators used “latitude equals elevation” as a tool in
navigating Lancaster heavy bombers in night attacks on Germany during World
War II as described in [9], written by a retired war-time RAF navigator.

Thus, many centuries of navigation by land and sea provides overwhelming
evidence for the following.

Empirical Law 1: Latitude is the Elevation Angle of Polaris _In the Northern
Hemisphere, the latitude of the observer is approximately equal to the
elevation angle measured by the observer, or in other words_

$\displaystyle\epsilon\approx\varphi$ (1)

.

We show below that this heuristic rule is in fact a provable theorem.

Practical navigators are obliged to make small corrections to the observed or
“apparent” elevation angle to account for errors intrinsic to the sextant and
aslo atmospheric refraction [11], observer height above the sea (“dip angle”)
and so on. The observation and corrections collectively generate the “true”
elevation angle denoted by $\epsilon$; it is this corrected elevation that
appears in Empirical Law 1.

##  3 The Leveled Observer

The surface of the planet is corrugated with small scale curvature in the form
of innumerable hills and valleys. A planar earth has zero large scale
curvature, but a spherical earth is curved on a planetary scale as well as on
a small scale. To obtain consistent star measurements, we assume a “leveled
observer” who employs a sextant with a spirit level (or equivalent surveying
tricks) applied so that a vector from the surface to his zenith is parallel to
the “effective gravity” (gravity plus the centrifugal force due to the
planet’s rotation).

The effective gravity is not constant; a mountain will alter the gravitational
field in its immediate vicinity. However, the variations in Earth’s surface
gravitational field are less than one part in 15,000 of its mean value. The
corrugations will then have no effect on our measurements or results except
perhaps in the fourth decimal place. Geodesists and oil prospectors need to be
concerned about these variations, but we need not.

##  4 Empirical Law 2: The Distance to the North Pole Is the Product of
Colatitude With 6360 Kilometers

By the late 16th century, Englishmen knew that the ratio of distances at sea
to degrees were constant along any great circle such as the equator or any
meridian, 111[ A meridian is a Great Circle of constant longitude; it passes
through both poles.] Sssuming that Earth was a sphere. Robert Hues wrote in
1594 that the distance along a great circle was 60 miles per degree, that is,
one nautical, mile per arcminute (pg. 374 of [20]).

Edmund Gunter wrote in 1623 that the distance along a great circle was 20
leagues per degree. (Pg. 374 of [20]), which is the same statement in a
different length unit.

Pierre-Louis Maupertuis (1698-1759) was a major figure in physics for stating
the Principle of Least-Action. In 1736, he was appointed chief of the French
Geodesic Mission to Lapland to measure the length of a degree of latitude.
Simultneously, a second team was sent to Ecuador. The measurements
triumphantly confirmed Newton’s prediction that Earth is a (slightly!) oblate
spheroid. In 1799, the meter was defined to be one-quarter of the length of a
meridian divided by 10,000,000.

Let “geodesic distance” denote the arclength of the curve that connects two
points by the shortest path that lies in the surface of the planet. This
definition of the meter implies that one degree of latitude shortens the
geodesic distance to the North Pole by (10,000 kilometers/90) $\sim$ 111
kilometers, consistent with Hues, Gunter and other navigators. Collectively,
these and other astronomers over many centuries arrived at the following.

Empirical Law 2: Nautical Miles to the Pole and Degrees of Latitude _Each
increase in latitude $\varphi$ by one degree decreases the geodesic distance
to the North Pole from the observer by 60 nautical miles or equivalently 111
kilometers._ In symbols, this is equivalently

$\displaystyle d(\varphi)=W\,\left(\dfrac{\pi}{2}-\,\varphi\right)$ (2)

where $W$ is a length scale determined in the next section. It is not possible
to _prove_ this empirical law in the same way that Empirical Law 1 will be
elevated to a theorem below. However, it is possible to prove that Emprical
Law 2 can be _measured_ as true.

##  5 A Trivial Derivation of the Radius of the Earth

The distance $d$ along the geodesic distance (great circle arc) from the
observer to the North Pole has no obvious length scale. Nevertheless, it is
still convenient to write $d$ as the product of a nondimensional distance and
an arbitrary length scale $W$. The latitude $\varphi$ is known for millions of
cities and geographical features and can be easily determined by a sextant
measurement of the elevation $\epsilon(\text{Pole Star})$ via
$\varphi=\epsilon(\text{Pole Star})$. It is therefore convenient to use
latitude as the nondimensional measure of distance from the observer to the
North Pole. The formula for $d$ becomes

$\displaystyle gin{eqnarray}~{}d(\varphi)$ $\displaystyle=$ $\displaystyle
W\,\left(\pi/2-\varphi\right).$ (3)

Empirical Law 2 requires that the pole-to-equator geodesic distance is the
length associated with one degree of latitude (111.1 km) multiplied by the
number of degrees of latitude (90) between the equator and pole, yielding
$d(\pi/2)=10,000\text{km}$. Eq. (3) becomes

$\displaystyle 10,000\text{km}$ $\displaystyle=$ $\displaystyle
W\dfrac{\pi}{2}$ (4)

which demands that

$\displaystyle W$ $\displaystyle=$ $\displaystyle 6360\,\text{km}$ (5)
$\displaystyle=$ $\displaystyle W^{earth}\,$ (6)

Note that $W=W^{earth}$ was obained without any presuppositions of planetary
shape or curvature; the radius of the spherical Earth has _emerged
spontaneously_.

##  6 Parallel Ray/Distant Source Approximation

> “The historical assumption of the sun rays parallelism is far from being a
> spontaneous assumption for the students.”

— by Nicolas Décamp and Cécile de Hosson on pg. 919 [6]

Whether light from a star falls upon the Earth or the Moon as a beam of
_parallel rays_ or as a beam of _divergent_ rays has a profound impact upon
the consequences. A perfectly parallel beam is impossible, but a beam of
_almost_ -parallel rays will nonetheless generate a much different situation
than a beam of strongly divergent rays. It is therefore essential to quantify
the meaning of “almost-parallell”.

Let $D^{star}$ denote the distance from the observer’s position to the target
star. Let $W^{star}$ and $W^{earth}$ be the widths (radii) of the light-
emission source (star) and absorber (Earth). A ray from the center of the star
to the center of the surface of the earth (“central ray”) is shown as the
dashed line in Fig. 1.

The ray which is the “least parallel” (of those hitting the planet) in the
diagram is the solid ray with the arrowhead. Its vertical range, as
represented by the right side of the triangle in the lower half of Fig. 1, is
the radius of the star plus the radius of the planet. Let $\eta$ denote the
angle in radians between the central ray and the least parallel ray. This is
an upper bound on the angles between rays from star to planet.

Figure 1: Top: Schematic of an emission source [star or the Sun] (yellow
disk), Earth (blue disk) and light rays. The dashed line with an arrowhead is
the path of a ray from the center of the star to the center of the Earth; we
shall call this the “central ray”. The degree of parallelness of other rays is
measured by the angle between the ray and the central ray. The least parallel
ray that hits the Earth is one which, while traveling a distance $D^{star}$ (
horizontal in the diagram) from the star to the Earth is moving in the
perpendicular direction (vertical in the figure) from one side of the star to
the opposite side of the planet. This implies that this least-parallel ray is
the hypotenuse of the triangle in the lower half of the diagram. The
horizontal top side of the triangle is the path of the central ray with length
$D^{star}$. The angle $\eta$ between these two rays is an upper bound on the
angle between the central ray and all other rays that pass from star to
planet.

######  Definition 4 (Nearly Parallel Rays)

All the rays from the star to the Earth will be “nearly parallel” if and only
if the bound-on-angles $\eta\ll 1$.

It will prove useful to define the “Eratosthenes parameter” as

$\displaystyle~{}\nu\equiv\frac{W^{star}+W^{earth}}{D^{star}}\qquad\text{[Eratosthenes
parameter]}$ (7)

######  Theorem 1 (Spread of Angles and the Eratosthenes Parameter $\nu$)

In terms of the Eratosthenes parameter $\nu$,

$\displaystyle~{}\eta$ $\displaystyle=$
$\displaystyle\arctan\left(\frac{W^{star}+W^{earth}}{D^{star}}\right)$ (8)
$\displaystyle=$ $\displaystyle\arctan(\nu).$ (9)

Proof:

The Law of Sines applied to the triangle in Fig. 1 yields [18]

$\displaystyle\eta$ $\displaystyle=$
$\displaystyle\arctan\left(\frac{W^{star}+W^{earth}}{D^{star}}\right),$ (10)

which is (8). The second line follows from recognizing that the argument of
the arctangent is just the Eratosthenes parameter $\nu$ as defined by (7).
$\blacksquare$

######  Theorem 2 (Angles for Nearly Parallel Rays)

If the distance to the star $D^{star}$ is large compared to the larger of the
radius of the radiating star and the radius of the target planet, then

  1. 1.

$\displaystyle\nu\ll 1$ (11)

  2. 2.

$\displaystyle\eta\approx\nu$ (12)

  3. 3.

$\displaystyle\eta\ll 1$ (13)

This implies that all rays emitted by the star and falling on the planet will
be “nearly parallel” as defined above.

  4. 4.

In the limit $\nu\rightarrow\,0$ or equivalently, $D^{star}\rightarrow\infty$,
the rays are parallel.

Proof: Recall that the definition of the Eratosthenes parameter is (7 )
$\nu\equiv(W^{star}+W^{earth})/D^{star}$. When $D^{star}$ is large compared to
the radii of the star and planet, it follows that $\nu\ll 1$, which is the
first proposition.

The Taylor expansion of the arctangent function, (4.4.42) of [1], transforms
(8), $\eta=\arctan\left(\nu\right)$, into

$\displaystyle\eta$ $\displaystyle=$
$\displaystyle\sum_{j=0}^{\infty}\,(-1)^{j}\,\dfrac{1}{1+2j}\,\nu^{1+2j}$ (14)

Retaining just the term of smallest degree gives, with a relative error of
$O(\nu^{2})$, $\eta\approx\nu$, which is the second proposition of the
theorem. This in turn requires that small $\nu$ implies small $\eta$, which is
the third proposition. The fourth part of the theorem follows from taking the
limit in the approximation $\eta=\nu+O(\nu^{3})$ and then recognizing that the
limit of “nearly parallel waves”, differing in angle from the central ray by
no more than $\eta$, is _parallel_ waves. $\blacksquare$

The nearly parallel approximation may be alternatively labeled the “Distant
Source” approximation because it becomes more and more accurate as $D^{star}$
increases.

For all stars, $\nu$ is tiny because the distance to the star is huge compared
to the diameter of the star.

For Polaris itself, $D^{star}=4\times 10^{15}$ km and $W^{star}=3\times
10^{6}$ km yield

$\displaystyle\nu\approx 0.75\,\times 10^{-9}\qquad$ (15)

For the Sun and earth combination, $D^{star}=1.496\times 10^{8}$ km, and
$W^{star}=6.96\times 10^{5}$ km,

$\displaystyle\nu\approx 0.0046\approx 1/215$ (16)

The smallness of $\nu$ vindicates Eratosthenes, who assumed parallel (or
nearly parallel) rays from the sun in his measurement of the radius of the
earth.[6, 21, 10]

Figure 2: Parallel rays falling on a Flat Earth (left) and a spherical planet
(right). All the observers on the Flat Earth measure exactly the same
elevation angle of a given star, regardless of observer location. In contrast,
observers at different points on the sphere measure _different_ elevation
angles. The elevation angle for observer C on the right is around 90 degrees
where the elevation angle for observer A is near 0; for him the rays are low,
just above the horizon.

######  Theorem 3 (Flat Earth Parallel Rays)

On a Flat Earth when the rays from the star are parallel, the elevation of
Polaris (and every other star) is _independent_ of the observer’s latitude.

Proof: Obvious from the left half of Fig. 2. $\blacksquare$

The Flat Earth-parallel ray prediction that every observer sees Polaris at 90
degrees of elevation is in gross contradiction to observations. Thus, the Flat
Earth is possible only if the parallel ray approximation _fails_ , which
requires that the distance to the star is of comparable magnitude to the
larger of the radii of the earth and the star.

In contrast, if the parallel ray approximation is accurate, observers on a
_spherical_ earth measure _different_ elevation angles depending on latitude.
Let us suppose that all observations are made on the same meridian as that of
the target star.

The stars appear to rotate about the Pole Star; a time-lapse photograph gives
a so-called “star trails” plot in which each star (except the Pole Star)
traces a circle. Each point on the star trail gives a different angle; the
desired elevation angle of a star is the angle of the point on its star trail
which is closest to the observer’s zenith, i. e., the highest point on the
star trail.

Figure 3: Parallel rays falling on a spherical planet. Two rays from the
center of the planet through the observers define an angle $\psi$. Observers B
and C will measure different elevation angles for the star. Because the rays
from the star (and indeed any star) are parallel, the difference
$\epsilon_{C}\,-\,\epsilon_{B}$ is due entirely to the fact that a line
segment [red] parallel to the horizontal direction of Observer B is the same
as for Observer C except for a rotation through the angle $\psi$.

Fig. 3 shows why different observers on a meridian of a spherical planet
obtain different angles for the beam of approximately parallel rays from a
target star. The horizontal direction for each observer can be represented by
a vector that is orthogonal to the local effective gravity vector. The crucial
point is that as we move from one observer to another, this vector must
rotate, and the rotation angle must be added to the elevation angle of the
previous observer to obtain the elevation angle at the location of the current
observer as shown in the figure. Thus, observer B measures an elevation angle
$90-\psi$ where $\psi$ is about 45 degrees while observer C has an elevation
angle of $\epsilon_{C}=90\,\text{degrees}$. The difference in elevation angle
arises because B and C are separated by about one-eighth of the circumference
of the planet or equivalently by an angle of 360/8 degrees. In general, the
measured elevation angles of two observers differ by an angle $\Psi$ which is
360 degrees multplied by the geodesic distance between the two observers
divided by the circumference of the earth.

##  7 Elevation Is Latitude

###  7.1 Conversion of an Empirical Law to a Theorem

######  Theorem 4 (Elevation Angle Is Latitude)

Let $\varphi$ denote latitude and $\epsilon$ denote the “elevation angle”,
also called the “altitude angle” or “altitude”. If the Pole Star is observed
from a spherical planet and the distance to the Pole Star is large compared to
the radius of the Pole Star (allowing the parallel rays approximation), then
after $\epsilon$ is corrected for the small errors of the sextant, the slight
displacement of Polaris from the northern Celestial Pole, and atmospheric
refraction, etc.

$\displaystyle\varphi=\epsilon$ (17)

Proof: To handle general positions of stars on the Celestial Sphere, rotate
the latitudinal coordinate so that C is at a latitude of 90 degrees N. in the
new coordinate $\breve{\varphi}$. (Longitude is not modified.) Observer B is
at a $\breve{\varphi}$ that differs from C’s by $\psi$ degrees. Therefore,

$\displaystyle\epsilon_{B}=90\,\,-\,\psi\qquad[\text{degrees}]$ (18)

But the latitude of observer B is also $90-\psi$. It follows that latitude and
elevation angle are always equal, which is the proposition of the theorem.
$\blacksquare$

This theorem is identical with Empirical Law 1, which is thus converted from
an empirical relationship — highly probable but vulnerable to disproof by as
yet undisovered counterexamples — to a mathematical theorem.

###  7.2 Sextant Corrections

As noted in the introduction, the accuracy of the sextant measurement can be
usefully improved by applying a number of corrections. The details can be
found in any handbook of celestial navigation [5, 4] or videos such as [12,
11]; it would take us too far afield to describe them here. Instead, the
correction for atmospheric refraction [next subsection] will have to
illustrate the idea of corrections for all species of corrections.

The important point is that most of these errors can be greatly reduced by
using the tables of the _Astronomical Almanac_ , a joint production of the
U.S. Nautical Almanac Office, U. S. Naval Observatory, Her Majesty’s Nautical
Almanac Office and the U. K. Hydrographic Office. The printed version contains
precise positions over time of astronomical objects both natural (stars,
planets, etc.) and man-made (artificial satellites) for a given year. This
ephemeris serves as a world-wide standard for such information. The online
version [http://asa.hmnao.com/] extends the printed version by providing data
in machine-readable form. A good source for the theory behind the _Almanac_ is
[15]; the history and pre-history of navigational almanacs is described in
[14]. The videos [12, 11] walk the reader through examples.

###  7.3 Atmospheric Refraction

When the material properties of the medium vary, this induces variations in
the speed of light, smoothly curving or discontinuously bending the
trajectories of light rays as discussed further in Subsection 11.1. In
particular, density variations in the atmosphere alter the speed of light.
Unfotunately, the density of the atmosphere varies with the seasons, local
weather phenomena and the exponential decrease of density required by
hydrostatic equilibrium. Without density measurements, the best one can do is
to apply a simple analytic approximation that represents a time-average of
refraction.

Let $\epsilon_{app}$ denote the apparent elevation angle. Let $\epsilon$
denote the true elevation (or synonymously,“ true altitude”). Let $\digamma$
denote the correction to $\epsilon_{app}$. The simplest analytical formula for
refraction is, in radians,

$\displaystyle~{}\digamma$ $\displaystyle=$
$\displaystyle-0.000293\cot(\epsilon_{app})$ (19)

which is given on the first page of [2]. It is not recommended when
$\epsilon_{app}$ is smaller than 15 degrees.

Bennett [2] compares many analytic approximations for the refraction angle
$\digamma(\epsilon_{app}$). We apply his method H, a more complicated formula
than (19), which is widely used, notably in the U. S. Naval Observatory’s
_Vector Astrometry_ software.

$\displaystyle\digamma\,=0.000291\,\cot\left(\epsilon_{app}+\dfrac{7.31}{\epsilon_{app}+4.4}\,\right)$
(20)

Note that the argument of the cotangent function must be computed in degrees
and then converted to radians.

Once the refraction $\digamma$ has been calculated, the best approximation to
the elevation is

$\displaystyle\epsilon\approx\epsilon_{app}+\digamma$ (21)

where $\epsilon_{app}$ is the “apparent” elevation and $\epsilon$ is the
“true” elevation, which is the sextant measurement modified by all non-
refractive corrections. This is also the best approximation to latitude,
$\varphi\approx\epsilon$.

Fig. 4 and Table 2 show that the refractive correction is less than an
arcminute for all $\epsilon_{app}$ greater than 45 degrees, which implies a
navigational error, if the refractive correction is omitted, of only one
nautical mile. The correction grows as the elevation angle decreases, which is
why navigators are advised to measure Polaris and other celestial objects only
when the target is at least 15 degrees above the horizon. However, $\epsilon$
is a good approximation to latitude even when Polaris is within a degree or
two of the horizon.

Figure 4:  Corrected and uncorrected elevation angles and their difference on
the full range of the elevation angle,
$\epsilon_{app}\in[0\,\text{degrees},90\,\text{degrees}]$ (top). Not that the
uncorrected and refraction-corrected curves are graphically indistinguishable
and the dotted error curve is invisible underneath the horizontal axis. Bottom
left: the same but on the reduced range
$\epsilon_{app}\in[0\,\text{degrees},10\,\text{degrees}]$. The lower right
graph plots $\epsilon-\epsilon_{app}$, also on the shortened range
$[0,10\,\text{degrees}]$. Table 2:  Refraction corrections using the formula
“H” of Bennett [2]

Note that refraction in nautical miles (rightmost column) is numerically
identical to the refraction in arc minutes.

$\epsilon_{app}$ [degrees] |  $\epsilon$ [degrees] |  $\digamma$ [radians] |  $\digamma$ [degrees] |  $\digamma$ [nautical miles]
---|---|---|---|---
0.5 | 0.021 | \- 0.008 | -0.48 | -28.8
1 | 0.59 | -0.01 | -0.41 | \- 24.3
2 | 1.70 | -0.0071 | -0.30 | -18.2
5 | 4.84 | \- 0.0029 | -0.165 | -9.89
10 | 9.91 | \- 0.0016 | -0.090 | \- 5.39
15 | 14.9 | -0.0011 | -0.061 | -3.63
30 | 29.97 | -0.00050 | -0.029 | -1.72
45 | 44.98 | -0.00029 | -0.017 | -1.00
60 | 59.99 | -0.00017 | -0.0096 | -0.57
75 | 75.0 | -0.000077 | -0.0044 | -0.27

##  8 Why Flat Earth Needs a Nearby-Polaris and Non-parallel Light From the
North Star

On a Flat Earth, the horizontal levels of all Flat Earth observers are
parallel. If the light from a star like Polaris is also approximately
parallel, then, as seen in Fig. 2 _all Flat Earth observers_ will see Polaris
_at the same point_ in the celestial sphere. This is in gross
contradistinction to what is actually observed, $\epsilon\approx\varphi$,
Theorem 4.

Therefore, Flat Earth requires that Polaris be so close to Earth — a few
thousand kilometers or less— so that its rays are _not_ approximately
parallel. Different observers will see different elevation angles.222This
means jettisoning several independent measurements of the distance to Polaris
by parallax and other methods that show the actual distance to Polaris is over
400 light-years. However, a central tenet of all Flat Earth sects and
denominations is that all space agencies and professional astronomers are part
of a vast conspiracy and hoax, and their publications and findings must be
discarded. Flat Earth must be approached as fruitful for what-if scenarios,
but is otherwise the exponential of the exponential of exponential of
delusion. The close-Polaris hypothesis allows the computation of $D^{star}$,
the distance to Polaris, by triangularization as described in the next
section.

If, on the contrary, Polaris is very distant, then the following corollary to
Theorem 2 shows that calculation of $D^{star}$ is impossible.

######  Corollary 1 (Triangularization Fails in the Parallel Limit)

In the limit of parallel rays, it is impossible to obtain an estimate of
$D^{star}$ from triangulation because the relevant triangles are infinitely
long and narow, degenerating into rays. This conclusion is _independent_ of
the shape of the planet.

Proof: Take the limit, for fixed $W^{star}$ and $W^{earth}$, of
$D^{star}\rightarrow\infty$ in Fig. 1. $\blacksquare$

##  9 Triangularization On the Flat Earth for a Local Polaris

Let $D$ denote the distance from the North Pole to the Pole Star.
Equivalently, $D$ is the length of a geodesic from the star to the North Pole.
Because the star is at the zenith of an observer standing on the North Pole (a
“leveled observer” as defined earlier), the line segment is perpendicular to
the surface of Flat Earth.

Let $d$ denote the distance from the observer’s position to the North Pole.
Let $\epsilon$ denote the elevation angle.

Fig. 5 shows a triangle whose vertices are (i) the observer’s position on the
surface of the planet (ii) the North Pole and (iii) the Pole Star. It is a
good first approximation to pretend Polaris is the Pole Star because the error
is only three-quarters of a degree.

Figure 5: Triangulating the Pole Star on on a Flat Earth.

The “Law of Triangles” [18] then shows that the ratio of the two shortest
sides is equal to the tangent of the angle $\epsilon$ between these two sides,
or in other words

$\displaystyle~{}D=d\tan(\epsilon)$ $\displaystyle\leftrightarrow$
$\displaystyle d=D\cot(\epsilon)$ (22)

We shall apply the left equation in Flat Earth Case I and the right equation
[$d(D,\epsilon)$] in Flat Earth Case II.

##  10 Flat Earth Case I : $d$ [distance to pole] is a linear function of
elevation angle $\epsilon$

$\displaystyle\arctan(D/d)\qquad\Leftrightarrow\qquad\epsilon=\operatorname{\arctan}\left(\dfrac{D}{d}\right)\qquad\text{Pole-
Star-on-the-Flat-Earth}$ (23)

where $D$ is the distance to Polaris and $d$ is the distance from the
observation point to the North Pole. This can be rewritten to give the unknown
$D$ in terms of the known quantities $d$ and $\epsilon(d)$:

$\displaystyle D=\,d\,\tan(\epsilon)\qquad\Leftrightarrow\qquad
D=d\,\tan\left(\epsilon(d)\right)\qquad\text{Pole-Star-on-the-Flat-Earth}$
(24)

Substituting

$\displaystyle d(\varphi)$ $\displaystyle=$ $\displaystyle
W\,\left(\dfrac{\pi}{2}\,-\,\varphi\,\right)$ (25)

where $W=6360$ km and replacing elevation angle $\epsilon$ by latitude
$\varphi$ yields

$\displaystyle~{}D(\varphi)=\dfrac{\pi}{2}W\,\left(\pi/2-\varphi\,\right)\,\tan\left(\varphi\,\right)$
(26)

This function is plotted in Fig. 6. Table 3 lists numerical values of the
distance to Polaris for different latitudes of the observer.

Figure 6: Flat Earth Case I: The plot shows the distance $D$ from the North
Pole to the Pole Star, scaled by $W=6360$ km, so that the nondimensional
quantity plotted is $D/W=(\pi/2-\varphi)\,\tan(\varphi)$. The figure assumes
that $d=W\,\left(\pi/2-\varphi\right)$.

For a Flat Earth evangelist, the very bad news is that Eq. 26 predicts a
_different_ location for Polaris for every latitude. An observer at 75 degrees
N. latitude sees Polaris more than 6220 km above the North Pole whereas an
observer in the tropics at 15 degrees N. sees Polaris only 2233 km above the
Flat Earth. One is reminded of the cartoon character Multiman, who could
almost instantly create an army of clones. Even his magic powers fail here,
however: placing clones of Polaris at both 2233 km and 6220 km above the
planet would mean the observer now sees _two_ Polarises. Since the number of
observers is not restricted to two but may be arbitrarily large, we are
immediately confronted with an _infinite_ number of yellow giant clones of
Polaris.

To amplify the absurdity, $\max(D)=W=6360$ km. This implies that Polaris is
always _near_ and _local_ , but can a star so close (with a distance of 10,000
km at most)

  1. 1.

be a point source

  2. 2.

be a star of the _first magnitude_ over the entire daylight sky

  3. 3.

while invisible at night?

Table 3:  Distance $D$ from the North Pole to the North Star, Polaris. The
star cannot be in multiple places at the same time, so this Flat Earth model
fails.

Distances $D$ have been independently checked and confirmed in [3]. Analytical
Form Fourier Latitude (degrees) $\epsilon,\varphi$ (radians) $d$ $D$ [Distance
to Polaris] 1 degree $(1/180)\pi$ 9889 173 km 15 degrees $\pi/12$ 8333 km 2233
km 45 degrees $\pi/4$ 5000 km 5000 km 75 degrees $(5/12)\pi$ 1667 km 6220 km
89 degrees $(89/180)\pi$ 111.11 km 6360 km

##  11 Flat Earth Case II: Polaris Location $D$ is Independent of Latitude

Polaris cannot simultaneously be at every possible distance $D$ above the
North Pole from 0 to 6360 km. Can this unreality be fixed?

The linear relationship between $d$, the distance from the observer to the
North Pole, and latitude $\varphi$, which is
$d(\varphi)=W\,\left(\dfrac{\pi}{2}\,-\,\varphi\,\right)$ as derived in the
previous section, is supported by at least half a millenia of experience at
sea by those whose lives literally depended on it! Nevertheless, in the spirit
that desperate illnesses often require drastic remedies, we now explore a
_nonlinear_ relationship $d(\varphi)$.

The right equation in Eq. (22) gives a relationship between $d$, $D$ and
$\epsilon$. If $D$ is to be a constant so that Polaris has a single, unique
distance from the planet, independent of the observer’s latitude, then the
constraint from (22) becomes, with $D$ as a constant,

$\displaystyle d(\varphi)=D\,\cot(\varphi).$ (27)

as illustrated in Fig. 7. However, as $\varphi\rightarrow 0$ (the equator),
$d$ diverges to infinity.

Since $d(0)=\infty$, the Northern Hemisphere, latitude $\varphi\in[0,\pi])$ is
mapped to the entire infinite plane. Where, then, to put the Southern
Hemisphere? One could perhaps put the Southern Hemisphere on the other side of
the infinite plane. Then a journey across the equator — easy in our reality —
would become a journey of infinite length. The Pan-American Highway from
Alaska to Chile would be stretched to a never-ending strip of pavement.

We are left with insoluble difficulties:

  1. 1.

The Northern Hemisphere is an _infinite_ disk.

  2. 2.

There is no place to put the Southern Hemisphere.

  3. 3.

$\lim_{\varphi\rightarrow 0}d(\varphi)=\infty$.

  4. 4.

The nonlinear relationship between $d$ and $\varphi$ is refuted by many
centuries of navigation on land, sea and air.

Figure 7: Flat Earth Case II. The distance $d$ from the North Pole to the
observer at latitude $\varphi$, scaled by $W$, as illustrated for three
different values of the distance $D$ (in kilometers) to Polaris.

###  11.1 Flat Earth Case III: Flat Earth with Refracted Light Rays

Another possibility is that atmospheric refraction of the light from Polaris
can remove the absurdities of the two previous Flat Earth models. An important
constraint is that refraction is zero in a vacuum. Therefore, refraction will
be confined to a layer of thickness $L$ where $L$ is at most a few kilometers.
Equivalently, refractive light-bending occurs only in the lower and middle
atmosphere. Since the distances $D$ from the North Pole to Polaris have been a
few hundred to a few thousand kilometers in the earlier sections, we shall
make the plausible assumption that

$\displaystyle L<<D$ (28)

We shall suppose for simplicity that the index of refraction (i) is equal to
one, its vacuum value, everywhere above the refractive layer and (ii) varies
discontinuously at $z=L$:

$\displaystyle n=\left\\{\begin{array}[]{c}1,\qquad\,\qquad z>L\\\
n_{bot}(\varphi),\qquad z\,\leq\,L\end{array}\right.$ (31)

where $z$ is height above the surface of the planet.

Snell’s Law is usually stated in terms of the angles $\theta_{i}$ and
$\theta_{o}$ for the ray incoming from the Pole Star and the outgoing,
refracted wave, respectively, as shown in Fig. 8. Denoting the indices of
refraction in the top and bottom layers by $n_{top}$ and $n_{bot}$,
respectively, Snell’s law is

$\displaystyle\dfrac{\sin(\theta_{o})}{\sin(\theta_{i})}=\dfrac{n_{top}}{n_{bot}}$
(32)

Figure 8: Snell’s Law of Refraction for Flat Earth Case III for the Pole Star.
The horizontal dotted line is the boundary [interface] between the two layers
of air; the dashed line is a vertical guideline perpendicular to the
interface. $\theta_{1}$ and $\theta_{o}$ are the angles of the incoming and
outgoing light rays. $\phi$ is colatitude

However, we previously employed angles that are the complements of
$\theta_{i}$ and $\theta_{o}$, $\delta=\pi/2-\theta_{i}$ and
$\epsilon=\pi/2-\theta_{o}$. Furthermore $n_{top}=1$. Making these
substitutions and invoking the identity that the sine of an angle is equal to
the cosine of its complementary angle, Snell’s Law can be rewritten as

$\displaystyle
n_{bot}(\varphi)=\dfrac{\cos\left(\delta(\varphi)\right)}{\cos(\epsilon(\varphi)}$
(33)

To obtain expressions for $\delta(\varphi)$ and $\epsilon(\varphi),$ we
execute the following steps:

  1. 1.

Empirical Law 1 is

$\displaystyle\epsilon=\varphi$ (34)

  2. 2.

Empirical Law 2 gives the Polaris-to-North-Pole distance as

$\displaystyle d=W\,\left(\dfrac{\pi}{2}\,-\,\varphi\right)$ (35)

  3. 3.

Trigonometric identities applied to the triangle with vertices at the North
Star, the North Pole and the observer’s position gives

$\displaystyle\delta=\operatorname{arctan}\left(\dfrac{W}{D}\left(\dfrac{\pi}{2}-\varphi)\right)\right)$
(36)

Then

$\displaystyle
n_{bot}(\varphi)=\dfrac{\cos\left(\operatorname{arctan}\left(\dfrac{W}{D}\left(\dfrac{\pi}{2}-\varphi\right)\right)\right)}{\cos(\varphi)}$
(37)

Figure 9: Flat Earth Case III. The downward-propagating light ray from Polaris
(broken line segment with chevrons) has an elevation angle $\delta$ in the
near-vacuum region $z\in[L,D]$. At $z=L$ [dashed line], the density jumps
discontinously with the density larger in the thin layer $z\in[0,L]$. The
light bends at the density jump so that the measured elevation angle is
altered by refraction from $\delta$ (for $z>L$) to give the smaller angle
$\epsilon$ seen by an observer on the ground. The observer is a distance $d$
from the North Pole. Note that although it is difficult to show these
triangles to scale [and these triagles are therefore omitted from the graph],
$D\gg L$ and $d\gg d_{L}$. These inequalities make it possible to accurately
apply triangle formulas with the approximation $L=0$, i. e., the refractive
layer were absent.

The index of refraction observed in the atmosphere is very close to one. The
National Institutes of Science and Technology offers an Engineering Metrology
Toolbox at emtoolbox.nist.gov. This includes a calculator to compute the index
of refraction for air as a function of temperature, humidity and pressure. One
finds that even for extreme pressure, density and temperature, the index of
refraction is less than 1.004. Furthermore, the atmosphere slows down the
incoming light; the index of refraction can never be less than 1. Thus, the
physically realizable values of the index of refraction are

$\displaystyle~{}1\,\leq\,n_{bot}\leq 1.004$ (38)

Fig. 10 shows the required index of refraction in the lower, refractive layer
for the Flat Earth to have a single, unique distance to Polaris. For almost
all observer latitudes, the needed index of refraction $n_{bot}$ is
unphysical, violating (38).

Figure 10: Flat Earth Case III. $n_{bot}(\varphi)$ for three different
geodesic distances $D$ of Polaris from the North Pole. Values with the index
of refraction less than one, which are not physically realizable, are marked
with disks.

A little reflection will show that the simplifying assumptions made earlier
have no effect on the conclusion: Refraction fails to salvage the Flat Earth.

##  12 Summary

For at least half a millenia, navigators have used the empirical rule that the
elevation of the Pole Star (Polaris) is equal to latitude. We have shown that
proposition can be elevated to a proved theorem.

A second ancient empirical law is that the distance in kilometers from the
observer to the North Pole, the geodesic distance measured along the spherical
surface of the planet, is the number of degrees of colatitude multiplied by
111.1 kilometers. This cannot be proved _a priori_ , but we show that this
relationship can be established by _measurments_.

However, once Empirical Law 2, that is, that distance to the pole is
colatitude multiplied by 111 km, is confirmed by measurement, we prove that
this rigorously and uniquely determines the circumference of the Earth to be
40,000 km.

We also extend these ideas and the parallel ray approximation to three
different ways of modeling a Flat Earth. All fail. As noted by many previous
videographers, too numerous to list, triangularization gives different
distances from Earth to Polaris for each latitude. Absurd!

Photographs from space, taken on a very expensive satellite or astronauts on
the moonm

are unnecessary. Simple mathematics and ancient empirical laws prove Earth a
sphere just as well.

## References

  * [1] M. Abramowitz and I. A. Stegun, Handbook of Mathematical Functions, Dover, New York, 1965.
  * [2] G. G. Bennett, The calculation of astronomical refraction in marine navigation, J. Navigation, 35 (1982), pp. 255–259.
  * [3] Blue Marble Science, Mitchell from Australia and Trigonometry Just Don’t Get Along.  video, September 27 2021.  Channel=Blue Marble Science , duration=15:50.
  * [4] D. Burch, Celestial Navigation: A Complete Home Study Course, Second Edition, Starpath Publications, 2d ed., 2015.  296 pp.
  * [5] T. Cunliffe, Celestial Navigation: Learn How to Master One of the Oldest Mariner’s Arts Paperback – Illustrated, Fernhurst Books, 3rd ed., 2010\.  80 pp.
  * [6] N. Décamp and C. de Hosson, Implementing Eratosthenes’ Discovery in the Classroom: Educational Difficulties Needing Attention, Science and Education, 21 (2012), pp. 911–920.
  * [7] H. C. Freiesleben, Early pole star tables, J. Navigation, 8 (1955), pp. 285–288.
  * [8] J. B. Hewson, A History of the Practice of Navigation, Brown, Son and Ferguson, Glasgow, 1951.
  * [9] P. Hoare, Navigating Lancasters in WW II, J. Navigation, 60 (2007), pp. 201–215.
  * [10] T. V. Papathomas, Celestial illusions and ancient astronomers: Aristarchus and Eratosthenes, in Proceedings of SPIE 5666, Human Vision and Electronic Imaging X, B. E. Rogowitz, T. N. Pappas, and S. J. Daily, eds., New York, 2005, SPIE, SPIE, pp. 4–8.
  * [11] Practical Navigator, Basic latitude by Polaris.  video, Sept 27 2017.  Channel=Practical Navigator , duration=11:59.
  * [12] Rigging Doctor , Can You Really FIND Your COORDINATES With a SEXTANT? — Sailing Wisdom.  video, August 6 2021.  channel=Rigging Doctor, duration= 12:28.
  * [13] P. K. Seidelmann, Before the nautical almanacs, in The History of Celestial Navigation: Rise of the Royal Observatory and Nautical Almanacs (Historical & Cultural Astronomy), P. K. Seidelmann and C. Y. Hohenkerk, eds., Springer, Cham, Switzerland, 2020, pp. 1–20.
  * [14] P. K. Seidelmann and C. Y. Hohenkerk, eds., The History of Celestial Navigation: Rise of the Royal Observatory and Nautical Almanacs (Historical & Cultural Astronomy, Springer, Cham, Switzerland, 2020.
  * [15] S. E. Urban and P. K. Seidelmann, eds., Explanatory Supplement to the Astronomical Almanac, University Science Books, New York, 3rd ed., 2012.
  * [16] G. Van Brummelen, The Mathematics of the Heavens and the Earth: The Early History of Trigonometry, Princeton University Press, Princeton, New Jersey, 2009.  352 pp.
  * [17] , Heavenly Mathematics: The Forgotten Art of Spherical Trigonometry, Princeton University Press, Princeton, New Jersey, 2013.  208 pp.
  * [18] , A Very Short Introduction to Trigonometry, vol. 626 of Very Short Introductions, Oxford University Press, Oxford, 2020.  192 pp.
  * [19] , The Doctrine of Triangles: A History of Modern Trigonometry, Princeton University Press, Princeton, New Jersey, 2021.  392 pp.
  * [20] D. W. Waters, The Art of Navigation in England in Elizabethan and Early Stuart Times, Yale University Press, New Haven, Connecticut, 1958.  696 pp.
  * [21] Yazdan and his class of tweens, Measure the circumference of the earth like Eratosthenes 240 BC.  video, Apr 3 2019.  channel=Mr Yazdan, duration=3:59.