A Geometric Origin for Quasi-periodic Oscillations in Black Hole X-Ray Binaries

Black hole X-ray binaries (BHXRBs) are systems with a primary black hole gravitationally bound to a nondegenerate companion star. These systems display transient behavior exhibiting high X-ray luminosities (L_X ∼ 10³⁸ erg s⁻¹) during the outburst state, lasting from a few days to many months, followed by a long quiescent state (L_X ∼ 10³⁰ erg s⁻¹; Remillard et al. 2006). The triggering of these X-ray outbursts has been modeled as an instability arising in the accretion disk when the accretion rate is not adequate for the continuous flow of matter to the black hole, and when a critical surface density is reached (Dubus et al. 2001). However, the disk instability model has not been able to explain the outbursts of much shorter or longer timescales; for example, BHXRB GRS 1915+105 has shown a high X-ray luminosity state for more than 10 yr (Fender & Belloni 2004). During the outburst phase, the X-ray intensity shows rapid variations with timescales ranging from milliseconds to a few seconds, which are most likely to arise in the proximity to the black hole (r ∼ r_I, where r_I represents the innermost stable spherical orbit (ISSO)). The power density spectrum (PDS) of the X-ray intensity, which is commonly used to probe this fast variability, exhibits distinct features called quasi-periodic oscillations (QPOs) during the outburst period with their peak frequency, ν₀, ranging from 0.01 to 450 Hz (Remillard et al. 2006; Belloni & Stella 2014). QPOs can be distinguished from other broad features of the PDS by their high-quality factor ν₀/FWHM ≳ 2. Hence, the study of properties and origin of QPOs in BHXRBs is crucial to understanding the properties of inner accretion flow close to the black hole, where general relativistic effects are ascendant.

QPOs in BHXRBs are categorized as low-frequency QPOs (LFQPOs) with ν₀ < 30 Hz, which are again classified as type A, B, or C based on their various properties, and high-frequency QPOs (HFQPOs) with ν₀ > 30 Hz (Motta 2016). These different types of QPOs are also known to show a remarkable association with various spectral states during the outburst phase (Fender et al. 2004; Remillard et al. 2006; Fender & Belloni 2012; Motta 2016). The launch of the Rossi X-ray Timing Explorer (RXTE) in 1995 with its high sensitivity significantly increased the detection of BHXRBs and made it possible to detect HFQPOs in their PDS in the late 1990s (Belloni & Stella 2014), for example, the detection of 300 and 450 Hz QPOs in GROJ 1655-40 (Remillard et al. 1999b; Strohmayer 2001a); QPOs in the range 102–284 Hz, at 188, 249–276 Hz and near 183, 283 Hz in XTEJ 1550-564 (Homan et al. 2001; Miller et al. 2001; Remillard et al. 2002); 67, 40, and 170 Hz in GRS 1915+105 (Morgan et al. 1997; Strohmayer 2001b; Belloni et al. 2006); 250 Hz in XTEJ 1650-500 (Homan et al. 2003); 240 and 160 Hz in H1743-322 (Homan et al. 2005; Remillard et al. 2006); and more. Some of these HFQPOs have been detected simultaneously along with their peak frequencies showing nearly 3:2 or 5:3 ratios, indicating a resonance phenomenon (Remillard et al. 2006; Belloni & Stella 2014). There is also an interesting case of BHXRB GROJ 1655-40, which showed three QPOs simultaneously—two HFQPOs and one type C LFQPO (Motta et al. 2014a). Understanding the origin of HFQPOs and their simultaneity has been the prime focus of observational studies as well as theoretical models.

The study of general relativistic effects is important for a theoretical understanding of the origin of QPOs and their connection with various spectral states during the X-ray outburst, as these signals appear to emanate very close to the black hole. Several existing models, based on the instabilities in the accretion disk and other geometrical effects, attempt to explain the origin of LFQPOs and HFQPOs. Most of these models assume that the disk inhomogeneities orbiting in the innermost regions of the accretion disk are the cause of high variability in the X-ray flux, resulting in QPOs in the PDS. A widely accepted model among them is the relativistic precession model (RPM; Stella & Vietri 1999; Stella et al. 1999), which ascribes two simultaneous HFQPOs to the azimuthal, ν_ϕ, and periastron precession frequencies, $\left({\nu }_{\phi }-{\nu }_{r}\right)$, and a third simultaneous type C LFQPO to the nodal precession frequency, $\left({\nu }_{\phi }-{\nu }_{\theta }\right)$, of a self-emitting blob of matter in the accretion disk. The RPM has been applied to the cases of BHXRBs GROJ 1655-40 (Motta et al. 2014a) and XTEJ 1550-564 (Motta et al. 2014b) to estimate the spin parameter and mass of the black hole, where they assumed the precession frequencies of nearly circular particle trajectories in the accretion disk around a Kerr black hole. Recently, in contrast with the localized assumption of the RPM, the most frequently detected type C QPOs in BHXRBs have been modeled as the Lense–Thirring frequency of a radially extended thick torus precessing as a rigid body (Ingram et al. 2009; Ingram & Done 2011, 2012). This model describes the increase in type C QPO frequency with the hard to soft spectral transition during outburst as coincident with the decrease in outer radius of the torus and also shows that the maximum type C QPO frequency should be close to 10–30 Hz (Motta et al. 2018). Other models concentrate on the 3:2 or 5:3 resonance phenomena of simultaneous HFQPOs under the regime of particle approach; for instance, the nonlinear resonance models (Kato 2004, 2008; Török et al. 2005, 2011) explain the phenomenon of simultaneous HFQPOs as an excitation due to nonlinear resonant coupling between the oscillations within the accretion disk. One such nonlinear resonance phenomenon is the parametric resonance between radial, ν_r, and vertical, ν_θ, oscillation frequencies of particles in the accretion disk (Abramowicz et al. 2003). Another explanation of HFQPOs is based on the Keplerian and radial frequencies of the deformation of the clumps of matter that is due to the simulated tidal interactions in the accretion disk (Germanà et al. 2009). A recent model involves the study of (magneto)hydrodynamic instabilities, in particular, to understand the 3:2 resonance of HFQPOs using the general relativistic and ray-tracing simulations (Tagger & Varnière 2006; Varniere et al. 2019).

The RPM takes into account the fundamental phenomenon of relativistic precession, which is dominant and inevitable in the strong-field regime around a black hole. Although the emission mechanism for the production of QPOs with a strong rms (∼20%) is hitherto unknown, it explains some important observational relations, for example, the Psaltis–Belloni–Klis relation (PBK; Psaltis et al. 1999), which is a positive correlation between the HFQPOs and the LFQPOs in different BHXRBs. In a few other BHXRBs, the characteristic frequency of a broad feature (not a QPO) in the PDS during the hard state shows the same correlation with the LFQPOs. This correlation has been explained using the RPM as a variation of the radius of origin around the Kerr black hole, tracing the QPO frequency.

In this paper, we expand the RPM from a restricted study of circular orbits and explore the fundamental frequency range of the nonequatorial eccentric, equatorial eccentric, and spherical particle trajectories around a Kerr black hole and associate them with the properties of QPOs. We call this the generalized RPM (GRPM). The general trajectory solutions around a Kerr black hole and their corresponding fundamental frequencies have been extensively studied before (Schmidt 2002; Fujita & Hikida 2009; Rana & Mangalam 2019a, 2019b). The existence of nonequatorial eccentric, equatorial eccentric, and spherical orbits near a rotating black hole is tangible, and hence the relativistic precession of these orbits can also be included in the model for the emission of QPOs. The quadrature form of the general trajectory solution {ϕ, θ, r, t} around a Kerr black hole (Carter 1968) and the corresponding fundamental frequencies {ν_ϕ, ν_r, ν_θ} (Schmidt 2002) are well known. Later, the complete analytic form for the trajectories and the fundamental frequencies was derived in terms of the Mino time (Mino 2003) and the standard elliptic integrals (Fujita & Hikida 2009). More recently, a more compact, analytic, and numerically faster form was derived, in terms of the standard elliptic integrals, for the particle trajectory solutions and their fundamental frequencies (Rana & Mangalam 2019a, 2019b). We use these analytic formulae for the fundamental frequencies via the GRPM for the periastron and nodal precession of nonequatorial eccentric, equatorial eccentric, and spherical trajectories around a Kerr black hole to associate them with the detected QPO frequencies. The RPM was previously predicted for circular {e = 0, Q = 0} orbits (Stella & Vietri 1999; Stella et al. 1999). We now include {$e\ne 0$, $Q\ne 0$} orbits in this paradigm and test the more general model in this paper. Finally, we show that the eccentric trajectory solution also satisfies the PBK correlation for the case of BHXRB GROJ 1655-40.

The paper is structured as follows. We first motivate the association of fundamental frequencies of the general eccentric and spherical trajectories with the QPOs in BHXRBs assuming the GRPM in Sections 2.1 and 2.2; see Figure 1 for the terminology used for eQ (general case), Q0 (spherical), e0 (eccentric equatorial), and 00 (circular orbits). We then take up the cases of BHXRBs M82 X-1, GROJ 1655-40, XTEJ 1550-564, 4U 1630-47, and GRS 1915+105, where HFQPOs have been discovered before. We discuss their observation history in Appendix D, and we discuss observations of each BHXRB that we use for our analysis in Section 3.1. Using the observed QPO frequencies in these BHXRBs, we calculate the corresponding orbital parameters. The method for the parameter estimation and its corresponding errors is discussed in Section 3.2 and in Appendix E. We discuss the results for general eccentric trajectories in Section 3.2.1, and those corresponding to the spherical orbit in Section 3.2.2. We also show in Section 4 that the PBK correlation is well explained by the eccentric trajectory solutions found in the case of BHXRB GROJ 1655-40. In Section 5, we compare our model with another model for the fluid flow in the general-relativistic thin accretion disk. We finally discuss and conclude our results in Section 6. A glossary of symbols used in this article is given in Table 1, and a concept flowchart of the paper is given in Figure 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 1.  Glossary of Symbols Used

Common Physical Parameters
c	Speed of light	G	Gravitational constant
M•	Mass of the black hole	a	Spin of the black hole scaled by $\left({{GM}}_{\bullet }^{2}/c\right)$
${ \mathcal M }$	M•/M⊙
Orbital Parameters
E	Energy per unit rest mass of the	Lz	z component of angular momentum
	particle, scaled by mc2		per unit rest mass of the particle,
			scaled by $\left({{GM}}_{\bullet }/c\right)$
L	Angular momentum per unit rest	Q	Carter's constant scaled by ${\left({{GM}}_{\bullet }^{2}/c\right)}^{2}$
	mass of the particle, scaled by $\left({{GM}}_{\bullet }/c\right)$
ra	Apastron distance of the orbit	rp	Periastron distance of the orbit
	scaled by $\left({{GM}}_{\bullet }/{c}^{2}\right)$		scaled by $\left({{GM}}_{\bullet }/{c}^{2}\right)$
e	Eccentricity parameter	μ	Inverse latus rectum parameter
rs	Radius of spherical orbit scaled	rI	ISSO radius scaled by $\left({{GM}}_{\bullet }/{c}^{2}\right)$
	by $\left({{GM}}_{\bullet }/{c}^{2}\right)$
Fundamental Frequencies
νϕ	Azimuthal frequency	νnp	Nodal precession frequency, $\left({\nu }_{\phi }-{\nu }_{\theta }\right)$
νθ	Vertical oscillation frequency	νpp	Periastron precession frequency, $\left({\nu }_{\phi }-{\nu }_{r}\right)$
νr	Radial frequency	ν0	Centroid frequency of the QPO
$\bar{\nu }$	Frequency scaled by the factor $\left({c}^{3}/{{GM}}_{\bullet }\right)$
Probability Analysis for Estimating Parameter Errors
P	Probability density (space)	${ \mathcal P }$	Normalized probability density (space)
${ \mathcal N }$	Normalization factor	${{ \mathcal J }}_{l}$	Jacobian of transformation from frequency
			to parameter space

Download table as:  ASCIITypeset image

The relativistic precession is a phenomenon that is due to strong gravity near a rotating black hole, and its consequence for QPOs originating very close to the black hole is studied. We motivate the association of QPOs in BHXRBs with the fundamental frequencies of general nonequatorial bound particle trajectories around a Kerr black hole through the GRPM. Figure 3 shows the periastron and nodal precession of an eccentric particle trajectory near the equatorial plane of a rotating black hole.

Zoom In Zoom Out Reset image size

We suggest that the instabilities in the inner region close to the rotating black hole might provide a radiating plasma cloud (it could be a blob or a torus with the collection of such trajectories degenerate in the parameter space) with enough energy and angular momentum to attain an eccentric ($e\ne 0$) trajectory, or a nonequatorial trajectory $(Q\ne 0$, Carter's constant, Carter 1968), or both simultaneously ($e\ne 0$, $Q\ne 0$). Carter's constant can be roughly interpreted as representative of the residual of the angular momentum in the x–y plane, $Q\propto {L}^{2}-{L}_{z}^{2}$, so we have Q = 0 for the equatorial orbits where L = L_z. We first try to find the suitable range for the parameters, {e, r_p, a, Q}, of these orbits that produce the fundamental frequencies to compare with the observed range of QPO frequencies in BHXRBs, where r_p represents the periastron point of the orbit and a represents the spin of the black hole. We divide our study of the trajectories into three categories (see Figure 1), where a particle follows one of these:

• 1.  
  A nonequatorial eccentric trajectory ($e\ne 0$, $Q\ne 0$) called eQ.
• 2.  
  An equatorial eccentric trajectory ($e\ne 0$, Q = 0) called e0.
• 3.  
  A nonequatorial and noneccentric, also called a spherical trajectory (e = 0, $Q\ne 0$), called Q0.

We are using dimensionless parameters (G = c = M_• = 1) as the convention in this article for simplicity, so ${r}_{p}\to {r}_{p}/\left({{GM}}_{\bullet }/{c}^{2}\right),{r}_{a}\to {r}_{a}/\left({{GM}}_{\bullet }/{c}^{2}\right),a\to J/\left({{GM}}_{\bullet }^{2}/c\right)$ and $Q\to Q/{\left({{GM}}_{\bullet }^{2}/c\right)}^{2}$, where J is the angular momentum and M_• is the mass of the black hole, and r_a is the apastron point of the bound orbit, while $e=\left({r}_{a}-{r}_{p}\right)/\left({r}_{a}+{r}_{p}\right)$, the eccentricity parameter, is dimensionless by definition (see Table 1). We also define another mass parameter ${ \mathcal M }={M}_{\bullet }/{M}_{\odot }$ scaled by solar mass for convenience. The most general nonequatorial trajectory (eQ) around a Kerr black hole comprises periastron precession in the orbital plane, superimposed on the precession of the orbital plane about the spin axis of the rotating black hole. Figure 4 shows one such trajectory around a Kerr black hole centered at the origin.

Zoom In Zoom Out Reset image size

There are a variety of bound Kerr orbits, for example, nonequatorial eccentric, separatrix, zoom-whirl, and spherical orbits, that have been systematically studied (e.g., Rana & Mangalam 2019a, 2019b and references within). Hence, here we first discuss the distribution of these orbits in the parameter space and then isolate the most plausible type of orbits, which should give us the observed range of QPO frequencies assuming the GRPM. A complete description of various types of trajectories is given in Table 2, where MBSO (MBCO) is the marginally bound spherical (circular) orbit, and ISCO is the innermost stable circular orbit. These bound orbits are distributed in particular regions in the parameter space and into different parameter ranges for different types of orbits. In Figure 5, we show how this distribution belongs in different regions in the (r, a) plane, where r = R/R_g represents distance from the black hole, and ${R}_{g}=\left({{GM}}_{\bullet }/{c}^{2}\right)$. These regions are separated by important radii, which are shown as various curves for the equatorial (Q = 0) and nonequatorial (Q = 4) trajectories in Figure 5, where we see that the (un)stable bound orbits are found in regions 1, 2, and 3. Region 4 is beyond the light radius, which extends down to the horizon radius $[{r}_{+}=\left(1+\sqrt{1-{a}^{2}}\right)$], where bound particle orbits are not present, which means any particle in this region would plunge into the black hole, and region 5 is inside the horizon surface. Hence, we restrict our exploration of suitable parameters for required QPO frequencies to regions 1 and 2, where stable circular (spherical), equatorial (nonequatorial) eccentric, zoom-whirl, and separatrix orbits are found.

Zoom In Zoom Out Reset image size

These bound orbits can also be shown as a region in the (e, μ) space, which is defined as

$$\begin{eqnarray}&&e=\displaystyle \frac{{r}_{a}-{r}_{p}}{{r}_{a}+{r}_{p}},\ \ \ \ \ \ \ \ \ \mu =\displaystyle \frac{{r}_{a}+{r}_{p}}{2{r}_{a}{r}_{p}},\end{eqnarray} \tag{ 1 }$$

where r_a is the apastron point of the orbit. This bound orbit region is shown as a shaded region in Figure 6. The condition for these bound orbits is given by (Rana & Mangalam 2019a, 2019b)

$$\begin{eqnarray}&&\begin{array}{l}\left[{\mu }^{3}{a}^{2}Q{\left(1+e\right)}^{2}+{\mu }^{2}\left(\mu {a}^{2}Q-{x}^{2}-Q\right)\left(3-e\right)\right.\left.\left(1+e\right)+1\right]\geqslant 0,\end{array}\end{eqnarray} \tag{ 2 }$$

where μ can also be written as $\mu =1/\left[{r}_{p}\left(1+e\right)\right]$, where the equality sign corresponds to the separatrix trajectories. This bound orbit region shown in Figure 6 only includes regions 1 and 2 of the (r_p, a) plane shown in Figure 5. The RPM has been applied to two cases of BHXRBs, assuming the precession of nearly circular orbits (negligible eccentricity¹ ) in the equatorial plane of a Kerr black hole (Motta et al. 2014a, 2014b). In general, the observed range of HFQPOs in BHXRBs is 40–500 Hz, whereas that of type C LFQPOs is 10 mHz to 30 Hz (Remillard et al. 2006; Belloni & Stella 2014). The formulae for fundamental particle frequencies of nearly circular and equatorial orbits are given by Bardeen et al. (1972) and Wilkins (1972); see Appendix C for the derivation of these formulae from the general frequency formulae of e0 (Equation (5)) and Q0 (Equation (7)) orbits:

$$\begin{eqnarray}\begin{array}{rclcl}{\nu }_{\phi }\left(r,a\right) & = & \displaystyle \frac{{c}^{3}}{2\pi {{GM}}_{\bullet }}\displaystyle \frac{1}{\left({r}^{3/2}+a\right)},{\bar{\nu }}_{\phi }\left(r,a\right) & = & \displaystyle \frac{{\nu }_{\phi }}{\left({c}^{3}/{{GM}}_{\bullet }\right)}=\displaystyle \frac{1}{2\pi \left({r}^{3/2}+a\right)},\end{array}\end{eqnarray} \tag{ 3a }$$

$$\begin{eqnarray}\begin{array}{rclcl}{\nu }_{r}\left(r,a\right) & = & {\nu }_{\phi }{\left(1-\displaystyle \frac{6}{r}-\displaystyle \frac{3{a}^{2}}{{r}^{2}}+\displaystyle \frac{8a}{{r}^{3/2}}\right)}^{1/2},{\bar{\nu }}_{r}\left(r,a\right) & = & \displaystyle \frac{{\nu }_{r}}{\left({c}^{3}/{{GM}}_{\bullet }\right)},\end{array}\end{eqnarray} \tag{ 3b }$$

$$\begin{eqnarray}\begin{array}{rclcl}{\nu }_{\theta }\left(r,a\right) & = & {\nu }_{\phi }{\left(1+\displaystyle \frac{3{a}^{2}}{{r}^{2}}-\displaystyle \frac{4a}{{r}^{3/2}}\right)}^{1/2},{\bar{\nu }}_{\theta }\left(r,a\right) & = & \displaystyle \frac{{\nu }_{\theta }}{\left({c}^{3}/{{GM}}_{\bullet }\right)},\end{array}\end{eqnarray} \tag{ 3c }$$

where $\{{\bar{\nu }}_{\phi },{\bar{\nu }}_{r},{\bar{\nu }}_{\theta }\}$ are the dimensionless frequencies, where we use the convention a > 0 for the prograde and a < 0 for the retrograde orbits in this article. Using these formulae and assuming the RPM, it was retrodicted for BHXRB GROJ 1655-40 and XTEJ 1550-564 that these signals originated very close to and outside the ISCO radius, at nearly r = 5.677 ± 0.035 for GROJ 1655-40 and r = 5.47 ± 0.12 for XTEJ 1550-564 (Motta et al. 2014a, 2014b). We show that the expected QPO frequency range associated with the 00 orbits in the RPM {ν_ϕ, ν_pp ≡ (ν_ϕ − ν_r), ν_np ≡ (ν_ϕ − ν_θ)} is valid for a wide range of r, where ν_ϕ, ν_pp, and ν_np correspond to the HFQPO-1, HFQPO-2, and type C LFQPO, respectively.² To illustrate this, we present a mass-independent model of these frequencies. In Table 3, we have shown the observed range of the HFQPO and LFQPO frequencies in BHXRBs along with a typical range in dimensionless values {${\bar{\nu }}_{\phi }$, ${\bar{\nu }}_{{pp}}$, ${\bar{\nu }}_{{np}}$}, obtained by scaling the observed frequencies of HFQPOs in BHXRBs using the corresponding known value of the black hole mass ${ \mathcal M }\sim 5\mbox{--}10$ (given in Table 5). For a BHXRB, the typical frequency range of the type C QPOs is 10 mHz to 30 Hz, and we have scaled this frequency range with ${ \mathcal M }=10$ (a typical mass value for BHXRB) to obtain the dimensionless frequency range. This provides an expected range of the geometrical orbital parameters independent of the black hole mass that implies largely a range of r_p. Figure 7 shows the contours of ${\bar{\nu }}_{\phi }$, ${\bar{\nu }}_{{pp}}$, and ${\bar{\nu }}_{{np}}$ for the 00 orbits, using Equations 3(a)–(c), in the (r, a) plane outside the ISCO radius (region 1 of Figure 5(a)). We find the following:

• 1.  
  The expected range of simultaneous QPO frequencies corresponds to a wide range of r ∼ 5–15 for the 00 orbits, which is typically the inner region of the accretion disk.
• 2.  
  The simultaneous QPOs, if associated with the 00 orbits, should originate very near the ISCO radius.
• 3.  
  We expect much higher QPO frequency values $\{{\bar{\nu }}_{\phi }\gtrsim 0.015,{\bar{\nu }}_{{pp}}\gtrsim 0.009,{\bar{\nu }}_{{np}}\gtrsim 0.001\}$ for the 00 orbits near the ISCO radius for a ≳ 0.5, as seen in Figure 7, which are outside the observed QPO frequency range.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 3.  Summary of the Observed QPO Frequency Range in BHXRBs and Their Corresponding Dimensionless Values Derived from Data Given in Table 5

Type of QPO	QPOs in the	Observed QPO Frequency	Dimensionless Frequency Range
	RPM and GRPM	Range in Hz	$\ \ \ \bar{\nu }=\nu \cdot {10}^{-3}/\left({c}^{3}/{{GM}}_{\bullet }\right)\ \$
HFQPO-1	νϕ	100–500	2–15
HFQPO-2	νpp	40–300	1–9
Type C LFQPO	νnp	10−2–30	10−5–1

Download table as:  ASCIITypeset image

Now, with this, we can explore the frequency range of the nonequatorial eccentric, equatorial eccentric, and spherical orbits using a similar approach assuming the GRPM in regions 1 and 2 of Figure 5 (shaded region of Figure 6).

2.1. Nonequatorial and Equatorial Eccentric Orbits: eQ and e0

We first discuss the useful formulae of the fundamental frequencies for the nonequatorial and equatorial eccentric particle trajectories derived in Rana & Mangalam (2019a, 2019b). Later, we use these formulae to explore the required frequency range for QPOs in BHXRBs, based on the GRPM, and determine the corresponding parameter range {e, r_p, a, Q} associated with these trajectories.

As shown in Figure 4, the orbital plane of a nonequatorial eccentric trajectory oscillates with respect to the spin axis of the black hole, along with the phenomenon of periastron precession taking place in the orbital plane. A complete analytic trajectory solution and the fundamental frequencies for such trajectories around a Kerr black hole were derived in terms of {e, μ, a, Q} parameters (Rana & Mangalam 2019a, 2019b), where μ is the inverse latus rectum of the orbit, and it can be written in terms of {e, r_p} as $\mu =1/\left[{r}_{p}\left(1+e\right)\right]$. The expressions of dimensionless fundamental frequencies for these trajectories are given by (Rana & Mangalam 2019a, 2019b)

$$\begin{eqnarray}&&{\bar{\nu }}_{\phi }\left(e,{r}_{p},a,Q\right)=\displaystyle \frac{\left[-{I}_{1}\left(e,{r}_{p},a,Q\right)-2{L}_{z}{I}_{8}\left(e,{r}_{p},a,Q\right)\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)+2{L}_{z}{I}_{8}\left(e,{r}_{p},a,Q\right){\rm{\Pi }}\left({z}_{-}^{2},\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)}{2\pi \left\{\left[{I}_{2}\left(e,{r}_{p},a,Q\right)+2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\left(e,{r}_{p},a,Q\right)\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\left(e,{r}_{p},a,Q\right)K\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ 4a }$$

$$\begin{eqnarray}&&{\bar{\nu }}_{r}\left(e,{r}_{p},a,Q\right)=\displaystyle \frac{F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)}{\left\{\left[{I}_{2}\left(e,{r}_{p},a,Q\right)+2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\left(e,{r}_{p},a,Q\right)\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\left(e,{r}_{p},a,Q\right)K\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ 4b }$$

$$\begin{eqnarray}&&{\bar{\nu }}_{\theta }\left(e,{r}_{p},a,Q\right)=\displaystyle \frac{a\sqrt{1-{E}^{2}}{z}_{+}{I}_{8}\left(e,{r}_{p},a,Q\right)}{2\left\{\left[{I}_{2}\left(e,{r}_{p},a,Q\right)+2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\left(e,{r}_{p},a,Q\right)\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\left(e,{r}_{p},a,Q\right)K\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ 4c }$$

where L_z is the z-component of a particle's angular momentum and E is its energy per unit rest mass, which can be explicitly expressed as the functions of {e, μ, a, Q} parameters (see Equations 5(a)–(e) in Rana & Mangalam 2019a). Here, ${I}_{1}\left(e,{r}_{p},a,Q\right)$, ${I}_{2}\left(e,{r}_{p},a,Q\right)$, and ${I}_{8}\left(e,{r}_{p},a,Q\right)$ are the radial integrals of motion given in their simplest analytic forms, along with the constants involved, by Equations 6(a)–(h), 7(a)–(l), 8(a)–(c), and 9(d) in Rana & Mangalam (2019a); $F\left(\varphi ,{p}^{2}\right)$, $K\left(\varphi ,{p}^{2}\right)$, and ${\rm{\Pi }}\left({q}^{2},\varphi ,{p}^{2}\right)$ used in Equations 4(a)–(c) are the standard elliptic integrals (Gradshteyn & Ryzhik 2007).

Next, in the case of equatorial eccentric orbits (e0), the expressions for the azimuthal and radial fundamental frequencies can be further reduced to a form simpler than Equations 4(a)–(c), which are given by (Rana & Mangalam 2019a, 2019b)

$$\begin{eqnarray}&&{\bar{\nu }}_{\phi }\left(e,{r}_{p},a\right)=\displaystyle \frac{{a}_{1}{\rm{\Pi }}\left(-{p}_{2}^{2},\tfrac{\pi }{2},{m}^{2}\right)+{b}_{1}{\rm{\Pi }}\left(-{p}_{3}^{2},\tfrac{\pi }{2},{m}^{2}\right)}{2\pi \left\{\tfrac{{a}_{2}}{2\left(1+{p}_{1}^{2}\right)}\left[-F\left(\tfrac{\pi }{2},{m}^{2}\right)+\tfrac{{p}_{1}^{2}\,K\left(\tfrac{\pi }{2},{m}^{2}\right)}{\left({m}^{2}+{p}_{1}^{2}\right)}\right]+{\rm{\Pi }}\left(-{p}_{1}^{2},\tfrac{\pi }{2},{m}^{2}\right)\left[{a}_{2}\tfrac{\left[{p}_{1}^{4}+2{p}_{1}^{2}\left(1+{m}^{2}\right)+3{m}^{2}\right]}{2\left(1+{p}_{1}^{2}\right)\left({m}^{2}+{p}_{1}^{2}\right)}+{b}_{2}\right]+{c}_{2}{\rm{\Pi }}\left(-{p}_{2}^{2},\tfrac{\pi }{2},{m}^{2}\right)+{d}_{2}{\rm{\Pi }}\left(-{p}_{3}^{2},\tfrac{\pi }{2},{m}^{2}\right)\right\}},\end{eqnarray} \tag{ 5a }$$

$$\begin{eqnarray}&&{\bar{\nu }}_{r}\left(e,{r}_{p},a\right)=\displaystyle \frac{1}{2\left\{\tfrac{{a}_{2}}{2\left(1+{p}_{1}^{2}\right)}\left[-F\left(\tfrac{\pi }{2},{m}^{2}\right)+\tfrac{{p}_{1}^{2}K\left(\tfrac{\pi }{2},{m}^{2}\right)}{\left({m}^{2}+{p}_{1}^{2}\right)}\right]+{\rm{\Pi }}\left(-{p}_{1}^{2},\tfrac{\pi }{2},{m}^{2}\right)\left[{a}_{2}\tfrac{\left[{p}_{1}^{4}+2{p}_{1}^{2}\left(1+{m}^{2}\right)+3{m}^{2}\right]}{2\left(1+{p}_{1}^{2}\right)\left({m}^{2}+{p}_{1}^{2}\right)}+{b}_{2}\right]+{c}_{2}{\rm{\Pi }}\left(-{p}_{2}^{2},\tfrac{\pi }{2},{m}^{2}\right)+{d}_{2}{\rm{\Pi }}\left(-{p}_{3}^{2},\tfrac{\pi }{2},{m}^{2}\right)\right\}},\end{eqnarray} \tag{ 5b }$$

$$\begin{eqnarray}&&\begin{array}{l}{\bar{\nu }}_{\theta }\left(e,{r}_{p},a\right)=\displaystyle \frac{2{\bar{\nu }}_{r}\left(e,\mu ,a\right){\mu }^{1/2}\sqrt{\left({x}^{2}+{a}^{2}+2{aEx}\right)}\cdot F\left(\tfrac{\pi }{2},{m}^{2}\right)}{\pi {\left[1-{\mu }^{2}{x}^{2}\left(3-{e}^{2}-2e\right)\right]}^{1/2}},\end{array}\end{eqnarray} \tag{ 5c }$$

where $x=\left({L}_{z}-{aE}\right)$, and {${p}_{1}^{2}$, ${p}_{2}^{2}$, ${p}_{3}^{2}$} are given by Equation 7(k) of Rana & Mangalam (2019a), while m² is given by Equation 13(c) of Rana & Mangalam (2019a), for the e0 orbits. See Appendix A for the derivation of Equation 5(c), which is a novel reduced form for ${\bar{\nu }}_{\theta }$.

Now, we use these frequency formulae, Equations 4(a)–(c), to deduce the suitable range of parameters {e, r_p, a, Q} for eQ and Equations 5(a)–(c) for e0 trajectories to find {e, r_p, a} to retrodict the observed range of QPOs in BHXRBs, which is provided in Table 3. In Figures 8–10, we have shown the variation of the quantities

$$\begin{eqnarray}&&\begin{array}{l}{\delta }_{\phi }\left(e,{r}_{p},a,Q\right)=\displaystyle \frac{\left[{\bar{\nu }}_{\phi }\left(e,{r}_{p},a,Q\right)-{\bar{\nu }}_{\phi }\left(e=0,{r}_{p},a,Q=0\right)\right]}{{\bar{\nu }}_{\phi }\left(e=0,{r}_{p},a,Q=0\right)},\end{array}\end{eqnarray} \tag{ 6a }$$

$$\begin{eqnarray}&&\begin{array}{l}{\delta }_{{pp}}\left(e,{r}_{p},a,Q\right)\\ \quad \ =\ \displaystyle \frac{\left[{\bar{\nu }}_{{pp}}\left(e,{r}_{p},a,Q\right)-{\bar{\nu }}_{{pp}}\left(e=0,{r}_{p},a,Q=0\right)\right]}{{\bar{\nu }}_{{pp}}\left(e=0,{r}_{p},a,Q=0\right)},\end{array}\end{eqnarray} \tag{ 6b }$$

$$\begin{eqnarray}&&\begin{array}{l}{\delta }_{{np}}\left(e,{r}_{p},a,Q\right)\\ \quad \ =\ \displaystyle \frac{\left[{\bar{\nu }}_{{np}}\left(e,{r}_{p},a,Q\right)-{\bar{\nu }}_{{np}}\left(e=0,{r}_{p},a,Q=0\right)\right]}{{\bar{\nu }}_{{np}}\left(e=0,{r}_{p},a,Q=0\right)},\end{array}\end{eqnarray} \tag{ 6c }$$

in the (r_p, a) plane for combinations of e = {0.25, 0.5} and Q = {0, 2, 4}. These quantities provide a fractional deviation between frequencies of general eccentric orbits and circular orbits for the same spin and periastron radius. For this comparison, we have calculated the frequency corresponding to a circular orbit at the same radius, r_p = r, for a fixed value of parameter a. Hence, the deviation, δ, between frequencies defined in this manner is dominated by the parameters e and Q. Also, these deviations are shown only in the region where ${\bar{\nu }}_{\phi }\left(e,{r}_{p},a,Q\right)$, ${\bar{\nu }}_{{pp}}\left(e,{r}_{p},a,Q\right)$, and ${\bar{\nu }}_{{np}}\left(e,{r}_{p},a,Q\right)$ are in the range of QPO frequencies allowed by the observations, as provided in Table 3. Hence, these plots together give us the information of deviation of frequencies from circularity, as the e and Q parameters are varied, along with information on the range of $\left(e,{r}_{p},a,Q\right)$ for general eccentric orbits allowed by the observed range of QPO frequency. The 3:2 and 5:3 ratios of the simultaneous HFQPOs, seen in a few BHXRBs, are also a remarkable phenomenon that we need to fathom; for example, 300 Hz and 450 Hz HFQPOs were seen in GROJ 1655-40 (Remillard et al. 1999b; Strohmayer 2001a) and 240 Hz and 160 Hz HFQPOs in H1743-322 (Homan et al. 2005; Remillard et al. 2006). Assuming the GRPM, this ratio is given by ${\nu }_{\phi }/{\nu }_{{pp}}={\bar{\nu }}_{\phi }/{\bar{\nu }}_{{pp}}$, which is a dimensionless quantity. The contours of this ratio are shown in Figure 11 for the six combinations from the set e = {0.2, 0.4}, Q = {0, 2, 4}. The blue contours in Figures 8–11 represent the ISSO radius, and black contours represent the MBSO radius as also indicated in Figure 5, whereas the magenta color contours represent the separatrix orbits, given by the equality in Equation (2), defining the innermost limit for r_p of an eccentric orbit with a given e.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

A summary of the results is given below:

• 1.  
  A novel and reduced form of ${\bar{\nu }}_{\theta }\left(e,{r}_{p},a\right)$ for e0 trajectories, given by Equation 5(c), is derived in Appendix A.
• 2.  
  Assuming the GRPM, (non)equatorial eccentric trajectories with small to moderate eccentricities, e ≲ 0.5, with Q ∼ 0–4 also generate the expected range of QPO frequencies, {${\bar{\nu }}_{\phi }$, ${\bar{\nu }}_{{pp}}$, ${\bar{\nu }}_{{np}}$}, in BHXRBs, as shown in Table 3. We have not taken very high values for the Q parameter, as the particle oscillation is expected to be close to the equatorial plane in typical BHXRB scenarios.
• 3.  
  The effective r_p ranges that produce the required QPO frequency ranges are Δr_p ∼ 2–15 for ${\bar{\nu }}_{\phi }$, Δr_p ∼ 2–10 for ${\bar{\nu }}_{{pp}}$, and Δr_p ∼ 4–20 for ${\bar{\nu }}_{{np}}$, where a varies from 0 to 1. While these Δr_p values are strongly dependent on e, they are only weakly dependent on the Q parameter. The frequency ${\bar{\nu }}_{{np}}$ (see Figure 10) increases with a, which implies that we expect to find high type C LFQPO values (nearly ${\bar{\nu }}_{{np}}\sim 0.001$) for the black holes with high spin.
• 4.  
  As e increases, the allowed region shifts close to the black hole. In other words, we expect (non)equatorial eccentric orbits close to the black hole to create the allowed frequency range, whereas circular orbits at comparatively larger radius cater to the same frequency range (see Figure 7). This is consistent with the finding that the GRPM favors slightly eccentric and strongly relativistic orbits. We also see that as e increases, the frequencies deviate and decrease from corresponding circular orbit frequencies; for example, ${\bar{\nu }}_{\phi }$ decreases by 30% for e = 0.25 to 60% for e = 0.5 (see Figure 8), ${\bar{\nu }}_{{pp}}$ decreases by 40% for e = 0.25 to 79% for e = 0.5 (see Figure 9), and ${\bar{\nu }}_{{np}}$ decreases by 40% for e = 0.25 to 80% for e = 0.5 (see Figure 10).
• 5.  
  The dependence of these frequencies on Q is very weak. Although the change is comparatively small, we see that these frequencies increase with Q. For example, the maximum increase in ${\bar{\nu }}_{\phi }$ is ∼3% (see Figure 8) and ∼10% for ${\bar{\nu }}_{{pp}}$ (see Figure 9), whereas it is ∼3% for ${\bar{\nu }}_{{np}}$ (see Figure 10) as Q changes from 0 to 4. Even for high Q values, say Q ∼ 10, the change in frequencies is of the same order.
• 6.  
  Expectedly, the associated frequencies increase as the r_p of a trajectory decreases for a given {e, a, Q}, where r_p of an eccentric trajectory is limited by the corresponding separatrix orbit, having the same {e, a, Q} values.
• 7.  
  As shown in Figure 11, the 3:2 or 5:3 ratios of HFQPOs originate in the region very close to the separatrix orbits, which is between MBSO and ISSO radii corresponding to typically Δr_p ∼ 2–6; this range is dependent on a since r_p decreases as a increases. The frequency ratio contours shift close to the black hole as e is increased, whereas these contours move toward large r_p as Q is increased. This indicates that nonequatorial orbits show a 3:2 or 5:3 ratio of HFQPO frequencies farther away from the black hole than the equatorial orbits, and eccentric orbits have such ratios comparatively closer to the black hole than the circular orbits. Therefore, eQ and 00 orbits close to the black hole can account for these ratios, as e and Q have canceling effects.

2.2. Spherical Orbits: Q0

Similar to the eQ trajectories, the spherical orbits (Q0) are also specific to the rotating black holes. They are the orbits with a constant radius, r_s, where the orbital plane precesses on a sphere about the spin axis of the black hole. Similar to the ISCO and MBCO radii for circular orbits, ISSO and MBSO radii exist for the spherical orbits that are functions of the a and Q parameters (Rana & Mangalam 2019a, 2019b). We explore the ranges of parameters, {r_s, a, Q}, for spherical orbits allowed by the observed frequency range of QPOs (see Table 3). The fundamental frequency formulae for the spherical orbits reduce to the form given by (see Appendix B for a derivation: Equations (B14), B15(b), and B16(c))

$$\begin{eqnarray}&&{\bar{\nu }}_{\phi }\left({r}_{s},a,Q\right)=\displaystyle \frac{\left\{\left[-\tfrac{\left(2{L}_{z}{r}_{s}-{L}_{z}{r}_{s}^{2}-2{r}_{s}{aE}\right)}{{\rm{\Delta }}}-{L}_{z}\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)+{L}_{z}\cdot {\rm{\Pi }}\left({z}_{-}^{2},\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}}{2\pi \left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}E\cdot K\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ 7a }$$

$$\begin{eqnarray}&&{\bar{\nu }}_{r}\left({r}_{s},a,Q\right)=\displaystyle \frac{\sqrt{{r}_{s}^{4}\left(1-{E}^{2}\right)+\left(3{{Qa}}^{2}-2{x}^{2}{r}_{s}-2{{Qr}}_{s}\right)}\cdot F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)}{2\pi {r}_{s}\left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}E\cdot K\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ 7b }$$

$$\begin{eqnarray}&&{\bar{\nu }}_{\theta }\left({r}_{s},a,Q\right)=\displaystyle \frac{a\sqrt{1-{E}^{2}}{z}_{+}}{4\left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}E\cdot K\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ 7c }$$

where ${\rm{\Delta }}={r}_{s}^{2}+{a}^{2}-2{r}_{s}$, and z_± are given by Equation 9(d) of Rana & Mangalam (2019a). In Figure 12, we show the contours of the quantities

$$\begin{eqnarray}&&{\delta }_{\phi }\left({r}_{s},a,Q\right)=\displaystyle \frac{\left[{\bar{\nu }}_{\phi }\left({r}_{s},a,Q\right)-{\bar{\nu }}_{\phi }\left({r}_{s},a,Q=0\right)\right]}{{\bar{\nu }}_{\phi }\left({r}_{s},a,Q=0\right)},\end{eqnarray} \tag{ 8a }$$

$$\begin{eqnarray}&&{\delta }_{{pp}}\left({r}_{s},a,Q\right)=\displaystyle \frac{\left[{\bar{\nu }}_{{pp}}\left({r}_{s},a,Q\right)-{\bar{\nu }}_{{pp}}\left({r}_{s},a,Q=0\right)\right]}{{\bar{\nu }}_{{pp}}\left({r}_{s},a,Q=0\right)},\end{eqnarray} \tag{ 8b }$$

$$\begin{eqnarray}&&{\delta }_{{np}}\left({r}_{s},a,Q\right)=\displaystyle \frac{\left[{\bar{\nu }}_{{np}}\left({r}_{s},a,Q\right)-{\bar{\nu }}_{{np}}\left({r}_{s},a,Q=0\right)\right]}{{\bar{\nu }}_{{np}}\left({r}_{s},a,Q=0\right)},\end{eqnarray} \tag{ 8c }$$

for QPOs in the (r_s, a) plane for spherical orbits with Q = {2, 4} assuming the GRPM, using Equations 7(a)–(c). The blue contours in Figures 12 and 13 represent the ISSO radii, and the black contours represent the MBSO radii. The results for spherical orbits are enumerated below:

• 1.  
  Novel and reduced forms for the equations of motion {$\phi \left({r}_{s},a,Q\right)$, $t\left({r}_{s},a,Q\right)$}, given by Equation (B12), and the fundamental frequencies {${\bar{\nu }}_{\phi }\left({r}_{s},a,Q\right)$, ${\bar{\nu }}_{r}\left({r}_{s},a,Q\right)$, ${\bar{\nu }}_{\theta }\left({r}_{s},a,Q\right)$}, given by Equation (7), for spherical trajectories are derived in Appendix B.
• 2.  
  Assuming the GRPM, we see that the spherical orbits with Q ∼ 0–4 are in the expected range of QPO frequencies for BHXRBs. The allowed range of r_s to source the QPOs is typically ∼3–18 for ${\bar{\nu }}_{\phi }$ (see Figures 12(a) and (b)), ∼3–12 for ${\bar{\nu }}_{{pp}}$ (see Figures 12(c) and (d)), and ∼3–20 for ${\bar{\nu }}_{{np}}$ (see Figures 12(e) and (f)), where a varies from 0 to 1.
• 3.  
  The frequencies change weakly with Q. The maximum changes in frequencies are ∼2%–3% for ${\bar{\nu }}_{\phi }$, ∼11%–23% for ${\bar{\nu }}_{{pp}}$, and ∼4%–8% for ${\bar{\nu }}_{{np}}$ as Q changes from 2 to 4 for the spherical orbits. The associated frequencies increase as r_s decreases for a given {a, Q}.
• 4.  
  We see from Figure 13 that the 3:2 or 5:3 ratio of HFQPOs, ${\bar{\nu }}_{\phi }/{\bar{\nu }}_{{pp}}$, for spherical orbits should emanate in the region r_s ∼ 3–7 for Q = 2 and r_s ∼ 3.5–7.5 for Q = 4. The ranges of r_s are also dependent on a, where r_s for a given ratio contour decreases as a increases.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Next, we take up a few cases of black hole systems that are known to have shown either two or three simultaneous QPOs in their PDS, and we extract the parameter values of the nonequatorial eccentric (eQ), equatorial eccentric (e0), and spherical orbits (Q0) corresponding to the observed QPO frequencies using our GRPM. The solution for a given GRPM class (eQ, Q0, e0) being attempted here is based on balancing the known (number of simultaneous frequencies, two or three) with the number of unknown parameters {e, r_p, a, Q} (see Table 4 illustrating this criterion). For the three frequency cases (M82 X-1 and GROJ 1655-40), we have to either supply a from available data or deduce this using a procedure that involves minimizing χ² in the unknown parameter volume. For the geometric study of orbits that is of importance here, we have taken the view that the best approximation to a is to be determined first, and then the solution vector {e, r_p, Q} (which is crucial for the orbital shape) for the peak probability is found. We have taken slightly different approaches for the two sources as exact solutions are found only in one of the two sources (M82 X-1), where we minimize χ² in the a dimension to isolate a. In the other case where no exact solution vector is found (GROJ 1655-40), and where it is computationally expensive to explore the full four-dimensional parameter volume of {e, r_p, Q, a} in a fine-grained manner, we have only done a primary coarse-grained search to find a sufficiently accurately and then proceeded to determine the unknown parameters {e, r_p, Q} by a fine-grid search. The two QPO frequency cases (XTEJ 1550-564, 4U 1630-47, and GRS 1915+105) are searched by direct fine-grid computations assuming a from available data (see Tables 4 and 5).

We describe our parameter search criteria below:

• 1.  
  For BHXRBs with three simultaneous QPOs, that is, M82 X-1 and GROJ 1655-40 (see Table 5), since a type C LFQPO is also present, which corresponds to the nodal oscillation frequency (ν_np), we search for all eQ, e0, and Q0 orbit solutions. We use Equations 4(a)–(c) and 5(a)–(c) to equate the QPO frequencies to {ν_ϕ, ν_pp, ν_np} and find the parameters {e, r_p, a} of eQ and e0 orbits for M82 X-1 and GROJ 1655-40. Next, we calculate the most probable spin of the black hole to estimate {e, r_p, Q} of the orbit. Similarly, we study the Q0 orbits as solutions to the QPOs using Equations 7(a)–(c) and find the parameters {r_s, a, Q} for these BHXRBs. Hence, the parameters searched for these cases are

  $$\begin{eqnarray}3\mathrm{QPOs}=\left\{\begin{array}{ll}{eQ}\ \mathrm{and}\ e0, & \{{ \mathcal M }=\mathrm{fixed}\ \mathrm{from}\ \mathrm{observations},e,{r}_{p},a,Q=\{0,1,2,3,4\}\},\\ Q0, & \{{ \mathcal M }=\mathrm{fixed}\ \mathrm{from}\ \mathrm{observations},e=0,{r}_{s},a,Q\}.\end{array}\right.\end{eqnarray} \tag{ 9a }$$
• 2.  
  For BHXRBs with two simultaneous QPOs, that is, XTEJ 1550-564, 4U 1630-47, and GRS 1915+105 (see Table 5), we expect that the solutions are likely to be equatorial as the LFQPO, or ν_np oscillation, is absent (this is consistent with no large-amplitude nodal oscillations and strictly equatorial orbits). Hence, we search for e0 solutions using Equations 5(a)–(b) for {ν_ϕ, ν_pp} to find {e, r_p} of the orbit. However, we also check for the Q0 orbital solution in these systems and estimate the parameters {r_s, Q} using {ν_ϕ, ν_pp}. Hence, the parameters searched for in these cases are

  $$\begin{eqnarray}2\mathrm{QPOs}=\left\{\begin{array}{ll}\ e0, & \{\left({ \mathcal M },a\right)=\mathrm{fixed}\ \mathrm{from}\ \mathrm{observations},e,{r}_{p},Q=0\},\\ Q0, & \{\left({ \mathcal M },a\right)=\mathrm{fixed}\ \mathrm{from}\ \mathrm{observations},e=0,{r}_{s},Q\}.\end{array}\right.\end{eqnarray} \tag{ 9b }$$

Table 5.  Summary of Existing BHXRBs That Exhibit Either Three or Two Simultaneous QPOs

S. No.	BHXRB	ν1 (Hz)	ν2 (Hz)	ν3 (Hz)	${\boldsymbol{ \mathcal M }}$	a	Model
							Classes
1.	M82 X-1	5.07 ± 0.06 (a)	3.32 ± 0.06 (a)	$\left(204.8\pm 6.3\times {10}^{-3}\right)$ (b)	428 ± 105 (a)	⋯	eQ, e0, Q0
2.	GROJ 1655-40	441 ± 2 (c)	298 ± 4 (c)	17.3 ± 0.1 (c)	5.4 ± 0.3 (d)	⋯	eQ, e0, Q0
3.	XTE J1550-564	268 ± 3 (e)	188 ± 3 (e)	⋯	9.1 ± 0.61 (f)	${0.34}_{-0.45}^{+0.37}$ (g)	e0, Q0
4.	4U 1630-47	179.3 ± 5.7 (h)	38.06 ± 7.3 (h)	⋯	10 ± 0.1 (i)	${0.985}_{-0.014}^{+0.005}$ (j)	e0, Q0
5.	GRS 1915+105	69.2 ± 0.15 (k)	41.5 ± 0.4 (k)	⋯	10.1 ± 0.6 (l)	0.98 ± 0.01 (m)	e0, Q0

Note. The first two rows represent the cases having twin HFQPOs with simultaneous type C QPOs. The remaining rows show the cases of BHXRBs having only twin HFQPOs. The columns show the source name, QPO frequencies, and previously measured mass through optical, infrared, or X-ray observations, previously known spin of the black hole measured by fit to the Fe Kα line or to the continuum spectrum (for 1 and 2 we calculate the parameter a from our method), and the class of GRPM applied to estimate the parameters. The references are indicated by lowercase letters (a–m).

References. (a) Pasham et al. (2014), (b) Pasham & Strohmayer (2013a), (c) Motta et al. (2014a), (d) Beer & Podsiadlowski (2002), (e) Miller et al. (2001), (f) Orosz et al. (2011), (g) Miller & Miller (2015), (h) Klein-Wolt et al. (2004), (i) Seifina et al. (2014), (j) King et al. (2014), (k) Strohmayer (2001b), (l) Steeghs et al. (2013), (m) Miller et al. (2013).

Download table as:  ASCIITypeset image

We have summarized the history of black hole systems considered here in Appendix D. In Section 3.1, we summarize the observations related to QPO detection, mass, and spin estimation and the parameters we estimated for each source. In Section 3.2, we explain the method used to estimate the parameters of these orbits and corresponding errors and then present the results for the (non)equatorial eccentric orbits in Section 3.2.1 and spherical orbits in Section 3.2.2.

3.1. Source Selection

3.2. Method Used and Results

We apply the GRPM to associate the fundamental frequencies of eQ, e0, and Q0 orbits with QPOs. In Appendix E, we describe a generic procedure that we have used to estimate errors in the orbital parameters. A flowchart of this method is provided in Figure 14. Next, we summarize the results corresponding to the eQ and e0 models in Section 3.2.1 and the Q0 model in Section 3.2.2.

Zoom In Zoom Out Reset image size

3.2.1. Nonequatorial and Equatorial Eccentric Orbits (eQ and e0)

We have taken the cases of five BHXRBs, known to have either three or two simultaneous detections of QPOs in their observations, to study the eccentric and nonequatorial trajectories as solutions to the QPOs assuming the GRPM. Here we summarize the results for the cases of three and two simultaneous QPOs separately, as discussed below:

• 1.  
  Three simultaneous QPOs: In our analysis, varying the dimensionless parameter $Q=\{0,1,2,3,4\}\propto \left({L}^{2}-{L}_{z}^{2}\right)$ gives us various trajectory solutions having different {e, r_p, a, Q} combinations. We first find the exact solutions for the parameters {e, r_p, a}, given in Table 6, by equating the centroid frequencies of three simultaneous QPOs (see Table 5) to {ν_ϕ, ν_pp, ν_np} for each value of Q = {0, 1, 2, 3, 4} using our analytic formulae (Equations 4(a)–(c)). We estimate errors for the parameters {e, r_p, a} using the method discussed in Appendix E (see Figure 14) for each value of Q. The results of fits to the integrated profiles {${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, ${{ \mathcal P }}_{1}\left(a\right)$} are summarized in Table 6. Since the spin of the black hole is not expected to vary, we estimate the most probable spin value for these black holes and then estimate the orbital parameters {e, r_p, Q} and their corresponding errors again using the same method discussed in Appendix E (see Figure 14). The results for each case are as follows:

  • (a)  
    M82 X-1: In this case, we find that the (non)equatorial trajectories with small to moderate eccentricities e ∼ 0.18–0.28 with r_p = 4.6–5.07 and a = 0.28–0.31 (see Table 6) are possible exact solutions for the observed QPO frequencies in M82 X-1, for Q between 0 and 4. Starting with these exact solutions, the most probable value of the spin is found first. In Figure 15(a), we show the spin variation in the parameter solutions for QPOs as a function of Q. Next, to estimate the most probable value of the spin, we minimize the function

    $$\begin{eqnarray}&&{\chi }_{a}^{2}=\displaystyle \sum _{i}\displaystyle \frac{{\left(a-{a}_{i}\right)}^{2}}{{\sigma }_{{a}_{i}}^{2}},\end{eqnarray} \tag{ 10a }$$

    which gives

    $$\begin{eqnarray}&&a=\displaystyle \frac{{\displaystyle \sum }_{i}\left({a}_{i}/{\sigma }_{{a}_{i}}^{2}\right)}{{\displaystyle \sum }_{i}\left(1/{\sigma }_{{a}_{i}}^{2}\right)},\end{eqnarray} \tag{ 10b }$$

    where i = 1–6 corresponds to six probable solutions for a, and the σ_i values are the corresponding 1σ errors, where five of these are given in Table 6 and the remaining one corresponds to the spherical orbit solution found for M82X-1 given in Table 9. By including these six solutions, we have spanned the complete (e, Q) parameter space, which is bounded by e0 and Q0 solutions. This gives us the most probable spin value of a = 0.2994. Hence, we fix the spin of the black hole to this most probable estimate and then calculate the remaining parameters {e, r_p, Q} and corresponding errors using the method given in Appendix E and Figure 14. We find the exact solution for QPOs at {e = 0.2302, r_p = 4.834, Q = 2.362} calculated by equating centroid QPO frequencies to {ν_ϕ, ν_pp, ν_np} while fixing a = 0.2994. The probability density distribution profiles {${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, ${{ \mathcal P }}_{1}\left(Q\right)$}, along with their model fit, and the probability contours in the parameter plane {e, Q}, {r_p, e}, and {Q, r_p} are shown in Figure 16. The results of the model fit to the integrated profiles are summarized in Table 7. The corresponding errors are quoted with respect to the exact solution of the parameters, which slightly differ from the peak of the integrated profiles {${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, ${{ \mathcal P }}_{1}\left(Q\right)$}, as expected (see Figure 16).
  • (b)  
    GROJ 1655-40: For this case, we did not find the exact solution for the parameters {e, r_p, a} when the centroid frequencies of QPOs, Table 5, are equated to {ν_ϕ, ν_pp, ν_np}. However, we generate the probability density profiles ${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, and ${{ \mathcal P }}_{1}\left(a\right)$ for each value of Q between 0 and 4. The results of fits for these profiles are summarized in Table 6. We found that the probability density peaks near very small eccentricities e ∼ 0.05–0.07 for various values of Q, whereas r_p ranges between 5.24 and 5.43 and a ranges between 0.282 and 0.291; see Table 6. The change in the value of the spin of the black hole as a function of Q is shown in Figure 15(b) for GROJ 1655-40. Next, we find the most probable value of the spin for this BHXRB. Since we did not find any exact solution for the parameters by equating centroid frequencies of QPOs to the frequency formulae, we calculated the χ² function given by

    $$\begin{eqnarray}\begin{array}{rcl}{\chi }^{2} & = & \displaystyle \frac{{\left({\nu }_{\phi }-{\nu }_{10}\right)}^{2}}{{\sigma }_{1}^{2}}+\displaystyle \frac{{\left({\nu }_{{pp}}-{\nu }_{20}\right)}^{2}}{{\sigma }_{2}^{2}}+\displaystyle \frac{{\left({\nu }_{{np}}-{\nu }_{30}\right)}^{2}}{{\sigma }_{3}^{2}}\end{array}\end{eqnarray} \tag{ 11 }$$

    in the four-dimensional parameter space {e, r_p, a, Q} using Equations 4(a)–(c) for {ν_ϕ, ν_pp, ν_np}, and we numerically found the minimum χ² = 2.814 for the parameter combination {e = 0.021, r_p = 5.51, a = 0.283, Q = 0}. This is a primary coarse-grained search to find a viable solution of a. Next, we assume a = 0.283 corresponding to the minimum χ² to calculate the final solution for the parameters {e, r_p, Q}, which are the key parameters for the geometric study, using the more accurate fine-grid method described in Appendix E and Figure 14. We find that the probability density peaks near {e = 0.071, r_p = 5.25, Q = 0}. The results of fitting to the {${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, ${{ \mathcal P }}_{1}\left(Q\right)$} profiles are summarized in Table 7, whereas these profiles with their model fit and the probability contours in the parameter plane {e, Q}, {r_p, e}, and {Q, r_p} are shown in Figure 17.

  Hence, we conclude for both M82 X-1 and GROJ 1655-40 that (non)equatorial trajectories (both eQ and e0) with small or moderate eccentricities in the region very close to the black hole are the solutions for the observed QPOs assuming our GRPM. A self-emitting blob of matter close to a Kerr black hole can have enough energy and angular momentum to attain an eccentric and nonequatorial trajectory. These results are also consistent with the conclusions made in Section 2.1 that the trajectories having small to moderate eccentricities with Q = 0–4 are also possible solutions for the observed range of QPO frequencies in BHXRBs.The errors in QPO frequencies cause a distribution in the solution space {e, r_p, Q} as solutions using our GRPM, as shown in Figures 16 and 17. We take various combinations of these parameters within the range of 1σ errors, as summarized in Table 7, as any such parameter combination is a probable solution for the frequencies within the width of QPOs observed in the power spectrum. In Figure 18, we have plotted together the trajectories for these parameter combinations for both BHXRBs M82 X-1 and GROJ 1655-40. Each trajectory has different parameter values {e, r_p, Q} and is indicated by a different color, where we fixed the spin of the black hole to a = 0.2994 for M82 X-1 and a = 0.283 for GROJ 1655-40. Hence, these trajectories, having fundamental frequencies very close to each other and within the width of the QPO, together simulate the strong rms of the observed QPOs. The trajectories together span a torus in the region 4.7–9.08 for M82 X-1 and 5.11–6.67 for GROJ 1655-40, which should be the emission region for QPOs, where we expect precession frequencies of both the eQ and e0 trajectories. The ISCO radius is ∼5 for both cases of BHXRB. We suggest that the simultaneous HFQPO and LFQPO emission should be from a region that is close to the inner edge of the accretion disk (r_in), where both eQ and e0 trajectories span a torus; the disk edge could be a source of blobs that are generating QPOs, as we will argue in Section 5. In contrast, a rigid body precession model is invoked by some authors (Ingram et al. 2009; Ingram & Done 2011, 2012), where Lense–Thirring precession of a rigid torus is suggested as the origin of the type C QPOs. Here, instead of the rigid precession of a solid torus, we propose that a collective precession of various trajectories, spanning a torus region, explains the origin of HFQPOs and LFQPOs simultaneously. We argue that HFQPOs originate when r_in comes in very close to the black hole at some point during the outburst (the soft state). In the hard state, r_in is farther out, and the type C QPO is more frequent and is more prone to vertical oscillations (ν_np). This scenario explains the increase in the frequency of type C QPOs with a decrease in r_in, while the spectrum transits from the hard to soft state.
• 2.  
  Two simultaneous QPOs: We have considered only equatorial eccentric trajectories, Q = 0, for these BHXRBs, as we can estimate only two parameters of the orbit corresponding to two simultaneous QPOs. First, we find the exact solutions for the parameters {e, r_p}, summarized in Table 8, by equating the centroid frequencies of two simultaneous QPOs (see Table 5) to {ν_ϕ, ν_pp} using our analytic formulae for Q = 0, Equations 5(a) and (b). Then we calculate the errors in the parameters {e, r_p} using the method discussed in Appendix E (see Figure 14). The results are summarized in Table 8. These results are described below:

  • (a)  
    XTEJ 1550-564: We find an equatorial trajectory with eccentricity e = 0.262 with r_p = 4.365 (see Table 8) as a solution for the observed QPO frequencies in XTEJ 1550-564. The calculated probability density profiles in e and r_p dimensions, ${{ \mathcal P }}_{1}\left(e\right)$ and ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, were found to be skew symmetric and were fit by an interpolating function. The corresponding errors were obtained by taking the integrated probability of 68.2% about the peak value of the probability density distributions. The quoted errors are calculated with respect to the exact solution of the parameters, which slightly deviates from the peak of the integrated profiles {${{ \mathcal P }}_{1}\left(e\right)$ and ${{ \mathcal P }}_{1}\left({r}_{p}\right)$}; see Figure 19 and Table 8. These profiles, corresponding model fit, and the probability contours in the (e, r_p) plane are shown in Figure 19.
  • (b)  
    4U 1630-47: We found an exact solution at {e = 0.734, r_p = 2.249} (see Table 8) by equating {ν_ϕ, ν_np} instead of {ν_ϕ, ν_pp} to the centroid QPO frequencies. This might be because the QPO with a lower frequency of ∼38 Hz (see Table 5) is too small to be an HFQPO. The calculated probability density profiles in the e and r_p dimensions, the corresponding model fit, and the probability contours in the (e, r_p) plane are shown in Figure 20. In this case, too, we see that ${{ \mathcal P }}_{1}\left(e\right)$ and ${{ \mathcal P }}_{1}\left({r}_{p}\right)$ profiles are skew, such that the integrated probability is 68.2% about the peak value of the probability density distributions, and the errors are quoted with respect to the exact solution of the parameters, which slightly deviates from the peak of the integrated profiles ${{ \mathcal P }}_{1}\left(e\right)$ and ${{ \mathcal P }}_{1}\left({r}_{p}\right)$ (see Figure 20 and Table 8). We see that a highly eccentric orbit is found as the most probable solution.
  • (c)  
    GRS 1915+105: We found an exact solution at {e = 0.918, r_p = 1.744}; see Table 8. We find a highly eccentric equatorial trajectory as the most probable solution that can give the observed QPO frequencies in GRS 1915+105. This result is similar to the case of 4U 1630-47, which leads us to observe that a black hole with a high spin value prefers a highly eccentric orbit solution to simultaneous QPOs. The calculated probability density profiles ${{ \mathcal P }}_{1}\left(e\right)$ and ${{ \mathcal P }}_{1}\left({r}_{p}\right)$ are well fit by the Gaussian. The corresponding model fit and the probability contours in the (e, r_p) plane are shown in Figure 21.

  Hence, we conclude that for XTEJ 1550-564, 4U 1630-47, and GRS 1915+105, the e0 model in the region r_p = 1.74–4.36 is the probable cause of the observed QPOs in the power spectrum. We found high eccentricity values for the orbits as solutions for QPOs in the cases of BHXRB 4U 1630-47 and GRS 1915+105, and this indicates that black holes with very high spin values prefer highly eccentric orbits in the QPO solutions.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 6.  Summary of Results Corresponding to (Non)equatorial Eccentric Solutions (eQ and e0) for BHXRBs M82 X-1 and GROJ 1655-40

BHXRB	Q	e Range	Resolution	Exact	Model Fit	rp Range	Resolution	Exact	Model Fit	a Range	Resolution	Exact	Model Fit
			Δe	Solution	to ${{ \mathcal P }}_{1}\left(e\right)$		Δrp	Solution	to ${{ \mathcal P }}_{1}\left({r}_{p}\right)$		Δa	Solution	to ${{ \mathcal P }}_{1}\left(a\right)$
				e0				rp0				a0
M82 X-1	0	0.23–0.32	0.001	0.277	${0.277}_{-0.045}^{+0.066}$	4.4–4.85	0.005	4.616	${4.616}_{-0.126}^{+0.069}$	0.26–0.32	0.001	0.290	0.290 ± 0.009
	1	0.21–0.31	0.001	0.259	${0.259}_{-0.045}^{+0.072}$	4.3–5	0.01	4.698	${4.698}_{-0.154}^{+0.068}$	0.265–0.315	0.001	0.294	0.294 ± 0.009
	2	0.19–0.29	0.001	0.239	${0.239}_{-0.046}^{+0.080}$	4.45–5.1	0.01	4.795	${4.795}_{-0.166}^{+0.066}$	0.27–0.32	0.001	0.298	0.298 ± 0.009
	3	0.16–0.26	0.001	0.214	${0.214}_{-0.045}^{+0.090}$	4.55–5.25	0.01	4.913	${4.913}_{-0.163}^{+0.081}$	0.28–0.33	0.001	0.302	0.302 ± 0.009
	4	0.12–0.24	0.001	0.187	${0.187}_{-0.047}^{+0.113}$	4.65–5.35	0.01	5.067	${5.067}_{-0.221}^{+0.076}$	0.285–0.335	0.001	0.308	0.308 ± 0.009
GROJ 1655-40	0	0–0.22	0.002	⋯	${0.07}_{-0.038}^{+0.042}$	4.6–5.7	0.01	⋯	${5.24}_{-0.186}^{+0.191}$	0.265–0.3	0.001	⋯	0.282 ± 0.003
	1	0–0.22	0.002	⋯	${0.062}_{-0.034}^{+0.040}$	4.6–5.8	0.015	⋯	${5.305}_{-0.185}^{+0.170}$	0.24–0.32	0.002	⋯	0.284 ± 0.003
	2	0–0.2	0.002	⋯	${0.056}_{-0.031}^{+0.038}$	4.7–5.85	0.015	⋯	${5.345}_{-0.169}^{+0.167}$	0.27–0.31	0.001	⋯	0.286 ± 0.003
	3	0–0.2	0.002	⋯	${0.052}_{-0.029}^{+0.036}$	4.75–5.9	0.015	⋯	${5.395}_{-0.170}^{+0.151}$	0.275–0.32	0.001	⋯	0.288 ± 0.003
	4	0–0.2	0.002	⋯	${0.05}_{-0.028}^{+0.034}$	4.8–5.95	0.015	⋯	${5.43}_{-0.162}^{+0.147}$	0.278–0.31	0.001	⋯	0.291 ± 0.003

Note. The columns describe the range of parameter volume considered for {e, r_p, a} with a chosen resolution to calculate the normalized probability density at each point inside the parameter volume using Equation E31(b), the exact solutions for {e, r_p, a} calculated using Equations 4(a)–(c), and the results of the model fit to ${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, and ${{ \mathcal P }}_{1}\left(a\right)$ for each value of Q between 0 and 4.

Download table as:  ASCIITypeset image

Table 7.  Summary of Results for {e, r_p, Q} Parameter Solutions and Corresponding Errors for QPOs in BHXRBs M82 X-1 and GROJ 1655-40

BHXRB	e Range	Resolution	Exact	Model Fit	rp Range	Resolution	Exact	Model Fit	Q Range	Resolution	Exact	Model Fit
		Δe	Solution	to ${{ \mathcal P }}_{1}\left(e\right)$		Δrp	Solution	to ${{ \mathcal P }}_{1}\left({r}_{p}\right)$		ΔQ	Solution	to ${{ \mathcal P }}_{1}\left(Q\right)$
			e0				rp0				Q0
M82 X-1	0.1–0.35	0.002	0.230	${0.230}_{-0.049}^{+0.067}$	4.2–5.4	0.02	4.834	${4.834}_{-0.268}^{+0.181}$	0–5	0.1	2.362	${2.362}_{-1.439}^{+1.519}$
GROJ 1655-40	0–0.18	0.001	⋯	${0.071}_{-0.035}^{+0.031}$	4.9–5.75	0.0125	⋯	${5.25}_{-0.142}^{+0.171}$	0–3	0.1	⋯	0+0.623

Note. The columns describe the range of parameter volume taken for {e, r_p, Q}, and the chosen resolution to calculate the normalized probability density at each point inside the parameter volume, the exact solutions, and the results of the model fit to the integrated profiles. The spin of the black hole is fixed to the most probable estimates, which are a = 0.2994 for M82 X-1 and a = 0.283 for GROJ 1655-40.

Download table as:  ASCIITypeset image

Table 8.  Summary of Results Corresponding to the Equatorial Eccentric Orbit Solutions for BHXRBs XTEJ 1550-564, 4U 1630-47, and GRS 1915+105

BHXRB	e Range	Resolution	Exact	Model Fit	rp Range	Resolution	Exact	Model Fit
		Δe	Solution	to ${{ \mathcal P }}_{1}\left(e\right)$		Δrp	Solution	to ${{ \mathcal P }}_{1}\left({r}_{p}\right)$
			e0				rp0
XTEJ 1550-564	0.01–0.7	0.0005	0.262	${0.262}_{-0.062}^{+0.090}$	2.5–6.5	0.005	4.365	${4.365}_{-0.279}^{+0.169}$
4U 1630-47	0.4–0.99	0.0005	0.734	${0.734}_{-0.048}^{+0.066}$	1–4.5	0.005	2.249	${2.249}_{-0.353}^{+0.249}$
GRS 1915+105	0.6–0.999	0.0005	0.918	0.918 ± 0.002	0.3–4	0.005	1.744	${1.744}_{-0.011}^{+0.025}$

Note. The columns describe the parameter range considered for {e, r_p}, its resolution, the exact solutions for {e, r_p} calculated using {ν_ϕ, ν_pp} for XTEJ 1550-564 and GRS 1915+105, and using {ν_ϕ, ν_np} for 4U 1630-47 using Equations 5(a), (b), and (c), and results of the model fit to ${{ \mathcal P }}_{1}\left(e\right)$ and ${{ \mathcal P }}_{1}\left({r}_{p}\right)$.

Download table as:  ASCIITypeset image

We show all of the eccentric trajectory solutions together for both Q = 0 and $Q\ne 0$ in Figure 22 in the (r_p, a) plane along with the radii ISCO (ISSO), MBCO (MBSO), light radius, and the horizon. We see that the calculated eccentric orbit solutions are found in region 1 of the (r_p, a) plane (as defined in Figure 5) and near ISCO for Q = 0 in the cases of BHXRB 4U 1630-47, GROJ 1655-40, and GRS 1915+105. The trajectory solutions are found in region 2 near ISCO for XTEJ 1550-564 (Q = 0) and near ISSO for M82 X-1 (Q = 2.362; as defined in Figure 5). These results are also consistent with the results discussed in Section 2.1, except that very high e values are found for trajectories in BHXRB 4U 1630-47 and GRS 1915+105. Hence, we conclude that the eccentric trajectory solutions with Q = 0 and $Q\ne 0$ for the observed QPOs in BHXRBs are found either in region 1 or region 2 of the (r_p, a) plane but close to the ISCO (ISSO) curve; we call this radius R₀. As all these orbit solutions are distributed near R₀, it is expected that this radius is very close to the inner edge radius, r_in, of the circular accretion disk, which could also be a source of blobs that are generating these QPOs. The torus region, shown in Figure 18, spans part of regions 1 and 2 near R₀, which can be represented as $\left({{R}_{0}}_{-{{\rm{\Delta }}}_{2}}^{+{{\rm{\Delta }}}_{1}}\right)$, where Δ_i represents a small deviation from R₀ (which need not be the center point of the torus in this scenario). This means that the orbits near R₀ are induced by the instabilities in the inner flow to be (non)equatorial and eccentric.

Zoom In Zoom Out Reset image size

3.2.2. Spherical Orbits

Here we summarize the results of associating the spherical orbits around a Kerr black hole with QPOs in BHXRBs. We limited this study to the cases of BHXRBs M82 X-1 and XTEJ 1550-564, as we found the exact solutions for the parameters {r_s, a, Q} or {r_s, a} for only these two BHXRBs when we solved {ν_ϕ = ν₁₀, ν_pp = ν₂₀, ν_np = ν₃₀} for M82 X-1 and {ν_ϕ = ν₁₀, ν_pp = ν₂₀} for XTEJ 1550-564 using Equations 7(a)–(c). We calculated errors for the parameters using the method discussed in Appendix E (also see Figure 14); these results are summarized in Table 9 and are presented below:

• 1.  
  M82 X-1: We found the exact solution for a spherical orbit at {r_s = 6.044, a = 0.321, Q = 6.113} for M82 X-1. The spherical trajectory with these parameter values is shown in Figure 23(a). The calculated probability density profiles and the model fit are shown in Figure 24. The ${{ \mathcal P }}_{1}\left({r}_{s}\right)$ and ${{ \mathcal P }}_{1}\left(Q\right)$ profiles were found to be skew symmetric, and the integrated probability is 68.2% about the peak of the probability density distribution between the error bars, while ${{ \mathcal P }}_{1}\left(a\right)$ is well fit by a Gaussian. We see that the spin of the black hole is also found very close to the spin solutions estimated in Section 3.2.1. We conclude that along with the eQ trajectories having moderate eccentricities, as discussed in Section 3.2.1, a spherical trajectory (Q0) at r_s = 6.044 with Q = 6.113 is also a viable solution that can produce the observed QPO frequencies in M82 X-1. The corresponding spin estimate a = 0.321 ± 0.0132 was utilized in Section 3.2.1 using Equation 10(b) to calculate the most probable value of the spin for M82 X-1.
• 2.  
  XTEJ 1550-564: A spherical trajectory solution was found at r_s = 5.538 and Q = 2.697 for BHXRB XTEJ 1550-564 that is shown in Figure 23(b), and the calculated probability density profiles, the Gaussian model fit, and the probability contours in the {r_s, Q} plane are shown in Figure 25. So, along with an e0 trajectory, as discussed in Section 3.2.1, a Q0 orbit is also a viable candidate for the observed QPOs in the temporal power spectrum of XTEJ 1550-564.

We found that the spherical trajectories are also possible solutions for QPOs in BHXRBs M82 X-1 (a = 0.321, Q = 6.113, r_s = 6.044, r_I = 5.258) and XTEJ 1550-564 (a = 0.34, Q = 2.697, r_s = 5.538, r_I = 4.988). This indicates that the spherical trajectory solutions are in region 1 of the (r, a) plane, as defined in Figure 5, for both BHXRBs, and they are very close to the ISSO radius, r_I. These results are also consistent with the results discussed in Section 2.2, where the QPO-generating region is close to the ISSO curve in the (r, a) plane. For the case of M82 X-1, the spherical trajectory solution has a different value of spin compared to the ones estimated in Section 3.2.1, but it is very close to the other estimates given in Table 6. This value of spin, together with other results in Table 6, is used to estimate the most probable value of spin of the black hole for M82 X-1, which is a = 0.2994. We also see that a low eccentric trajectory prefers a high Q value and vice versa, as seen from the results shown in Table 6. As the Q value of the orbit is increased, the eccentricity of the trajectory solution decreases for both BHXRBs M82 X-1 and GROJ 1655-40. This trend is also followed here: for the spherical orbit (e = 0), Q ∼ 6 is found as a solution for M82 X-1 and Q ∼ 2.7 for XTEJ 1550-564, whereas a moderately eccentric trajectory solution was found with Q = 0 for XTEJ 1550-564; see Table 8.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 9.  Summary of Results Corresponding to the Spherical Orbit Solutions for BHXRBs M82 X-1 and XTEJ 1550-564

BHXRB	rs Range	Resolution	Exact	Model Fit	Q Range	Resolution	Exact	Model Fit	a Range	Resolution	Exact	Model Fit
		Δrs	Solution	to ${{ \mathcal P }}_{1}\left({r}_{s}\right)$		ΔQ	Solution	to ${{ \mathcal P }}_{1}\left(Q\right)$		Δa	Solution	to ${{ \mathcal P }}_{1}\left(a\right)$
			rs0				Q0				a0
M82 X-1	5.75–6.35	0.005	6.044	${6.044}_{-0.072}^{+0.071}$	2–10	0.03	6.113	${6.113}_{-1.645}^{+2.124}$	0.29–0.36	0.001	0.321	0.321 ± 0.013
XTEJ 1550-564	3–8	0.005	5.538	5.538 ± 0.054	0.01–5	0.01	2.697	${2.697}_{-1.627}^{+1.738}$	⋯	⋯	⋯	⋯

Note. The columns describe the range of parameter volume considered for {r_s, a, Q} and its resolution to calculate the normalized probability density using Equation E31(b), the exact solutions for {r_s, a, Q} calculated using Equations 7(a)–(c), the value of parameters corresponding to the peak of the integrated profiles in {r_s, a, Q}, and results of the model fit to ${{ \mathcal P }}_{1}\left({r}_{s}\right)$, ${{ \mathcal P }}_{1}\left(Q\right)$, and ${{ \mathcal P }}_{1}\left(a\right)$.

Download table as:  ASCIITypeset image

We conclude that various kinds of Kerr orbits, for example, spherical {e = 0, $Q\ne 0$}, equatorial eccentric {$e\ne 0$, Q = 0}, and nonequatorial eccentric {$e\ne 0$, $Q\ne 0$}, are also viable solutions for QPOs in BHXRBs. Hence, such trajectories with similar fundamental frequencies can together give a strong QPO signal in the temporal power spectrum.

A tight correlation between the frequencies of two components in the PDS of various sources, including black hole and neutron star X-ray binaries, was discovered (Psaltis et al. 1999). Such a correlation among various variability components of the PDS in both types of sources suggests a common and important emission mechanism for these signals. This correlation is either between two QPOs, an LFQPO, and either of the two HFQPOs, or it is between an LFQPO and high-frequency broadband noise components. We adopt the definition of Belloni et al. (2002) for these variability components: L_LF for LFQPO, and L_l and L_u for lower and upper HFQPOs or broad noise components. A systematic study of 571 RXTE observations was carried out for BHXRB GROJ 1655-40 between 1996 March and 2005 October (Motta et al. 2014a), and they also found such correlation between the type C QPOs and high-frequency QPOs and broadband components (either L_l or L_u; see Tables 1 and 2 and Figure 5 of Motta et al. 2014a). In this study, they calculated mass, spin of the black hole, and the radius at which QPOs originated {${ \mathcal M }=5.31$, a = 0.29, r = 5.68} (Motta et al. 2014a) using {L_u = ν_ϕ, L_l = ν_pp, L_LF = ν_np}, assuming that circular equatorial orbits are the origin of three simultaneous QPOs in the RPM (00 model as defined in Figure 1). Using the estimated values of ${ \mathcal M }$ and a, they fit the PBK correlation of variability components in GROJ 1655-40 by varying r.

Here we apply the e0 model solution calculated in Section 3.2.1 assuming {L_u = ν_ϕ, L_l = ν_pp, L_LF = ν_np}, using the observation ID having three simultaneous QPOs detected in GROJ 1655-40 (shown in Table 5), to fit the PBK correlation. We fix the mass of the black hole to ${ \mathcal M }=5.4$ (Beer & Podsiadlowski 2002) and the spin of the black hole to the most probable value, a = 0.283, estimated by minimizing the χ² function, given by Equation (11). We fix e and Q to the values estimated by the fine-grid method {e = 0.071, Q = 0} and vary r_p to calculate the frequencies. In Figure 26, we show the correlations of the frequencies corresponding to the parameters {e = 0.071, a = 0.283, Q = 0}, which are in good agreement with the PBK correlation. In Figure 27, these frequencies are shown as functions of r_p. We see that the data points for L_u components fit very well (see Figure 26(a)), whereas the L_l components show a good fit in the high-frequency range (see Figures 26(b), (c)). The L_LF components also show good agreement with the eccentric orbit solution (see Figure 26(d)).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Thirty-four L_l and L_LF components were detected simultaneously in the same observation ID (see Table 1 of Motta et al. 2014a). To calculate r_p, we first solve for L_LF = ν_np for the solution vector {e = 0.071, Q = 0, a = 0.283, ${ \mathcal M }=5.4$}; this locates the r_p, where oscillations are present, to a good approximation. Using these r_p values, we simultaneously solve {ν_pp = L_l, ν_np = L_LF} using the centroid frequencies of these components and estimate the exact solutions for parameters {e, Q} with {a = 0.283, ${ \mathcal M }=5.4$}. In 10 out of 34 cases, we found low-eccentricity eQ solutions for these PDS components, where the calculated parameters are shown in Table 10. We find orbits with high Q values at large r_p (this is expected as $Q\propto {L}^{2}-{L}_{z}^{2}$) as solutions for these PDS components. This exercise confirms the existence of eQ in addition to e0 solutions for QPOs.

In this section, we discuss our torus picture of eccentric trajectories, and we examine the model of fluid flow in the general-relativistic thin disk around a Kerr black hole (Penna et al. 2012; Mohan & Mangalam 2014) with the aim of finding a source of the e0, eQ, and Q0 trajectories. In this model, the region around the rotating black hole was divided into various regimes: (1) the plunge region between the ISCO radius and black hole horizon dominated by gas pressure and electron scattering based opacity, (2) the edge region at and very near the ISCO radius dominated by gas pressure and electron scattering based opacity, (3) the inner region outside the edge region with small radii comparable to ISCO dominated by radiation pressure and electron scattering based opacity, (4) the middle region outside the inner region where gas pressure again dominates over the radiation pressure and electron scattering based opacity, (5) the outer region far from the black hole horizon and outside the middle region dominated by gas pressure and electron scattering based opacity. The analytic forms for the important quantities like flux of radiant energy, F, temperature, T, and radial velocity in the locally nonrotating frame, β_r, were given for these different regions (as functions of r, a, viscosity, α, accretion rate, $\dot{m}=\dot{{M}_{\bullet }}/{\dot{M}}_{\mathrm{Edd}}$, and M_•), where nonzero stresses were incorporated at the inner edge of the disk in this model (Penna et al. 2012). Also, the expression for quality factor ${Q}_{\phi }\left(r,a,{\beta }_{r}\right)$ was derived for ν_ϕ QPO frequencies in the equatorial plane, which is given by (Mohan & Mangalam 2014, typo fixed in Equation (10))

$$\begin{eqnarray}&&{Q}_{\phi }\left(r,a,{\beta }_{r}\right)=\displaystyle \frac{-\sqrt{A}}{3\pi {\beta }_{r}{\rm{\Delta }}{r}^{1/2}}{\left[1-\displaystyle \frac{{\left(A{\rm{\Omega }}-2{ar}\right)}^{2}}{{{\rm{\Sigma }}}^{2}{\rm{\Delta }}}\right]}^{-1/2},\end{eqnarray} \tag{ 12 }$$

where $A={\left({r}^{2}+{a}^{2}\right)}^{2}-{a}^{2}{\rm{\Delta }}{\sin }^{2}\theta$, Δ = r² + a² − 2r, ${\rm{\Sigma }}={r}^{2}+{a}^{2}{\cos }^{2}\theta$, and ${\rm{\Omega }}=1/\left({r}^{3/2}+a\right)$, and where θ = π/2 is assumed in Equation (12). Using this formula, one can obtain the quality factor of the QPO in various regions close to the black hole by substituting the β_r of the corresponding region as defined above. The expressions for β_r in the edge and inner regions are given by (Equations (12), (13) of Mohan & Mangalam 2014)

$$\begin{eqnarray}&&\begin{array}{l}{\beta }_{r,\mathrm{edge}}=-7.1\times {10}^{-5}{\alpha }^{4/5}{m}_{1}^{-1/5}{\dot{m}}^{2/5}{r}^{-2/5}{{ \mathcal B }}^{4/5}{{ \mathcal C }}^{-1/2}{{ \mathcal D }}^{3/10}{{\rm{\Phi }}}^{-3/5},\end{array}\end{eqnarray} \tag{ 13a }$$

$$\begin{eqnarray}&&{\beta }_{r,\mathrm{inner}}=-124.416\ \alpha \ {\dot{m}}^{2}{r}^{-5/2}{{ \mathcal A }}^{2}{{ \mathcal B }}^{-3}{{ \mathcal C }}^{-1/2}{{ \mathcal D }}^{-1/2}{{ \mathcal S }}^{-1}{\rm{\Phi }},\end{eqnarray} \tag{ 13b }$$

where m₁ = M_•/10M_⊙; ${ \mathcal C }=1-3{r}^{-1}+2{{ar}}^{-3/2}$ (there is a typo in the expression of ${ \mathcal C }$, Equation A4(c), in Penna et al. 2012); and ${ \mathcal A }$, ${ \mathcal B }$, ${ \mathcal D }$, ${ \mathcal S }$, and Φ are given in Penna et al. (2012; Equations A4(a), (b), (d), (o) and (3.6)).

In Figures 28(a) and (b), we have shown the contours for β_r and Q_ϕ for the edge region in the (r, a) plane, and the p^gas/p^rad ratio as a function of r in Figure 28(c). One can discern the transition from the inner to edge region by the sudden increase of the p^gas/p^rad ratio, as seen in Figure 28(c), which is given by (Penna et al. 2012, Equation 3.7(g))

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{p}^{\mathrm{gas}}}{{p}^{\mathrm{rad}}}=1.983\times {10}^{-8}{m}_{1}^{-1/4}{\alpha }^{-1/4}{\dot{m}}^{-2}{r}^{21/8}{{ \mathcal A }}^{-5/2}{{ \mathcal B }}^{9/2}{ \mathcal D }{{ \mathcal S }}^{5/4}{{\rm{\Phi }}}^{-2}.\end{array}\end{eqnarray} \tag{ 14 }$$

Zoom In Zoom Out Reset image size

In Table 11, we give the range of {r, Q_ϕ, β_r, p^gas/p^rad} for the edge and inner regions for different combinations of a and $\dot{m}$, fixing {m₁ = 1, α = 0.1} for BHXRBs, with a low accretion rate ($\dot{m}\simeq 0.1$) corresponding to the hard spectral state and a high accretion rate ($\dot{m}\simeq 0.3$) corresponding to the soft spectral state of BHXRBs. We see a sharp rise in p^gas/p^rad values in the edge region in Figure 28(c). The ranges of Q_ϕ in both the edge and inner regions are very high compared to those observed in BHXRBs (Q_ϕ = 5–40). We, therefore, suggest that the QPOs are coming from a region very close to and inside ISCO; we identify this with the torus region, consisting of geodesics (Penna et al. 2012), and hence Q_ϕ is different. This is also supported by the observation that the edge-flow-sourced geodesics populate the torus region obtained here for M82 X-1 (r = 4.7–9.08) and GROJ 1655-40 (r = 5.1–6.67); see Figure 18. Specifically, the sharp pressure ratio gradient suggests that the edge region can be a launchpad for the instabilities that then oscillate with fundamental frequencies, causing geodesic flows in the torus region inside ISCO (r < r_ISCO), where the fluid motion is close to Hamiltonian flow. A further understanding of this proposal (or conjecture) can be gained by carrying out a detailed model or simulation of the GRMHD flow in the edge region.

Table 11.  Ranges of r, Pressure Ratio, p^gas/p^rad, Quality Factor, Q_ϕ, and Radial Velocity, β_r, in the Edge and Inner Regions of Fluid Flow in the Relativistic Thin Accretion Disk around a Kerr Black Hole (Penna et al. 2012; Mohan & Mangalam 2014), Where We Have Fixed {m₁ = 1, α = 0.1} for BHXRBs

Region	$\left(a=0.3,\dot{m}=0.1\right)$	$\left(a=0.5,\dot{m}=0.1\right)$	$\left(a=0.3,\dot{m}=0.3\right)$	$\left(a=0.5,\dot{m}=0.3\right)$
	$\left(r,{p}^{\mathrm{gas}}/{p}^{\mathrm{rad}},{\beta }_{r},{Q}_{\phi }\right)$	$\left(r,{p}^{\mathrm{gas}}/{p}^{\mathrm{rad}},{\beta }_{r},{Q}_{\phi }\right)$	$\left(r,{p}^{\mathrm{gas}}/{p}^{\mathrm{rad}},{\beta }_{r},{Q}_{\phi }\right)$	$\left(r,{p}^{\mathrm{gas}}/{p}^{\mathrm{rad}},{\beta }_{r},{Q}_{\phi }\right)$
Edge	4.98–5.93	4.23–4.87	4.98–5.25	4.23–4.35
	1.002–29.84	1.003–18.41	1.026–1.921	1.003–1.186
	$-\left(2.84\mbox{--}10.2\right)\times {10}^{-5}$	$-\left(3.81\mbox{--}11.37\right)\times {10}^{-5}$	$-\left(1.03\mbox{--}1.34\right)\times {10}^{-4}$	$-\left(1.36\mbox{--}1.49\right)\times {10}^{-4}$
	914.46–2624.12	1019.24–2473.07	694.29–844.59	773.62–814.27
Inner	5.93–85.22	4.87–87.81	5.25–226.2	4.35–229.45
	0.0589–0.998	0.0373–0.999	0.0065–0.998	0.0041–0.998
	$-\left(1.2052\mbox{--}69.28\right)\times {10}^{-6}$	$-\left(1.1601\mbox{--}110.71\right)\times {10}^{-6}$	$-\left(1.127\mbox{--}626.95\right)\times {10}^{-6}$	$-\left(1.107\mbox{--}1001.39\right)\times {10}^{-6}$
	688.09–9830.85	501.112–10046.3	76.21–6333.26	55.52–6397.35

Download table as:  ASCIITypeset image

The QPOs in BHXRBs have been an important probe for comprehending the inner accretion flow close to the rotating black hole. Many theoretical models have been proposed to explain its origin and in particular LFQPOs and HFQPOs (Kato 2004; Török et al. 2005; Tagger & Varnière 2006; Germanà et al. 2009; Ingram et al. 2009; Ingram & Done 2011, 2012). These various models have been able to explain different properties of QPOs. For example, one of these models attributes the HFQPOs to the Rossby instability under the general relativistic regime (Tagger & Varnière 2006), whereas another model attributes type C QPOs to the Lense–Thirring precession of a rigid torus of matter around a Kerr black hole (Ingram et al. 2009; Ingram & Done 2011, 2012). Although these models can explain either LFQPOs or HFQPOs, they do not explain the simultaneity of these QPOs, as previously observed in BHXRB GROJ 1655-40 (Motta et al. 2014a). The RPM, which is based on the geometric phenomenon of the relativistic precession of particle trajectories, explains these simultaneous QPOs as {ν_ϕ, $\left({\nu }_{\phi }-{\nu }_{r}\right)$, $\left({\nu }_{\phi }-{\nu }_{\theta }\right)$} of a self-emitting blob of matter (or instability) in a bound orbit near a Kerr black hole. We have extended the RPM for QPOs in BHXRBs to study and associate the fundamental frequencies of the bound particle trajectories near a Kerr black hole, which are eQ, e0, and Q0 solutions with the frequencies of QPOs. We call this the generalized RPM (GRPM). Recently, novel and compact analytic forms for the trajectories around a Kerr black hole and their fundamental frequencies were derived (Rana & Mangalam 2019a, 2019b). We applied these formulae to the GRPM to extract the QPO frequencies. Graphical examples of these trajectories around a Kerr black hole are shown in Figures 18, 23, and 29. A summary of these results is given in Table 12.

Zoom In Zoom Out Reset image size

Table 12.  Summary of Orbital Solutions Found for QPOs Observed in Five BHXRBs Using the GRPM in This Article, and the Corresponding Region of the (r_p, a) Plane Where QPOs Originate

BHXRB	Number	Model	e	rp	a	Q	MBSO	ISCO	ISSO	Region in
	of QPOs	Class								(rp, a) Plane
M82 X-1	3	eQ	${0.230}_{-0.049}^{+0.067}$	${4.834}_{-0.268}^{+0.181}$	0.299	${2.362}_{-1.439}^{+1.519}$	3.424	4.981	5.096	2
		Q0	0	${6.044}_{-0.072}^{+0.071}$	0.321 ± 0.013	${6.113}_{-1.645}^{+2.124}$	3.475	4.903	5.258	1
GROJ 1655-40	3	eQ	${0.071}_{-0.035}^{+0.031}$	${5.25}_{-0.142}^{+0.171}$	0.283	0+0.623	⋯	5.039	⋯	1
XTE J1550-564	2	e0	${0.262}_{-0.062}^{+0.090}$	${4.365}_{-0.279}^{+0.169}$	0.34	0	⋯	4.835	⋯	2
		Q0	0	5.538 ± 0.054	0.34	${2.697}_{-1.627}^{+1.738}$	3.35	4.835	4.988	1
4U 1630-47	2	e0	${0.734}_{-0.048}^{+0.066}$	${2.249}_{-0.353}^{+0.249}$	0.985	0	⋯	1.541	⋯	1
										1
GRS 1915+105	2	e0	0.918 ± 0.002	${1.744}_{-0.011}^{+0.025}$	0.98	0	⋯	1.614	⋯	1

Download table as:  ASCIITypeset image

We add the following caveats and conclusions:

• 1.  
  Novel and useful formulae: We have derived novel forms for the spherical trajectory solutions {$\phi \left({r}_{s},a,Q\right)$, $t\left({r}_{s},a,Q\right)$}, given by Equation (B12), and their fundamental frequencies {${\bar{\nu }}_{\phi }\left({r}_{s},a,Q\right)$, ${\bar{\nu }}_{r}\left({r}_{s},a,Q\right)$, ${\bar{\nu }}_{\theta }\left({r}_{s},a,Q\right)$}, given by Equation (7). A reduced form of the vertical oscillation frequency, ${\bar{\nu }}_{\theta }\left(e,{r}_{p},a\right)$ given by Equation 5(c), for equatorial eccentric orbits is also derived in Appendix A. These new and compact formulae are useful for various theoretical studies of Kerr orbits, besides other astrophysical applications (e.g., Rana & Mangalam 2020).
• 2.  
  Orbital solutions: The fundamental frequencies of the eQ, e0, and Q0 trajectories are in the range of QPO signals observed in BHXRBs, so these are viable solutions for explaining the observed QPOs in BHXRBs M82 X-1, GROJ 1655-40, XTEJ 1550-564, 4U 1630-47, and GRS 1915+105 in the GRPM paradigm. We see that these trajectory solutions are found in either region 1 or 2 of the (r, a) plane, as defined in Figure 5 and shown in Figure 22. The values of the black hole spin for BHXRBs M82 X-1 and GROJ 1655-40 were fixed to their most probable values calculated in Section 3.2.1 and to the previously observed values for the other BHXRBs for eccentric orbit solutions. For BHXRBs with two QPOs, fixing the spin to previously known values increases the uncertainty in the estimated orbital parameters, because the spin values assumed have uncertainties associated with the X-ray spectroscopic methods that are influenced by systematics, with the general finding that the solution lies near ISCO. However, our exercise still supports the GRPM. A spin value was also calculated for M82 X-1 for a Q0 solution. A summary of these parameter solutions and corresponding MBSO, ISCO, and ISSO radii for all BHXRBs is given in Table 12.
• 3.  
  Trajectories in the torus: We found trajectories, having different parameter combinations within the estimated range of errors in the orbital parameters and having fundamental frequencies within the width of the observed QPOs, as solutions for QPOs in BHXRBs M82 X-1 and GROJ 1655-40. We also found that the distinct parameter solutions found for these cases follow a trend that, as the eccentricity of the orbit decreases, the Q value increases for a given QPO frequency pair. This behavior can also be understood from Figures 8–10, where the frequencies increase as Q increases, but decrease as e increases for a given r_p. This implies that to obtain the degenerate parameter solutions for the same set of frequencies, a low eccentricity $\ \Longleftrightarrow \$ high Q trend is expected. We also found that these trajectories span a torus region near the Kerr black hole, as shown in Figure 18, which together give rise to the same peaks in the power spectrum. This should also explain the strong rms seen for the HFQPOs and type C LFQPOs. Another possibility of a rigidly precessing torus was suggested (Ingram et al. 2009; Ingram & Done 2011, 2012); our proposal consists of a nonprecessing torus, which includes all viable solutions of the GRPM: eQ, e0, and Q0 trajectories.
• 4.  
  Torus region: The emission of simultaneous QPOs is expected from a region where different trajectories having similar fundamental frequencies span a torus, as shown in Figure 18, and they can together show a strong peak in the power spectrum. The inner radius of the circular accretion disk is expected to be close to this torus region in such a scenario. In Figure 29, we depict this geometric model where the emission region of the simultaneous QPOs is shown as a torus region close to the inner edge of the accretion disk. This torus region is expected to be outside the MBSO radius, and the ISSO radius is expected to be in between the torus region for the eccentric orbit solutions, as observed in the case of M82 X-1. The torus region can be represented as ${{R}_{0}}_{-{{\rm{\Delta }}}_{2}}^{+{{\rm{\Delta }}}_{1}}$, where R₀ is an e = 0 orbit (ISCO or ISSO) and Δ_i represents the region very close to R₀. The width of the torus region in this model is given by ${\rm{\Delta }}r=\left({{\rm{\Delta }}}_{1}+{{\rm{\Delta }}}_{2}\right)$. All of the orbit solutions are found to be distributed near R₀; hence, it is expected that this radius corresponds to the inner edge radius, r_in, of the circular accretion disk. This torus region exists in region 1 and(or) 2 near the R₀ radius. Due to the instabilities in the inner flow, we argue that the nearly e0 orbits near the R₀ radius transcend to eQ orbits. Based on the geometry of the orbits and the emission region, we plan to build a detailed GRMHD-based model to expand on the GRPM paradigm. More cases of BHXRBs with three simultaneous QPOs, if detected in the future, will help us test our models.
• 5.  
  Highly eccentric solutions: For QPOs in BHXRBs 4U 1630-47 and GRS 1915+105, we found highly eccentric e0 solutions. This indicates that black holes with high spin values prefer highly eccentric trajectories as solutions to the QPOs. This behavior can also be understood from Figures 8–10, where we see that for black holes with very high spins, the QPOs originate very close to the black hole, and the solution contours move close to the black hole as e increases. This implies that more eccentric orbits are preferred for a given frequency pair of QPOs for a black hole with very high spin. We do not find any spherical orbit solution for QPOs in these two BHXRBs, which confirms that the orbital solution is purely equatorial, but such highly eccentric solutions are unlikely. We expect more and better estimates of the orbital solutions in the future if a more precise estimate of the spin is available, or if three simultaneous QPOs are discovered in BHXRBs 4U 1630-47 and GRS 1915+105. For the case that we studied in this paper of 4U 1630-47, the lower frequency of the QPO pair probably has a different origin than the high-frequency feature suggested by Klein-Wolt et al. (2004). However, even in such a scenario, the frequency range of this QPO still implies an origin near the torus region in our model. There was also another pair of QPOs observed in 4U 1630-47 (Klein-Wolt et al. 2004), for which there was no exact solution found in the orbital parameter space.
• 6.  
  Nonequatorial solutions: In the case of BHXRBs M82 X-1 and XTE 1550-564, we found both eQ (e0 for XTE 1550-564) and Q0 solutions, and the spin determinations are slightly different for the two different types of trajectory solutions. These solutions were found close to and outside their corresponding ISSO radii. The mass of the black hole in M82 X-1 was fixed to the intermediate-mass black hole (IMBH) range, ${ \mathcal M }=428$, because the QPOs observed in the low-frequency range (3–5 Hz) were found to be very stable, unlike LFQPOs, implying that they are HFQPO counterparts of BHXRBs, and hence indicating an IMBH (Pasham et al. 2014). Although this mass estimation stems from the mass-scaling relation of QPOs, which is not very reliable, a more accurate estimate of ${ \mathcal M }$, if found in the IMBH range, will not significantly change the result. However, if, in the future, a more reliable and precise estimate places it in the stellar-mass range, then the outcome from the GRPM will be dramatically different. The QPOs observed in XTE 1550-564 by Miller et al. (2001) were later shown to be the result of the data averaging by Motta et al. (2014b), where the same QPO moved up in frequency, appearing as a distinct QPO. As in the case of 4U 1630-47, the range of this QPO frequency still implies an origin near the torus region.
• 7.  
  Spectral states: We suggest that HFQPOs originate when r_in comes very close (near ISCO/ISSO) to the black hole during the soft spectral state of the outburst. When r_in is farther out as in the hard state, the resulting type C QPO frequency is of the order of millihertz. As a type C QPO occurs more frequently and is prone to vertical oscillations, the increase in its frequency is explained as an increase in ν_np when r_in decreases, with the spectral transition from the hard to soft state.
• 8.  
  Circularity: The RPM was previously applied to understand the QPOs observed in BHXRBs GROJ 1655-40 and XTEJ 1550-564 (Motta et al. 2014a, 2014b) using the fundamental frequencies of 00 orbits. We have found an eQ solution for GROJ 1655-40 very close to an equatorial orbit having a very small eccentricity e ∼ 0.071 (see Table 12), which is in very close agreement with the solution found by Motta et al. (2014a), where their estimated mass of the black hole, ${ \mathcal M }=5.307$, is also very close to our assumption, ${ \mathcal M }=5.4$ (see Table 5). Our most probable spin estimated for GROJ 1655-40, a = 0.283, is almost the same as found by Motta et al. (2014a), a ∼ 0.286, but our solution provides a more precise estimation of e and Q values while confirming a near 00 orbit solution as assumed by Motta et al. (2014a). For the case of XTEJ 1550-564, the mass of the black hole was assumed to be ${ \mathcal M }=9.1$ by Motta et al. (2014b), as also in our model. Our assumption for the spin, $a={0.34}_{-0.45}^{+0.37}$ (Orosz et al. 2011), is also nearly the same as the value estimated by Motta et al. (2014b), but our model gives the e0 and Q0 solutions for XTEJ 1550-564, having moderate $e={0.262}_{-0.062}^{+0.090}$ and $Q={2.697}_{-1.627}^{+1.738}$ values, respectively (see Table 12). This indicates that the assumption of circularity is not always valid.
• 9.  
  Solution degeneracy: To study the impact of the GRPM (with nonzero e and Q), we have explored the behavior of {δ_ϕ, δ_pp, δ_np}(e, r_p, a, Q) as defined in Equation (6) as deviations from the 00 behavior (circularity). We find that the frequencies are strongly dependent on e but not so much on Q (see Figures 8–10). This is elaborated upon in points 3 and 4 in Section 2.1, and in points 2 and 3 in Section 2.2 for spherical orbits. The GRPM has a built-in degeneracy in the parameter space {e, r_p, a, Q}, called the isofrequency pairs, for a given combination of QPO frequencies. This degeneracy is a known behavior of trajectories around a Kerr black hole (Warburton et al. 2013), where different combinations of {E, L_z, Q} can have the same set {ν_ϕ, ν_r, ν_θ} for a fixed a. An evidence of this degeneracy is also seen in Figures 8–10, where the contours of {δ_ϕ, δ_pp, δ_np}(e, r_p, a, Q) have multiple solutions; that is, for a given δ value, there are different combinations of {e, Q} that have distinct contours on the (r_p, a) plane. Unlike RPM, the mass of the black hole is always assumed from the previous estimates in the GRPM, which is a valid assumption because the underlying physics or behavior of the Kerr orbits is independent of M_•. The GRPM, along with the statistical method (Figure 14, Appendix E) that is applied, provides a more precise estimation of the spin of the black hole.
• 10.  
  Frequency ratio: The 3:2 and 5:3 ratios of the simultaneous HFQPOs are a phenomenon observed in a few cases of BHXRBs: 300 and 450 Hz HFQPOs in GROJ 1655-40 (Remillard et al. 1999b; Strohmayer 2001a), and 240 and 160 Hz HFQPOs in H1743-322 (Homan et al. 2005; Remillard et al. 2006). Such claims, other than the case of GROJ 1655-40, are probably not real (Belloni et al. 2012). Hence, the possibility of such ratios still causes skepticism. However, if true, the GRPM suggests that the origin of these ratios is very close to the torus region and r_in.
• 11.  
  The PBK correlation: In Section 4, we show that the e0 solution {e = 0.071, a = 0.283, Q = 0, ${ \mathcal M }=5.4$}, estimated using a fine-grid method in Section 3.2.1, fits the PBK correlation that was previously observed in BHXRB GROJ 1655-40 (Motta et al. 2014a). This fit is shown in Figure 26. We also found that 10 observation IDs, where L_l and L_lf (broad frequency components) were detected simultaneously (Motta et al. 2014a), show low-eccentricity eQ solutions, where the calculated parameters are shown in Table 10. The calculated Q values are consistent with large r_p and small e values. This exercise suggests that eQ solutions for QPOs are viable.
• 12.  
  Probing the disk edge with a GR fluid model: We study a model of fluid flow in the general-relativistic thin accretion disk (Penna et al. 2012; Mohan & Mangalam 2014). We find that the disk edge flows into a torus region containing Hamiltonian geodesics that was obtained for M82 X-1 (r = 4.7–9.08) and GROJ 1655-40 (r = 5.11–6.67). Specifically, the sharp gradient of the p^gas/p^rad pressure ratio, seen in Figure 28(c), suggests that the edge region is a launch pad for the instabilities that orbit with fundamental frequencies of the geodesics in the edge and geodesic regions, which then follow the geodesics inside the torus region and also close to the edge region, where Hamiltonian dynamics is applicable, that is built into the GRPM. The range of {r, Q_ϕ, β_r, p^gas/p^rad} for the edge and inner regions for different combinations of a and $\dot{m}$, fixing {m₁ = 1, α = 0.1}, is given in Table 11, and the contours of β_r and Q_ϕ in the (r, a) plane for the edge region are shown in Figures 28(a) and (b). The ranges of Q_ϕ (tuned to Δν, the width of the observed QPO), which is defined by orbits in the torus which was provided by observed frequency centroids, in both the edge and inner regions are very high compared to those observed in BHXRBs (Q_ϕ = 5–40). We are suggesting that the QPOs originate in the geodesic region. We also see that the edge is adjacent to the torus region (consisting of geodesics) found for M82 X-1 and GROJ 1655-40, implying that the QPOs are originating from geodesics close to the edge region. Hence, the particle and gas dynamics models together justify the scenario sketched in Figure 29, of a unified fluid-particle picture that is the following: the source of the particles in the torus are dynamical instabilities of plasma blobs ejected from the edge region. These blobs have zero α and therefore obey the Hamiltonian dynamics. The clue that the torus region physically overlaps the edge and geodesic regions is a subject of future detailed GRMHD models (and simulations).
• 13.  
  Isofrequency combinations: In the cases with three simultaneous QPOs, once a is fixed (to the most probable value or the previously estimated value), it is easy to predict the remaining parameters {e, r_p, Q} using three QPO frequencies. In the case of M82 X-1 and GROJ 1655-40, when a was fixed to the most probable value (Table 7), we obtained a single solution for {e, r_p, Q} and their errors {Δe, Δr_p, ΔQ}, where this range of parameters spans the torus region based on the GRPM. However, there is a finite possibility (Warburton et al. 2013) that distinct solutions for the {e, r_p, Q} triad are obtained for the same triple QPO frequency set, subject to the bound orbit conditions: 0 ≤ e < 1, Q ≥ 0, and Equation (2). This completely depends on the values of the QPO frequency set that are further subject to the constraints of bound orbit conditions. In the cases where only two simultaneous QPOs exist, it is difficult to predict whether an e = 0 orbit will be preferred over an e > 0 orbit, or a Q = 0 orbit will be preferred over a Q > 0 orbit, or vice versa. This will be clear when more cases of three simultaneous QPOs are found and whether they yield distinct solution sets for {e, r_p, a, Q}, thereby indicating if the torus region at the disk edge is indeed the geometric origin of QPOs. From our numerical experiment, we find a distinct exact solution for {e, r_p, Q} for the three simultaneous QPOs case, where a was fixed to the most probable value. The RPM restricts the search to {e = 0, Q = 0} orbital solutions, while the GRPM expands it to more general but astrophysically possible {$e\ne 0,Q\ne 0$} solutions and thereby subsumes the RPM within its framework. Hence, the GRPM provides more realistic orbit solutions around a Kerr black hole that are outside the scope of the RPM, thus giving more impetus to probes of physical models of the origin of QPOs.
• 14.  
  Caveats: The results predicted by the GRPM are subject to the veracity of the observed data that are inputs to our model. For example, in the case of 4U 1630-47 and GRS 1915+105, very highly eccentric orbit solutions obtained by the GRPM are unlikely; this implies that very high spin values in these cases are probably unreliable. Similarly, if M82 X-1 does not host an IMBH but a stellar-mass black hole or a neutron star, then the results predicted by the GRPM will change drastically. Also, for 4U 1630-47 and XTEJ 1550-564, where the input frequencies of QPOs are not very reliable (Klein-Wolt et al. 2004; Motta et al. 2014b), as discussed before, the results obtained by the GRPM might not be physically meaningful. As most of the measured frequencies do exist in a similar range, then their geometric origin in the torus region (as predicted by the GRPM) is valid.
• 15.  
  Future work: In the near future, we expect suitable observational results from the currently operative Indian X-ray satellite, AstroSat, and from future missions, such as eXTP, which is proposed to have instruments with much higher sensitivity for fast variations and X-ray timing. If simultaneous QPO signals are observed from these missions, we expect to test our GRPM further.

We would like to thank the anonymous referee for detailed and insightful suggestions that have improved our paper significantly. We acknowledge the DST SERB grant No. CRG 2018/003415 for financial support. We would like to thank Dr. Prashanth Mohan for his useful suggestions and comments. We would like to thank Saikat Das for helping us with Figures 3 and 29. We acknowledge the use and support of the IIA-HPC facility.

Here, we derive the θ oscillation frequency for the equatorial eccentric orbits about the equatorial plane. Using Equations 4(b) and (c), we can write

$$\begin{eqnarray}&&\displaystyle \frac{{\bar{\nu }}_{\theta }}{{\bar{\nu }}_{r}}=\displaystyle \frac{a\sqrt{1-{E}^{2}}{z}_{+}{I}_{8}\left(e,\mu ,a,Q\right)}{2F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)},\end{eqnarray} \tag{ A1a }$$

where the substitution of ${I}_{8}\left(e,\mu ,a,Q\right)$ from Equation 6(h) of Rana & Mangalam (2019a) into the above equation yields

$$\begin{eqnarray}&&\displaystyle \frac{{\bar{\nu }}_{\theta }}{{\bar{\nu }}_{r}}=\displaystyle \frac{\mu \left(1-{e}^{2}\right)a\sqrt{1-{E}^{2}}{z}_{+}F\left(\tfrac{\pi }{2},{k}^{2}\right)}{\sqrt{C-A+\sqrt{{B}^{2}-4{AC}}}F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)}.\end{eqnarray} \tag{ A1b }$$

By the substitution of A, B, and C using Equations 7(f)–(h) of Rana & Mangalam (2019a), and using Q = 0 for the equatorial orbits, we find

$$\begin{eqnarray}&&\begin{array}{l}\sqrt{C-A+\sqrt{{B}^{2}-4{AC}}}={\mu }^{1/2}\left(1-{e}^{2}\right)\left[1\right.-{\left.{\mu }^{2}{x}^{2}\left(3-{e}^{2}-2e\right)\right]}^{1/2}.\end{array}\end{eqnarray} \tag{ A2a }$$

Also, from Equation 9(d) of Rana & Mangalam (2019a), we see that

$$\begin{eqnarray}\begin{array}{rclcl}{z}_{-} & = & 0,{z}_{+} & = & \displaystyle \frac{\sqrt{{L}_{z}^{2}+{a}^{2}\left(1-{E}^{2}\right)}}{a\sqrt{\left(1-{E}^{2}\right)}}=\displaystyle \frac{\sqrt{{x}^{2}+{a}^{2}+2{aEx}}}{a\sqrt{\left(1-{E}^{2}\right)}},\end{array}\end{eqnarray} \tag{ A2b }$$

for Q = 0, which implies that

$$\begin{eqnarray}&&F\left(\displaystyle \frac{\pi }{2},\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)=\displaystyle \frac{\pi }{2}.\end{eqnarray} \tag{ A2c }$$

Hence, Equations A1(b)–A2(c) together reduce ${\bar{\nu }}_{\theta }/{\bar{\nu }}_{r}$ for equatorial orbits to

$$\begin{eqnarray}&&\displaystyle \frac{{\bar{\nu }}_{\theta }}{{\bar{\nu }}_{r}}=\displaystyle \frac{2{\mu }^{1/2}\sqrt{{x}^{2}+{a}^{2}+2{aEx}}\cdot F\left(\tfrac{\pi }{2},{k}^{2}\right)}{\pi {\left[1-{\mu }^{2}{x}^{2}\left(3-{e}^{2}-2e\right)\right]}^{1/2}}.\end{eqnarray} \tag{ A3 }$$

We see from Equations 7(f)–(j) of Rana & Mangalam (2019a) that ${k}^{2}=\left({n}^{2}-{m}^{2}\right)/\left(1-{m}^{2}\right)$ can be written in terms of A, B, and C as

$$\begin{eqnarray}&&{k}^{2}=\displaystyle \frac{2\sqrt{{B}^{2}-2{AC}}}{\left(-A+C+\sqrt{{B}^{2}-2{AC}}\right)},\end{eqnarray} \tag{ A4 }$$

where the substitution of A, B, and C for Q = 0 gives

$$\begin{eqnarray}&&{k}^{2}={m}^{2}=\displaystyle \frac{4{{ex}}^{2}{\mu }^{2}}{\left[1-{\mu }^{2}{x}^{2}\left(3-{e}^{2}-2e\right)\right]}.\end{eqnarray} \tag{ A5 }$$

Hence, we can write ${\bar{\nu }}_{\theta }$ for the equatorial orbits as

$$\begin{eqnarray}&&\begin{array}{l}{\bar{\nu }}_{\theta }\left(e,\mu ,a\right)=\ \displaystyle \frac{2{\bar{\nu }}_{r}\left(e,\mu ,a\right){\mu }^{1/2}\sqrt{\left({x}^{2}+{a}^{2}+2{aEx}\right)}\cdot F\left(\tfrac{\pi }{2},{k}^{2}\right)}{\pi {\left[1-{\mu }^{2}{x}^{2}\left(3-{e}^{2}-2e\right)\right]}^{1/2}},\end{array}\end{eqnarray} \tag{ A6 }$$

where ${\bar{\nu }}_{r}\left(e,\mu ,a\right)$ is given by Equation 5(b) and k² is given by Equation (A5).

• 1.  
  Azimuthal angle and coordinate time: The integrals of motion for a general nonequatorial trajectory of a particle with rest mass m₀ around a Kerr black hole have been derived using the Hamilton–Jacobi method, in terms of the Boyer–Lindquist coordinates (r, ϕ, θ, t) (Carter 1968; Schmidt 2002):

  $$\begin{eqnarray}\begin{array}{rcl}\phi -{\phi }_{0} & = & -\displaystyle \frac{1}{2}{\displaystyle \int }_{{r}_{0}}^{r}\displaystyle \frac{1}{{\rm{\Delta }}\sqrt{R}}\displaystyle \frac{\partial R}{\partial {L}_{z}}{{dr}}^{{\prime} }-\displaystyle \frac{1}{2}{\displaystyle \int }_{{\theta }_{0}}^{\theta }\displaystyle \frac{1}{\sqrt{{\rm{\Theta }}}}\displaystyle \frac{\partial {\rm{\Theta }}}{\partial {L}_{z}}d{\theta }^{{\prime} }\,=\ -\displaystyle \frac{1}{2}{I}_{1}-\displaystyle \frac{1}{2}{H}_{1},\end{array}\end{eqnarray} \tag{ B7a }$$

  $$\begin{eqnarray}\begin{array}{rcl}t-{t}_{0} & = & \displaystyle \frac{1}{2}{\displaystyle \int }_{{r}_{0}}^{r}\displaystyle \frac{1}{{\rm{\Delta }}\sqrt{R}}\displaystyle \frac{\partial R}{\partial E}{{dr}}^{{\prime} }+\displaystyle \frac{1}{2}{\displaystyle \int }_{{\theta }_{0}}^{\theta }\displaystyle \frac{1}{\sqrt{{\rm{\Theta }}}}\displaystyle \frac{\partial {\rm{\Theta }}}{\partial E}d{\theta }^{{\prime} }\,=\,\displaystyle \frac{1}{2}{I}_{2}+\displaystyle \frac{1}{2}{H}_{2},\end{array}\end{eqnarray} \tag{ B7b }$$

  $$\begin{eqnarray}&&{\int }_{{r}_{0}}^{r}\displaystyle \frac{{{dr}}^{{\prime} }}{\sqrt{R}}={\int }_{{\theta }_{0}}^{\theta }\displaystyle \frac{d{\theta }^{{\prime} }}{\sqrt{{\rm{\Theta }}}}\Rightarrow {I}_{8}={H}_{3},\end{eqnarray} \tag{ B7c }$$

  where R and Θ are given by

  $$\begin{eqnarray}&&R={\left[\left({{r}^{{\prime} }}^{2}+{a}^{2}\right)E-{{aL}}_{z}\right]}^{2}-{\rm{\Delta }}\left[{{r}^{{\prime} }}^{2}+{\left({L}_{z}-{aE}\right)}^{2}+Q\right],\end{eqnarray} \tag{ B7d }$$

  $$\begin{eqnarray}&&{\rm{\Theta }}=Q-\left[\left(1-{E}^{2}\right){a}^{2}+\displaystyle \frac{{L}_{z}^{2}}{{\sin }^{2}{\theta }^{{\prime} }}\right]{\cos }^{2}{\theta }^{{\prime} }.\end{eqnarray} \tag{ B7e }$$

  We have from Equation B7(c) that

  $$\begin{eqnarray}&&\displaystyle \frac{{{dr}}^{{\prime} }}{\sqrt{R}}=\displaystyle \frac{d{\theta }^{{\prime} }}{\sqrt{{\rm{\Theta }}}};\end{eqnarray} \tag{ B8 }$$

  the substitution of the above equation into Equations B7(a), (b) for the spherical orbits reduces the expressions of the azimuthal angle and coordinate time to

  $$\begin{eqnarray}\begin{array}{rclcl}\phi -{\phi }_{0} & = & -\displaystyle \frac{1}{2}\left[\displaystyle \frac{1}{{\rm{\Delta }}}\displaystyle \frac{\partial R}{\partial {L}_{z}}{H}_{3}+{H}_{1}\right],t-{t}_{0} & = & \displaystyle \frac{1}{2}\left[\displaystyle \frac{1}{{\rm{\Delta }}}\displaystyle \frac{\partial R}{\partial E}{H}_{3}+{H}_{2}\right].\end{array}\end{eqnarray} \tag{ B9 }$$

  Since r = r_s is constant for the spherical orbits, the expressions $\tfrac{1}{{\rm{\Delta }}}\tfrac{\partial R}{\partial {L}_{z}}$ and $\tfrac{1}{{\rm{\Delta }}}\tfrac{\partial R}{\partial E}$ can be written as

  $$\begin{eqnarray}\begin{array}{rclcl}\displaystyle \frac{1}{{\rm{\Delta }}}\displaystyle \frac{\partial R}{\partial {L}_{z}} & = & \displaystyle \frac{2\left(2{L}_{z}{r}_{s}-{L}_{z}{r}_{s}^{2}-2{r}_{s}{aE}\right)}{{\rm{\Delta }}},\displaystyle \frac{1}{{\rm{\Delta }}}\displaystyle \frac{\partial R}{\partial E} & = & \displaystyle \frac{2\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}},\end{array}\end{eqnarray} \tag{ B10 }$$

  and the integrals H₁, H₂, and H₃ have been previously derived to be (Fujita & Hikida 2009; Rana & Mangalam 2019a)

  $$\begin{eqnarray}&&\begin{array}{l}{H}_{1}\left(\theta ,{\theta }_{0},e,\mu ,a,Q\right)=\displaystyle \frac{2{L}_{z}}{{z}_{+}a\sqrt{1-{E}^{2}}}\left\{F\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-\ F\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.\\ \quad \ +\,{\rm{\Pi }}\left({z}_{-}^{2},\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\left.\,-\,{\rm{\Pi }}\left({z}_{-}^{2},\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\},\end{array}\end{eqnarray} \tag{ B11a }$$

  $$\begin{eqnarray}&&\begin{array}{l}{H}_{2}\left(\theta ,{\theta }_{0},e,\mu ,a,Q\right)=\displaystyle \frac{2{{Eaz}}_{+}}{\sqrt{1-{E}^{2}}}\left\{K\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.-F\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\\ \qquad \ -\ K\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)+\left.F\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\},\end{array}\end{eqnarray} \tag{ B11b }$$

  $$\begin{eqnarray}&&\begin{array}{l}{H}_{3}\left(\theta ,{\theta }_{0},e,\mu ,a,Q\right)=\displaystyle \frac{1}{a\sqrt{1-{E}^{2}}{z}_{+}}\left\{F\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.\\ \,\,-\,\left.F\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\},\end{array}\end{eqnarray} \tag{ B11c }$$

  where z_± are given by Equation 9(d) of Rana & Mangalam (2019a). Hence, the substitution of Equations (B10) and (B11) into Equation (B9) yields the expressions of (ϕ − ϕ₀, t − t₀) for the spherical orbits, given by

  $$\begin{eqnarray}\begin{array}{rcl}\phi -{\phi }_{0} & = & \displaystyle \frac{1}{a\sqrt{1-{E}^{2}}{z}_{+}}\left\{\displaystyle \frac{\left({a}^{2}{L}_{z}-2{{aEr}}_{s}\right)}{{\rm{\Delta }}}\right.\left[F\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-F\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right]\\ & & -\ {L}_{z}\left[{\rm{\Pi }}\left({z}_{-}^{2},\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.-{\rm{\Pi }}\left.\left.\left({z}_{-}^{2},\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right]\right\},\end{array}\end{eqnarray} \tag{ B12a }$$

  $$\begin{eqnarray}\begin{array}{rcl}t-{t}_{0} & = & \displaystyle \frac{1}{a\sqrt{1-{E}^{2}}{z}_{+}}\left\{{{Ea}}^{2}{z}_{+}^{2}\left[K\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.\right.-\left.K\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right]\\ & & +\ \left[F\left(\arcsin \left(\displaystyle \frac{\cos {\theta }_{0}}{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.-\left.F\left(\arcsin \left(\displaystyle \frac{\cos \theta }{{z}_{-}}\right),\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right]\\ & & \cdot \ \left.\left[{{Ea}}^{2}{z}_{+}^{2}+\displaystyle \frac{E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}}{{\rm{\Delta }}}\right]\right\}.\end{array}\end{eqnarray} \tag{ B12b }$$
• 2.  
  Fundamental frequencies: The closed forms for fundamental frequencies associated with the nonequatorial eccentric bound trajectories have been previously derived (Schmidt 2002; Rana & Mangalam 2019a) and are given by Equations 4(a)–(c). We first reduce the common denominator of these expressions to the case of spherical orbits. If we take ${I}_{8}\left(e,\mu ,a,Q\right)$ common from the denominator, it gives

  $$\begin{eqnarray}&&\begin{array}{l}\quad \ \left[\left({I}_{2}+2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\right)F\left(\displaystyle \frac{\pi }{2},\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-2{a}^{2}{z}_{+}^{2}{{EI}}_{8}\,K\left(\displaystyle \frac{\pi }{2},\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right]\\ \quad \ =\ {I}_{8}\left[\left(\displaystyle \frac{{I}_{2}}{{I}_{8}}+2{a}^{2}{z}_{+}^{2}E\right)F\left(\displaystyle \frac{\pi }{2},\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right.-\left.2{a}^{2}{z}_{+}^{2}{EK}\left(\displaystyle \frac{\pi }{2},\displaystyle \frac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right],\end{array}\end{eqnarray} \tag{ B13 }$$

  where by definition ${I}_{2}/{I}_{8}=\tfrac{1}{{\rm{\Delta }}}\tfrac{\partial R}{\partial E}$ for spherical orbits, which is given by Equation (B10). Hence, Equations (B13), (B10), and 4(c) combine to give the vertical oscillation frequency for the spherical orbits:

  $$\begin{eqnarray}&&\bar{{\nu }_{\theta }}\left({r}_{s},a,Q\right)=\displaystyle \frac{a\sqrt{1-{E}^{2}}{z}_{+}}{4\left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}{EK}\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}}.\end{eqnarray} \tag{ B14 }$$

  Next, using Equation (B13), the azimuthal frequency, Equation 4(a), can be written as

  $$\begin{eqnarray}&&\bar{{\nu }_{\phi }}\left({r}_{s},a,Q\right)=\displaystyle \frac{\left\{\left[-\tfrac{{I}_{1}}{{I}_{8}}-2{L}_{z}\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)+2{L}_{z}\cdot {\rm{\Pi }}\left({z}_{-}^{2},\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}}{4\pi \left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}{EK}\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ B15a }$$

  where ${I}_{1}/{I}_{8}=\tfrac{1}{{\rm{\Delta }}}\tfrac{\partial R}{\partial {L}_{z}}$, which is given by Equation (B10). Hence, the azimuthal frequency for the spherical orbits is given by

  $$\begin{eqnarray}&&\bar{{\nu }_{\phi }}\left({r}_{s},a,Q\right)=\displaystyle \frac{\left\{\left[-\tfrac{\left(2{L}_{z}{r}_{s}-{L}_{z}{r}_{s}^{2}-2{r}_{s}{aE}\right)}{{\rm{\Delta }}}-{L}_{z}\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)+{L}_{z}\cdot {\rm{\Pi }}\left({z}_{-}^{2},\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}}{2\pi \left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}{EK}\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}}.\end{eqnarray} \tag{ B15b }$$

  Similarly, the radial oscillation frequency, Equation 4(b), can be written for the spherical orbits by using Equation (B13) as

  $$\begin{eqnarray}&&\bar{{\nu }_{r}}\left({r}_{s},a,Q\right)=\displaystyle \frac{F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)}{2{I}_{8}\left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}{EK}\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}},\end{eqnarray} \tag{ B16a }$$

  where, for spherical orbits the integral I₈ reduces to a constant, as shown below.We see that the expression for k², Equation (A4), reduces to zero because A = B = 0 (Equations 7(f), (g) of Rana & Mangalam 2019a) for spherical orbits (e = 0). Hence, ${I}_{8}\left(e=0,\mu ,a,Q\right)$ (Equation 6(h) of Rana & Mangalam 2019a) reduces to

  $$\begin{eqnarray}\begin{array}{rcl}{I}_{8} & = & \displaystyle \frac{2\mu }{\sqrt{C}}F\left(\displaystyle \frac{\pi }{2},{k}^{2}=0\right)\,=\,\displaystyle \frac{\pi {r}_{s}}{\sqrt{{r}_{s}^{4}\left(1-{E}^{2}\right)+\left(3{{Qa}}^{2}-2{x}^{2}{r}_{s}-2{{Qr}}_{s}\right)}}.\end{array}\end{eqnarray} \tag{ B16b }$$

  Hence, the radial oscillation frequency for spherical orbits reduces to

  $$\begin{eqnarray}&&\bar{{\nu }_{r}}\left({r}_{s},a,Q\right)=\displaystyle \frac{\sqrt{{r}_{s}^{4}\left(1-{E}^{2}\right)+\left(3{{Qa}}^{2}-2{x}^{2}{r}_{s}-2{{Qr}}_{s}\right)}\cdot F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)}{2\pi {r}_{s}\left\{\left[\tfrac{\left[E\left({a}^{2}{r}_{s}^{2}+{r}_{s}^{4}+2{a}^{2}{r}_{s}\right)-2{L}_{z}{{ar}}_{s}\right]}{{\rm{\Delta }}}+{a}^{2}{z}_{+}^{2}E\right]F\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)-{a}^{2}{z}_{+}^{2}{EK}\left(\tfrac{\pi }{2},\tfrac{{z}_{-}^{2}}{{z}_{+}^{2}}\right)\right\}}.\end{eqnarray} \tag{ B16c }$$

Here, we reduce the fundamental frequency formulae to the known case of equatorial circular orbits (00). We show this reduction from both the equatorial eccentric (e0) and the spherical (Q0) orbits below:

• 1.  
  Reduction from e0 orbits: We see that for circular orbits (e = 0), the expressions of m², ${{p}_{1}}^{2}$, ${{p}_{2}}^{2}$, and ${{p}_{3}}^{2}$ (Equations 7(i), (k) of Rana & Mangalam 2019a) reduce to

  $$\begin{eqnarray}&&{m}^{2}={{p}_{1}}^{2}={{p}_{2}}^{2}={{p}_{3}}^{2}=0.\end{eqnarray} \tag{ C17 }$$

  We first make the subtitution m² = 0 in Equations 5(a)–(c), which gives

  $$\begin{eqnarray}&&{\bar{\nu }}_{\phi }=\displaystyle \frac{{a}_{1}{\rm{\Pi }}\left(-{p}_{2}^{2},\tfrac{\pi }{2},0\right)+{b}_{1}{\rm{\Pi }}\left(-{p}_{3}^{2},\tfrac{\pi }{2},0\right)}{2\pi \left\{{\rm{\Pi }}\left(-{p}_{1}^{2},\tfrac{\pi }{2},0\right)\left[{a}_{2}\tfrac{\left({p}_{1}^{2}+2\right)}{2\left(1+{p}_{1}^{2}\right)}+{b}_{2}\right]+{c}_{2}{\rm{\Pi }}\left(-{p}_{2}^{2},\tfrac{\pi }{2},0\right)+{d}_{2}{\rm{\Pi }}\left(-{p}_{3}^{2},\tfrac{\pi }{2},0\right)\right\}},\end{eqnarray} \tag{ C18a }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{r}=\displaystyle \frac{1}{2\left\{{\rm{\Pi }}\left(-{p}_{1}^{2},\tfrac{\pi }{2},0\right)\left[{a}_{2}\tfrac{\left({p}_{1}^{2}+2\right)}{2\left(1+{p}_{1}^{2}\right)}+{b}_{2}\right]+{c}_{2}{\rm{\Pi }}\left(-{p}_{2}^{2},\tfrac{\pi }{2},0\right)+{d}_{2}{\rm{\Pi }}\left(-{p}_{3}^{2},\tfrac{\pi }{2},0\right)\right\}},\end{eqnarray} \tag{ C18b }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{\theta }=\displaystyle \frac{{\bar{\nu }}_{r}{\mu }^{1/2}\sqrt{\left({x}^{2}+{a}^{2}+2{aEx}\right)}}{\sqrt{1-3{\mu }^{2}{x}^{2}}}.\end{eqnarray} \tag{ C18c }$$

  Next, the substitution of ${{p}_{1}}^{2}={{p}_{2}}^{2}={{p}_{3}}^{2}=0$ in Equation (C18) yields

  $$\begin{eqnarray}&&{\bar{\nu }}_{\phi }=\displaystyle \frac{{a}_{1}+{b}_{1}}{2\pi \left({a}_{2}+{b}_{2}+{c}_{2}+{d}_{2}\right)},\end{eqnarray} \tag{ C19a }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{r}=\displaystyle \frac{1}{\pi \left({a}_{2}+{b}_{2}+{c}_{2}+{d}_{2}\right)},\end{eqnarray} \tag{ C19b }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{\theta }=\displaystyle \frac{{\bar{\nu }}_{r}{\mu }^{1/2}\sqrt{\left({x}^{2}+{a}^{2}+2{aEx}\right)}}{\sqrt{1-3{\mu }^{2}{x}^{2}}}.\end{eqnarray} \tag{ C19c }$$

  By substituting e = 0 in Equation (16) of Rana & Mangalam (2019b), we find that

  $$\begin{eqnarray}&&{a}_{1}+{b}_{1}=\displaystyle \frac{2{\mu }^{1/2}\left({L}_{z}-2x\mu \right)}{\sqrt{1-3{\mu }^{2}{x}^{2}}\left(1-2\mu +{a}^{2}{\mu }^{2}\right)},\end{eqnarray} \tag{ C20a }$$

  $$\begin{eqnarray}&&\begin{array}{l}{a}_{2}+{b}_{2}+{c}_{2}+{d}_{2}=\displaystyle \frac{2\left(E+{{Ea}}^{2}{\mu }^{2}-2{ax}{\mu }^{3}\right)}{{\mu }^{3/2}\sqrt{1-3{\mu }^{2}{x}^{2}}\left(1-2\mu +{a}^{2}{\mu }^{2}\right)}.\end{array}\end{eqnarray} \tag{ C20b }$$

  Now, by substituting Equation (C20) in Equation (C19), we get

  $$\begin{eqnarray}&&{\bar{\nu }}_{\phi }=\displaystyle \frac{{\mu }^{2}\left({L}_{z}-2x\mu \right)}{2\pi \left(E+{{Ea}}^{2}{\mu }^{2}-2{ax}{\mu }^{3}\right)},\end{eqnarray} \tag{ C21a }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{r}=\displaystyle \frac{{\mu }^{3/2}\sqrt{1-3{\mu }^{2}{x}^{2}}\left(1-2\mu +{a}^{2}{\mu }^{2}\right)}{2\pi \left(E+{{Ea}}^{2}{\mu }^{2}-2{ax}{\mu }^{3}\right)},\end{eqnarray} \tag{ C21b }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{\theta }=\displaystyle \frac{{\mu }^{2}\left(1-2\mu +{a}^{2}{\mu }^{2}\right)\sqrt{\left({x}^{2}+{a}^{2}+2{aEx}\right)}}{2\pi \left(E+{{Ea}}^{2}{\mu }^{2}-2{ax}{\mu }^{3}\right)}.\end{eqnarray} \tag{ C21c }$$

  The expressions of E, L_z, and x for 00 orbits are given by (Bardeen et al. 1972)

  $$\begin{eqnarray}&&E=\displaystyle \frac{\left({r}_{c}^{2}-2{r}_{c}+a\sqrt{{r}_{c}}\right)}{{r}_{c}{\left({r}_{c}^{2}-3{r}_{c}+2a\sqrt{{r}_{c}}\right)}^{1/2}},\end{eqnarray} \tag{ C22a }$$

  $$\begin{eqnarray}&&{L}_{z}=\displaystyle \frac{\sqrt{{r}_{c}}\left({r}_{c}^{2}+{a}^{2}-2a\sqrt{{r}_{c}}\right)}{{r}_{c}{\left({r}_{c}^{2}-3{r}_{c}+2a\sqrt{{r}_{c}}\right)}^{1/2}},\end{eqnarray} \tag{ C22b }$$

  $$\begin{eqnarray}&&x=\displaystyle \frac{{r}_{c}\left({r}_{c}^{1/2}-a\right)}{{\left({r}_{c}^{2}-3{r}_{c}+2a\sqrt{{r}_{c}}\right)}^{1/2}},\end{eqnarray} \tag{ C22c }$$

  where r_c is the radius of the circular orbit. These expressions can also be obtained by substituting {e = 0, Q = 0, μ = 1/r_c} in the more general expressions given by Equation (5) of Rana & Mangalam (2019a). Finally, by substituting E, L_z, x, and μ = 1/r_c from Equation (C22) into Equation (C21), we recover the frequency formulae for 00 orbits:

  $$\begin{eqnarray}&&{\bar{\nu }}_{\phi }=\displaystyle \frac{1}{2\pi \left({r}_{c}^{3/2}+a\right)},\end{eqnarray} \tag{ C23a }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{r}={\bar{\nu }}_{\phi }{\left(1-\displaystyle \frac{6}{{r}_{c}}-\displaystyle \frac{3{a}^{2}}{{r}_{c}^{2}}+\displaystyle \frac{8a}{{r}_{c}^{3/2}}\right)}^{1/2},\end{eqnarray} \tag{ C23b }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{\theta }={\bar{\nu }}_{\phi }{\left(1+\displaystyle \frac{3{a}^{2}}{{r}_{c}^{2}}-\displaystyle \frac{4a}{{r}_{c}^{3/2}}\right)}^{1/2},\end{eqnarray} \tag{ C23c }$$

  as given by Equation (3).
• 2.  
  Reduction from Q0 orbits: We find that for circular orbits (Q = 0), the expressions of z_± (Equation 9(d) of Rana & Mangalam 2019a) reduce to

  $$\begin{eqnarray}&&{z}_{-}=0,\ \ \ \ \ {z}_{+}=\displaystyle \frac{\sqrt{{L}_{z}^{2}+{a}^{2}\left(1-{E}^{2}\right)}}{a\sqrt{1-{E}^{2}}}.\end{eqnarray} \tag{ C24 }$$

  The substitution of Equation (C24) in the frequency formulae of Q0 orbits, Equation (7), yields

  $$\begin{eqnarray}&&{\bar{\nu }}_{\phi }=\displaystyle \frac{\left(-2{L}_{z}{r}_{c}+{L}_{z}{r}_{c}^{2}+2{r}_{c}{aE}\right)}{2\pi \left[E\left({a}^{2}{r}_{c}^{2}+{r}_{c}^{4}+2{a}^{2}{r}_{c}\right)-2{L}_{z}{{ar}}_{c}\right]},\end{eqnarray} \tag{ C25a }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{r}=\displaystyle \frac{\sqrt{{r}_{c}^{4}\left(1-{E}^{2}\right)-2{x}^{2}{r}_{c}}{\rm{\Delta }}}{2\pi {r}_{c}\left[E\left({a}^{2}{r}_{c}^{2}+{r}_{c}^{4}+2{a}^{2}{r}_{c}\right)-2{L}_{z}{{ar}}_{c}\right]},\end{eqnarray} \tag{ C25b }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{\theta }=\displaystyle \frac{\sqrt{{L}_{z}^{2}+{a}^{2}\left(1-{E}^{2}\right)}{\rm{\Delta }}}{2\pi \left[E\left({a}^{2}{r}_{c}^{2}+{r}_{c}^{4}+2{a}^{2}{r}_{c}\right)-2{L}_{z}{{ar}}_{c}\right]}.\end{eqnarray} \tag{ C25c }$$

  Using the expressions of E, L_z, and x from Equation (C22), we find that

  $$\begin{eqnarray}&&\begin{array}{l}\left[E\left({a}^{2}{r}_{c}^{2}+{r}_{c}^{4}+2{a}^{2}{r}_{c}\right)-2{L}_{z}{{ar}}_{c}\right]=\displaystyle \frac{{r}_{c}^{3/2}{\rm{\Delta }}\left({r}_{c}^{3/2}+a\right)}{{\left({r}_{c}^{2}-3{r}_{c}+2{{ar}}_{c}^{1/2}\right)}^{1/2}},\end{array}\end{eqnarray} \tag{ C26a }$$

  $$\begin{eqnarray}&&\left(-2{L}_{z}{r}_{c}+{L}_{z}{r}_{c}^{2}+2{r}_{c}{aE}\right)=\displaystyle \frac{{r}_{c}^{3/2}{\rm{\Delta }}}{{\left({r}_{c}^{2}-3{r}_{c}+2{{ar}}_{c}^{1/2}\right)}^{1/2}},\end{eqnarray} \tag{ C26b }$$

  $$\begin{eqnarray}&&\begin{array}{l}\sqrt{{r}_{c}^{4}\left(1-{E}^{2}\right)-2{x}^{2}{r}_{c}}=\displaystyle \frac{{r}_{c}^{3/2}{\left({r}_{c}^{2}-6{r}_{c}-3{a}^{2}+8{{ar}}_{c}^{1/2}\right)}^{1/2}}{{\left({r}_{c}^{2}-3{r}_{c}+2{{ar}}_{c}^{1/2}\right)}^{1/2}},\end{array}\end{eqnarray} \tag{ C26c }$$

  $$\begin{eqnarray}&&\sqrt{{L}_{z}^{2}+{a}^{2}\left(1-{E}^{2}\right)}=\displaystyle \frac{\sqrt{{r}_{c}^{3}+3{a}^{2}{r}_{c}-4{{ar}}_{c}^{3/2}}}{{\left({r}_{c}^{2}-3{r}_{c}+2{{ar}}_{c}^{1/2}\right)}^{1/2}}.\end{eqnarray} \tag{ C26d }$$

  Finally, substituting these factors, given by Equation (C26), in Equation (C25), we recover the expressions for 00 orbits, which are given by

  $$\begin{eqnarray}&&{\bar{\nu }}_{\phi }=\displaystyle \frac{1}{2\pi \left({r}_{c}^{3/2}+a\right)},\end{eqnarray} \tag{ C27a }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{r}={\bar{\nu }}_{\phi }{\left(1-\displaystyle \frac{6}{{r}_{c}}-\displaystyle \frac{3{a}^{2}}{{r}_{c}^{2}}+\displaystyle \frac{8a}{{r}_{c}^{3/2}}\right)}^{1/2},\end{eqnarray} \tag{ C27b }$$

  $$\begin{eqnarray}&&{\bar{\nu }}_{\theta }={\bar{\nu }}_{\phi }{\left(1+\displaystyle \frac{3{a}^{2}}{{r}_{c}^{2}}-\displaystyle \frac{4a}{{r}_{c}^{3/2}}\right)}^{1/2},\end{eqnarray} \tag{ C27c }$$

  as given in Equation (3).

We summarize the history of each BHXRB below:

• 1.  
  M82 X-1: This is the brightest X-ray source in the M82 galaxy. This source is thought to harbor an intermediate-mass black hole because of its very high X-ray luminosity, average 2–10 keV luminosity ∼5 × 10⁴⁰ erg s⁻¹, and variability characteristics (Patruno et al. 2006; Casella et al. 2008; Pasham & Strohmayer 2013b), although other models claim that it might contain a black hole of mass ∼20 M_⊙ (Okajima et al. 2006). However, the discovery of twin-peak and stable QPOs at 3.32 ± 0.06 Hz and 5.07 ± 0.06 Hz in M82 X-1, which are nearly in 3:2 ratio, gave a shred of affirmative evidence that these QPOs are analogs of HFQPOs in stellar BHXRBs (Pasham et al. 2014). Following and extrapolating the inverse-mass scaling that holds for HFQPOs in stellar-mass BHXRBs (McClintock & Remillard 2006), it was found that the mass of the black hole in M82 X-1 could be 428 ± 105 M_⊙ (Pasham et al. 2014), making it an intermediate-mass black hole system.
• 2.  
  GROJ 1655-40: This is one among the few BHXRBs in the Milky Way galaxy whose BH mass is known with good precision through dynamical studies of infrared and optical observations during the quiescent state (Beer & Podsiadlowski 2002). GROJ 1655-40 is also one of the BHXRBs known to produce relativistic radio jets having a double-lobed radio structure (Mirabel & Rodríguez 1994). The first detection of two simultaneous HFQPOs near ∼450 and 300 Hz in GROJ 1655-40 was reported by Strohmayer (2001a). The detection of 300 Hz QPO was reported in BHXRB GROJ 1655-40 (Remillard et al. 1999b), and later the detection of a simultaneous 450 Hz QPO along with 300 Hz in the same observations was confirmed (Strohmayer 2001a). A systematic study of the LFQPOs and HFQPOs in 571 RXTE observations taken between the years 1996 and 2005 was carried out by Motta et al. (2014a), who detected three simultaneous QPOs (two HFQPOs and one LFQPO) at 441 ± 2 Hz, 298 ± 4 Hz, and 17.3 ± 0.1 Hz in one of these observations. Using these QPO frequencies, the mass, the spin of the black hole, and the radius of the equatorial circular orbit where these QPOs originated were estimated using Equations 3(a)–(c) assuming the RPM (Motta et al. 2014a).
• 3.  
  XTEJ 1550-564: This BHXRB was first detected by ASM/RXTE on 1998 September 7. Since then, it has undergone four X-ray outbursts between the years 1998 and 2002 as observed by RXTE, among which the 1998 September to 1999 May outburst was the most luminous one. XTEJ 1550-564 is also among the few BHXRBs that have shown HFQPOs; for example, QPOs with frequencies in the range 185–237 Hz were detected during the 1998–1999 outburst (Remillard et al. 1999a; Homan et al. 2001). After a quiescent period of a few months, XTEJ 1550-564 again underwent a short X-ray outburst in the period 2000 April to May following a fast rise and an exponential decay of the X-ray luminosity. The simultaneous occurrence of two HFQPOs at 268 ± 3 and 188 ± 3 Hz frequencies during the 2000 outburst was reported (Miller et al. 2001), indicating a resonance phenomenon. However, no LFQPOs were detected simultaneously with these two HFQPOs. A systematic study of all archival RXTE observations of XTEJ 1550-564 was carried out by Motta et al. (2014b), who reported the detection of an HFQPO at ∼183 Hz along with a simultaneous type C LFQPO at ∼13 Hz and type B LFQPO at ∼5 Hz, but no second peak of HFQPO was detected during this observation.
• 4.  
  4U 1630-47: This soft X-ray transient was discovered by the Uhuru satellite (Jones et al. 1976), which is known to have an inclination of ∼60°–75° (Kuulkers et al. 1998). This source is one among the few BHXRBs to show HFQPOs during its 1998 outburst in the frequency range ∼100–300 Hz, and also twin simultaneous HFQPOs with frequency ratio 1:4 (Klein-Wolt et al. 2004). It shows a regular X-ray outburst cycle after every ∼600–690 days (Jones et al. 1976; Priedhorsky 1986). The QPO frequencies in this system during the 1998 X-ray outburst were observed to increase during the rising phase, followed by a phase where the frequencies were found to be stable near ∼180 Hz, and then a decrease in QPO frequencies was observed during the decay of the outburst.
• 5.  
  GRS 1915+105: This BHXRB is known to be a very bright system during the whole RXTE period, showing its peculiar behavior and superluminal radio outflows (Mirabel & Rodríguez 1994). This is also the first BHXRB to show an HFQPO at a characteristic constant frequency of ∼67 Hz (Morgan et al. 1997) in the RXTE observations taken during 1996 April to May. Later, simultaneous ∼67 and ∼40 Hz QPOs were discovered in the RXTE observations taken during 1997 July and November (Strohmayer 2001b). A systematic study of all RXTE observations of GRS 1915+105 discovered 51 observations that showed detection of HFQPOs, out of which 48 observations showed the centroid frequency of QPOs in the range 63–71 Hz (Belloni & Altamirano 2013a). Another pair of simultaneous HFQPOs was also discovered at ∼34 and ∼68 Hz (Belloni & Altamirano 2013b).

Here, we describe a generic procedure that we have used to estimate errors in the orbital parameters. A flowchart of this method is provided in Figure 14.

• 1.  
  We assume that the QPO frequencies, ν₁, ν₂, and ν₃, are Gaussian distributed with mean values at ν₁₀, ν₂₀, and ν₃₀ (with ν₁₀ > ν₂₀ > ν₃₀), which are equal to the observed QPO centroid frequencies (see Table 5). For BHXRBs with two simultaneous QPOs, we only have ν₁ and ν₂. The joint probability density distribution of these frequencies will be given by

  $$\begin{eqnarray}&&P\left(\nu \right)=\displaystyle \prod _{i=1}^{l}{P}_{i}\left({\nu }_{i}\right),\end{eqnarray} \tag{ E28a }$$

  where l = 3 and l = 2 for BHXRBs with three and two simultaneous QPOs, respectively. Here, ${P}_{i}\left({\nu }_{i}\right)$ represents the Gaussian distribution of the ith QPO frequency, given by

  $$\begin{eqnarray}&&{P}_{i}\left({\nu }_{i}\right)=\displaystyle \frac{1}{\sqrt{2\pi {\sigma }_{i}^{2}}}\exp \left[-\displaystyle \frac{{\left({\nu }_{i}-{\nu }_{i0}\right)}^{2}}{2{\sigma }_{i}^{2}}\right].\end{eqnarray} \tag{ E28b }$$
• 2.  
  We find the Jacobian of the transformation from frequency to orbital parameter space using the formulae of fundamental frequencies, which are given by

  $$\begin{eqnarray}{{ \mathcal J }}_{l}=\displaystyle \frac{\partial {\nu }_{i}}{\partial {x}_{j}};\ {{ \mathcal J }}_{l}=\left\{\begin{array}{ll}{{ \mathcal J }}_{2}, & 2\ \mathrm{simultaneous}\ \mathrm{QPOs},\\ {{ \mathcal J }}_{3}, & 3\ \mathrm{simultaneous}\ \mathrm{QPOs},\end{array}\right.\end{eqnarray} \tag{ E29a }$$

  where {i, j} = 1 to l, x_j represent the orbital parameters, and ${ \mathcal J }$ is given by

  $$\begin{eqnarray}{{ \mathcal J }}_{3}=\left[\begin{array}{ccc}\tfrac{\partial {\nu }_{1}}{\partial {x}_{1}} & \tfrac{\partial {\nu }_{1}}{\partial {x}_{2}} & \tfrac{\partial {\nu }_{1}}{\partial {x}_{3}}\\ \tfrac{\partial {\nu }_{2}}{\partial {x}_{1}} & \tfrac{\partial {\nu }_{2}}{\partial {x}_{2}} & \tfrac{\partial {\nu }_{2}}{\partial {x}_{3}}\\ \tfrac{\partial {\nu }_{3}}{\partial {x}_{1}} & \tfrac{\partial {\nu }_{3}}{\partial {x}_{2}} & \tfrac{\partial {\nu }_{3}}{\partial {x}_{3}}\end{array}\right]\ \ \ \ \ \mathrm{and}\ \ \ \ \ {{ \mathcal J }}_{2}=\left[\begin{array}{cc}\tfrac{\partial {\nu }_{1}}{\partial {x}_{1}} & \tfrac{\partial {\nu }_{1}}{\partial {x}_{2}}\\ \tfrac{\partial {\nu }_{2}}{\partial {x}_{1}} & \tfrac{\partial {\nu }_{2}}{\partial {x}_{2}}\end{array}\right].\end{eqnarray} \tag{ E29b }$$

  For general eccentric trajectories ($Q\ne 0$), which are implemented for BHXRBs with three QPOs, we have {x₁, x₂, x₃} = {e, r_p, a}, whereas for equatorial eccentric trajectories (Q = 0), implemented for BHXRBs with two QPOs, we have {x₁, x₂} = {e, r_p}. Similarly, for the spherical orbit case, these parameters are {x₁, x₂, x₃} = {r_s, Q, a} or {x₁, x₂} = {r_s, Q}. The Jacobian is completely expressible in terms of the standard elliptic integrals and can be easily evaluated from Equation (E29) and using the frequency formulae Equations (4), (5), and (7), where ν₁ = ν_ϕ, ${\nu }_{2}=\left({\nu }_{\phi }-{\nu }_{r}\right)$, and ${\nu }_{3}=\left({\nu }_{\phi }-{\nu }_{\theta }\right)$ according to the RPM and GRPM. The analytic expressions for the elements of the Jacobian are too long to reproduce here, but they are used to make our computations faster.
• 3.  
  Next, we write the probability density distribution in the parameter space given by

  $$\begin{eqnarray}&&P\left([x]\right)=P\left(\nu \right)| {{ \mathcal J }}_{l}| ,\end{eqnarray} \tag{ E30 }$$

  where [x] represent the set of parameters {x₁, x₂, x₃} for l = 3 and {x₁, x₂} for l = 2, and {ν₁, ν₂, ν₃} or {ν₁, ν₂} are substituted in terms of parameters using our analytic formulae.We take Q = {0, 1, 2, 3, 4} for the general {e, Q} trajectory solutions that are implemented for the sources M82 X-1 and GROJ 1655-40. For each fixed value of Q, we find the corresponding probability density distribution in the parameter space using Equation (E30).
• 4.  
  We calculate the exact solutions for parameters by solving ν_ϕ = ν₁₀, ν_pp = ν₂₀, and ν_np = ν₃₀ using Equations 4(a)–(c) for nonequatorial eccentric trajectories, Equations 5(a)–(c) for equatorial eccentric, and Equations 7(a)–(c) for the spherical trajectories. We fix ${ \mathcal M }$ for l = 3, and both ${ \mathcal M }$ and a for l = 2 to the previous values; see Table 5. We find 1σ errors in the parameters by taking an appropriate parameter volume around the exact solution, and we generate sets of parameter combinations with resolution Δx_j in this volume. The chosen parameter range, exact solutions, and corresponding resolutions are summarized in Tables 6, 8, and 9. We then calculate the probability density using Equation (E30) for all of the generated parameter combinations and normalize the probability density by the normalization factor

  $$\begin{eqnarray}\begin{array}{rclcl}{ \mathcal N } & = & \displaystyle \frac{{\displaystyle \sum }_{k}P\left({\left[x\right]}_{k}\right){\rm{\Delta }}{V}_{k}}{V},{\rm{\Delta }}{V}_{k} & = & \displaystyle \prod _{j=1}^{l}{\rm{\Delta }}{x}_{j,k},\ V=\displaystyle \sum _{k}{\rm{\Delta }}{V}_{k},\end{array}\end{eqnarray} \tag{ E31a }$$

  where k varies from 1 to the number of total parameter combinations taken in the parameter volume, and ${\left[x\right]}_{k}$ is the kth combination of the parameters in the parameter volume. Hence, the normalized probability density is given by

  $$\begin{eqnarray}&&{ \mathcal P }\left([x]\right)=\displaystyle \frac{P\left([x]\right)}{{ \mathcal N }}.\end{eqnarray} \tag{ E31b }$$
• 5.  
  The allowed parameter combinations for the bound orbits are governed by the condition Equation (2). For a spherical orbit, we have e = 0. Hence, we ensure that the parameters (e, r_p, a, Q) for eccentric and (r_s, a, Q) for spherical trajectories follow the bound orbit condition. If any parameter combination does not obey the bound orbit condition, then ${ \mathcal P }\left([x]\right)$ is taken to be zero at that point in the parameter volume.
• 6.  
  Next, we integrate the normalized probability density, ${ \mathcal P }\left([x]\right)$, Equation E31(b), in two dimensions to obtain the profile in the remaining third dimension of the parameters for BHXRBs with three simultaneous QPOs, and similarly by integrating in one dimension for the two QPO cases, we obtain the profile in the other dimension. So we finally obtain the one-dimensional distributions ${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, and ${{ \mathcal P }}_{1}\left(a\right)$.
• 7.  
  Finally, we fit the normalized probability density profiles in each of the parameter dimensions to find the corresponding mean values, and quoted errors are obtained such that they contain a probability of 68.2% about the peak value of the probability density. The results of these fits are given in Tables 6, 8, and 9.
• 8.  
  For BHXRBs M82 X-1 and GROJ 1655-40, we find various orbital solutions showing varying {a, Q} values. As the spin of the black hole should be fixed, we choose the most probable value of a, and then we estimate the remaining parameters {e, r_p, Q}, their profiles {${{ \mathcal P }}_{1}\left(e\right)$, ${{ \mathcal P }}_{1}\left({r}_{p}\right)$, ${{ \mathcal P }}_{1}\left(Q\right)$}, and the corresponding errors using the same procedure given above in steps 1 to 6, where the orbital parameters are now given by {x₁, x₂, x₃} = {e, r_p, Q}.
• 9.  
  Although we have made accurate calculations described above, to obtain a rough and quick estimate of the errors, we may use the following procedure. Assuming that the probability density is Gaussian distributed independently in e, r_p, and a parameters, the normalized joint probability density distribution is given by

  $$\begin{eqnarray}&&\begin{array}{l}{ \mathcal P }\left(e,{r}_{p},a\right)=\displaystyle \frac{1}{{\left(2\pi \right)}^{3/2}{\sigma }_{e}{\sigma }_{{r}_{p}}{\sigma }_{a}}\ \exp \left\{-\displaystyle \frac{1}{2}\left[{\left(\displaystyle \frac{e-{e}_{0}}{{\sigma }_{e}}\right)}^{2}+{\left(\displaystyle \frac{{r}_{p}-{r}_{p0}}{{\sigma }_{{r}_{p}}}\right)}^{2}+{\left(\displaystyle \frac{a-{a}_{0}}{{\sigma }_{a}}\right)}^{2}\right]\right\},\end{array}\end{eqnarray} \tag{ E32a }$$

  where the distribution is centered at the exact solution (e₀, r_p0, a₀), and σ_e, ${\sigma }_{{r}_{p}}$, and σ_a are the corresponding 1σ errors, derived using the method described above. The total probability contained in a volume ${ \mathcal V }$ in (e, r_p, a) space is given by

  $$\begin{eqnarray}\begin{array}{c}\begin{array}{c}\begin{array}{ccc}p & = & \displaystyle \frac{1}{{\left(2\pi \right)}^{3/2}{\sigma }_{e}{\sigma }_{{r}_{p}}{\sigma }_{a}}\displaystyle \int \displaystyle \int {\displaystyle \int }_{{ \mathcal V }}\exp \left\{-\displaystyle \frac{1}{2}\left[{\left(\displaystyle \frac{e-{e}_{0}}{{\sigma }_{e}}\right)}^{2}\right.\right.+\,\left.\left.{\left(\displaystyle \frac{{r}_{p}-{r}_{p0}}{{\sigma }_{{r}_{p}}}\right)}^{2}+{\left(\displaystyle \frac{a-{a}_{0}}{{\sigma }_{a}}\right)}^{2}\right]\right\}{\rm{d}}e\cdot {\rm{d}}{r}_{p}\cdot {\rm{d}}a\end{array}\end{array}\end{array}\end{eqnarray} \tag{ E32b }$$

  so that the total probability p inside an ellipsoid in (e, r_p, a) space specified by

  $$\begin{eqnarray}&&\left[{\left(\displaystyle \frac{e-{e}_{0}}{{\sigma }_{e}}\right)}^{2}+{\left(\displaystyle \frac{{r}_{p}-{r}_{p0}}{{\sigma }_{{r}_{p}}}\right)}^{2}+{\left(\displaystyle \frac{a-{a}_{0}}{{\sigma }_{a}}\right)}^{2}\right]={s}_{3}^{2}\end{eqnarray} \tag{ E32c }$$

  is given by

  $$\begin{eqnarray}&&p=\sqrt{\displaystyle \frac{2}{\pi }}{\int }_{0}^{{s}_{3}}\exp \left(\displaystyle \frac{-{s}^{2}}{2}\right){s}^{2}\,{\rm{d}}s=\displaystyle \frac{2}{\sqrt{\pi }}\gamma \left(\displaystyle \frac{3}{2},\displaystyle \frac{{s}_{3}^{2}}{2}\right),\end{eqnarray} \tag{ E32d }$$

  where $\gamma \left(\displaystyle \frac{3}{2},\displaystyle \frac{{s}_{3}^{2}}{2}\right)$ is the incomplete gamma function.Similarly, for two QPO cases, the joint probability density distribution can be written as

  $$\begin{eqnarray}\begin{array}{rcl}{ \mathcal P }\left(e,{r}_{p}\right) & = & \displaystyle \frac{1}{2\pi {\sigma }_{e}{\sigma }_{{r}_{p}}}\exp \left\{-\displaystyle \frac{1}{2}\left[{\left(\displaystyle \frac{e-{e}_{0}}{{\sigma }_{e}}\right)}^{2}\right.\right.+\ \left.\left.{\left(\displaystyle \frac{{r}_{p}-{r}_{p0}}{{\sigma }_{{r}_{p}}}\right)}^{2}\right]\right\}.\end{array}\end{eqnarray} \tag{ E33a }$$

  The total probability contained in a surface ${ \mathcal S }$ in (e, r_p) space is given by

  $$\begin{eqnarray}\begin{array}{c}\begin{array}{rcl}p & = & \displaystyle \frac{1}{2\pi {\sigma }_{e}{\sigma }_{{r}_{p}}}\displaystyle \int {\displaystyle \int }_{{ \mathcal S }}\,\exp \left\{-\displaystyle \frac{1}{2}\left[{\left(\displaystyle \frac{e-{e}_{0}}{{\sigma }_{e}}\right)}^{2}+{\left(\displaystyle \frac{{r}_{p}-{r}_{p0}}{{\sigma }_{{r}_{p}}}\right)}^{2}\right]\right\}{\rm{d}}e\cdot {{\rm{d}}r}_{p}.\end{array}\end{array}\end{eqnarray} \tag{ E33b }$$

  The total probability inside an ellipse, specified by

  $$\begin{eqnarray}&&\left[{\left(\displaystyle \frac{e-{e}_{0}}{{\sigma }_{e}}\right)}^{2}+{\left(\displaystyle \frac{{r}_{p}-{r}_{p0}}{{\sigma }_{{r}_{p}}}\right)}^{2}\right]={s}_{2}^{2},\end{eqnarray} \tag{ E33c }$$

  is given by

  $$\begin{eqnarray}&&p={\int }_{0}^{{s}_{2}}\exp \left(\displaystyle \frac{-{s}^{2}}{2}\right)s\,{\rm{d}}s=1-\exp \left(\displaystyle \frac{-{s}_{2}^{2}}{2}\right).\end{eqnarray} \tag{ E33d }$$

  For a given p, we can calculate ${s}_{3}^{2}$ and ${s}_{2}^{2}$, and hence evaluate the error ellipsoid corresponding to p. ${s}_{3}^{2}$ and ${s}_{2}^{2}$ are shown as a functions of p in Figure 30. This can be used to get rough estimates of the error distribution of the parameters. However, we calculate them exactly in Section 3.2.

Zoom In Zoom Out Reset image size