# The prototype X-ray binary GX 339–4: using TeV $\gamma$-rays to assess LMXBs
as Galactic cosmic ray accelerators

D. Kantzas1,2, S. Markoff1,2, M. Lucchini3,1, C. Ceccobello4, V. Grinberg5,6,
R. M. T. Connors7 & P. Uttley1
1Anton Pannekoek Institute for Astronomy (API), University of Amsterdam,
Science Park 904, 1098 XH Amsterdam, the Netherlands
2GRavitational AstroParticle Physics Amsterdam (GRAPPA), University of
Amsterdam, Science Park 904, 1098 XH Amsterdam, the Netherlands
3MIT Kavli Institute for Astrophysics and Space Research, Massachusetts
Institute of Technology, Cambridge, MA 02139, USA
4Department of Space, Earth and Environment, Chalmers University of
Technology, Onsala Space Observatory, 439 92 Onsala, Sweden
5Institute for Astronomy und Astrophysics, University of Tübingen, Sand 1,
72076 Tübingen, Germany
6European Space Agency (ESA), European Space Research and Technology Centre
(ESTEC), Keplerlaan 1, 2201 AZ Noordwijk, the Netherlands
7Cahill Center for Astronomy and Astrophysics, California Institute of
Technology, 1200 California Boulevard, Pasadena, CA 91125, USA E-mail:
<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Since the discovery of cosmic rays (CRs) over a century ago, their origin
remains an open question. Galactic CRs with energy up to the knee ($10^{15}$
eV) are considered to originate from supernova remnants, but this scenario has
recently been questioned due to lack of TeV $\gamma$-ray counterparts in many
cases. Extragalactic CRs on the other hand, are thought to be associated with
accelerated particles in the relativistic jets launched by supermassive
accreting black holes at the center of galaxies. Scaled down versions of such
jets have been detected in X-ray binaries hosting a stellar black hole
(BHXBs). In this work, we investigate the possibility that the smaller-scale
jets in transient outbursts of low-mass BHXBs could be sources of Galactic
CRs. To better test this scenario, we model the entire electromagnetic
spectrum of such sources focusing on the potential TeV regime, using the
’canonical’ low-mass BHXB GX 339–4as a benchmark. Taking into account both the
leptonic radiative processes and the $\gamma$-rays produced via neutral pion
decay from inelastic hadronic interactions, we predict the GeV and TeV
$\gamma$-ray spectrum of GX 339–4 using lower-frequency emission as
constraints. Based on this test-case of GX 339–4 we investigate whether other,
nearby low-mass BHXBs could be detected by the next-generation very-high-
energy $\gamma$-ray facility the Cherenkov Telescope Array, which would
establish them as additional and numerous potential sources of CRs in the
Galaxy.

###### keywords:

acceleration of particles – radiation mechanisms: non-thermal – X-rays:
individual: GX 339–4

††pubyear: 2021††pagerange: The prototype X-ray binary GX 339–4: using TeV
$\gamma$-rays to assess LMXBs as Galactic cosmic ray accelerators–A

##  1 Introduction

Accreting supermassive black holes located at the centres of galaxies are the
most powerful engines in the Universe, and some of the most interesting
laboratories to investigate the physics of extreme gravity. Of particular
importance are those active galactic nuclei (AGN) that exhibit relativistic
and collimated jets. The underlying physics that unites the accretion of black
holes with the large scale jets is still an unanswered problem. These
relativistic jets are considered powerful enough to accelerate particles to
very high energy, making them likely a source of extragalactic cosmic rays
(CRs) that reach energies of at least $10^{19}\,$eV (Hillas, 1984; Abbasi et
al., 2018; Perrone, 2020).

CRs are elementary particles and/or atoms of extraterrestrial origin. The
resulting CR spectrum covers ten orders of magnitudes in particle energy and
shows two very well known characteristic spectral features where the slope
changes. The first one is the ‘knee’ that is located around $10^{15}\,$eV, and
the second feature is the ‘ankle’ that is located around $10^{17}\,$eV.
Current models assume that CRs up to the knee are produced within the Milky
Way, while CRs from above the ankle are of extragalactic origin (Hillas, 1984;
Drury, 2012; Blasi, 2013). Supernova remnants have long been considered the
dominant source of Galactic CRs based on their size and measured magnetic
fields (Hillas, 1984; Völk et al., 2003; Vink, 2012; Ackermann et al., 2013),
but due to the lack of TeV $\gamma$-ray counterparts the debate is still open
(Aharonian et al., 2019). Given the ability of AGN jets to accelerate cosmic
rays, another promising alternative source could be the Galactic jets launched
in X-ray binaries comprised of a stellar accreting black hole and a companion
star (BHXBs; Mirabel & Rodriguez 1994; Fender 2001; McClintock et al. 2006).
Such Galactic jets share the physical properties of AGN jets but on much
smaller scales (Heinz & Sunyaev, 2002; Romero et al., 2003; Romero & Orellana,
2005; Fender et al., 2005; Romero & Vila, 2008; Vila & Romero, 2010; Vila et
al., 2012; Pepe et al., 2015; Cooper et al., 2020; Kantzas et al., 2020).

The presence of jets in low mass BHXBs is transient, and is tightly connected
to the properties of the accretion flow. In the so-called hard state, BHXBs
display a flat or inverted radio-to-IR spectrum associated with jet
synchrotron emission analogous to AGN jets (Blandford & Königl, 1979;
Hjellming & Johnston, 1988; Falcke & Biermann, 1995; Markoff et al., 2001;
Fender et al., 2006; Corbel et al., 2003, 2012). BHXBs transit from quiescent
to hard and soft states within ‘human-like’ timescales, hence we can observe
the jet launching and jet quenching in real-time (see e.g., Russell et al.,
2020). The dynamical timescales are roughly proportional to the mass of the
black hole, so it would take typically millions of times longer to detect
similar state transitions in AGN.

Accelerated particles in AGN jets are the source of the non-thermal radiation
detected over the entire electromagnetic spectrum, from radio to TeV
$\gamma$-rays (see e.g. Tavecchio et al., 1998; Celotti et al., 2001;
Aharonian, 2004; Georganopoulos et al., 2006; Marscher et al., 2008;
Ghisellini et al., 2009). However, the exact radiative mechanism has been
under debate for a long time because it is tightly connected to the jet
composition and the exact particle acceleration mechanism, which remain
debated. Two scenarios are generally considered depending on the jet launching
mechanism. First, a purely leptonic jet powered by the black hole spin
(Blandford & Znajek, 1977) may accelerate electrons/positrons that are
responsible for the entire multi-wavelength spectrum (Maraschi et al., 1992;
Dermer & Schlickeiser, 1993; Levinson & Blandford, 1995; Blandford & Levinson,
1995; Marcowith et al., 1995; Böttcher & Schlickeiser, 1997; Georganopoulos et
al., 2002; Ghisellini et al., 2010). Second, a lepto-hadronic jet powered by
the accretion disc (Blandford & Payne, 1982) (or which starts out leptonic and
entrains hadronic mass) may accelerate leptons and baryons that contribute in
different energy bands via different mechanisms (Mannheim, 1993; Rachen &
Biermann, 1993; Mücke et al., 2003; Böttcher et al., 2013; Liodakis &
Petropoulou, 2020).

Recent GeV observations of the high-mass BHXBs Cygnus X–3 (Tavani et al.,
2009) and Cygnus X–1 (Tavani et al., 2009; Malyshev et al., 2013; Zanin et
al., 2016), and TeV observations of SS 433 (Abeysekara et al., 2018) suggest
that some Galactic jets can accelerate particles to high energy. However, it
is not known whether all BHXBs, especially the more abundant population of
low-mass BHXBs can routinely produce $\gamma$-rays. Until now, only the high-
mass BHXBs that are characterised by the presence of a strong stellar wind
that interacts with the jet have been detected in the GeV and TeV bands (see
e.g., Bodaghee et al., 2013).It is thus important to investigate whether the
far more populous low-mass BHXBs can also produce $\gamma$-rays. In this paper
we approach this question by studying the ‘canonical’ low-mass BHXB source GX
339–4, extending our previous work on the ‘canonical’ high-mass BHXB Cygnus
X–1 (Kantzas et al., 2020). Similar to AGN jets, the emitting mechanism
responsible for any $\gamma$-rays remains unclear, with both leptonic and
hadronic processes considered feasible. We are also interested in exploring
how the different composition scenarios may affect the jet dynamics and the
interpretation of the jet properties.

In this work, we employ a multi-zone jet model to study the hadronic
interactions within the jets, as well as the effect on the dynamics and the
electromagnetic signature of low-mass BHXB jets. We examine the bright
outburst of GX 339–4 in 2010 to model the radio-to-X-ray spectrum with the
goal of predicting the TeV radiation originating in the jets. Using the case
of GX 339–4 as a model, we assess the likelihood of other, closer low-mass
BHXBs to be potential sources for the next generation $\gamma$-ray facilities,
particularly the Cherenkov Telescope Array (CTA). Such TeV emission may be the
signature of efficient CR acceleration inside the BHXB jets, and hence the
entire Galactic population of BHXB jets may contribute to the Galactic CR
spectrum.

In Section 2 we discuss the physical properties of GX 339–4 and its spectral
behaviour. In Section 3 we describe the model we use to study the spectrum of
GX 339–4. We present our results in Section 4, discuss their implication in
Section 5 and come to our final conclusions in Section 6.

##  2 GX 339–4

GX 339–4 is a ’canonical’ low-mass BHXB discovered in 1973 (Markert et al.,
1973). It undergoes outbursts every two-to-three years that last from a few
weeks to months (Belloni et al., 1999; Corbel & Fender, 2002; Corbel et al.,
2003; Zdziarski et al., 2004; Homan et al., 2005; Belloni et al., 2006; Motta
et al., 2009; Corbel et al., 2012). During outbursts, GX 339–4 rises out of
quiescence and launches compact jets that contribute to the radio-to-optical
spectrum as the source continues into the hard state (Fender, 2001; Fender et
al., 2004; Corbel et al., 2000, 2003, 2012; Corbel & Fender, 2002; Casella et
al., 2010; Homan et al., 2005; Gandhi et al., 2011). Such consistent,
repetitive behavior along with extensive and often simultaneous
multiwavelength monitoring makes GX 339–4 a perfect target to better
understand the properties of relativistic jets.

Although GX 339–4 is a well-studied source, its physical parameters are not
well constrained because of the weakness of its companion star. Based on
optical photometry the orbital period is estimated to be between 14.8 hours
and 16.8 hours (Callanan et al., 1992; Cowley et al., 2002, respectively). The
inclination angle is still unknown but is constrained to $<60$ degrees because
of the lack of eclipsing (Cowley et al., 2002), and the lack of a detection of
the companion star means the mass of the black hole is also uncertain. Various
current estimates put the mass (in M⊙) between $4-16$ (Shidatsu et al., 2011),
5.8 $\pm$ 0.8 for an orbital period of 1.75 days (Hynes et al., 2003), $>7$
(Muñoz Darias et al., 2008) or 9.8 for a mass function of 1.91 $\pm$ 0.08 M⊙
(Heida et al., 2017). We adopt the most recent value of
$M_{\rm{bh}}=9.8\,\rm{M_{\odot}}$ of Heida et al. (2017). Hynes et al. (2004)
set the distance of GX 339–4 higher than 6 kpc, and Zdziarski et al. (2004)
derived a value of 8 kpc while Parker et al. (2016) found a distance of $8\pm
0.9\,$kpc, which is the distance we adopt here.

###  2.1 Observational constraints in the hard state

GX 339–4 has been detected in the optical bands, but the origin of this
emission is still not clear. Tetarenko et al. (2020) recently studied its
multiwavelength emission and concluded that the optical emission in bright
outbursts like the one of 2010 cannot originate exclusively from irradiation
of the accretion disc, because unreasonable amounts of energy would be
required. Thermal synchrotron emission from the base of the jets could then be
considered a good candidate for the optical emission. On the other hand, GX
339–4 shows a flat spectrum in the radio with a spectral break in the IR band
that corresponds to the transition of optically thick to optically thin
synchrotron emission (Corbel & Fender, 2002; Gandhi et al., 2011).
Extrapolating the optically thin IR emission to the X-ray band, significantly
underpredicts the optical flux (Maitra et al., 2009; Gandhi et al., 2011;
Tetarenko et al., 2019). Hence, if the optical emission originates in the
jets, it must come from a different region compared to the IR (Markoff et al.,
2003; Corbel et al., 2013).

Reflection features, including a broad iron emission line, are also evident in
the X-ray spectrum of GX 339–4 (Nowak et al., 2002; García et al., 2015; Fürst
et al., 2015; Parker et al., 2016; García et al., 2019; Dziełak et al., 2019).
A jet synchrotron component that is beamed perpendicularly away from the
accretion disc is very unlikely to produce significant relativistic reflection
(Markoff & Nowak, 2004; Reig & Kylafis, 2021). Furthermore, Uttley et al.
(2011) studied the energy-dependent time lags and found that the instabilities
in the accretion disc may be responsible for driving the continuum variability
on short and longer-than-second timescales. The large time-lags are due to the
travel-time between the illuminating region and the disc where the X-rays are
reprocessed, and can be only tens of gravitational radii at most. That
indicates that the X-ray continuum should be governed by a single component,
and a thermal corona close to the black hole could sufficiently explain it
(but also see Mahmoud et al., 2019, for a two-component corona).

Based on these results, we approach the modelling assuming the most
conservative case for the jet power: that the radio through IR up to the break
is self-absorbed synchrotron from the extended jets, the optical emission is
synchrotron emission from thermal particles at the base of the jets, and that
the X-ray reflecting power-law is from a separate coronal region.

###  2.2 Observational data

In this work we use archival quasi-simultaneous data to model the
multiwavelength spectrum of GX 339–4 from radio to X-rays during the hard
state of the 2010 outburst. We use the radio data obtained by the Australian
Telescope Compact Array (ATCA) on MJD 55263 (Corbel et al., 2012), IR data
obtained by the Wide-field Infrared Survey Explorer (WISE) on MJD 55266
(Gandhi et al., 2011), optical data obtained by the Small & Moderate Aperture
Research Telescope System (SMARTS) on MJD 55263, and X-rays from the Neil
Gehrels Swift Observatory/X-ray Telescope (Swift/XRT) on MJD 55262 (Corbel et
al., 2012) and Rossi X-ray Timing Explorer/Proportional Counter Array
(RXTE/PCA) on MJD 55263 (Corbel et al., 2012). We use the 0.5–4.0 keV XRT and
the 3–45 keV PCA X-ray data. The IR data are not simultaneous and were
obtained 3 days later, but we use them because they show a spectral break
crucial for our interpretation of the whole spectrum (see below). There was no
significant variability in this time, hence this is a decent assumption to
combine these data (Corbel et al., 2012; Corbel et al., 2013; Connors et al.,
2019). We also use the upper-limits in the GeV band set by the Fermi/Large
Area Telescope (LAT) $\gamma$-ray telescope during the 2010 outburst to
further constrain the highest energy regime of the spectrum (Bodaghee et al.,
2013). We provide the energy/frequency ranges and the corresponding flux
density of all the data we use in Table 1.

##  3 Modelling

In this section we briefly discuss our model, focusing on the interpretation
of the free parameters we fit for. A more detailed description of the model
can be found in Kantzas et al. (2020) and in Lucchini et al. (2022).

###  3.1 Jet properties

We assume that two compact jets are launched by the accreting black hole with
jet base radius $R_{0}$. The power injected into the jets in the comoving
frame $L_{\rm{jet}}$ defines the number density of the cold (non-relativistic)
protons in the plasma at the base of the jets as,

$n_{0}=\frac{L_{\rm{jet}}}{2\beta_{0,s}\Gamma_{0,s}c\,\pi
R_{0}^{2}{\rm{(m_{p}c^{2}+\langle\gamma_{e}\rangle m_{e}c^{2}(1+1/\beta))}}},$
(1)

where $\beta_{0,s}\Gamma_{0,s}c$ is the comoving velocity of the plasma in the
jet base assumed to be equal to the speed of sound in a relativistic fluid
(Falcke & Biermann, 1995; Markoff et al., 2008; Crumley et al., 2017; Lucchini
et al., 2022), $\beta={U_{\rm{e}}}/{U_{\rm{B}}}$ is the plasma beta where
$U_{\rm e}$ is the energy density of the electrons and $U_{\rm B}$ is the
magnetic field energy density.For simplicity, we assume equal number density
of electrons and protons, but we discuss the implication of this assumption in
Section 5. We further assume that the electron population at the jet base is
injected in a thermal Maxwell-Jüttner (MJ) distribution with a peak-energy of
$2.23\,k_{\rm{B}}T_{\rm{e}}$.

We vary the plasma beta at the jet base to define the strength of the magnetic
field, which scales inversely with distance along the jet $z$. Assuming the
electron enthalpy is not significant, we define the magnetisation of the jet
as

$\sigma=\frac{B_{0}^{2}}{4\pi n_{0}{\rm{m_{p}c^{2}}}},$ (2)

where $B_{0}$ is the strength of the magnetic field at the jet base. We do not
consider any particular magnetic field configuration (toroidal or poloidal)
but merely describe the magnetic field by its total strength $B$.

###  3.2 Particle acceleration

At some distance $z_{\rm{diss}}$ along the jet axis, energy is dissipated into
accelerating a fraction of the thermal particles into a non-thermal power-law.
We assume that the accelerated particles carry a fixed fraction of the jet
power, and in particular, we conservatively fix the power of the non-thermal
leptons to be 0.02$L_{\rm jet}$ and of the protons to be 0.05$L_{\rm jet}$.

We allow $z_{\rm{diss}}$ to vary as a fitted parameter and, for the case of
the leptonic populations, we assume constant re-acceleration along the jet,
but we constrain the proton acceleration to occur only between $z_{\rm diss}$
and $10\,z_{\rm diss}$ in order to limit the required power. Because the most
compact part of the jet produces the non-thermal particles, this dissipation
region also corresponds to the region where the synchrotron radiation breaks
from flat/inverted due to self-absorption, to steep/optically thin. After
predictions by Markoff et al. (2001), Corbel & Fender (2002) confirmed that
this break typically falls in the NIR band during hard states, and we chose
the epoch here because of high-quality observations by Gandhi et al. (2011)
that could pinpoint the synchrotron break frequency to be
$\rm{4.6^{+3.5}_{-2.0}\times 10^{13}\,Hz}$. To match this frequency, we fix
the particle acceleration region at 2600 $r_{\rm{g}}$ from the black hole (see
also Connors et al., 2019).

The accelerated particles follow a power-law in energy of the form

${\rm{d}}n\left(E\right)\propto
E^{-p}\times\exp{\left(-E/E_{\rm{max}}\right)}.$ (3)

In principle, the power-law index $p$ depends on the acceleration mechanism
and may differ between electrons and protons, but we choose to use the same
for both populations for simplicity.

In equation 3, $E_{{\rm{max}}}$ is the maximum particle energy constrained by
energy losses and/or escape. In this work, the maximum electron energy is
limited by synchrotron losses and the maximum proton energy is limited by the
lateral escape from the jet region. The maximum attainable energy is self-
consistently calculated along the jet by equating the characteristic
timescales of the losses to the acceleration timescale. The characteristic
acceleration timescale $t_{\rm{acc}}=4E/(3f_{\rm{sc}}\rm{ec}B)$ depends on the
acceleration efficiency parameter $f_{\rm{sc}}$ that we take to be close to
maximum, namely $f_{\rm{sc}}=0.01$ (Jokipii, 1987; Aharonian, 2004). We plot
the characteristic timescales versus the particle kinetic energy for the
population of the accelerated protons in appendix A.

The fractional number of accelerated particles with respect to the total
number of particles $f_{\rm{nth}}$ depends on the acceleration mechanism as
well. This number may not be constant along the jet. We parametrize the
density of the accelerated particles following:

$n_{\rm{nth}}=n_{\rm{th}}\,f_{\rm{nth}}\,\left(\frac{\log_{10}(z_{\rm{diss}})}{\log_{10}(z)}\right)^{f_{\rm{pl}}},$
(4)

where $f_{\rm{pl}}>0$ is a free parameter accounting for our ignorance about
the exact nature of the dissipation. ($n_{\rm{nth}}$) $n_{\rm{th}}$ is the
number density of the (non-)thermal particles. The physical motivation behind
such an assumption is the fact that it leads to the characteristic inverted
spectrum between radio and optical wavelengths detected in BHXBs (see
discussion in Lucchini et al. 2021).

The minimum energy of the accelerated particles depends on the injected
distributions in the base. We assume that the minimum energy for accelerated
protons is the rest mass energy ($\rm{m_{p}c^{2}}$). This choice is intended
purely to limit the number of free parameters; we discuss its implication
below. We take the peak of the MJ distribution $2.23k_{\rm{B}}T_{\rm{e}}$ to
be the minimum energy of the accelerated electrons. We further define a
heating parameter $f_{\rm{heat}}$

$E_{\rm{e,min}}=2.23f_{\rm{heat}}k_{\rm{B}}T_{\rm{e}}.$ (5)

The physical motivation behind this assumption is that along with the electron
acceleration, some extra heating has been reported by numerical simulations
(Sironi & Spitkovsky, 2009; Gedalin et al., 2012; Plotnikov et al., 2013;
Sironi et al., 2013; Sironi & Spitkovsky, 2014; Melzani et al., 2014; Crumley
et al., 2017). The value of this parameter is not well constrained, but we set
it to be $f_{\rm heat}$< 10 (Sironi & Spitkovsky, 2009, 2011; Crumley et al.,
2019).

###  3.3 Radiative Processes

####  3.3.1 Leptonic processes

Following Kantzas et al. (2020), the leptonic radiative processes we take into
account are cyclo-synchrotron radiation and inverse Compton scattering (ICS),
where the cyclo-synchrotron photons are further upscattered via the
synchrotron-self Compton mechanism (SSC) along the jets. Further photon
targets for the ICS are the disc photons. We also take into account a precise
treatment of pair production due to photon annihilation and pair annihilation
to electron-positron pairs (Coppi & Blandford, 1990; Böttcher & Schlickeiser,
1997).

####  3.3.2 Hadronic Processes

Accelerated protons interact with the bulk cold protons of the jet and, via
proton-proton (pp) interactions, lead to pion production. Neutral pions decay
into $\gamma$-rays and charged pions into secondary electrons and neutrinos
via the muon decay channel (Mannheim & Schlickeiser, 1994). Photomeson
interactions between the accelerated protons and target photons (p$\gamma$)
lead to similar distributions of secondary particles. The target photons we
consider here are: the thermal radiation of the accretion disc and the non-
thermal radiation originating in the compact jet. Finally, we also account for
photopair interactions that lead to the formation of pairs, after the
inelastic collision between protons and photons. We use the semi-analytical
formalism of Kelner et al. (2006) and Kelner & Aharonian (2008) for pp and
p$\gamma$ interactions, respectively. For the full description of the
treatment of the cascades, see Kantzas et al. (2020). For the case of GX
339–4, no photon field is significant enough to attenuate the GeV and TeV
emission (also see the discussion below).

###  3.4 Accretion disc and thermal corona

We assume a standard geometrically thin, optically thick accretion disc
truncated at some innermost radius $R_{\rm{in}}$ with temperature
$T_{\rm{in}}$ (Shakura & Sunyaev, 1973; Frank et al., 2002). We describe the
disc luminosity $L_{d}$ in terms of Eddington luminosity
$L_{\rm{Edd}}=4\pi{\rm{G}}M_{\rm{bh}}{\rm{m_{p}c}}/\sigma_{\rm{T}}$. We
further assume the existence of a hot electron plasma of temperature
$T_{\rm{cor}}$, in a spherical region centered on the black hole, normalized
by a radius $R_{\rm{cor}}$, and of optical depth $\tau_{\rm{cor}}=n_{\rm
e}R_{\rm cor}\sigma_{\rm T}$. These hot electrons upscatter the disc photons
to higher energies. We require the existence of such a plasma to be able to
model both the X-ray spectrum and properly account for the measured hard
timing lags as mentioned in Section 1 (and see e.g., Connors et al., 2019, and
discussion below).

Observatory | log Frequency (Hz) | log Energy (eV) | Flux Density (mJya) | Reference
---|---|---|---|---
ATCA | $\begin{array}[]{cc}9.74&9.94\end{array}$ | $\begin{array}[]{cc}-4.64&-4.44\end{array}$ | $\begin{array}[]{cc}10.2\pm 0.1&11.3\pm 0.1\end{array}$ | Corbel et al. 2012
WISE | $\begin{array}[]{cc}13.13&13.41\\\ 13.81&13.95\end{array}$ | $\begin{array}[]{cc}-1.25&-0.97\\\ -0.57&-0.43\end{array}$ | $\begin{array}[]{cc}87\pm 8&80\pm 7\\\ 64\pm 5&55\pm 4\end{array}$ | Gandhi et al. 2011
SMARTS | $\begin{array}[]{cc}14.25&14.40\\\ 14.57&14.73\end{array}$ | $\begin{array}[]{cc}-0.13&0.01\\\ 0.18&0.35\end{array}$ | $\begin{array}[]{cc}47\pm 5&50\pm 5\\\ 54\pm 5&92\pm 29\end{array}$ | Buxton et al. 2012
SWIFT/RXT | 17.08–18.0 | 2.7–3.7 |  $0.2$ at $3\,$keV | Corbel et al. 2012
RXTE/PCA | 17.9–18.9 | 3.5–4.5 |  $0.2$ at $3\,$keV | Corbel et al. 2012
Table 1:  The observational multiwavelength data we use in this work.
amJy$\,=10^{-26}\,\rm{erg\,cm^{-2}\,s^{-1}\,Hz^{-1}}$ parameter | value | description
---|---|---
$M_{\rm{BH}}\,\left(\rm{M_{\odot}}\right)$ | 9.8 | mass of the black hole†
$\theta_{\rm{incl}}$ | 40∘ | inclination angle†
$D\,\rm{\left(kpc\right)}$ | 8 | distance of the source⋆
$h=z_{0}/R_{0}$ | 2 | initial jet height to radius ratio
$z_{\rm{diss}}\,\left(r_{\rm g}\right)$ | $2600$ | particle acceleration region∗
$z_{\rm{diss,max}}$ | $10z_{\rm diss}$ | maximum proton acceleration region
$z_{\rm{max}}\,\left(r_{\rm g}\right)$ | $10^{8}$ | maximum jet height
$f_{\rm sc}$ | 0.01 | particle acceleration efficiency parameter
$P_{\rm e}/L_{\rm jet}$ | 0.02 | power of non-thermal electrons
$P_{\rm p}/L_{\rm jet}$ | 0.05 | power of non-thermal protons
$R_{\rm{in,disc}}$ | $\rm{R_{0}}$ | disc innermost radius (see Table 3)
$R_{\rm{out,disc}}\,\left(r_{\rm g}\right)$ | $10^{5}$ | disc outermost radius
$\rm{N_{H}\,(10^{22}\,cm^{-2})}$ | $0.6$ | absorption coefficient‡
refl | 0.29 | reflection fraction⋄
Table 2: The fixed parameters of our models, see text for further discussion.
†Heida et al. (2017), ⋆Parker et al. (2016), ∗Gandhi et al. (2011), ‡García et
al. (2019), ⋄Magdziarz & Zdziarski (1995)

##  4 Results

In this section, we present the results for the best fits of our model to the
multiwavelength spectrum of GX 339–4. We explore three different model
scenarios: one purely leptonic, and two lepto-hadronic models. For the purely
leptonic model, we assume that the non-thermal electrons follow a power-law
with $p=2.2$ (Corbel & Fender, 2002; Gandhi et al., 2011). For the two
hadronic models, we explore both a soft ($p=2.2$) and a hard ($p=1.7$)
particle power-law, respectively. For all models we fix some common parameters
as shown in Table 2. We choose the ratio between the height of the jet base
and its radius to be constant and equal to 2 (Maitra et al., 2011; Crumley et
al., 2017). The maximum height of the jet is fixed at a large enough value, so
it does not influence the spectrum in the radio band via the self-absorption
cutoff, and we choose the maximum reasonable particle acceleration efficiency
parameter $f_{\rm{sc}}=0.1$, which results in maximum proton energies of the
order of tens of TeV in the hadronic models. We tie the truncation radius of
the thin accretion disc to the jet base radius to reduce model degeneracy
because the disc does not contribute to the electromagnetic spectrum at all.
We use the tbabs model to account for the neutral photoelectric absorption in
the intergalactic medium, using the cross-sections by Verner et al. (1996) and
the cosmic abundances by Wilms et al. (2000), where the absorption coefficient
$N_{\rm{H}}$ sets the X-ray absorption column. We use the non-relativistic
reflect function to treat in a simplified way the reflection detected in GX
339–4, parametrised primarily via the reflection fraction
$\textit{refl}=\Omega/2\pi$, which indicates the amplitude of the reflected
spectrum (Magdziarz & Zdziarski, 1995). We choose this simple model in order
to minimize the free parameters used to describe the X-ray spectrum, which is
well-fit by a power-law. Our focus is on constraining the jet physics that
drives the $\gamma$-ray band, thus we retain most of the free parameters for
that model.

We use the Interactive Spectral Interpretation System (ISIS; Houck & Denicola
2000) to forward fold the model into X-ray detector space, and to find the
statistical best fit to the data presented in Section 2.2. We use the emcee
function to explore the parameter space using a Markov Chain Monte Carlo
(MCMC) method (Foreman-Mackey et al., 2013). We initiate 20 walkers per free
parameter and perform $10^{4}$ loops. We reject the first 50% of the run as
the “burn-in” period. We provide the $1\,\sigma$ uncertainties in Table 3,
along with the results of the best fit for each model.

In Figures 3–3 we show the best fits of the multi-wavelength spectrum of GX
339–4 for the three different models we explore. In Fig. 3, we show the purely
leptonic model, whereas in Figures 3 and 3 the results of the lepto-hadronic
models.

parameter \model | leptonic | hadronic soft | hadronic hard
---|---|---|---
$p_{\rm{e}}$ | 2.2 | 2.2 | 1.7
$p_{\rm{p}}$ | - | 2.2 | 1.7
$L_{\rm{jet}}\begin{array}[]{c}\left(10^{-3}\,L_{\rm{Edd}}\right)\\\ \rm{\left(\times 10^{36}\,erg\,s^{-1}\right)}\end{array}$ | $\begin{array}[]{c}2.5_{-2}^{+5}\\\ 3_{-2}^{+6}\end{array}$ | $\begin{array}[]{c}70^{+100}\\\ 90_{-83}^{+117}\end{array}$ | $\begin{array}[]{c}50^{+60}\\\ 70_{-60}^{+70}\end{array}$
$R_{0}\left(r_{\rm{g}}\right)$ | $100^{+100}$ | $110^{+100}$ | $90^{+90}$
$T_{\rm{e}}\,(\rm{keV})$ | $1600_{-600}^{+2400}$ | $2100_{-2000}^{+2300}$ | $2000_{-1900}^{+2000}$
$f_{\rm{pl}}$ | $4_{-3}^{+5}$ | $4_{-3}^{+5}$ | $4_{-3}^{+5}$
$f_{\rm{heat}}$ | $5_{-5}^{+7}$ | $8_{-6}^{+9}$ | $7_{-6}^{+9}$
$\beta$ | $0.2_{-0.1}^{+1.5}$ | $0.2_{-0.1}^{+0.4}$ | $0.04^{+0.04}$
$L_{d}\,(10^{-3}\,L_{\rm Edd})$ | $5_{-4}^{+9}$ | $23_{-22}^{+24}$ | $2^{+3}$
$T_{\rm{cor}}\,\rm{\left(keV\right)}$ | $170_{-100}^{+200}$ | $55^{+350}$ | $60_{-50}^{+70}$
$R_{\rm{cor}}\,\left(r_{\rm{g}}\right)$ | $300_{-200}^{+300}$ | $160_{-155}^{+460}$ | $460_{-440}^{+470}$
$\tau_{\rm{cor}}$ | $0.6_{-0.4}^{+0.6}$ | $0.7_{-0.6}^{+0.8}$ | $0.7^{+0.7}$
$\chi^{2}$/DoF | 250/233 | 240.8/233 | 190/233
$\sigma$ | 1.7 | 0.03 | 0.1
$B_{0}\,\left(\rm{G}\right)$ | $2\times 10^{5}$ | $1\times 10^{6}$ | $2\times 10^{5}$
$B\,\left(\rm{G}\right)@\,z_{\rm{diss}}$ | $1\times 10^{4}$ | $6\times 10^{3}$ | $1\times 10^{4}$
$E_{\rm{p,max}}\,\rm{\left(eV\right)}$ | - | $2.8\times 10^{13}$ | $2.7\times 10^{13}$
$E_{\rm{e,max}}\,\rm{\left(eV\right)}$ | $1\times 10^{8}$ | $5.2\times 10^{10}$ | $5.3\times 10^{10}$
Table 3: Parameters for the three fitted models, distinguished via the power-
law index of the accelerated electrons $p_{\rm{e}}$ and protons $p_{\rm{p}}$.
We show the free parameters and the 1 $\sigma$ uncertainties as discussed in
Section 3 before the double line. Below the double line are indicative
evaluated quantities of the plasma magnetisation, the magnetic field, the
total luminosity of the accelerated proton/electron population and the maximum
energy of the protons/electrons at the particle acceleration region.

Figure 1: The best fit with the $\chi$-residuals of the multiwavelength
spectrum of the 2010 outburst of GX 339–4 assuming a purely leptonic model.
The solid black line shows the total intrinsic emission, the red line shows
the X-ray absorbed emission, and the rest of the components are explained in
the legend.

Figure 2: Similar to Fig. 3 but for the lepto-hadronic model with
$p_{\rm{e}}$= 2.2 and $p_{\rm{p}}$= 2.2 power-law index.

Figure 3: Similar to Fig. 3 but for the lepto-hadronic model with
$p_{\rm{e}}$= 1.7 and $p_{\rm{p}}$= 1.7 power-law index.

The unique contribution of the hadronic processes can only be seen in the TeV
$\gamma$-ray band, because the purely leptonic model cannot produce
significant emission at GeV and above. In Figures 5 and 5 we show the
predicted GeV to 100 TeV $\gamma$-ray spectrum of GX 339–4. The primary-
accelerated electrons dominate in the GeV regime via SSC. The hadronic
processes dominate in the TeV energy band, in particular, the neutral pion
decay from both pp and p$\gamma$ collisions as well as the synchrotron
radiation of secondary pairs from the latter. Because we set the acceleration
efficiency parameter $f_{\rm{sc}}$ to a high value, the protons are able to
achieve high energies of the order of $\sim 10^{13}\,$eV, producing
$\gamma$-rays of the order of TeV.

Figure 4: The $\gamma$-ray spectrum for the lepto-hadronic model with
$p_{\rm{e}}$= 2.2 and $p_{\rm{p}}$= 2.2 power-law index. We compare the
predicted spectrum to the Fermi/LAT upper limits of the 2010 outburst
(Bodaghee et al., 2013), the sensitivity of LHAASO after one year of operation
(Bai et al., 2019), and the predicted 50-hour sensitivity of the North and
South site of CTA (from www.cta-observatory.org). The solid red line shows the
pion bump from p$\gamma$, the dashed red line shows the synchrotron radiation
from secondary leptons from p$\gamma$, the dash-dotted green line shows the
pion bump from pp interactions, and the dash-double dotted line shows the
synchrotron radiation from secondary leptons from pp.

Figure 5: Similar to Fig. 5 but for the lepto-hadronic model with
$p_{\rm{e}}$= 1.7 and $p_{\rm{p}}$= 1.7 power-law index.

##  5 Discussion

###  5.1 Multi-wavelength spectrum and jet dynamics

Our results with our new lepto-hadronic multi-zone jet model confirm earlier
results that stratified jets can self-consistently reproduce the radio-to-X-
ray spectrum, together with a thin accretion disc including reflection
(Markoff et al., 2001; Zhang et al., 2010; Kylafis & Reig, 2018; Connors et
al., 2019; Lucchini et al., 2022). However, compared to earlier works (e.g.,
Markoff et al., 2005), we can also better reproduce the significantly inverted
radio-to-IR spectrum by introducing a decreasing particle acceleration
efficiency along the jets (Lucchini et al., 2021). We see however in Table 3
that the parameter $f_{pl}$ controlling this effect cannot be well-constrained
by the data, and we can only set an upper-limit.

Apart from particle acceleration, we require significant electron heating of
the thermal population(Sironi & Spitkovsky, 2009; Gedalin et al., 2012;
Plotnikov et al., 2013; Sironi et al., 2013; Sironi & Spitkovsky, 2014;
Melzani et al., 2014; Crumley et al., 2019) to reproduce both the optical and
IR bands as jet synchrotron emission. In particular, we find that the scenario
where optical emission originates from the jet base and the IR emission
originates from the particle acceleration region $z_{\rm diss}$ is consistent
with the data. An alternative scenario is that both the IR and the optical
emission originate in a hot flow that consists of thermal and non-thermal
electrons (Poutanen & Veledina, 2014; Kosenkov et al., 2020), a scenario that
better describes the soft states (Kosenkov & Veledina, 2018). Further
simultaneous IR-to-optical observations in the hard state would be able to
test this scenario, as well as simultaneous polarisation measurements across
the entire optical/IR band (although see e.g., Russell & Fender 2008 for
measurements prior to the 2010 outburst).

In both the leptonic and lepto-hadronic scenarios the shape of the radio-to-X-
ray spectrum of GX 339–4 looks identical and the radiative mechanisms are also
the same. The spectral shape is determined primarily by the jet geometry and
dynamics, which are similar between the scenarios. However, for the case of
the lepto-hadronic models, where we assume equal number density of accelerated
electrons and protons, we require much more power injected into the jet base
than for the purely leptonic model, which is a well-known issue with hadronic
models (see e.g., Pepe et al., 2015; Abeysekara et al., 2018; Kantzas et al.,
2020).

To fit the optical emission with thermal synchrotron emission from the base of
the jets while the accelerated particles fit the radio-to-IR, we require high
electron temperature. This radiation leads to a curved IC spectrum in the soft
X-rays, so another component is required to explain the hard power-law. If it
can be confirmed that the optical emission is jet synchrotron (via
polarisation for instance), then the need for a second component to fit the
X-rays will be more robust. For this reason we have added a simple thermal
corona model, which together with reflection, can well account for the X-ray
spectrum, but is otherwise independent of the jet parameters. In reality,
these components should be linked, but it is well known that spectral
information alone is often not enough to probe the detailed geometry of the
corona, which is the case in our work here as well (see e.g., Del Santo et
al., 2008; Droulans et al., 2010; Reig & Kylafis, 2015, 2021; Kylafis & Reig,
2018; Connors et al., 2019; Cao et al., 2021).

When protons are accelerated, the hadronic interactions contribute with
additional flux in the $\gamma$-ray regime of the spectrum. For the scenario
with a hard proton power-law index of $p_{\rm p}$= 1.7, producing significant
TeV flux detectable by CTA would require a non-physical amount of power
dissipated into proton acceleration. By constraining the non-thermal proton
power to 5% of the jet power, we see that the TeV flux does not exceed the CTA
sensitivity (see Fig. 5). A more typical power-law index of $p_{\rm p}$= 2.2
produces even less GeV and TeV flux. In addition, both of these models require
strongly matter-dominated outflows even at their launching point
($\sigma\lesssim 0.1$). Such a low magnetisation raises issues of physicality
for these models, since the final bulk Lorentz factor of the flow is expected
to be on the order of the initial magnetisation $\sigma$ (Komissarov et al.,
2007, 2009; Tchekhovskoy et al., 2008, 2009; Chatterjee et al., 2019).
Specifically, BHXB jets consistently show at least mildly relativistic
velocities of $\Gamma\sim 2-3$ in several systems (Mirabel & Rodriguez, 1994;
Fender, 2001; Fender et al., 2004; Casella et al., 2010; Miller-Jones et al.,
2012). Such a low initial magnetisation would struggle to explain the bulk
acceleration of the flow unless further energy is available by, e.g., thermal
pressure. However, numerical simulations show that a jet ’sheath’ forms where
the originally Poynting-flux dominated ’spine’ interacts and entrains the
surrounding disc wind, resulting in a region with much lower magnetisation
(McKinney, 2006; Móscibrodzka et al., 2016; Nakamura et al., 2018; Chatterjee
et al., 2019). The instabilities that form along this boundary are expected to
be sites of reconnection and particle acceleration (Rieger & Duffy, 2004;
Faganello et al., 2010; Rieger, 2019; Sironi et al., 2020). Thus, although our
approach is quite simplistic, it would be consistent with the emission
occurring along this boundary as suggested by recent radio observations of AGN
jets, such as M87 (Hada et al., 2016) or Cen A (Janssen et al., 2021), and
GRMHD simulations (e.g. Móscibrodzka & Falcke, 2013; Davelaar et al., 2018).
Although BHXB jets cannot be resolved by current facilities, similar scenarios
may apply to them since the systems are likely to be governed by the same
physical laws (Heinz & Sunyaev, 2003; Merloni et al., 2003; Falcke et al.,
2004).

###  5.2 Particle distributions

In Fig. 8 and 8 we plot the total distribution of the primary electrons and
protons, respectively, integrated along the jets. The MJ-only distribution at
the jet base dominates the lower energy regime, with its peak defined by the
free parameter $T_{\rm{e}}$ (see Table 3), while the higher energy electrons
originate mostly at the first particle acceleration region $z_{\rm diss}$. The
shifting of the thermal peak between the two shows the effect of the
$f_{\rm{heat}}$ parameter. The fact that the slope is steeper than $p_{\rm
e}=2$ indicates that the synchrotron cooling break occurs below $\sim 10^{9}$
eV.

In Fig. 8 we plot the differential number density of the secondary pairs from
pp and p$\gamma$ for the lepto-hadronic model with $p_{\rm p}$= 2.2. We also
include for comparison the total distribution of the primary pairs of the
jets. We note that the secondary pairs from p$\gamma$ are synchrotron cooled
and hence their spectrum is flat. The excess of particles around $\sim
10^{12}\,\rm eV$ is responsible for the TeV flux of Fig. 5.

Assuming a maximum value of $f_{\rm{sc}}$= 0.01, we see that the compact jets
of GX 339–4 can accelerate CRs up to 100 TeV. Consequently, if this is true
and moreover the entire population of BHXBs can accelerate CRs up to 100 TeV,
then BHXBs may contribute to the Galactic CR spectrum up to the knee depending
on their total number (see also Cooper et al., 2020).

Figure 6: The total electron number density of the jets multiplied by the
square of the electron energy for the lepto-hadronic model (solid black line)
with $p_{\rm{e}}=1.7$ and $p_{\rm p}$= 1.7 power-law index. We also show the
contribution of the Maxwell-Jüttner distribution at the jet base (MJ; dashed
blue line), and the contribution of the MJ plus the power-law tail of
accelerated electrons (MJ+PL; dash-dotted green line) at the particle
acceleration region $z_{\rm{diss}}$.

Figure 7: Similar to Fig. 8 but for the population of protons for the lepto-
hadronic model with $p_{\rm{e}}$= 2.2 and $p_{\rm{p}}$= 2.2 power-law index.

Figure 8: Similar to Fig. 8 but for the secondary pairs of the lepto-hadronic
model with $p_{\rm{e}}$= 2.2 and $p_{\rm{p}}$= 2.2. We also plot the
distribution of the primary electrons for comparison

###  5.3 Non-thermal proton power

The uncomfortably high proton powers needed for lepto-hadronic jet models has
been a topic of discussion for many years (see e.g., Böttcher et al., 2013;
Zdziarski & Böttcher, 2015; Liodakis & Petropoulou, 2020; Kantzas et al.,
2020). As discussed in Section 5.1, given what we see in AGN jet observations
and simulations, we would expect proton acceleration to happen primarily at
the interface between the spine and the sheath of the jet, a region of limited
volume (Rieger & Duffy, 2019). In our current setup, as a first approximation,
we can limit the volume where proton acceleration occurs by reducing the
extent of this region with respect to the total jet length. In particular,
similar to previous studies (Romero & Vila, 2008; Vila & Romero, 2010; Zhang
et al., 2010; Pepe et al., 2015; Hoerbe et al., 2020), we terminate the proton
acceleration at a distance 10 $z_{\rm diss}$ from the region where
acceleration initiates. As a consequence, we see that even for a hard power
law index of $p_{\rm p}$=1.7, the TeV emission of GX 339–4 due to hadronic
processes will not be detectable by CTA, but the energy budget remains within
reasonable values.

A further way to constrain the total power of the accelerated protons is by
increasing the minimum energy of the accelerated particles (Zdziarski &
Böttcher, 2015; Pepe et al., 2015). We nevertheless decide to use as the
minimum energy for the accelerated leptons the peak of the MJ distribution and
for the accelerated protons the rest mass energy (see Section 3), but will
explore this in more detailed future work.

Recent high resolution magneto-hydrodynamic simulations have shown that jets
can be significantly mass-loaded via instabilities at distances well beyond
the launching point (Chatterjee et al., 2019). This progressive mass-loading
could significantly reduce the total proton power and make the hadronic models
more viable, but this is a project we will pursue in the future.

###  5.4 $\gamma$-ray attenuation on the optical/IR emission

In both lepto-hadronic models, the optical emission is produced in the jet
base due to synchrotron emission from the thermal leptons. The GeV-to-TeV
$\gamma$-ray emission on the other hand, is produced in the particle
acceleration region and above, which is located at some distance of 3000
$r_{\rm g}$ from the black hole, two orders of magnitude further away from the
jet base. Moreover, the $\gamma$-ray is beamed away making it difficult for
any attenuation on this optical emission. The IR emission of GX 339–4 is
produced in the particle acceleration region where the $\gamma$-ray emission
originates as well. We therefore examine any $\gamma$-ray attenuation on the
IR emission.

We calculate the optical depth of a 3 TeV $\gamma$-ray that has the maximum
likelihood to interact with the $\sim$0.08 eV IR emission using equation 16 of
Mastichiadis (2002):

$\tau_{\gamma\gamma}=\frac{R_{\rm diss}}{4\pi}\int\epsilon_{\rm ph}n_{\rm
ph}(\epsilon_{\rm ph})\int d\Omega(1-\cos\theta)\,\sigma_{\gamma\gamma}\approx
10^{-8},$ (6)

where $\epsilon_{\rm ph}$ is the target photon energy and $n_{\rm ph}$ is the
target photon number density of the particle acceleration region. Such a small
values indicates that the particle acceleration region is optically thin to
TeV $\gamma$-rays.

###  5.5 $\gamma$-rays from BHXBs

Despite the fact that GX 339–4 is considered a ’canonical’ low-mass BHXB, it
is also amongst the most distant ones. There are Galactic low-mass BHXBs that
are as close as approximately 1–3 kpc, e.g., GRO J0422+32 (Webb et al., 2000;
Gelino & Harrison, 2003; Hynes, 2005), XTE J1118+480 (Gelino et al., 2006;
Hernández et al., 2008), XTE J1650–500 (Homan et al., 2006; Orosz et al.,
2004), GRO J1655–40 (Hjellming & Rupen, 1995; Shahbaz et al., 1999; Beer &
Podsiadlowski, 2002), GRS 1716–249 (Remillard & McClintock, 2006), GS 2000+251
(Casares et al., 1995; Barret et al., 1996; Harlaftis et al., 1996), V404 Cyg
(Miller-Jones et al., 2009), VLA J2130+12 (Kirsten et al., 2014; Tetarenko et
al., 2016b), Swift J1357.2–0933 (Shahbaz et al., 2013; Torres et al., 2015),
MAXI J1348–630 (Chauhan et al., 2020) and many more at unknown distances that
might also be as low as 2–3 kpc (see Liu et al., 2007; Kreidberg et al., 2012;
Tetarenko et al., 2016a).

For this reason we also check whether some BHXBs at a distance of 3 kpc with
the same $\gamma$-ray luminosity and spectrum as GX 339–4 could be detected by
CTA. Assuming that the jets in this putative source have identical properties
to GX 339–4, the $\gamma$-ray flux of a nearer source scales as $\left(d_{\rm
GX~{}339-4}/d_{\rm source}\right)^{2}\,F_{\gamma}$, where $d_{\rm GX~{}339-4}$
and $d_{\rm source}$ are the distances of GX 339–4 and the source,
respectively, and $F_{\gamma}$ is the $\gamma$-ray flux of GX 339–4. We plot
this $\gamma$-ray flux in Fig. 9 and compare it to the simulated sensitivity
of CTA for various energies, as a function of observation time
111http://www.cta-observatory.org/science/cta-performance/. The energy range
we study here coincides with the energy range of Fermi/LAT which as we can see
in Fig. 9 is orders of magnitude less sensitive than CTA for short integration
times. We assume that the $\gamma$-ray flux remains constant for up to one day
and its uncertainty is of the order of 30 percent. We see CTA is sensitive
enough to detect the 100 GeV emission of a GX 339–4-like source at 3 kpc
distance, with an exposure of approximately one hour, assuming the emission
remains persistent for that long. Consequently, CTA should be able to detect
GeV $\gamma$-rays from several future bright outbursts of nearby Galactic
BHXBs assuming that the accelerated particles form hard spectra within the
relativistic jets produced at peak hard/hard-intermediate states.

Figure 9: The $\gamma$-ray light curves in two energy bins as indicated in the
legend. The horizontal green lines indicate the predicted flux of a BHXB with
the same luminosity as GX 339–4 but located at a distance of 3 kpc instead. We
assume that the accelerated particles follow a power-law with index $p_{\rm
e}$= 1.7 and $p_{\rm p}$= 1.7, and the emitted flux remains constant for one
day. CTA can detect such a GeV emission within the first hour of the outburst,
but Fermi/LAT is not sensitive enough to detect such an outburst.

We finally examine a more specific example, in particular that of MAXI
J1820+070, which is at 2.96 kpc (Gandhi et al., 2019; Atri et al., 2020).
During its outburst in 2018, the source was monitored across the multi-
wavelength spectrum, from radio to X-rays (Tucker et al., 2018). Here, we
merely benchmark the spectral energy distribution instead of optimising to
determine the best fit, with the goal of illustrating the similarities and
differences with our results on GX 339–4. We use the radio-to-X-ray spectrum,
as presented by Tetarenko et al. (2021). We set the black hole mass at 8.5
$\rm{M_{\odot}}$ (Torres et al., 2020), the inclination angle at 63∘ and the
injected jet power at 15% of the Eddington luminosity (Atri et al., 2020). We
take the same model parameters we found for the best fit of GX 339–4 for the
case of $p_{\rm e}=p_{\rm p}$= 1.7 and present the spectral energy
distribution of the 2018 outburst in Fig. 10 and 11. We see that the radio-to-
X-ray spectrum is similar to the one of GX 339–4, namely the radio spectrum is
due to non-thermal synchrotron radiation, the optical band is due to thermal
synchrotron in agreement with Tetarenko et al. (2021) (although see Veledina
et al. (2019) for further contributors), and the X-ray spectrum is due to a
thermal corona. In contrast to GX 339–4, the p$\gamma$ emission exceeds the
CTA sensitivity in the sub-TeV regime. We further compare our predicted
spectrum in Fig. 11 to the upper limits set by Fermi/LAT and the Cherenkov
telescopes MAGIC, VERITAS and HESS (Hoang et al., 2019). We see that the
predicted emission exceeds the upper limits of HESS and marginally those of
VERITAS, but it is worth mentioning that these upper limits are derived after
26.9 and 12.2 hours, respectively (Hoang et al., 2019). We are unable to
capture the timing signature of the TeV emission with the current version of
our model, but we moreover do not know yet whether the high-energy emission of
these sources is persistent for up to 20-30 hours (Bodaghee et al., 2013). If
MAXI J1820+070’s TeV emission persists for at least a couple of hours during
its next outburst, it could then be a possible target-of-opportunity for CTA.
Moreover, based on the population-synthesis results of Olejak, A. et al.
(2020) and on the recent X-ray observations of Hailey et al. (2018) and Mori
et al. (2021), Cooper et al. (2020) estimated that a few thousands BHXBs may
reside in the Galactic disc capable of accelerating protons to high energy
(also see Fender et al. 2005). If these sources spend approximately 1% of
their outburst in the hard to hard-intermediate state (Tetarenko et al.,
2016a), then CTA might be able to detect a few tens of BHXBs in its first
years of operation.

In our current analysis, we assume equal number density of electrons and
protons in the jets, similar to previous studies (Vila & Romero, 2010; Connors
et al., 2019). Following this assumption, we derive the jet kinetic power to
be $2\times 10^{37}\rm erg\,s^{-1}$. Tetarenko et al. (2021) though suggest
that the jets of MAXI J1820+070 cannot be proton dominated and constrain the
ratio of protons to positrons to be $\sim 0.6$ otherwise the jet kinetic
power, which they estimate to be $6\times 10^{37}\rm erg\,s^{-1}$, may reach
18 times the accretion power. We aim to further study the impact of the pair-
to-proton ratio to jet evolution and emission in a forthcoming work.

Figure 10:  The predicted flux density based on the lepto-hadronic scenario
with $p_{\rm p}$=1.7 for the 2018 outburst of MAXI J1820+070. We compare the
total emitted spectrum to the data of Tetarenko et al. (2021). The rest of the
components are the same as in Fig. 3. Figure 11:  The spectral energy
distribution based on the lepto-hadronic scenario with $p_{\rm p}$=1.7 for the
2018 outburst of MAXI J1820+070. We compare the total emitted spectrum to the
upper limits of Fermi/LAT (F), MAGIC (M), VERITAS (V) and HESS (H) from Hoang
et al. (2019), and the predicted 50-hour sensitivity of the North and South
site of CTA (from www.cta-observatory.org). The rest of the components are the
same as in Fig. 5.

##  6 Summary and conclusions

Astrophysical jets are ideal laboratories to understand the underlying physics
of particle acceleration and the physical processes responsible for the non-
thermal emission. It is still unclear whether BHXB jets can accelerate
particles to high enough energy to shine in the $\gamma$-ray regime of the
electromagnetic spectrum. Such emission strongly depends on the composition of
the jets, which remains poorly constrained for either Galactic or
extragalactic jets. A possible hadronic composition would support BHXB jets as
candidate sources of Galactic CRs and shed light on this long-standing open
question. Understanding the jet composition is clearly crucial not only for a
better understanding of the non-thermal radiation and total power
requirements, but also for our understanding of the jet launching and bulk
acceleration properties.

To further understand the properties of Galactic jets and predict any TeV
signature, we studied the ‘canonical’ low-mass BHXB GX 339–4 during the bright
outburst of 2010. We presented the best fit of our jet model to the
multiwavelength emission and found that the whole radio-to-GeV electromagnetic
spectrum can be due to primary leptonic processes. To explain both the radio
and the IR/Optical bands, we require a heating mechanism similar to what we
see in PIC simulations (Sironi & Spitkovsky, 2009, 2011; Crumley et al.,
2019). We further found that the jets of GX 339–4 can accelerate protons to a
non-thermal power law up to a few hundreds of TeV. Depending on the power-law
index, we saw that the accelerated protons can produce a strong TeV emission
via neutral pion decay and synchrotron radiation of secondary pairs. In the
case of a hard power law of protons in particular, we found that the
photomeson processes dominate the pp interactions and the synchrotron emission
of secondary pairs dominates the sub-TeV band.

GX 339–4 is however a distant source, located at 8 kpc and the predicted TeV
flux will not be strong enough to be detected by future $\gamma$-ray
facilities, such as CTA. We rescaled the emitted spectrum to a distance of 3
kpc and compared it to the predicted timing sensitivity of CTA. We find that
CTA would be able to detect such emission with an hour of integrated
observations in the energy range above 100 GeV, which would be an indication
that protons are accelerated into a hard power law. We further tested this
scenario by bench-marking the electromagnetic spectrum of a nearby source,
such as the newly discovered BHXB MAXI J1820+070. We found that this source
might be a potential target-of-opportunity for future CTA observations to hint
BHXBs as TeV sources and CR accelerators.

## Acknowledgements

We would like to thank K. Chatterjee for many fruitful conversations on jet
physics and A. López-Oramas for the insightful discussion on $\gamma$-ray
Astronomy. DK, SM and ML are grateful for support from the Dutch Research
Council (NWO) VICI grant (no. 639.043.513). CC acknowledges support from the
Swedish Research Council (VR). R.M.T.C acknowledges support from NASA grant
NNG08FD60C. This research made use of `ASTROPY` (http://www.astropy.org), a
community-developed core `PYTHON` package for Astronomy (Astropy Collaboration
et al., 2013; Price-Whelan et al., 2018), `MATPLOTLIB` (Hunter, 2007), `NUMPY`
(Oliphant, 2006), `SCIPY` (Virtanen et al., 2020), ISIS functions
(ISISscripts) provided by ECAP/Remeis observatory and MIT
(http://www.sternwarte.unierlangen.de/isis/), and the CTA instrument response
functions provided by the CTA Consortium and Observatory (see http://www.cta-
observatory.org/science/cta-performance/ for more details).

## Data availability

All observational data in this paper are publicly available (see Table 1). The
output of our model and the plotting scripts are available in Zenodo, at
https://dx.doi.org/

## References

  * Abbasi et al. (2018) Abbasi R. U., et al., 2018, ApJ, 867, L27
  * Abeysekara et al. (2018) Abeysekara A., et al., 2018, Nature, 562, 82
  * Ackermann et al. (2013) Ackermann M., et al., 2013, Science, 339, 807
  * Aharonian (2004) Aharonian F. A., 2004, Very high energy cosmic gamma radiation: a crucial window on the extreme Universe. World Scientific
  * Aharonian et al. (2019) Aharonian F., Yang R., de Oña Wilhelmi E., 2019, Nature astronomy, 3, 561
  * Astropy Collaboration et al. (2013) Astropy Collaboration et al., 2013, A&A, 558, A33
  * Atri et al. (2020) Atri P., et al., 2020, MNRAS: Letters, 493, L81
  * Bai et al. (2019) Bai X., et al., 2019, arXiv preprint arXiv:1905.02773
  * Barret et al. (1996) Barret D., McClintock J. E., Grindlay J. E., 1996, ApJ, 473, 963
  * Beer & Podsiadlowski (2002) Beer M. E., Podsiadlowski P., 2002, MNRAS, 331, 351
  * Belloni et al. (1999) Belloni T., Méndez M., van der Klis M., Lewin W. H. G., Dieters S., 1999, ApJ, 519, L159
  * Belloni et al. (2006) Belloni T., et al., 2006, MNRAS, 367, 1113
  * Blandford & Königl (1979) Blandford R., Königl A., 1979, ApJ, 232, 34
  * Blandford & Levinson (1995) Blandford R. D., Levinson A., 1995, ApJ, 441, 79
  * Blandford & Payne (1982) Blandford R. D., Payne D. G., 1982, MNRAS, 199, 883
  * Blandford & Znajek (1977) Blandford R. D., Znajek R. L., 1977, MNRAS, 179, 433
  * Blasi (2013) Blasi P., 2013, Nuclear Physics B Proceedings Supplements, 239, 140
  * Bodaghee et al. (2013) Bodaghee A., Tomsick J. A., Pottschmidt K., Rodriguez J., Wilms J., Pooley G. G., 2013, ApJ, 775, 98
  * Böttcher & Schlickeiser (1997) Böttcher M., Schlickeiser R., 1997, arXiv preprint astro-ph/9703069
  * Böttcher et al. (2013) Böttcher M., Reimer A., Sweeney K., Prakash A., 2013, ApJ, 768, 54
  * Buxton et al. (2012) Buxton M. M., Bailyn C. D., Capelo H. L., Chatterjee R., Dinçer T., Kalemci E., Tomsick J. A., 2012, The Astronomical Journal, 143, 130
  * Callanan et al. (1992) Callanan P. J., Charles P. A., Honey W. B., Thorstensen J. R., 1992, MNRAS, 259, 395
  * Cao et al. (2021) Cao Z., Lucchini M., Markoff S., Connors R. M. T., Grinberg V., 2021, Evidence for an expanding corona based on spectral-timing modelling of multiple black hole X-ray binaries (arXiv:2110.10547)
  * Casares et al. (1995) Casares J., Charles P. A., Marsh T. R., 1995, MNRAS, 277, L45
  * Casella et al. (2010) Casella P., et al., 2010, MNRAS: Letters, 404, L21
  * Celotti et al. (2001) Celotti A., Ghisellini G., Chiaberge M., 2001, MNRAS, 321, L1
  * Chatterjee et al. (2019) Chatterjee K., Liska M., Tchekhovskoy A., Markoff S. B., 2019, MNRAS, 490, 2200
  * Chauhan et al. (2020) Chauhan J., et al., 2020, arXiv e-prints, p. arXiv:2009.14419
  * Connors et al. (2019) Connors R. M. T., et al., 2019, MNRAS, 485, 3696
  * Cooper et al. (2020) Cooper A. J., Gaggero D., Markoff S., Zhang S., 2020, MNRAS, 493, 3212
  * Coppi & Blandford (1990) Coppi P., Blandford R., 1990, MNRAS, 245, 453
  * Corbel & Fender (2002) Corbel S., Fender R. P., 2002, ApJ, 573, L35
  * Corbel et al. (2000) Corbel S., Fender R. P., Tzioumis A., Nowak M., McIntyre V., Durouchoux P., Sood R., 2000, arXiv preprint astro-ph/0003460
  * Corbel et al. (2003) Corbel S., Nowak M. A., Fender R. P., Tzioumis A. K., Markoff S., 2003, A&A, 400, 1007
  * Corbel et al. (2012) Corbel S., Coriat M., Brocksopp C., Tzioumis A. K., Fender R. P., Tomsick J. A., Buxton M. M., Bailyn C. D., 2012, MNRAS, 428, 2500
  * Corbel et al. (2013) Corbel S., et al., 2013, MNRAS: Letters, 431, L107
  * Cowley et al. (2002) Cowley A. P., Schmidtke P. C., Hutchings J. B., Crampton D., 2002, The Astronomical Journal, 123, 1741
  * Crumley et al. (2017) Crumley P., Ceccobello C., Connors R. M. T., Cavecchi Y., 2017, A&A, 601, A87
  * Crumley et al. (2019) Crumley P., Caprioli D., Markoff S., Spitkovsky A., 2019, MNRAS, 485, 5105
  * Davelaar et al. (2018) Davelaar J., Móscibrodzka M., Bronzwaer T., Falcke H., 2018, A&A, 612, A34
  * Del Santo et al. (2008) Del Santo M., Malzac J., Jourdain E., Belloni T., Ubertini P., 2008, MNRAS, 390, 227
  * Dermer & Schlickeiser (1993) Dermer C. D., Schlickeiser R., 1993, ApJ, 416, 458
  * Droulans et al. (2010) Droulans R., Belmont R., Malzac J., Jourdain E., 2010, ApJ, 717, 1022
  * Drury (2012) Drury L. O., 2012, Astroparticle Physics, 39-40, 52
  * Dziełak et al. (2019) Dziełak M. A., Zdziarski A. A., Szanecki M., De Marco B., Niedźwiecki A., Markowitz A., 2019, MNRAS, 485, 3845
  * Faganello et al. (2010) Faganello M., Pegoraro F., Califano F., Marradi L., 2010, Physics of Plasmas, 17, 062102
  * Falcke & Biermann (1995) Falcke H., Biermann P. L., 1995, A&A, 293, 665
  * Falcke et al. (2004) Falcke H., Körding E., Markoff S., 2004, A&A, 414, 895
  * Fender (2001) Fender R. P., 2001, MNRAS, 322, 31
  * Fender et al. (2004) Fender R. P., Belloni T. M., Gallo E., 2004, MNRAS, 355, 1105
  * Fender et al. (2005) Fender R. P., Maccarone T. J., van Kesteren Z., 2005, MNRAS, 360, 1085
  * Fender et al. (2006) Fender R. P., Stirling A., Spencer R., Brown I., Pooley G., Muxlow T., Miller-Jones J., 2006, MNRAS, 369, 603
  * Foreman-Mackey et al. (2013) Foreman-Mackey D., Hogg D. W., Lang D., Goodman J., 2013, Publications of the Astronomical Society of the Pacific, 125, 306
  * Frank et al. (2002) Frank J., King A., King B., Raine D., 2002, Accretion Power in Astrophysics. Accretion Power in Astrophysics, Cambridge University Press, https://books.google.gr/books?id=GGM_t-xn8ocC
  * Fürst et al. (2015) Fürst F., et al., 2015, ApJ, 808, 122
  * Gandhi et al. (2011) Gandhi P., et al., 2011, ApJ, 740, L13
  * Gandhi et al. (2019) Gandhi P., Rao A., Johnson M. A. C., Paice J. A., Maccarone T. J., 2019, MNRAS, 485, 2642
  * García et al. (2015) García J. A., Steiner J. F., McClintock J. E., Remillard R. A., Grinberg V., Dauser T., 2015, ApJ, 813, 84
  * García et al. (2019) García J. A., et al., 2019, ApJ, 885, 48
  * Gedalin et al. (2012) Gedalin M., Smolik E., Spitkovsky A., Balikhin M., 2012, EPL (Europhysics Letters), 97, 35002
  * Gelino & Harrison (2003) Gelino D. M., Harrison T. E., 2003, ApJ, 599, 1254
  * Gelino et al. (2006) Gelino D. M., Balman Ş., Kızıloğlu U., Yılmaz A., Kalemci E., Tomsick J. A., 2006, ApJ, 642, 438
  * Georganopoulos et al. (2002) Georganopoulos M., Aharonian F. A., Kirk J. G., 2002, A&A, 388, L25
  * Georganopoulos et al. (2006) Georganopoulos M., Perlman E. S., Kazanas D., McEnery J., 2006, ApJ, 653, L5
  * Ghisellini et al. (2009) Ghisellini G., Tavecchio F., Ghirlanda G., 2009, MNRAS, 399, 2041
  * Ghisellini et al. (2010) Ghisellini G., Tavecchio F., Foschini L., Ghirlanda G., Maraschi L., Celotti A., 2010, MNRAS, 402, 497
  * Hada et al. (2016) Hada K., et al., 2016, ApJ, 817, 131
  * Hailey et al. (2018) Hailey C. J., Mori K., Bauer F. E., Berkowitz M. E., Hong J., Hord B. J., 2018, Nature, 556, 70
  * Harlaftis et al. (1996) Harlaftis E. T., Horne K., Filippenko A. V., 1996, Publications of the Astronomical Society of the Pacific, 108, 762
  * Heida et al. (2017) Heida M., Jonker P. G., Torres M. A. P., Chiavassa A., 2017, ApJ, 846, 132
  * Heinz & Sunyaev (2002) Heinz S., Sunyaev R., 2002, A&A, 390, 751
  * Heinz & Sunyaev (2003) Heinz S., Sunyaev R. A., 2003, MNRAS, 343, L59
  * Hernández et al. (2008) Hernández J. I. G., Rebolo R., Israelian G., Filippenko A. V., Chornock R., Tominaga N., Umeda H., Nomoto K., 2008, ApJ, 679, 732
  * Hillas (1984) Hillas A. M., 1984, Annual review of A&A, 22, 425
  * Hjellming & Johnston (1988) Hjellming R., Johnston K., 1988, ApJ, 328, 600
  * Hjellming & Rupen (1995) Hjellming R., Rupen M., 1995, Nature, 375, 464
  * Hoang et al. (2019) Hoang J., et al., 2019, arXiv preprint arXiv:1908.06958
  * Hoerbe et al. (2020) Hoerbe M. R., Morris P. J., Cotter G., Becker Tjus J., 2020, Monthly Notices of the Royal Astronomical Society, 496, 2885
  * Homan et al. (2005) Homan J., Buxton M., Markoff S., Bailyn C. D., Nespoli E., Belloni T., 2005, ApJ, 624, 295
  * Homan et al. (2006) Homan J., Wijnands R., Kong A., Miller J. M., Rossi S., Belloni T., Lewin W. H. G., 2006, MNRAS, 366, 235
  * Houck & Denicola (2000) Houck J. C., Denicola L. A., 2000, in Manset N., Veillet C., Crabtree D., eds, ASP Conf. Ser. Vol. 216, Astronomical Data Analysis Software and Systems IX. Astron. Soc. Pac., San Francisco, p. 591
  * Hunter (2007) Hunter J. D., 2007, Computing in Science & Engineering, 9, 90
  * Hynes (2005) Hynes R. I., 2005, ApJ, 623, 1026
  * Hynes et al. (2003) Hynes R. I., Steeghs D., Casares J., Charles P. A., O'Brien K., 2003, ApJ, 583, L95
  * Hynes et al. (2004) Hynes R. I., Steeghs D., Casares J., Charles P., O’Brien K., 2004, ApJ, 609, 317
  * Janssen et al. (2021) Janssen M., et al., 2021, Nature Astronomy, 5, 1017
  * Jokipii (1987) Jokipii J., 1987, ApJ, 313, 842
  * Kantzas et al. (2020) Kantzas D., et al., 2020, MNRAS, 500, 2112
  * Kelner & Aharonian (2008) Kelner S., Aharonian F., 2008, Physical Review D, 78, 034013
  * Kelner et al. (2006) Kelner S., Aharonian F. A., Bugayov V., 2006, Physical Review D, 74, 034018
  * Kirsten et al. (2014) Kirsten F., Vlemmings W., Freire P., Kramer, Michael Rottmann, Helge Campbell, Robert M. 2014, A&A, 565, A43
  * Komissarov et al. (2007) Komissarov S. S., Barkov M. V., Vlahakis N., Königl A., 2007, MNRAS, 380, 51
  * Komissarov et al. (2009) Komissarov S. S., Vlahakis N., Königl A., Barkov M. V., 2009, MNRAS, 394, 1182
  * Kosenkov & Veledina (2018) Kosenkov I. A., Veledina A., 2018, MNRAS, 478, 4710
  * Kosenkov et al. (2020) Kosenkov I. A., Veledina A., Suleimanov V. F., Poutanen J., 2020, A&A, 638, A127
  * Kreidberg et al. (2012) Kreidberg L., Bailyn C. D., Farr W. M., Kalogera V., 2012, ApJ, 757, 36
  * Kylafis & Reig (2018) Kylafis N. D., Reig P., 2018, A&A, 614, L5
  * Levinson & Blandford (1995) Levinson A., Blandford R., 1995, ApJ, 449, 86
  * Liodakis & Petropoulou (2020) Liodakis I., Petropoulou M., 2020, ApJ, 893, L20
  * Liu et al. (2007) Liu Q. Z., van Paradijs J., van den Heuvel E. P. J., 2007, A&A, 469, 807
  * Lucchini et al. (2021) Lucchini M., Russell T. D., Markoff S. B., Vincentelli F., Gardenier D., Ceccobello C., Uttley P., 2021, MNRAS, 501, 5910
  * Lucchini et al. (2022) Lucchini M., et al., 2022, MNRAS, subm. (arXiv:2108:12011),
  * Magdziarz & Zdziarski (1995) Magdziarz P., Zdziarski A. A., 1995, MNRAS, 273, 837
  * Mahmoud et al. (2019) Mahmoud R. D., Done C., De Marco B., 2019, MNRAS, 486, 2137
  * Maitra et al. (2009) Maitra D., Markoff S., Brocksopp C., Noble M., Nowak M., Wilms J., 2009, MNRAS, 398, 1638
  * Maitra et al. (2011) Maitra D., Miller J. M., Markoff S., King A., 2011, ApJ, 735, 107
  * Malyshev et al. (2013) Malyshev D., Zdziarski A. A., Chernyakova M., 2013, MNRAS, 434, 2380
  * Mannheim (1993) Mannheim K., 1993, A&A, 269, 67
  * Mannheim & Schlickeiser (1994) Mannheim K., Schlickeiser R., 1994, A&A, 286, 983
  * Maraschi et al. (1992) Maraschi L., Ghisellini G., Celotti A., 1992, ApJ, 397, L5
  * Marcowith et al. (1995) Marcowith A., Henri G., Pelletier G., 1995, MNRAS, 277, 681
  * Markert et al. (1973) Markert T. H., Canizares C. R., Clark G. W., Lewin W. H. G., Schnopper H. W., Sprott G. F., 1973, ApJ, 184, L67
  * Markoff & Nowak (2004) Markoff S., Nowak M. A., 2004, ApJ, 609, 972
  * Markoff et al. (2001) Markoff S., Falcke H., Fender R., 2001, A&A, 372, L25
  * Markoff et al. (2003) Markoff S., Nowak M., Corbel S., Fender R., Falcke H., 2003, A&A, 397, 645
  * Markoff et al. (2005) Markoff S., Nowak M. A., Wilms J., 2005, ApJ, 635, 1203
  * Markoff et al. (2008) Markoff S., et al., 2008, ApJ, 681, 905
  * Marscher et al. (2008) Marscher A. P., et al., 2008, Nature, 452, 966
  * Mastichiadis (2002) Mastichiadis A., 2002, Radiative Processes in Relativistic Outflows. Springer Berlin Heidelberg, Berlin, Heidelberg, pp 1–23, doi:10.1007/3-540-46025-X_1
  * McClintock et al. (2006) McClintock J., Remillard R., Lewin W., Van Der Klis M., 2006, ed. WHG Lewin, M van der Klis, Cambridge Univ, p. 71
  * McKinney (2006) McKinney J. C., 2006, MNRAS, 368, 1561
  * Melzani et al. (2014) Melzani M., Walder R., Folini D., Winisdoerffer C., Favre J. M., 2014, A&A, 570, A112
  * Merloni et al. (2003) Merloni A., Heinz S., Di Matteo T., 2003, MNRAS, 345, 1057
  * Miller-Jones et al. (2009) Miller-Jones J. C. A., Jonker P. G., Dhawan V., Brisken W., Rupen M. P., Nelemans G., Gallo E., 2009, ApJ, 706, L230
  * Miller-Jones et al. (2012) Miller-Jones J. C. A., et al., 2012, MNRAS: Letters, 419, L49
  * Mirabel & Rodriguez (1994) Mirabel I., Rodriguez L., 1994, Nature, 371, 46
  * Mori et al. (2021) Mori K., et al., 2021, The Astrophysical Journal, 921, 148
  * Móscibrodzka & Falcke (2013) Móscibrodzka M., Falcke H., 2013, A&A, 559, L3
  * Móscibrodzka et al. (2016) Móscibrodzka M., Falcke H., Shiokawa H., 2016, A&A, 586, A38
  * Motta et al. (2009) Motta S., Belloni T., Homan J., 2009, MNRAS, 400, 1603
  * Muñoz Darias et al. (2008) Muñoz Darias T., Casares J., Martínez-Pais I. G., 2008, MNRAS, 385, 2205
  * Mücke et al. (2003) Mücke A., Protheroe R., Engel R., Rachen J., Stanev T., 2003, Astroparticle Physics, 18, 593
  * Nakamura et al. (2018) Nakamura M., et al., 2018, ApJ, 868, 146
  * Nowak et al. (2002) Nowak M. A., Wilms J., Dove J. B., 2002, MNRAS, 332, 856
  * Olejak, A. et al. (2020) Olejak, A. Belczynski, K. Bulik, T. Sobolewska, M. 2020, A&A, 638, A94
  * Oliphant (2006) Oliphant T. E., 2006, A guide to NumPy. Vol. 1, Trelgol Publishing USA
  * Orosz et al. (2004) Orosz J. A., McClintock J. E., Remillard R. A., Corbel S., 2004, ApJ, 616, 376
  * Parker et al. (2016) Parker M. L., et al., 2016, ApJ, 821, L6
  * Pepe et al. (2015) Pepe C., Vila G. S., Romero G. E., 2015, A&A, 584, A95
  * Perrone (2020) Perrone L., 2020, Journal of Physics: Conference Series, 1342, 012018
  * Plotnikov et al. (2013) Plotnikov I., Pelletier G., Lemoine M., 2013, MNRAS, 430, 1280
  * Poutanen & Veledina (2014) Poutanen J., Veledina A., 2014, Space Science Reviews, 183, 61–85
  * Price-Whelan et al. (2018) Price-Whelan A. M., et al., 2018, AJ, 156, 123
  * Rachen & Biermann (1993) Rachen J. P., Biermann P. L., 1993, Astronomy and Astrophysics, 272, 161
  * Reig & Kylafis (2015) Reig P., Kylafis N. D., 2015, A&A, 584, A109
  * Reig & Kylafis (2021) Reig P., Kylafis N. D., 2021, A&A, 646, A112
  * Remillard & McClintock (2006) Remillard R. A., McClintock J. E., 2006, Annual Review of Astronomy and Astrophysics, 44, 49
  * Rieger (2019) Rieger F. M., 2019, Galaxies, 7
  * Rieger & Duffy (2004) Rieger F. M., Duffy P., 2004, ApJ, 617, 155
  * Rieger & Duffy (2019) Rieger F. M., Duffy P., 2019, ApJ, 886, L26
  * Romero & Orellana (2005) Romero G. E., Orellana M., 2005, A&A, 439, 237
  * Romero & Vila (2008) Romero G. E., Vila G. S., 2008, A&A, 485, 623
  * Romero et al. (2003) Romero G. E., Torres D. F., Bernadó M. K., Mirabel I., 2003, A&A, 410, L1
  * Russell & Fender (2008) Russell D. M., Fender R. P., 2008, MNRAS, 387, 713
  * Russell et al. (2020) Russell T. D., et al., 2020, MNRAS, 498, 5772
  * Shahbaz et al. (1999) Shahbaz T., van der Hooft F., Casares J., Charles P. A., van Paradijs J., 1999, MNRAS, 306, 89
  * Shahbaz et al. (2013) Shahbaz T., Russell D. M., Zurita C., Casares J., Corral-Santana J. M., Dhillon V. S., Marsh T. R., 2013, MNRAS, 434, 2696
  * Shakura & Sunyaev (1973) Shakura N. I., Sunyaev R. A., 1973, A&A, 24, 337
  * Shidatsu et al. (2011) Shidatsu M., et al., 2011, Publications of the Astronomical Society of Japan, 63, S785
  * Sironi & Spitkovsky (2009) Sironi L., Spitkovsky A., 2009, ApJ, 698, 1523
  * Sironi & Spitkovsky (2011) Sironi L., Spitkovsky A., 2011, ApJ, 726, 75
  * Sironi & Spitkovsky (2014) Sironi L., Spitkovsky A., 2014, ApJ, 783, L21
  * Sironi et al. (2013) Sironi L., Spitkovsky A., Arons J., 2013, ApJ, 771, 54
  * Sironi et al. (2020) Sironi L., Rowan M. E., Narayan R., 2020, arXiv e-prints, p. arXiv:2009.11877
  * Tavani et al. (2009) Tavani M., et al., 2009, Nature, 462, 620
  * Tavecchio et al. (1998) Tavecchio F., Maraschi L., Ghisellini G., 1998, ApJ, 509, 608
  * Tchekhovskoy et al. (2008) Tchekhovskoy A., McKinney J. C., Narayan R., 2008, MNRAS, 388, 551
  * Tchekhovskoy et al. (2009) Tchekhovskoy A., McKinney J. C., Narayan R., 2009, ApJ, 699, 1789
  * Tetarenko et al. (2016a) Tetarenko B. E., Sivakoff G. R., Heinke C. O., Gladstone J. C., 2016a, ApJ Supplement Series, 222, 15
  * Tetarenko et al. (2016b) Tetarenko B. E., et al., 2016b, ApJ, 825, 10
  * Tetarenko et al. (2019) Tetarenko A., Casella P., Miller-Jones J., Sivakoff G., Tetarenko B., Maccarone T., Gandhi P., Eikenberry S., 2019, MNRAS, 484, 2987
  * Tetarenko et al. (2020) Tetarenko B. E., Dubus G., Marcel G., Done C., Clavel M., 2020, MNRAS, 495, 3666
  * Tetarenko et al. (2021) Tetarenko A., et al., 2021, MNRAS
  * Torres et al. (2015) Torres M. A. P., Jonker P. G., Miller-Jones J. C. A., Steeghs D., Repetto S., Wu J., 2015, MNRAS, 450, 4292
  * Torres et al. (2020) Torres M. A. P., Casares J., Jiménez-Ibarra F., Álvarez-Hernández A., Muñoz-Darias T., Padilla M. A., Jonker P. G., Heida M., 2020, ApJ, 893, L37
  * Tucker et al. (2018) Tucker M. A., et al., 2018, ApJ, 867, L9
  * Uttley et al. (2011) Uttley P., Wilkinson T., Cassatella P., Wilms J., Pottschmidt K., Hanke M., Böck M., 2011, MNRAS: Letters, 414, L60
  * Veledina et al. (2019) Veledina A., et al., 2019, A&A, 623, A75
  * Verner et al. (1996) Verner D. A., Ferland G. J., Korista K. T., Yakovlev D. G., 1996, ApJ, 465, 487
  * Vila & Romero (2010) Vila G. S., Romero G. E., 2010, MNRAS, 403, 1457
  * Vila et al. (2012) Vila Romero, G. E. Casco, N. A. 2012, A&A, 538, A97
  * Vink (2012) Vink J., 2012, A&A Review, 20, 1
  * Virtanen et al. (2020) Virtanen P., et al., 2020, Nature Methods, 17, 261
  * Völk et al. (2003) Völk H. J., Berezhko E. G., Ksenofontov L. T., 2003, A&A, 409, 563
  * Webb et al. (2000) Webb N. A., Naylor T., Ioannou Z., Charles P. A., Shahbaz T., 2000, MNRAS, 317, 528
  * Wilms et al. (2000) Wilms J., Allen A., McCray R., 2000, ApJ, 542, 914
  * Zanin et al. (2016) Zanin R., Fernández-Barral A., de Oña Wilhelmi E., Aharonian F., Blanch O., Bosch-Ramon V., Galindo D., 2016, A&A, 596, A55
  * Zdziarski & Böttcher (2015) Zdziarski A. A., Böttcher M., 2015, MNRAS: Letters, 450, L21
  * Zdziarski et al. (2004) Zdziarski A. A., Gierlinski M., Mikołajewska J., Wardziński G., Smith D. M., Harmon B. A., Kitamoto S., 2004, MNRAS, 351, 791
  * Zhang et al. (2010) Zhang J. F., Feng Y. G., Lei M. C., Tang Y. Y., Tian Y. P., 2010, MNRAS, 407, 2468

Figure 12: The inverse of the characteristic timescales for various physical
processes in the jets as indicated by the legend versus the proton kinetic
energy. The top plots correspond to a power-law index of $p_{\rm p}$= 2.2 and
the bottom plots correspond to a power-law index of $p_{\rm p}$= 1.7. The left
plots correspond to the particle acceleration region and the right plots to
the final jet segment. The vertical line in each plot shows the maximum energy

##  Appendix A Proton characteristic timescales

In Fig. 12 we show the characteristic cooling timescales of synchrotron,
escape, pp and p$\gamma$ for the accelerated protons inside the jets, in
comparison to the acceleration timescale. When any of the energy-loss
timescales intersects the acceleration timescale we derive the maximum proton
energy for every jet segment. The residence timescale is the one that defines
the maximum proton energy inside the jets.