Bethe-Heitler signature in proton synchrotron models for gamma-ray bursts

We study the effect of Bethe-Heitler (BeHe) pair production on a proton synchrotron model for the prompt emission in gamma-ray bursts (GRBs). The possible parameter space of the model is constrained by consideration of the synchrotron radiation from the secondary BeHe pairs. We find two regimes of interest. 1) At high bulk Lorentz factor, large radius and low luminosity, proton synchrotron emission dominates and produces a spectrum in agreement with observations. For part of this parameter space, a subdominant (in the MeV band) power-law is created by the synchrotron emission of the BeHe pairs. This power-law extends up to few tens or hundreds of MeV. Such a signature is a natural expectation in a proton synchrotron model, and it is seen in some GRBs, including GRB 190114C recently observed by the MAGIC observatory. 2) At low bulk Lorentz factor, small radius and high luminosity, BeHe cooling dominates. The spectrum achieves the shape of a single power-law with spectral index $\alpha = -3/2$ extending across the entire GBM/Swift energy window, incompatible with observations. Our theoretical results can be used to further constrain the spectral analysis of GRBs in the guise of proton synchrotron models.


INTRODUCTION
The emission mechanism at the origin of the observed signal during the prompt phase of GRBs remains unknown. Among the prime contenders are photospheric emission, released when the plasma becomes optically thin (Goodman 1986;Paczynski 1986;Mészáros & Rees 2000;Drenkhahn & Spruit 2002), and synchrotron emission produced by relativistic particles accelerated by shocks or magnetic reconnection once the flow is optically thin (Rees & Meszaros 1994;Sari et al. 1996;Daigne & Mochkovitch 1998;Zhang & Yan 2011). In addition, protons may also contribute, either directly by synchrotron emission, or indirectly by emission from the secondaries produced in photo-hadronic and photo-pair processes (Asano et al. 2009;Crumley & Kumar 2013;Florou et al. 2021).
When comparing models to spectral data, the most crucial difference between the aforementioned models are the prediction for the low energy spectral slope α, usually associated with the low-energy slope of the Band model (Band et al. 1993). For photospheric emission models, the slope is expected to be around α = 0.4 (Beloborodov 2010;Pe'er & Ryde 2011;Bégué et al. 2013;Parsotan & Lazzati 2018), unless the ejecta becomes transparent during the acceleration phase (Goodman 1986;Paczynski 1986;Bégué & Vereshchagin 2014;Ryde et al. 2017), in which case a steeper slope up to α = 1 can be achieved. Slopes shallower than α = 0.4 can be obtained when considering geometrical effects such as emission from a structured jet (Lundman et al. 2013), or subphotospheric dissipation (Pe'er & Waxman 2005;Giannios 2006;Vurm & Beloborodov 2016). Observationally, the footprint of photospheric emission is seen in many GRB spectra, see e.g. Ryde & Pe'er (2009) ;Acuner et al. (2020); Dereli- . Moreover, analysis of GRB 090902B strongly supports a model where the emission is produced at the photosphere of a highly relativistic outflow (Ryde et al. 2010;Pe'Er et al. 2012). In the past years, photospheric models have been directly fitted to data achieving good agreement (Ahlgren et al. 2015;Vianello et al. 2018;Samuelsson et al. 2021).
Synchrotron models predict a low energy slope to be α = −2/3 in the slow cooling regime and α = −3/2 in the fast cooling regime. Slightly steeper slopes could be obtained when considering the effects of inverse Compton cooling in the Klein-Nishina regime (Bošnjak et al. 2009;Nakar et al. 2009;Daigne et al. 2011). Yet, most GRB spectra fitted with the Band function (Band et al. 1993) are found incompatible with a synchrotron model: this is known as the "synchrotron line-of-death" (Preece et al. 1998). Recently, it was found that fitting a synchrotron model directly to GRB spectra alleviates this problem 1 (Burgess et al. 2020;Acuner et al. 2020). The main reason for the disagreement is that the Band function is a poor approximation of the synchrotron emission around the peak, giving poor constraints when comparing the fitted results to the expectations from synchrotron models in a limiting energy window. In addition, two independent analysis using Swift X-ray (Oganesyan et al. 2018) and optical data (Oganesyan et al. 2019) showed that the spectra of several GRBs require an additional break around observed energy 1 keV, leading to the straightforward identification of the injection and cooling frequencies of a synchrotron model. The spectral slopes below and between the breaks were also found compatible with the expectation from synchrotron models with power-law distributed charged particles.
The closeness of the two identified breaks requires, within the framework of synchrotron models, that the emitting particles be in the marginally fast cooling regime, with their cooling Lorentz factor γ c nearly equal to their injection Lorentz factor γ m . This requirement is difficult to account for if the radiating particles are electrons (Beniamini et al. 2018). Possible solutions include the jet in jet model (Narayan & Kumar 2009;Beniamini et al. 2018) or emission in a time dependent magnetic field (Uhm & Zhang 2014).
Alternatively, it was proposed by Ghisellini et al. (2020) that protons could be the particles radiating synchrotron and producing the main prompt MeV-peak. The observed requirement of marginally fast cooling is then naturally fulfilled for emission radius in the order of 10 13 − 10 14 cm and bulk Lorentz factor of a few hundreds (Ghisellini et al. 2020), as expected for optically thin emission models of GRBs (Rees & Meszaros 1994;Daigne & Mochkovitch 1998). On the other hand, proton synchrotron models do not explain the observed spectral peak energy clustering (von Kienlin et al. 2020) and require a large magnetic luminosity L B 10 55 erg s −1 (Florou et al. 2021).
Synchrotron emission from protons and from the secondaries produced by photo-pion ( pγ → p + π 0 , pγ → n + π + and other channels producing two or more pions (Mücke et al. 2000;Lipari et al. 2007;Hümmer et al. 2010)) and BeHe (pγ → pe + e − ) interactions was thoroughly studied in connection with ultra-high energy cosmic ray acceleration (Böttcher & Dermer 1998;Totani 1998;Razzaque et al. 2010), high-energy PeV neutrino production (Petropoulou 2014) and high energy photon component observed by LAT and more recently by HESS and MAGIC (Gupta & Zhang 2007;Asano et al. 2009;Crumley & Kumar 2013;Sahu & Fortín 2020). However, all those studies have in common the leptonic origin of the main MeV peak component, when it is not set to be a fiducial Band model. Florou et al. (2021) numerically studied a proton synchrotron model as the source of the main MeV peak similar to the one proposed by Ghisellini et al. (2020). They concluded that emission from the secondaries produced by either the BeHe-process or photo-pion interactions would be too bright to account for the optical constraints. This is especially important since this result is inconsistent with the claim that optical observations support synchrotron emission models (Oganesyan et al. 2019). However, the bursts used by Oganesyan et al. (2019) and afterwards by Florou et al. (2021) are very long duration bursts (T 90 > 70 − 80s), which is required to have simultaneous optical observations. Thus, this small subset of bursts is not necessarily representative of the full GRB population, nor of the emission mechanism producing the early episodes of a GRB. Indeed, it was suggested that the emission mechanism might change throughout the burst episodes (e.g. Zhang et al. (2018); Li (2019)). It is therefore interesting to find predictions of proton synchrotron models that do not rely on optical data, to be able to test the model on shorter duration bursts.
In this paper, we assume that the prompt emission is due to proton synchrotron and derive constraints on the model. We present analytical estimates of the effect of BeHe pair production and the pairs subsequent radiation. We identify two emission regimes: 1) proton synchrotron dominated emission regime with little contribution from other processes, 2) BeHe pair dominated emission regime, leading to a spectrum incompatible with observations. We further describe the transition between these two extreme regimes in which a subdominant power-law from the BeHe pair synchrotron radiation appears across the MeV band, as observed in some GRBs (e.g. Vianello et al. (2018); Chand et al. (2020)).
The paper is organized as follows. In Section 2, we identify the parameter space for each of the three regimes mentioned above by comparing the timescales of synchrotron emission to that of BeHe pair production. Section 3 details the modification to the spectrum due to synchrotron radiation from the pairs produced by the BeHe process. Discussion with an emphasis on GRB 190114C is given in Section 4.

CONSTRAINING PROTON SYNCHROTRON EMISSION MODELS BY BETHE-HEITLER COOLING
In this section, we compare the timescale of proton synchrotron emission to that of BeHe pair production. If the protons cool too quickly from BeHe pair creation, the secondary emission from the pairs can greatly affect the observed spectrum. The valid parameter space for proton synchrotron models can thus be constrained.
Consider an emission region expanding relativistically with Lorentz factor Γ, emitting radiation at a distance r from a central engine, and threaded by a magnetic field of comoving strength B. The comoving dynamical time is given by t dyn = r/(Γc), where c is the speed of light. In a marginally fast cooling scenario, relativistic particles, here protons, are assumed to be steadily injected into a power-law with index −p above some injection Lorentz factor γ p,m . We present our results for p = 2.5 and p = 3.5. On the one hand, the value of p ∼ 2.5 is expected in many dissipation and acceleration scenarios (e.g. Bednarz & Ostrowski (1998); Kirk et al. (2000)), albeit softer values can also be obtained from simulations (Sironi et al. 2013;Crumley et al. 2019;Comisso et al. 2020). On the other hand, synchrotron fits to the GRB spectra require an average value of p = 3.5 (Burgess et al. 2020). Marginally fast cooling implies that γ p,m ∼ γ p,c , where γ p,c is the characteristics proton cooling Lorentz factor. We write γ p,m = ξγ p,c . In this paper, we assume ξ 1, i.e., the protons are fast cooling albeit marginally. This implies that protons efficiently radiate most of their energy, while satisfying the observed spectral constraints.
The comoving cooling time for protons with Lorentz factor γ p emitting synchrotron radiation is and the frequency of the synchrotron spectral peak is ν peak = (4/3)Γqγ 2 p,m B/(πcm p ), where we used ν peak as the frequency without redshift correction, i.e. in the frame of the burst. The observed frequency is ν obs = ν peak /(1 + z).
In those equations, m p and m e are the proton and electron masses, q is the elementary charge and σ T is the Thompson cross section.
Setting the dynamical timescale and the cooling timescale equal, t dyn ∼ t synch , gives the magnetic field and the proton cooling Lorentz factor B = 2 3 √ 6πqcm where we have replaced γ p,m by ξγ p,c and have used the notation Q x = Q/10 x . Here, ν MeV = hν peak /1 MeV is the peak energy from the proton synchrotron, normalised to the value 1 MeV in agreement with observations, and h is Planck's constant. The comoving photon peak energy is hν m = hν peak /(2Γ) = 5.0 keV ν MeV Γ −1 2 . Most of the accelerated protons 2 have Lorentz factor γ p,c . Therefore, for the bulk number of accelerated protons interacting with photons at the peak, one gets which satisfies the threshold requirement for the BeHe-process (γ p,c hν m > 2m e c 2 ) unless r MeV < 0.015 × Γ 2 2 ξ 4 3 . We note that the lowest energy protons cannot satisfy the energy threshold for photo-pion interaction γ p,c hν m < 135 MeV, and it is therefore expected that neutrino production in this model be small. We comment further on the relevant cooling times in the discussion section. The model can be further constrained by the tight constraints from the IceCube (Aartsen et al. 2017) and Antares (Albert et al. 2017) experiments, as shown by Florou et al. (2021) and Pitik et al. (2021).
Having verified that all accelerated protons are energetic enough to satisfy the threshold of BeHe pair creation, we now estimate the cooling of proton by the BeHe. This timescale is a function of the comoving photon spectrum near the peak of their distribution, which itself depends on the comoving proton density 3 . Let L obs be the observed isotropic photon luminosity of the burst. Assuming the main emission mechanism is proton synchrotron, the number of radiating protons N p is where P obs synch (γ p,m ) is the observed synchrotron power emitted by a single proton with Lorentz factor γ p,m and u B = B 2 /(8π) is the comoving magnetic energy density.
For the comoving volume, we use V = 4πr 2 (r/Γ) (e.g. Pe'er (2015)), and therefore, the comoving density of radiating protons is given by To normalise the photon spectrum, it is assumed that the whole power radiated by protons with Lorentz factor γ p,m is emitted at ν m . Thus, the peak spectral energy density is where P synch = Γ −2 P obs synch . Since protons with Lorentz factor γ p,m mostly interact via BeHe process with photons close to the peak, only the shape of the photon spectrum around the peak affects the cooling rate by the BeHe process. The photon distribution close to the peak is well approximated by the synchrotron radiation of the protons even when t BeHe ∼ t synch , where t BeHe is the BeHe cooling timescale. This can be understood because the photons produced by the BeHe pairs are at different energies (see section 3), where the cross-section is smaller. Furthermore, when t BeHe ∼ t synch , the proton distribution function is not strongly changed below γ p,m . Therefore, in the marginally fast cooling scenario considered in this paper (ξ ≥ 1), the comoving photon spectrum around the peak is obtained as (Sari et al. 1998) where p is the index of the proton spectrum.
In Appendix A, we obtain the cooling rate of protons by BeHe pair production (pγ → pe + e − ) following the prescription of Chodorowski et al. (1992). The cooling is a function of photon energy in the proton rest frame, thus, one has to integrate the photon distribution over energy and angle. This is done when generating the figures, which therefore show exact results in the case of an isotropic photon distribution. Here, we present approximate analytical estimates to demonstrate how the cooling varies with the parameters. Using Equations (A5) and (A6), and assuming ξ = 1 for simplicity, the BeHe timescale is given by where x m = hν m /(m e c 2 ) and κ 0 = 40 is found to provide an adequate approximation for the cooling rate, see Appendix A. The parameters κ 0 is introduced to simplify the expression for the BeHe cooling. It is roughly the photon energy (in the proton rest frame) corresponding to the maximum cooling rate by the BeHe process. To obtain the numerical value in the bottom expression, p was set to 2.5. The expression and numerical value for a different value of p can be obtained by using Equations (A5) and (A6). The first line of Equation (9) is for protons that mostly interact with photons below ν m , while the second line describes protons interacting with photons with frequency higher than ν m . An estimate of the ratio between the BeHe cooling time and the synchrotron cooling time at the injection Lorentz factor γ p,m is given by where we used the fact that t synch = t dyn for γ p,m . This results therefore indicates similar timescale for fiducial parameters. As noted above, the similarities of the time scales implies that the peak of the proton synchrotron is not substantially modified by the BeHe process. Figure 1 shows the ratio of cooling times at γ p,m assuming ξ = 1 for different parameter choices. It is obtained by direct integration of Equation (A1). In computing this figure, we have assumed that the proton distribution function is only modified by synchrotron losses. This assumption breaks when BeHe cooling strongly dominates in the lowest domain of each panel in Figure 1. From top to bottom, the emitted luminosity is 10 54 , 10 53 , and 10 52 erg s −1 and from left to right the observed spectral peak energy is 100 keV, 300 keV, and 1 MeV, respectively. Figure 1 is made with p = 2.5. A softer value of p increases the timescale ratio for the high-energy branch, i.e., it only affects the left-most part in the panels in Figure 1. For p = 3.5 as compared to 2.5, the timescale ratio increases by a factor ∼ 2 for an order of magnitude decrease in radius.
From this figure, as well as from Equation (10), one can identify two extreme regimes. For low luminosities L obs , high Lorentz factor Γ and large radius r, the protons are largely unaffected by BeHe pair creation (yellow region in Figure  1). In this scenario, the observed spectrum is due to the synchrotron emission from the marginally fast cooling protons as described in Ghisellini et al. (2020), without any modification by BeHe. This shows that explaining GRB prompt spectra with proton synchrotron requires high bulk Lorentz factor Γ 300, in agreement with the analysis of Florou et al. (2021) who used optical constraints. On the other end, when Γ and r are small and/or L is high, the BeHe process dominates the cooling (dark region in Figure 1). In this regime, the synchrotron photons from the very fast cooling pairs quickly outnumber the proton synchrotron photons, leading to even more rapid BeHe pair creation. Therefore, most of the available proton energy is extracted by the BeHe pairs. The cooling Lorentz factor of the pairs is γ ±,c ∼ 1, corresponding to a cooling break in the observed spectrum at ∼ 10 eV, whereas the νF ν -peak energy associated to the synchrotron from the created BeHe pairs is at ∼ 100 MeV (see Equation (14)). Thus, the observed spectrum consists of a single power-law with F ν ∝ ν −1/2 between these two energies, clearly incompatible with observed GRB spectra. In between the two extreme regimes when the cooling timescales are comparable, signatures from both processes can be seen in the spectrum, and we explore this scenario in Section 3.

SPECTRAL SIGNATURE OF BETHE-HEITLER PAIRS
In this section, we obtain predictions for the comoving pair distribution and their emission spectrum in the case where the cooling via synchrotron and BeHe are comparable. In this situation, the proton distribution at γ p,m is only marginally affected by BeHe cooling. This implies that the photon spectrum at the peak energy around 1 MeV (which is the optimal photon energy for BeHe pair creation; see Equation (4) and Appendix A) is not strongly modified by synchrotron radiation from the secondaries. If the secondary emission from the pairs do substantially contribute to the BeHe cooling of the protons, the BeHe pair creation becomes exponential in time and we are instead in the regime where BeHe dominates. Here, we use the proton synchrotron photons as targets to compute the rate at which pairs are created, namely in the case t BeHe > ∼ t synch . A proton with Lorentz factor γ p produces electrons and positrons with typical Lorentz factor γ ± = κ e (m p /m e )γ p , where κ e is the inelasticity. The dependence of the inelasticity on γ p x, where x is the target photon energy in units of electron rest mass, can be found in Mastichiadis et al. (2005). For photons at the peak energy, γ p x is given by Equation (4) and is of the order of a few to a few hundreds. Looking at Figure 1 of Mastichiadis et al. (2005) for those values of γ p x, the inelasticity is found to vary between 10 −3 and 10 −4 , giving an average pair Lorentz factor between Figure 1. Comparison of cooling rate between the BeHe and synchrotron processes at γp,m. From top to bottom the emitted luminosity is 10 54 , 10 53 and 10 52 erg s −1 , while from left to right the observed spectral peak frequency is 100 keV, 300 keV and 1 MeV. The purple, red and black thick lines correspond to tBeHe = t synch , tBeHe = 0.1t synch , and tBeHe = 10t synch respectively. The thin black lines show the variability time expected from the Lorentz factor and radius tvar ∼ r/(Γ 2 c) for selected variability time 10 −2 s, 10 −1 s, 1s and 10s. The figure is made with p = 2.5. A value of p = 3.5 slightly increases the valid parameter space for proton synchrotron models by increasing the timescale ratio tBeHe/t synch at small radii. γ ± ∼ 2γ p and γ ± ∼ γ p /5. Considering that most of the protons have Lorentz factor γ p,m ∼ γ p,c , we consider that all pairs are created with Lorentz factor γ ±,m ≡ κ e (m p /m e )γ p,m . In other words, we neglect the contribution of higher energy protons in the creation of pairs with higher energies. This effect only changes the very high energy photon spectrum, which is likely to be absorbed by pair creation. In addition, since the pairs are fast cooling, which follows as the protons are marginally fast cooling, they obtain a Lorentz factor much smaller than their initial Lorentz factor in one dynamical timescale γ ±,c γ ±,m , and therefore the exact details of their injection is lost via their cooling. The pair production rate by BeHe is given by Chodorowski et al. (1992), but cannot be analytically integrated for a general photon spectrum. Analytical estimates in some specific cases were provided by Petropoulou & Mastichiadis (2015). We provide the integral expression used in our numerical computation in Appendix B. In order to get analytical estimates of the number of pairs, we write the ratio between the synchrotron power and BeHe power to be equal to the ratio of their timescales: where P tot BeHe it the power emitted by all BeHe pairs. Using P tot BeHe = γ ±,m m e c 2ṅ ± V , which is valid since the pairs are fast cooled, one obtainṡ We note the very strong dependence on the parameters, specifically the radius and the peak energy. This means that in principle both scenarios with high and low pair yield are possible.
Assuming that all pairs are produced at Lorentz factor γ ±,m and that pair annihilation is negligible (see Section 4), the continuity equation for the pairs can be solved to obtain the pair distribution. This computation is done in Appendix C, and gives namely a single power-law with index −2 extending from γ ±,m down to γ ∼ 1. Here, P e is the synchrotron power emitted by an electron and H is the Heaviside function. We note that at low energies the pair distribution should substantially deviate from this power-law because of strong synchrotron self-absorption heating and pair annihilation. Our analysis also neglects the electrons originally present in the flow (see the discussion in Section 4). We now estimate the emerging spectrum. The emitted synchrotron spectrum from the pair distribution in Equation (13) is a fast cooling power-law with F ν ∝ ν −1/2 . It extends from an observed peak frequency of down to sub-keV energies. The synchrotron frequency of the pairs linearly depends on the observed peak frequency. It also indirectly depends on the other model parameters via the value of the inelasticity. Larger values of κ e results from smaller values of γ p,m ν m , i.e., larger Lorentz factor and peak frequency, and/or smaller radius and luminosity, see Equation (4). The shape of the spectrum above ν ±,m depends on the shape of the proton distribution function and of the pair injection details. Since we are only giving analytical estimates, it is out of the scope of this paper to account for a detailed analysis at those energies. We note however that emission at GeV and eventually TeV energies are constrained by the LAT instrument on-board Fermi (e.g. Guetta et al. (2011)). The normalization of the spectrum is either obtained from the pair distribution function in Equation (13), or by considering the ratio of the synchrotron power to the BeHe power in Equation (11). Indeed, the value and parameter dependence of the ratio between the proton synchrotron peak and the BeHe peak, ν m F νm /ν ±,m F ν±,m , are well described by the ratio of the timescales given in Equation (10). This is true as long as t synch is not much larger than t BeHe , so that the target photons for the BeHe process are those produced by proton synchrotron. Figure 2 shows an example spectrum when the two timescales are comparable. A subdominant power-law extends from ∼ 100 MeV all across the observation window. This extra component is solely due to synchrotron radiation from the BeHe pairs and it is simultaneous to the main MeV emission. It is potentially detectable at low energy (few tens of keV) and in the LLE data. In making this figure, we assumed r 14 = 1, Γ 2 = 1, ν MeV = 1, L 52 = 10, ξ = 1 and p = 2.5 (left) or p = 3.5 (right). They both correspond to t BeHe /t synch ∼ 1.9 (the approximate expression in Equation (10) gives t BeHe /t synch ∼ 0.85). The figure was made with BeHe inelasticity κ e ∼ 10 −4 to determine ν m,± . This value of κ e is appropriate when γ p,m hν m /(m e c 2 ) ∼ 100, as obtained for this choice of parameters. The total number of pairs was calculated using the integral formulation of Chodorowski et al. (1992) (see Appendix B).
In this example, the overall synchrotron spectrum is strongly modified at low and high energies by an extra power-law produced by the synchrotron emission of the BeHe pairs. This is a clear spectral signature of the proton synchrotron emission model, which can help to differentiate proton synchrotron models from electron synchrotron models. Several properties of this power-law are well determined and weakly sensitive to the parameter of the model. First, its slope is set to be α = −1.5 since it is produced by electrons and positrons in the fast cooling regime. Second, it extends from low (sub-keV) energy to high energy with a peak at few tens or hundreds of MeV. Therefore, this component crosses the entire GBM energy window. Third, the maximum energy of this component is linearly correlated with the peak frequency of the MeV component, as shown by Equation (14). Finally, the intensity of this component does not have a It is clear that the emission from the secondary BeHe pairs greatly affect the overall shape of the spectrum at low ( 30 keV) and high ( 10 MeV) energies . The frequency of the second peak is independent on the power-law index. Above the photon peak at ∼ 100 MeV, the photon spectrum is not specified by our analysis as it depends on the details of the pair injection spectrum, and as such on the exact shape of the proton distribution function above γp,m. This is represented by a dashed line.
strong dependence on the uncertain value of the proton index p, but strongly depends on the luminosity, Lorentz factor and emission radius. The combination of all those characteristics provides a clear smoking-gun of proton synchrotron model.

DISCUSSION AND CONCLUSION
Proton synchrotron models are an attractive solution to explain marginally fast cooling spectra from GRBs (Ghisellini et al. 2020). We have performed a detailed investigation of the effect of BeHe cooling on the protons and of the subsequent radiation of the BeHe pairs. For high bulk Lorentz factor Γ, large radius r and low luminosity L obs , proton synchrotron emission dominates and no BeHe pair signature is expected. Conversely, synchrotron emission from the BeHe pairs dominates when the luminosity is high, and Γ and r are relatively low. This constitutes an additional test for the model: if no pair signature is observed in high luminosity bursts, large radii and high Lorentz factors are necessary. The magnetic luminosity is given by where we have replaced the Lorentz factor Γ by the variability time t v = r/(cΓ 2 ) normalised to 1s. The large dependence on the radius and variability time implies a very large magnetic luminosity if no BeHe pair signature is observed. For instance, for a variability timescale of the order a few seconds as observed in the burst sample of Burgess et al. (2020), an observed luminosity of L obs = 10 54 erg s −1 requires that r ∼ 10 15 cm to suppress the BeHe pair creation (see Figure  1). Using Equation (15), this implies a magnetic luminosity of a few 10 56 erg s −1 . Such a high magnetic luminosity suggests that the jet is magnetically dominated, leading to an acceleration rate with Γ ∝ r 1/3 (Drenkhahn & Spruit 2002;Bégué et al. 2017), and possibly magnetic reconnection as the energy dissipation mechanism, see e.g. Lyutikov & Blackman (2001); Giannios & Uzdensky (2019).
In between the two regimes outlined above, we expect a parameter space where both signatures can be observed simultaneously in the MeV band (see Figure 2). This component might have already been observed in several GRBs. Indeed, it is reminiscent of the population of bursts whose spectra seems to have two components: a main emission peak together with a subdominant component well approximated by a power-law observed from a few keV to several tens of MeV (Ackermann et al. 2010;Guiriec et al. 2015). This sub-dominant power-law can be interpreted as the radiation from BeHe pairs, identified in our work as the spectral signature of proton synchrotron models.
Most notably, the spectra of GRB 190114C detected by MAGIC (MAGIC Collaboration et al. 2019) are composed of two components in the energy range 1keV -1 GeV 4 , with low energy slopes seemingly compatible with slow α = −2/3 and fast α = −3/2 cooling synchrotron radiation in some time bins (Chand et al. (2020), however see Ajello et al. (2020)). In addition, a spectral cut-off component at energies between 50 to 100 MeV was reported in Chand et al. (2020). We speculate that this burst might possess the clear signature of proton synchrotron with BeHe cooling discussed in this paper: 1) one main peak produced by proton synchrotron in the slow cooling regime, 2) a second component compatible with fast cooling synchrotron radiation from BeHe pairs and 3) a cutoff between 50-100 MeV corresponding to the injection limit of BeHe pairs. Furthermore, analysis of panel c) and d) in Figure 2 of Chand et al. (2020) shows that the ratio between the low and the high peak frequency is about 100, in rough agreement with Equation (14).
In our analysis, there are additional emission processes not dealt with that might modify the spectrum. We briefly discuss some of them here. First, we neglected the modification of the high energy peak by γγ absorption. Indeed, Equation (14) shows that the synchrotron peak of the pairs is marginally below the energy threshold for pair creation. Therefore, only the spectrum at energies larger than the peak can be affected by this process. Pair recombination could produce an observable characteristic around observed frequency ν obs,γγ ∼ 2Γ m e c 2 /h = 2.4 × 10 22 Hz Γ 2 . Yet, the time for pair recombination approximated by t e ± →γγ = 2/(σ T n ± v ) is such that t e ± →γγ t dyn = 2.0 × 10 2 Γ MeV L −2 52 , where we have approximated the density of pairs as n ± =ṅ ± t dyn using the upper branch of Equation (12) and use for the average electron velocity v = c. Thus, pairs do not substantially recombine for our fiducial parameters. However, we note the large dependence on the parameters, therefore, such a signature could in some parameter space regions be present. Furthermore, we did not treat the radiation from the initial population of electrons. This was discussed in Ghisellini et al. (2020), who argued that under the assumption that the same number of protons and electrons are accelerated, the electrons would not be seen if they achieve the same injection spectrum as the protons since their luminosity would be a factor m e /m p lower. If instead of similar injection spectrum, both electrons and protons carry the same energy, electrons would radiate their energy in the TeV band, which should trigger a leptonic cascade. In addition, thermalization of the background electrons via synchrotron self-absorption could also change the low energy spectrum. We expect this process to change the spectrum mostly in the optical band, therefore our conclusions would be largely unaffected as we have focused on the keV to GeV band.
We assumed that photo-pion interaction is inefficient in our model since the threshold for this process is not reached for the fiducial parameters, see Equation (4). However, it was argued by Florou et al. (2021) that cooling by photo-pion dominates cooling by the BeHe process. We show in Appendix D the additional conditions that the parameters should satisfy for BeHe cooling to dominate. We find that this holds true for a large set of parameters. Yet, if a substantial amount of energy is transferred from protons to charged pions by photo-hadronic interaction, the spectrum above the peak would be modified.
In a future work, we will compute a table model with SOPRANO, a code designed to simulate lepto-hadronic processes in optically thin environment (Gasparyan et al. 2022). There, we perform direct fits of a proton synchrotron model with BeHe pair production to GRB 190114C, in order to understand if the additional power-law and its cutoff are in agreement with the model presented here. The sample of bursts in Burgess et al. (2020) will also be studied, in order to constrain the parameters and further test the model.
To conclude, proton synchrotron models have a clear smoking-gun signature: the synchrotron radiation produced by the BeHe pairs. The emerging spectrum is composed of a main synchrotron peak produced by protons, and a power-law with index α = −1.5 extending from sub-keV energies to few tens or hundreds of MeV. The peak frequency of this component is linearly linked to the frequency of the MeV peak frequency. Identification of this extra power-law in spectra will help to constrain the emission mechanism of GRB jets and their parameters. The cooling rate for a proton of Lorentz factor γ p is given in a simple form by Chodorowski et al. (1992): where α is the fine structure constant, r 0 the classical electron radius, κ = 2γ p x is the maximum energy for a photon with energy x = hν/m e c 2 in the proton rest frame, and φ(κ)/κ 2 is the energy loss rate of a single-energy proton with Lorentz factor γ p in an isotropic photon background. It is approximately given by Equations (3.14) and (3.18) of Chodorowski et al. (1992). The photon distribution n x is the photon number per energy in units electron rest mass. In Equation (A1), the photon distribution function is used at the value κ/(2γ p ). It is given by (compare with Equation (8)) where x c = hν c /(m e c 2 ) and x m = hν m /m e c 2 are the cooling and injection frequencies in units of the electron rest mass energy, ξ = ν c /ν m and u νm is the peak spectral energy density as given by Equation (7). In order to obtain analytical results that still factorize the photon spectrum, we approximate φ(κ)/κ 2 by a δ-function that is normalized to the integral between 2 < κ < 500 (those numbers are arbitrarily, but we numerically checked that they encompass most of the cooling contribution): where we find A ± = 420 and κ 0 = 40, roughly corresponding to the position of the maximum of φ(κ)/κ 2 . Using this approximation in Equation (A1) gives We now use the expression for the photon spectrum given in Equation (8) to obtain in marginally fast (ξ > 1) cooling The cooling time by the BeHe process is finally obtained by In this appendix, we seek an approximate expression for the pair production under the assumption of marginally fast cooling ξ ≡ 1. The pair production rate is given by Chodorowski et al. (1992) ∂N ± ∂t = c 3 8π where ψ(κ) is given by Equations 2.3 of Chodorowski et al. (1992). Multiplying this equation by an electron energy γ ± m e c 2 and integrating over electron energies corresponds to the cooling rate of a proton whose expression is given by Equation (A1). Replacing with the expression of the synchrotron spectrum in the marginally fast cooling scenario with γ m = γ c gives: where κ * ≥ 2 is such that is the energy of a photon at the spectral peak in the rest frame of the proton. Using the expression for the photon spectrum given by Equation (8), it becomes: We assume that electron-positron pairs are produced with a single Lorentz factor γ ±,m and we proceed by solving the kinetic equation describing the evolution of the pair distribution function ∂n ± ∂t + 1 m e c 2 ∂ ∂γ e (P e (γ)n ± ) =ṅ ± δ(γ − γ ±,m ) where the pair distribution function n ± (γ e , t) is a function of both time and electron energy, andṅ ± is the pair production rate by the BeHe process, given by Equation (13). Let be a small positive variable and let us integrate the above equation between γ ±,m − and γ ±,m + . Keeping only the zeroth order term in , we obtain 1 m e c 2 [P e n ± ] γ±,m+ γ±,m− =ṅ ± (C12) describing a jump in the solution at γ ±,m . Since particles are cooling, n ± is null for γ > γ ±,m . We now solve the equation for γ < γ ±,m . Further assuming that electron are fast cool, far from γ ±,c γ ±,m , the equation can be further simplified to ∂ ∂γ (P e n ± ) = 0 (C13) Therefore, the solution is of the form where we introduced γ ±,m for convenience. We use Equation (C12) to obtain Q = m e c 2 P e (γ ±,m )ṅ ± (C15)