© The Authors 2023

1. Introduction

After about five decades of research, some insights into gamma-ray bursts (GRBs) have been obtained. Generally speaking, there are two different phases in GRBs: the prompt emission, and the afterglow. The afterglow can be relatively well interpreted by the external shock model (Mészáros & Rees 1997; Sari et al. 1998; Huang et al. 2006; Geng et al. 2016; Lazzati et al. 2018; Xu et al. 2022). However, the radiation process of the prompt emission is still under debate. Different models have been proposed to explain the complicated prompt emission, such as the internal shock model (Rees & Meszaros 1994; Kobayashi et al. 1997; Bošnjak et al. 2009), the dissipative photosphere model (Rees & Mészáros 2005; Pe’er & Ryde 2011), and the internal-collision-induced magnetic reconnection and turbulence (ICMART) model (Zhang & Yan 2011; Zhang & Zhang 2014). The most commonly discussed model is the internal shock model, which is naturally expected for a highly variable central engine. In this model, the collision and merger of shells create relativistic shocks to accelerate particles. Then the accelerated particles will cause the observed GRB prompt emission.

The spectrum of GRB prompt emission has traditionally been described by the empirical Band function (Band et al. 1993). It has been suggested that synchrotron radiation may be responsible for the non-thermal component (Rees & Meszaros 1994; Sari et al. 1998). However, for a standard fast-cooling spectrum, the theoretical low-energy spectral index is too soft (Sari et al. 1998). Several studies have shown that this problem can be mitigated when a decaying magnetic field and the detailed cooling process are considered (Uhm & Zhang 2014; Zhang et al. 2016; Geng et al. 2018). Some authors argued that the low-energy index of the empirical Band function may be misleading (Burgess et al. 2020). They proposed that the synchrotron emission mechanism can interpret most spectra and can satisfactorily fit the observational data (Zhang et al. 2020).

In addition to the spectrum, the diversity of GRB energy is another mystery. For a typical long GRB (duration longer than 2 s), the isotropic energy is around 10⁵² − 10⁵⁴ erg. However, there also exist some low-luminosity GRBs (LLGRBs) with energies down to 10⁴⁸ erg (GRB 980425) for long GRBs and 10⁴⁶ erg (GRB 170817A) for short GRBs (duration shorter than 2 s). Some authors claimed that LLGRBs may form a distinct population of GRBs (Liang et al. 2007). However, a later study reported no clear separation between LLGRBs and standard high-luminosity GRBs in a larger GRB sample (Sun et al. 2015). On the other hand, LLGRBs can naturally be explained by off-axis jets. Due to the relativistic beaming effect, a GRB will become dimmer when the viewing angle (θ_obs) is larger than the jet opening angle (θ_jet) (Granot et al. 2002; Huang et al. 2002; Yamazaki et al. 2003a). The interesting short GRB 170817A clearly proves that at least some LLGRBs are viewed off-axis (Granot et al. 2017).

The GRB empirical relations can help us understand the physical nature of GRBs. The most famous relation is the so-called Amati relation (Amati et al. 2002). This relation connects the isotropic energy (E_iso) and the peak photon energy (E_p). The index of the Amati relation (E_p − E_iso) is about 0.5 (Amati 2006; Nava et al. 2012; Demianski et al. 2017; Minaev & Pozanenko 2020). However, it was found that this relation is not followed by some LLGRBs. These outliers of the Amati relation appear to follow a flatter track on the E_p − E_iso plane (Farinelli et al. 2021). Previous studies showed that the off-axis effect may play a role in this phenomenon (Ramirez-Ruiz et al. 2005; Dado & Dar 2012).

Analytical derivation of the on-axis and off-axis Amati relation indices has been attempted in several researches (Granot et al. 2002; Eichler & Levinson 2004; Ramirez-Ruiz et al. 2005; Dado & Dar 2012). Eichler & Levinson (2004) derived the index of the Amati relation as 0.5 by considering a uniform, axisymmetric jet with a hole cut out of it (i.e., a ring-shaped fireball). Dado & Dar (2012) argued that the Amati relation should have an index of 1/2 ± 1/6 in the framework of the so-called cannonball model. Here, we show that the standard fireball model will also lead to an index of 0.5 for the on-axis Amati relation. For the off-axis Amati relation, previous studies showed that the index should be 1/3 based on the effect of the viewing angle alone (Granot et al. 2002; Ramirez-Ruiz et al. 2005). However, the effect of the Lorentz factor was usually ignored or not fully considered for the off-axis cases. Here we argue that the Lorentz factor is equally responsible for determining the index of the off-axis Amati relation. We suggest that this index should be between 1/4 and 4/13 after fully considering both the effect of the viewing angle and the Lorentz factor.

Researchers have also tried to reproduce the Amati relation by means of numerical calculations (Yamazaki et al. 2004; Kocevski 2012; Mochkovitch & Nava 2015). Recently, Farinelli et al. (2021) used an empirical comoving-frame spectrum (a spectrum described by a smoothly broken power-law function) to simulate the prompt emission. Their radiation flux was calculated by averaging over the pulse duration. They found that the Amati relation should be $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.5}$ for on-axis and $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.25}$ for off-axis cases. Here, we consider the detailed synchrotron spectrum instead of an empirical one in our study, so that we can gain more insights into the detailed physics. Furthermore, the isotropic peak luminosity L_p can be precisely calculated in our model. As a result, in addition to deriving the Amati relation, we can also examine the Yonetoku relation between E_p and L_p (Yonetoku et al. 2004), which was previously discussed only for the on-axis cases (Zhang et al. 2009; Ito et al. 2019).

Our paper is organized as follows. In Sect. 2 we briefly describe our model. Then, in Sect. 3, we present an analytic derivation for the quantities of E_p, E_iso, and L_p for both on-axis and off-axis cases. The relations between these parameters are also derived. Numerical results are presented in Sect. 4. Next, in Sect. 5, Monte Carlo simulations are performed to confirm the theoretically derived correlations. The theoretical results are compared with the observational data in Sect. 6. Conclusions and discussion are presented in Sect. 7.

2. Jet model

We mainly focus on long GRBs. We used a simple top-hat jet model similar to that of Farinelli et al. (2021). For simplicity, the radiation was assumed to only last for a very short time interval (δt) in the local burst frame. In other words, photons are emitted from the shell at almost the same time. Hence the photon arrival time is largely decided by the curvature effect. A photon emitted at a radius r with a polar angle θ will reach the observer at

$$\begin{array}{c}\hfill t_{\mathrm{obs}}(\theta )=\frac{r}{c}(cos{\theta }_1-cos\theta )(1+z),\end{array}$$(1)

where c is the speed of light, and z is the redshift of the source. θ₁ refers to the angle from which the photons first arrive. We have θ₁ = 0 for an on-axis jet (θ_jet ≥ θ_obs), and θ₁ = θ_obs − θ_jet for an off-axis jet (θ_jet < θ_obs). The emission is produced through internal shocks, which can naturally dissipate the kinetic energy of a baryonic fireball (Zhang 2018). As a result, we have r ∼ γ²d, where γ is the bulk Lorentz factor of the shell, and d is the initial separation between the clumps ejected by the central engine (Kobayashi et al. 1997). A schematic illustration of our model is presented in Fig. 1.

Fig. 1. Geometry of our jet model. The upper panel shows the case for an on-axis jet, where θjet ≥ θobs. The lower panel shows an off-axis jet with θjet < θobs.

Usually, the light curve of a single GRB pulse is composed of a fast rising phase and an exponential decay phase (Li & Zhang 2021). Complicated processes may be involved in the prompt GRB phase, such as the hydrodynamics of the outflow and the cooling of electrons (Zhang et al. 2007). Analytical solutions will be too difficult to derive if these ingredients are included. Especially in the context of the internal shock model, our assumption of a small δt implies that the pulse is simply produced by the collision of two very thin shells. In realistic cases, a reverse shock may form and propagate backward during the collision. The pulse profile could then be affected by the hydrodynamics of the flow, not only by the curvature of the emitting shell. Depending on the initial distribution of the Lorentz factor, the effect can remain moderate, however, and the two-shell model we considered is approximately valid.

Electrons accelerated by internal shocks should be in the fast-cooling regime. They will follow a distribution of (Geng et al. 2018)

$$\begin{array}{c}\hfill \frac{\mathrm{d}{N}_{\mathrm{e}}}{\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }}=\{\begin{array}{cc}{C}_1{\gamma }_{\mathrm{e}}^{\prime -2},\hfill & {\gamma }_{\mathrm{m}}^{\prime }<{\gamma }_{\mathrm{e}}^{\prime }<{\gamma }_{\mathrm{inj}}^{\prime },\hfill \\ C_2{\gamma }_{\mathrm{e}}^{\prime -p-1},\hfill & {\gamma }_{\mathrm{inj}}^{\prime }<{\gamma }_{\mathrm{e}}^{\prime }<{\gamma }_{\max }^{\prime },\hfill \end{array}\end{array}$$(2)

with

$$\begin{array}{c}\hfill \begin{array}{c}{C}_1=N_{\mathrm{e}}{[{\gamma }_{\mathrm{m}}^{\prime -1}-(1-\frac{1}{p}){\gamma }_{\mathrm{inj}}^{\prime -1}-\frac{1}{p}{\left(\frac{{\gamma }_{\mathrm{inj}}^{\prime }}{{\gamma }_{\max }^{\prime }}\right)}^p{\gamma }_{\mathrm{inj}}^{\prime -1}]}^{-1},\\ C_2=C_1{\gamma }_{\mathrm{inj}}^{\prime -p-1},\end{array}\end{array}$$(3)

where ${\gamma }_{\mathrm{e}}^{\prime }$ is the Lorentz factor of electrons, and N_e is their total number. The spectral index p usually ranges from 2–3 (Huang et al. 2000). ${\gamma }_{\mathrm{m}}^{\prime }$ and ${\gamma }_{\max }^{\prime }$ are the minimum and maximum electron Lorentz factor, respectively. ${\gamma }_{\mathrm{m}}^{\prime }$ is less than the initial Lorentz factor of the injected electrons (${\gamma }_{\mathrm{inj}}^{\prime }$) in the fast-cooling regime (Burgess et al. 2020), while ${\gamma }_{\max }^{\prime }$ can be calculated approximately as ${\gamma }_{\max }^{\prime }\sim {10}^8{(\frac{B\prime }{1G})}^{-0.5}$ (Dai & Lu 1999; Huang et al. 2000), where B′ is the magnetic field strength in the comoving frame.

In the comoving frame, the synchrotron radiation power of these electrons at frequency ν′ is (Rybicki & Lightman 1979)

$$\begin{array}{c}\hfill P_{{\nu }^{\prime }}^{\prime }\left(N_{\mathrm{e}}\right)=\frac{\sqrt{3}{e}^3{B}^{\prime }}{m_{\mathrm{e}}{c}^2}{\displaystyle {\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\max }^{\prime }}\left(\frac{\mathrm{d}{N}_{\mathrm{e}}}{\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }}\right)F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)d{\gamma }_{\mathrm{e}}^{\prime },}\end{array}$$(4)

where e is the electron charge, and m_e is the electron mass. ${\nu }_{\mathrm{c}}^{\prime }=3{\gamma }_{\mathrm{e}}^{\prime 2}e{B}^{\prime }/4\pi m_{\mathrm{e}}c$ and $F(x)=x{\int }_x^{+\infty }{K}_{5/3}(k)\mathrm{d}k$, where $x=\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}\prime }$ and K_5/3(k) is the Bessel function.

In the local burst frame, the angular distribution of the radiation power is (Rybicki & Lightman 1979)

$$\begin{array}{c}\hfill \frac{\mathrm{d}{P}_{\nu }}{\mathrm{d}\mathrm{\Omega }}={\mathcal{D}}^3\frac{\mathrm{d}{P}_{{\nu }^{\prime }}^{\prime }}{\mathrm{d}{\mathrm{\Omega }}^{\prime }}={\mathcal{D}}^3\frac{P_{{\nu }^{\prime }}^{\prime }}{4\pi },\end{array}$$(5)

where 𝒟 = 1/[γ(1 − β cos θ)] is the Doppler factor, and ν′=(1 + z)ν_obs/𝒟. The radiation is assumed to be isotropic in the comoving frame.

Due to the different light-traveling time, photons emitted at different angles (θ) will be received by the observer at different times. This is called the curvature effect. Summing up the contribution from the whole jet, we can obtain the observed γ-ray light curve. We first consider a time interval of Δt in the local burst frame. It corresponds to a thin ring of [θ, θ + Δθ]. When the electrons are uniformly distributed in the shell, the number of electrons in the ring is

$$\begin{array}{c}\hfill N_{\theta }=\frac{\varphi (\theta )sin\theta \mathrm{\Delta }\theta }{2\pi f_{\mathrm{b}}}{N}_{\mathrm{tot}},\end{array}$$(6)

where f_b = 1 − cos θ_jet is the beaming factor, N_tot is the total number of electrons in the shell, and ϕ(θ) is the annular angle of the ring. N_tot can be derived from the isotropic-equivalent mass of the shell m_sh as N_tot = 2πf_bm_sh/m_p, where m_p is the mass of the proton. The annular angle ϕ(θ) takes the form of

$$\begin{array}{c}\hfill \varphi (\theta )=\{\begin{array}{cc}2\pi ,\hfill & \theta \le |{\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}}|,\hfill \\ 2arccos(\frac{cos{\theta }_{\mathrm{jet}}-cos{\theta }_{\mathrm{obs}}cos\theta }{sin{\theta }_{\mathrm{obs}}sin\theta }),\hfill & \theta >|{\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}}|.\hfill \end{array}\end{array}$$(7)

Since the radiation process only lasts for a short interval of δt, the energy emitted into one unit solid angle is

$$\begin{array}{c}\hfill \frac{\mathrm{\Delta }{E}_{\nu }}{\mathrm{d}\mathrm{\Omega }}={\mathcal{D}}^3\frac{P_{{\nu }^{\prime }}^{\prime }\delta t}{4\pi }.\end{array}$$(8)

This energy corresponds to a duration of Δt = r sin θΔθ/c in the local burster frame. The luminosity is then L_ν = ΔE_ν/Δt. For an observer at distance D_L, the observed flux at t_obs(θ) should be $F_{\nu }(\theta )=(1+z)L_{\nu }/4\pi D_{\mathrm{L}}^2$, that is,

$$\begin{array}{c}\hfill F_{\nu }(\theta )=\frac{(1+z)}{4\pi D_{\mathrm{L}}^2}\frac{c}{r}\frac{{\mathcal{D}}^3{P}_{{\nu }^{\prime }}^{\prime }(N_{\theta })\delta t}{sin\theta \mathrm{\Delta }\theta }.\end{array}$$(9)

The total isotropic energy can be calculated as

$$\begin{array}{c}\hfill E_{\mathrm{iso}}=\frac{4\pi D_{\mathrm{L}}^2}{(1+z)}{\displaystyle {\int }_0^{t_{\mathrm{end}}}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}{F}_{\nu }(\theta )\mathrm{d}\nu \mathrm{d}{t}_{\mathrm{obs}}(\theta ),}\end{array}$$(10)

where t_end = r(cos θ₁ − cos(θ_jet + θ_obs))(1 + z)/c is the end time of the pulse in the observer frame. The isotropic peak luminosity is

$$\begin{array}{c}\hfill L_{\mathrm{p}}=4\pi D_{\mathrm{L}}^2{\displaystyle {\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}{F}_{\nu }({\theta }_{\mathrm{p}})\mathrm{d}\nu ,}\end{array}$$(11)

where θ_p is the angle corresponding to the peak time of the light curve. The peak photon energy E_p can be derived from the time-integrated spectrum. In the next section, we present an analytic derivation for the prompt γ-ray emission, paying special attention to the three parameters E_iso, L_p, and E_p.

3. Analytic derivation of E_p, E_iso, and L_p

We first consider the on-axis cases. From Eqs. (1) and (10), we can derive E_iso as

$$\begin{array}{cc}\hfill E_{\mathrm{iso}}& =c\delta t{\displaystyle {\int }_0^{t_{\mathrm{end}}}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}\frac{{\mathcal{D}}^3{P}_{{\nu }^{\prime }}^{\prime }(N_{\theta })}{rsin\theta \mathrm{\Delta }\theta }\mathrm{d}\nu \mathrm{d}{t}_{\mathrm{obs}}(\theta )}\hfill \\ \hfill & \propto {\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}\frac{(1+z)\frac{r}{c}sin\theta \mathrm{d}\theta P_{{\nu }^{\prime }}^{\prime }(N_{\theta })}{{\gamma }^3{(1-\beta cos\theta )}^{3}rsin\theta \mathrm{\Delta }\theta }\mathrm{d}\nu }\hfill \\ \hfill & \propto \frac{1}{{\gamma }^3}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{1}{{(1-\beta cos\theta )}^3}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}(1+z)\frac{P_{{\nu }^{\prime }}^{\prime }(N_{\theta })}{\mathrm{\Delta }\theta }\mathrm{d}\nu \mathrm{d}\theta ,}\hfill \end{array}$$(12)

where θ₂ = θ_jet + θ_obs. Further combining Eq. (4), we have

$$\begin{array}{cc}\hfill E_{\mathrm{iso}}\propto \frac{(1+z)B^{\prime }}{{\gamma }^3}\frac{\sqrt{3}{e}^3}{m_{\mathrm{e}}{c}^2}& {\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{1}{{(1-\beta cos\theta )}^3}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}\frac{1}{\mathrm{\Delta }\theta }}\hfill \\ \hfill & {\displaystyle {\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\max }^{\prime }}\left(\frac{d{N}_{\theta }}{\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }}\right)F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)d{\gamma }_{\mathrm{e}}^{\prime }\mathrm{d}\nu \mathrm{d}\theta .}\hfill \end{array}$$(13)

We note that C₁ in Eq. (3) is largely dependent on ${\gamma }_{\mathrm{m}}^{\prime }$ and ${\gamma }_{\mathrm{inj}}^{\prime }$. In this study, a simple case of ${\gamma }_{\mathrm{m}}^{\prime }={\gamma }_{\mathrm{inj}}^{\prime }/10$ is assumed, which naturally leads to $C_1\propto {\gamma }_{\mathrm{inj}}^{\prime }$. Then with Eqs. (2), (3), and (6), we can further obtain

$$\begin{array}{cc}\hfill E_{\mathrm{iso}}& \propto \frac{(1+z)B^{\prime }}{{\gamma }^3}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{1}{{(1-\beta cos\theta )}^3}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}\frac{1}{\mathrm{\Delta }\theta }}\hfill \\ \hfill & {\displaystyle [{\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\mathrm{inj}}^{\prime }}{N}_{\theta }{\gamma }_{\mathrm{inj}}^{\prime }{\gamma }_{\mathrm{e}}^{\prime -2}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }+}\hfill \\ \hfill & {\displaystyle {\int }_{{\gamma }_{\mathrm{inj}}^{\prime }}^{{\gamma }_{\max }^{\prime }}{N}_{\theta }{\gamma }_{\mathrm{inj}}^{\prime p}{\gamma }_{\mathrm{e}}^{\prime -p-1}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }]\mathrm{d}\nu \mathrm{d}\theta }\hfill \\ \hfill & \propto \frac{(1+z)B^{\prime }{m}_{\mathrm{sh}}}{{\gamma }^3}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{\varphi (\theta )sin\theta }{{(1-\beta cos\theta )}^3}}\hfill \\ \hfill & {\displaystyle [{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}{\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\mathrm{inj}}^{\prime }}{\gamma }_{\mathrm{inj}}^{\prime }{\gamma }_{\mathrm{e}}^{\prime -2}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }\mathrm{d}\nu +}\hfill \\ \hfill & {\displaystyle {\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}{\int }_{{\gamma }_{\mathrm{inj}}^{\prime }}^{{\gamma }_{\max }^{\prime }}{\gamma }_{\mathrm{inj}}^{\prime p}{\gamma }_{\mathrm{e}}^{\prime -p-1}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }\mathrm{d}\nu ]\mathrm{d}\theta .}\hfill \end{array}$$(14)

For an on-axis observer, most of the observed photons are emitted by electrons within a small angle around the line of sight. According to Eq. (7), it is safe for us to simply take ϕ(θ) as 2π. Changing the integral order in Eq. (14), we obtain

$$\begin{array}{cc}\hfill E_{\mathrm{iso}}& \propto \frac{(1+z)B^{\prime }{m}_{\mathrm{sh}}}{{\gamma }^3}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{sin\theta }{{(1-\beta cos\theta )}^3}}\hfill \\ \hfill & {\displaystyle [{\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\mathrm{inj}}^{\prime }}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}{\gamma }_{\mathrm{inj}}^{\prime }{\gamma }_{\mathrm{e}}^{\prime -2}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }\mathrm{d}\nu +}\hfill \\ \hfill & {\displaystyle {\int }_{{\gamma }_{\mathrm{inj}}^{\prime }}^{{\gamma }_{\max }^{\prime }}{\int }_{{\nu }_1/(1+z)}^{{\nu }_2/(1+z)}{\gamma }_{\mathrm{inj}}^{\prime p}{\gamma }_{\mathrm{e}}^{\prime -p-1}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }\mathrm{d}\nu ]\mathrm{d}\theta }\hfill \\ \hfill & \propto \frac{(1+z)B^{\prime }{m}_{\mathrm{sh}}}{{\gamma }^3}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{sin\theta }{{(1-\beta cos\theta )}^3}}\hfill \\ \hfill & {\displaystyle [{\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\mathrm{inj}}^{\prime }}{\gamma }_{\mathrm{inj}}^{\prime }{\gamma }_{\mathrm{e}}^{\prime -2}\frac{{\gamma }_{\mathrm{e}}^{\prime 2}{B}^{\prime }}{(1+z)\gamma (1-\beta cos\theta )}}\hfill \\ \hfill & {\displaystyle {\int }_{\frac{{\nu }_1^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}/(1+z)}^{\frac{{\nu }_2^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}/(1+z)}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}(\frac{{\nu }_1^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }})\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }+}\hfill \\ \hfill & {\displaystyle {\int }_{{\gamma }_{\mathrm{inj}}^{\prime }}^{{\gamma }_{\max }^{\prime }}{\gamma }_{\mathrm{inj}}^{\prime p}{\gamma }_{\mathrm{e}}^{\prime -p-1}\frac{{\gamma }_{\mathrm{e}}^{\prime 2}{B}^{\prime }}{(1+z)\gamma (1-\beta cos\theta )}}\hfill \\ \hfill & {\displaystyle {\int }_{\frac{{\nu }_1^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}/(1+z)}^{\frac{{\nu }_2^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}/(1+z)}F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }}\right)\mathrm{d}(\frac{{\nu }_1^{\prime }}{{\nu }_{\mathrm{c}}^{\prime }})\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }]\mathrm{d}\theta }\hfill \\ \hfill & \propto \frac{B^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{sin\theta }{{(1-\beta cos\theta )}^4}[{\int }_{{\gamma }_{\mathrm{m}}^{\prime }}^{{\gamma }_{\mathrm{inj}}^{\prime }}{k}_0{\gamma }_{\mathrm{inj}}^{\prime }\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }+}\hfill \\ \hfill & {\displaystyle {\int }_{{\gamma }_{\mathrm{inj}}^{\prime }}^{{\gamma }_{\max }^{\prime }}{k}_0{\gamma }_{\mathrm{inj}}^{\prime p}{\gamma }_{\mathrm{e}}^{\prime -p+1}\mathrm{d}{\gamma }_{\mathrm{e}}^{\prime }]\mathrm{d}\theta .}\hfill \end{array}$$(15)

In the last step above, the integral of $F\left(\frac{{\nu }^{\prime }}{{\nu }_{\mathrm{c}}\prime }\right)$ has been simplified as a constant (k₀). Equation (15) can be further reduced as

$$\begin{array}{cc}\hfill E_{\mathrm{iso}}& \propto \frac{B^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{sin\theta }{{(1-\beta cos\theta )}^4}[k_1{\gamma }_{\mathrm{inj}}^{\prime }({\gamma }_{\mathrm{inj}}^{\prime }-{\gamma }_{\mathrm{m}}^{\prime })}\hfill \\ \hfill & +k_2{\gamma }_{\mathrm{inj}}^{\prime p}({\gamma }_{\mathrm{inj}}^{\prime 2-p}-{\gamma }_{\max }^{\prime 2-p})]\mathrm{d}\theta \hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{sin\theta }{{(1-\beta cos\theta )}^4}\mathrm{d}\theta ,}\hfill \end{array}$$(16)

where k₁ and k₂ are integration constants. When θ ≪ 1 and β ∼ 1, we have 1/(1 − β cos θ)∼2γ²/(θ²γ² + 1). The integral in Eq. (16) can therefore be approximated as

$$\begin{array}{cc}\hfill {\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}}\frac{sin\theta }{{(1-\beta cos\theta )}^4}\mathrm{d}\theta & \approx {\displaystyle {\int }_{{\theta }_1}^{{\theta }_2}\frac{\theta }{{(\frac{{\theta }^2{\gamma }^2+1}{2{\gamma }^2})}^4}\mathrm{d}\theta }\hfill \\ \hfill & =16{\gamma }^6{\int }_{{\theta }_1}^{{\theta }_2}\frac{\theta \gamma }{{({\theta }^2{\gamma }^2+1)}^4}\mathrm{d}(\theta \gamma )\hfill \\ \hfill & =16{\gamma }^6{-\frac{1}{6{({(\theta \gamma )}^2+1)}^3}|}_{{\theta }_1=0}^{{\theta }_2={\theta }_{\mathrm{obs}}+{\theta }_{\mathrm{jet}}\gg \frac{1}{\gamma }}\hfill \\ \hfill & \approx \frac{8}{3}{\gamma }^6.\hfill \end{array}$$(17)

Finally, combining Eqs. (16) and (17), it is easy to obtain $E_{\mathrm{iso}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{B}^{\prime 2}{\gamma }^2$.

Now we derive the peak luminosity of L_p. Since the main difference between Eqs. (10) and (11) is the integral of time, we can first consider the luminosity corresponding to a particular angle θ as

$$\begin{array}{cc}\hfill L(\theta )& \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{r{\gamma }^4}\frac{1}{{(1-\beta cos\theta )}^4}\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{d{\gamma }^6}\frac{1}{{(\frac{{\theta }^2{\gamma }^2+1}{2{\gamma }^2})}^4}\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{\gamma }^2}{d}\frac{1}{{({\theta }^2{\gamma }^2+1)}^4}.\hfill \end{array}$$(18)

The flux peaks at θ_p ∼ 0, thus the peak luminosity is $L_{\mathrm{p}}=L({\theta }_{\mathrm{p}})\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{\mathrm{d}}^{-1}{B}^{\prime 2}{\gamma }^2$.

The peak photon energy in the observer frame can be derived by considering the standard synchrotron emission mechanism,

$$\begin{array}{cc}\hfill E_{\mathrm{p},\mathrm{obs}}& \approx \frac{1}{(1+z)}\frac{1}{\gamma (1-\beta cos\theta )}\frac{3he{B}^{\prime }}{4\pi m_{\mathrm{e}}c}{\gamma }_{\mathrm{inj}}^{\prime 2}\hfill \\ \hfill & \propto \frac{B^{\prime }{\gamma }_{\mathrm{inj}}^{\prime 2}}{(1+z)}\frac{2\gamma }{{\theta }^2{\gamma }^2+1}\hfill \\ \hfill & \approx \frac{B^{\prime }{\gamma }_{\mathrm{inj}}^{\prime 2}}{(1+z)}\frac{2\gamma }{{\theta }_{\mathrm{p}}^2{\gamma }^2+1}\hfill \\ \hfill & \propto \frac{1}{(1+z)}\gamma B^{\prime }{\gamma }_{\mathrm{inj}}^{\prime 2}.\hfill \end{array}$$(19)

Therefore, the peak photon energy in the comoving frame is $E_{\mathrm{p}}=(1+z)E_{\mathrm{p},\mathrm{obs}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}B\prime \gamma$.

Next, we consider the off-axis cases. Again, we first focus on L_p. For an off-axis jet, the main difference is that ϕ(θ) can no longer be taken as 2π. Instead, it should be calculated as $\varphi (\theta )\equiv 2arccos(\frac{cos{\theta }_{\mathrm{jet}}-cos{\theta }_{\mathrm{obs}}cos\theta }{sin{\theta }_{\mathrm{obs}}sin\theta })$. Equation (18) now becomes

$$\begin{array}{c}\hfill L(\theta )\propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{\mathrm{d}{\gamma }^6}\frac{1}{{(1-\beta cos\theta )}^4}2arccos\left(\frac{cos{\theta }_{\mathrm{jet}}-cos{\theta }_{\mathrm{obs}}cos\theta }{sin{\theta }_{\mathrm{obs}}sin\theta }\right).\end{array}$$(20)

Here θ ∈ (θ_obs − θ_jet, θ_obs + θ_jet). We define $R(\theta )=\frac{cos{\theta }_{\mathrm{jet}}-cos{\theta }_{\mathrm{obs}}cos\theta }{sin{\theta }_{\mathrm{obs}}sin\theta }$. Then it is obvious that R(θ)∈(0, 1), which leads to the approximation $arccos(R(\theta ))\sim \sqrt{2}\sqrt{1-R(\theta )}$. We note that R(θ) can be simplified as

$$\begin{array}{cc}\hfill R(\theta )& =\frac{cos{\theta }_{\mathrm{jet}}-cos{\theta }_{\mathrm{obs}}cos\theta }{sin{\theta }_{\mathrm{obs}}sin\theta }\hfill \\ \hfill & \approx \frac{1-{\theta }_{\mathrm{jet}}^2/2-(1-{\theta }_{\mathrm{obs}}^2/2)(1-{\theta }^2/2)}{{\theta }_{\mathrm{obs}}\theta }\hfill \\ \hfill & \approx \frac{{\theta }_{\mathrm{obs}}^2-{\theta }_{\mathrm{jet}}^2+{\theta }^2}{2{\theta }_{\mathrm{obs}}\theta }.\hfill \end{array}$$(21)

When $\theta >\frac{1}{\gamma }$, we have $\frac{1}{{(1-\beta cos\theta )}^4}\sim \frac{16{\gamma }^8}{{({\theta }^2{\gamma }^2+1)}^4}\sim \frac{16}{{\theta }^8}$. Combining Eqs. (20) and (21), we obtain

$$\begin{array}{cc}\hfill L(\theta )& \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{\mathrm{d}{\gamma }^6}\frac{1}{{\theta }^8}\sqrt{1-R(\theta )}\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{\mathrm{d}{\gamma }^6}\frac{1}{{\theta }^8}\sqrt{\frac{{\theta }_{\mathrm{jet}}^2-{({\theta }_{\mathrm{obs}}-\theta )}^2}{2{\theta }_{\mathrm{obs}}\theta }}.\hfill \end{array}$$(22)

In this case, L_p will peak at ${\theta }_{\mathrm{p}}=\frac{1}{15}(16{\theta }_{\mathrm{obs}}-\sqrt{255{\theta }_{\mathrm{jet}}^2+{\theta }_{\mathrm{obs}}^2})\sim \frac{16}{15}({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})$. We note that the condition of ${\theta }_{\mathrm{p}}>\frac{1}{\gamma }$ leads to $\frac{16}{15}({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})>\frac{1}{\gamma }$. Under this condition, L_p can be further written as

$$\begin{array}{cc}\hfill L_{\mathrm{p}}& \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{\mathrm{d}{\gamma }^6}\frac{1}{{\theta }_{\mathrm{p}}^8}\sqrt{\frac{{\theta }_{\mathrm{jet}}^2-{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{p}})}^2}{2{\theta }_{\mathrm{obs}}{\theta }_{\mathrm{p}}}}\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{\mathrm{d}{\gamma }^6}\frac{1}{{\theta }_{\mathrm{p}}^{8.5}}\sqrt{\frac{{\theta }_{\mathrm{jet}}^2-{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{p}})}^2}{2{\theta }_{\mathrm{obs}}}}\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{\mathrm{d}{\gamma }^6}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-8.5}I({\theta }_{\mathrm{obs}},{\theta }_{\mathrm{jet}},\theta ={\theta }_{\mathrm{p}}).\hfill \end{array}$$(23)

The last factor is $I({\theta }_{\mathrm{obs}},{\theta }_{\mathrm{jet}},\theta ={\theta }_{\mathrm{p}})=\sqrt{\frac{{\theta }_{\mathrm{jet}}^2-{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{p}})}^2}{2{\theta }_{\mathrm{obs}}}}$. Its effect is not important compared with (θ_obs − θ_jet)^−8.5, and therefore, we ignore this factor for simplicity. Finally, we obtain $L_{\mathrm{p}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{d}^{-1}{B}^{\prime 2}{\gamma }^{-6}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-8.5}$.

Now we continue to calculate E_iso. Since most of the γ-ray energy is released in the decay stage of the light curve, we only need to integrate the energy over a θ ranging from θ_p to θ_obs + θ_jet, that is,

$$\begin{array}{cc}\hfill E_{\mathrm{iso}}& \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_{\mathrm{p}}}^{{\theta }_{\mathrm{obs}}+{\theta }_{\mathrm{jet}}}\frac{sin\theta }{{(1-\beta cos\theta )}^4}}\hfill \\ \hfill & 2arccos\left(\frac{cos{\theta }_{\mathrm{jet}}-cos{\theta }_{\mathrm{obs}}cos\theta }{sin{\theta }_{\mathrm{obs}}sin\theta }\right)\mathrm{d}\theta \hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_{\mathrm{p}}}^{{\theta }_{\mathrm{obs}}+{\theta }_{\mathrm{jet}}}\frac{1}{{\theta }^7}\sqrt{\frac{{\theta }_{\mathrm{jet}}^2-{({\theta }_{\mathrm{obs}}-\theta )}^2}{2{\theta }_{\mathrm{obs}}\theta }}\mathrm{d}\theta }\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_{\mathrm{p}}}^{{\theta }_{\mathrm{obs}}+{\theta }_{\mathrm{jet}}}\frac{1}{{\theta }^{7.5}}I({\theta }_{\mathrm{obs}},{\theta }_{\mathrm{jet}},\theta )\mathrm{d}\theta }\hfill \\ \hfill & \approx \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\displaystyle {\int }_{{\theta }_{\mathrm{p}}}^{{\theta }_{\mathrm{obs}}+{\theta }_{\mathrm{jet}}}\frac{1}{{\theta }^{7.5}}I({\theta }_{\mathrm{obs}},{\theta }_{\mathrm{jet}},\theta ={\theta }_p)\mathrm{d}\theta }\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}I({\theta }_{\mathrm{obs}},{\theta }_{\mathrm{jet}},\theta ={\theta }_{\mathrm{p}}){{\theta }^{-6.5}|}_{{\theta }_{\mathrm{p}}}^{{\theta }_{\mathrm{obs}}+{\theta }_{\mathrm{jet}}}\hfill \\ \hfill & \propto \frac{B^{\prime 2}{\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}}{{\gamma }^4}{\theta }_{\mathrm{p}}^{-6.5}.\hfill \end{array}$$(24)

Here the value of I(θ_obs, θ_jet, θ) will not change significantly in the range from θ_p to θ_obs + θ_jet. Hence it can be taken as a constant for simplicity. Since ${\theta }_{\mathrm{p}}\sim \frac{16}{15}({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})$, we obtain the result as $E_{\mathrm{iso}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{B}^{\prime 2}{\gamma }^{-4}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-6.5}$.

Finally, we calculate E_p. Similar to Eq. (19), the peak energy in the rest frame is

$$\begin{array}{cc}\hfill E_{\mathrm{p}}& \approx B^{\prime }{\gamma }_{\mathrm{inj}}^{\prime 2}\frac{2\gamma }{{\theta }_{\mathrm{p}}^2{\gamma }^2+1}\hfill \\ \hfill & \approx B^{\prime }{\gamma }_{\mathrm{inj}}^{\prime 2}\frac{2\gamma }{{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^2{\gamma }^2}\hfill \\ \hfill & \propto \frac{B^{\prime }{\gamma }_{\mathrm{inj}}^{\prime 2}}{\gamma }{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-2}.\hfill \end{array}$$(25)

We see that $E_{\mathrm{p}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}B\prime {\gamma }^{-1}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-2}$.

To summarize, in the on-axis cases, we have

$$\begin{array}{cc}& E_{\mathrm{p}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{B}^{\prime }\gamma ,\hfill \\ \hfill & E_{\mathrm{iso}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{B}^{\prime 2}{\gamma }^2,\hfill \\ \hfill & L_{\mathrm{p}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{d}^{-1}{B}^{\prime 2}{\gamma }^2,\hfill \end{array}$$(26)

while in the off-axis cases (θ_obs − θ_jet > 1/γ), we obtain

$$\begin{array}{cc}& E_{\mathrm{p}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{B}^{\prime }{\gamma }^{-1}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-2},\hfill \\ \hfill & E_{\mathrm{iso}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{B}^{\prime 2}{\gamma }^{-4}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-6.5},\hfill \\ \hfill & L_{\mathrm{p}}\propto {\gamma }_{\mathrm{inj}}^{\prime 2}{m}_{\mathrm{sh}}{d}^{-1}{B}^{\prime 2}{\gamma }^{-6}{({\theta }_{\mathrm{obs}}-{\theta }_{\mathrm{jet}})}^{-8.5}.\hfill \end{array}$$(27)

The equations derived above are very simple but intriguing. There are eight input parameters in our model: γ, B′, ${\gamma }_{\mathrm{inj}}^{\prime }$, θ_jet, θ_obs, m_shc², d, and the electron spectral index p. Seven of them are closely related to the prompt emission parameters of E_p, E_iso, and L_p, while p is largely irrelevant. Of all the input parameters, γ plays the most important role in affecting the observed spectral peak energies and fluxes for on-axis cases (Ghirlanda et al. 2018). Usually, γ may vary in a wide range for different GRBs, from lower than 100 to more than 1000. In contrast, ${\gamma }_{\mathrm{inj}}^{\prime }$ and m_sh only vary in much narrower ranges. Taking ${\gamma }_{\mathrm{inj}}^{\prime }$ and m_sh approximately as constants, from Eq. (26), we can easily derive the on-axis Amati relation as $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.5}$ and the on-axis Yonetoku relation as $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.5}$. In both relations, the power-law indices are 0.5, which is largely consistent with observations (Nava et al. 2012; Demianski et al. 2017).

In the off-axis cases, previous studies mainly focused on the effect of the variation of the viewing angle among different GRBs. Consequently, the power-law index of the off-axis Amati relation is derived as 1/3 based on the sharp-edge homogeneous jet geometry (Granot et al. 2002; Ramirez-Ruiz et al. 2005). Here, from our Eq. (27), we see that E_iso and L_p sensitively depend on both the viewing angle and the Lorentz factor. This indicates that both γ and θ_obs are important parameters that will affect the slope of the Amati relation and the Yonetoku relation. In Eq. (27), if only the variation of the viewing angle is considered (i.e., the indices of (θ_obs − θ_jet) in Eq. (27) are taken into account, but the item of γ is omitted), an index of 4/13 will be derived for the off-axis Amati relation, which is very close to the previous result of 1/3. When the combined effect of varying γ and θ_obs is included, the off-axis Amati relation is derived as $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{1/4\sim 4/13}$, and the corresponding off-axis Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{1/6\sim 4/17}$. The power-law indices become smaller than in the on-axis cases.

4. Numerical results

Some approximations have been made in deriving Eqs. (26) and (27). In this section, we carry out numerical simulations to confirm whether the above analytical derivations are correct. For convenience, a set of standard values were taken for the eight input parameters in our simulations, as shown in Table 1. When studying the effect of one particular parameter, we only changed this parameter, but fixed all other parameters at the standard values.

To examine the correctness of our analytical derivations, we chose one input parameter as a variable, and standard values were taken for all other input parameters so that the dependence of E_p, E_iso, and L_p on that particular parameter can be illustrated.

We first consider the effect of the bulk Lorentz factor γ for both on-axis and off-axis cases. Here we let γ vary between 50 and 1000. Figure 2 shows the time-integrated spectrum obtained for different γ. The peak photon energy is correspondingly marked with a black star on the curve.

Fig. 2. Time-integrated spectrum for different γ values (marked in the figure). Panel a: Spectrum of on-axis cases. Panel b: Spectrum of off-axis cases. The black star on each curve marks the position of the peak energy (Ep).

For on-axis cases (θ_jet ≥ θ_obs), the numerical results are plotted in the left panel of Fig. 2. When the value of γ increases, the peak photon energy and the flux increase simultaneously. However, a different trend is found for the off-axis cases (θ_jet < θ_obs), as shown in the right panel of Fig. 2. Both the flux and E_p decrease as γ increases. The flux is significantly lower than that of the on-axis cases. The results are consistent with those of Farinelli et al. (2021).

Figure 2 clearly shows that both E_p and the flux are positively correlated with γ in the on-axis cases, but they are negatively correlated with γ in the off-axis cases. Next, we wish to determine the precise relation between the prompt emission parameters and the input parameters. Hence, we plot E_p, E_iso, and L_p as functions of γ, B′, ${\gamma }_{\mathrm{inj}}^{\prime }$, |θ_obs − θ_jet|, and p in Fig. 3.

Fig. 3. Numerical results showing the dependence of the prompt emission features on the input parameters. The calculated prompt emission parameters (Ep, Eiso, and Lp) are plotted against the input parameters (γ, B′, ${\gamma }_{\mathrm{inj}}^{\prime }$, |θobs − θjet|, and p) for both the on-axis cases (with cross symbols) and the off-axis cases (with circle symbols). For comparison, the power-law indices derived from our analytical solutions are marked correspondingly in each panel.

Generally, the simulation results are well consistent with our analytical results of Eqs. (26) and (27). Especially E_p, E_iso, and L_p have a positive dependence on the Lorentz factor γ in the on-axis cases, but they have a negative dependence on γ in the off-axis cases. It is also interesting to note that the prompt emission parameters are nearly independent of θ_obs − θ_jet in the on-axis cases, since the line of sight is always within the jet cone. However, in the off-axis cases, the dependence of the prompt emission parameters on θ_obs − θ_jet can be described as a two-phase behavior. When θ_obs − θ_jet ≪ 1/γ, the prompt emission parameters are independent of θ_obs − θ_jet. When θ_obs − θ_jet > 1/γ, however, they decrease sharply with the increase of θ_obs − θ_jet, as indicated in Eq. (27). Finally, the last column of Fig. 3 shows that the parameter p has little impact on the prompt emission parameters.

5. Monte Carlo simulations

In this section, we perform Monte Carlo simulations to generate a large number of mock GRBs to further test the existence of the Amati relation and the Yonetoku relation. For this purpose, we performed Monte Carlo simulations as follows. First, a group of eight input parameters is generated randomly, assuming that each parameter follows a particular distribution. This group of parameters defines a mock GRB. Second, our model is applied to numerically calculate the corresponding values of E_p, E_iso, and L_p for this mock GRB. The above two steps are repeated until we obtain a large sample of mock GRBs. Finally, the mock bursts are plotted on the E_p − E_iso plane and the E_p − L_p plane, and a best-fit correlation is obtained for them.

For convenience, we designate the distribution function of a particular parameter a as D(a). Here a refers to the input parameters (i.e., γ, B′, ${\gamma }_{\mathrm{inj}}^{\prime }$, θ_jet, θ_obs, p, m_shc², and d). Three kinds of distributions are adopted in this study: (i) A normal Gaussian distribution, which is noted as N(μ, σ), where μ and σ refer to the mean value and standard deviation, respectively. (ii) A uniform distribution, which is noted as U(min,max), where min and max define the range of the parameter. Finally, (iii) a log-normal distribution, noted as D(log(a)), which means that log(a) follows a normal Gaussian distribution of N(μ, σ). The detailed distribution for each of the input parameters adopted in our Monte Carlo simulations is described below.

The ranges of some input parameters, such as γ, θ_jet, and θ_obs, can be inferred from observations. For the bulk Lorentz factor, a distribution of D(log(γ)) = N(2.2, 0.8) was assumed, which is consistent with observational constraints (Liang et al. 2010; Ghirlanda et al. 2018). As for the angle parameters, we set θ_jet as a uniform distribution of D(θ_jet) = U(0.04, 0.2) for the on- and off-axis cases (Frail et al. 2001; Wang et al. 2018). The distribution of θ_obs was taken as D(θ_obs) = U(0.0, 0.2) in the on-axis cases, together with the restriction of θ_jet ≥ θ_obs. On the other hand, in the off-axis cases, θ_obs was taken as D(θ_obs) = U(0.1, 0.5), together with θ_jet < θ_obs.

Other parameters cannot be directly measured from observations. Their distributions are then taken based on some theoretical assumptions. The isotropic ejecta mass is usually thought to range from 10⁻³ M_⊙ to 10⁻¹ M_⊙ for a typical long GRB. Hence we assumed that m_shc² has a log-normal distribution. We set its mean value as 10⁻² M_⊙ c² and its 2σ range as 10⁻³ M_⊙ c² to 10⁻¹ M_⊙ c², that is, D(log(m_shc²/erg)) = N(52.2, 0.5). The distribution of d was taken as D(log(d/cm)) = N(10, 0.2), which, combined with the above distribution of γ, gives the range of the internal shock radius as r ∼ γ²d ∼ 10¹³ − 10¹⁶ cm. It is largely consistent with the theoretical expectation of 10¹¹ − 10¹⁷ cm (Zhang 2018).

Some parameters (B′, ${\gamma }_{\mathrm{inj}}^{\prime }$, p) concern the microphysics of relativistic shocks. The comoving magnetic field strength (B′) is hard to estimate from observations (Zhang & Yan 2011; Burgess et al. 2020). Here, a log-normal distribution with a mean value of 10 G and a standard deviation of 0.4 G was assumed for it, that is, D(log(B′/G)) = N(1, 0.4). The mean value of ${\gamma }_{\mathrm{inj}}^{\prime }$ was taken as 10⁵, with a distribution of $D(log({\gamma }_{\mathrm{inj}}^{\prime }))=N(5,0.2)$ (Burgess et al. 2020). For the parameter p, we simply assumed a normal distribution of D(p) = N(2.5, 0.4).

With the above distribution functions of the input parameters, two samples of mock GRBs were generated through Monte Carlo simulations. One sample included 200 on-axis events, and the other sample included 200 off-axis cases. E_p, E_iso, and L_p were calculated for each mock event.

The distribution of the simulated bursts on the E_p − E_iso and E_p − L_p planes is shown in Fig. 4. In the on-axis cases, the best-fit slope is 0.495 ± 0.015 for the E_p − E_iso correlation, and it is 0.511 ± 0.016 for the E_p − L_p correlation. These two indices are very close to our theoretical results of 0.5, proving credibility for our analytical derivations. The intrinsic scatters of both relations are about 0.21, indicating that the correlations are rather tight.

Fig. 4. Monte Carlo simulation results. (a,b) Distribution of simulated bursts on the Ep − Eiso plane and the Ep − Lp plane for the on-axis cases. (c,d) Distribution of simulated bursts on the Ep − Eiso plane and the Ep − Lp plane for off-axis cases. The solid lines show the best-fitting results for the simulated data points, and the dashed lines represent the 3σ confidence level. To conclude, for on-axis bursts, the best-fit Amati relation is $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.495\pm 0.015}$ and the best-fit Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.511\pm 0.016}$. For off-axis bursts, the best-fit Amati relation is $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.282\pm 0.006}$, while the best-fit Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.206\pm 0.006}$.

In the off-axis cases, the best-fit slope is 0.282 ± 0.006 for the E_p − E_iso correlation, which is consistent with our expected range of 1/4 ∼ 4/13 for the off-axis Amati relation. Similarly, on the E_p − L_p plane, the best-fit slope is 0.206 ± 0.006, which also agrees well with our theoretical range of 1/6 ∼ 4/17 for the off-axis Yonetoku relation. Another interesting point in Figs. 4c and d is that most of the simulated off-axis GRBs have an E_iso lower than 10⁴⁸ erg and an L_p lower than 10⁴⁶ erg s⁻¹. Only about 5% of our simulated off-axis GRBs have greater E_iso and L_p. We thus argue that off-axis GRBs are usually LLGRBs. It also indicates that the majority of currently observed normal GRBs should be on-axis GRBs. In Fig. 5, all the simulated on-axis and off-axis bursts are plotted in the same panel. Only a very small number of off-axis events are strong enough to be detected. Their distribution obviously deviates from that of the on-axis events on both the E_p − E_iso and E_p − L_p plane.

Fig. 5. Distribution of simulated off-axis bursts as compared with that of on-axis ones. a) All mock bursts on the Ep − Eiso plane. b) All mock bursts on the Ep − Lp plane. The dash-dotted line and the dotted line correspond to the best-fit result for on- and off-axis samples, respectively. The dashed lines represent the corresponding 3σ confidence level. Most of the off-axis bursts are too weak to appear in the two panels here.

Using the mock bursts, we can also test the Ghirlanda relation (Ghirlanda et al. 2004), which links E_p with the collimation-corrected energy of E_γ = (1 − cos θ_jet)E_iso. The distribution of the simulated bursts on the E_p − E_γ plane is shown in Fig. 6. The data points are best fit by a power-law function of $E_{\mathrm{p}}\propto E_{\gamma }^{0.459\pm 0.018}$ for the on-axis case (Fig. 6a), with an intrinsic scatter of 0.263 ± 0.031. This agrees well with the recently updated Ghirlanda relation of $E_{\mathrm{p}}\propto E_{\gamma }^{0.44\pm 0.07}$ for a sample of 55 observed GRBs (Wang et al. 2018). For the off-axis case (Fig. 6b), we find $E_{\mathrm{p}}\propto E_{\gamma }^{0.269\pm 0.008}$, with a larger scatter of 0.272 ± 0.029. This indicates that the Ghirlanda relation may be largely connected with the Amati relation, but compared with the latter, the slope of the Ghirlanda relation is slightly flatter, and the data points are also obviously more scattered. It is difficult task to compare the off-axis Ghirlanda relation with observations because it needs the precise determination of both the beaming angle and the viewing angle for a number of off-axis GRBs, which themselves are rather dim.

Fig. 6. Distribution of the mock GRBs on the Ep − Eγ plane for the on-axis case. Panel a: on-axis case. Panel b: off-axis case. The solid line shows the best-fit result for the simulated data points, and the dashed lines represent the 3σ confidence level. The on-axis Ghirlanda relation is $E_{\mathrm{p}}\propto E_{\gamma }^{0.459\pm 0.018}$, and the off-axis Ghirlanda relation reads $E_{\mathrm{p}}\propto E_{\gamma }^{0.269\pm 0.008}$.

6. Comparison with observations

In Sect. 5 we illustrated the simulated E_p − E_iso and E_p − L_p relations for both on-axis and off-axis cases. We now proceed to compare our results with observational data.

A sample containing 172 long GRBs was used for this purpose, 162 of which were collected from the previous study by Demianski et al. (2017) and Nava et al. (2012). E_p and E_iso parameters are available for all these 162 events. The L_p parameter is available for only 45 GRBs (Nava et al. 2012). Additionally, we collected another 10 bursts that were argued to be possible outliers of the normal Amati relation in previous studies. It has been suggested that they might be off-axis GRBs (Yamazaki et al. 2003b; Ramirez-Ruiz et al. 2005). The observational data of these 10 events are listed in Table 2.

All the observed GRB samples are plotted on the E_p − E_iso and E_p − L_p planes in Fig. 7. Most of the GRBs are well consistent with our on-axis results, thus they both follow the Amati relation and the Yonetoku relation.

Fig. 7. Distribution of all the observed GRBs. Panel a: Ep vs. Eiso. Panel b: Ep vs. Lp. The dash-dotted lines represent our theoretical results for on-axis bursts, and the thick dotted line shows our theoretical results for off-axis bursts. The 3σ confidence levels are marked correspondingly. The solid lines are the best-fit results for normal GRBs.

On the E_p − E_iso plane (Fig. 7a), four of the ten possible outliers are located between our on- and off-axis lines. It is thus difficult to judge whether these events are on-axis or off-axis GRBs simply from the Amati relation. However, the remaining six bursts can only be matched by our off-axis line. They should indeed be off-axis events. Figures 4c and d clearly show that most off-axis GRBs should be LLGRBs with E_iso < 10⁴⁸ erg, with only a small portion falling in the range of 10⁴⁸ − 10⁵² erg (see also Fig. 5a). This is easy to understand. The input parameters should take some extreme values for an off-axis GRB to be strong enough. Especially these off-axis events should generally have a low θ_obs − θ_jet value (see Eq. (27)), which means that they are still slightly off-axis. In Fig. 7a, only a small number of GRBs are outliers of the on-axis Amati relation. This low event rate is also consistent with our simulation results.

For the E_p − L_p correlation (Fig. 7b), we note that almost all the ten possible outliers clearly deviate from the on-axis Yonetoku relation. From combining the Amati relation and the Yonetoku relation, we therefore argue that the ambiguous sample should also be off-axis GRBs. Furthermore, Fig. 7 seems to indicate that the Yonetoku relation may be a better tool for distinguishing between on- and off-axis GRBs.

7. Conclusions and discussion

The prompt emission of GRBs was studied for both on- and off-axis cases. Especially the three prompt emission parameters E_p, E_iso, and L_p were considered. Their dependence on the input model parameters was obtained via both analytical derivations and numerical simulations, which are consistent with each other. We confirmed that E_p, E_iso, and L_p are independent of θ_obs and θ_jet as long as our line of sight is within the homogeneous jet cone (i.e., the on-axis cases). However, they strongly depend on the value of θ_obs − θ_jet when the jet is viewed off-axis. Additionally, their dependence on γ is very different for on- and off-axis cases, as shown in Eqs. (26) and (27). Through Monte Carlo simulations, we found that the Amati relation is $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.495\pm 0.015}$ in on-axis cases. Correspondingly, the Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.511\pm 0.016}$. On the other hand, in off-axis cases, the Amati relation is $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.282\pm 0.006}$ and the Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.206\pm 0.006}$. The simulated samples were directly compared with the observational samples. We found that they are well consistent with each other in the slopes and intrinsic scatters. The E_p − E_γ relation was also tested and was found to be $E_{\mathrm{p}}\propto E_{\gamma }^{0.459\pm 0.018}$ for on-axis GRBs and $E_{\mathrm{p}}\propto E_{\gamma }^{0.269\pm 0.008}$ for off-axis GRBs.

We have focused on the empirical correlations of long GRBs in the above analysis. Here we present some discussion of short GRBs. Generally speaking, long GRBs usually contain many pulses (tens or even up to hundreds) in their light curves, while short GRBs contain much fewer pulses. For short GRBs, we can therefore perform a similar analysis, but only when the number of pulses considered is reduced significantly. In our model, the total energy release (E_iso) is proportional to the number of pulses, while the peak luminosity (L_p) and the spectral parameter E_p are nearly irrelevant. In other words, a short GRB will have a much smaller E_iso, but its L_p and E_p are largely unchanged. As a result, short GRBs will follow a parallel track on the E_p − E_iso plane as compared with long GRBs, the only difference is that their E_iso are significantly smaller. Alternatively speaking, the power-law index of the Amati relation should be the same for short GRBs and long GRBs, consistent with currently available observations (Zhang et al. 2009, 2018). As for the Yonetoku relation on the E_p − L_p plane, we argue that it should be almost identical for both long and short GRBs. However, the number of currently available short GRBs that is suitable for exploring the Yonetoku relation (i.e., for which E_p and L_p are measured) in detail is still too small.

A simple top-hat jet geometry was adopted in this study. It should be noted that the GRB outflow might be a structured jet (Ioka & Nakamura 2018). A recent study suggested that E_p and the flux will not change much, regardless of the choice of a top-hat jet or a not-too-complicated structured jet (Farinelli et al. 2021). However, in some more complicated scenarios, structured jets may further have an angle-dependent energy density (Lamb et al. 2021). They may be choked jets or even in jet-cocoon systems (Ioka & Nakamura 2018; Mooley et al. 2018; Troja et al. 2019). The empirical correlations of GRBs in these situations are beyond the scope of this study and could be further considered in the future.

In our calculations, the prompt emission was assumed to be due to nonthermal radiation from internal shocks. Other mechanisms such as photosphere radiation and magnetic reconnection above the photosphere may also contribute to the observed flux of GRBs. Generally, there are two types of photosphere models: nondissipative photosphere models, and dissipative photosphere models. The former can produce a narrow quasi-Planckian component and could account for the thermal component observed in some GRBs (Pe’er 2008; Beloborodov 2011), although this component does not play a dominant role in most GRBs (Guiriec et al. 2011; Axelsson et al. 2012). In the dissipative photosphere models, various subphotosphere dissipation processes (Rees & Mészáros 2005), especially Comptonization, may affect the spectrum (Thompson 1994; Pe’er et al. 2006; Pe’er & Ryde 2011; Veres & Mészáros 2012), which may lead to a nonthermal spectrum. Some authors have even used this model to explain the Band spectrum of GRBs (Thompson 1994; Beloborodov 2013; Lundman et al. 2013). Reconnection above the photosphere has also been considered, especially in the view of the ICMART model, focusing on Poynting-flux-dominated outflows (Zhang & Yan 2011; Zhang & Zhang 2014). It is still unclear how commonly photosphere radiation and magnetic reconnection are involved in the prompt phase of GRBs. If these radiation mechanisms are included, the derivation will become much more complicated, but it deserves a trial in the future.

The Amati and Yonetoku relations are essential for us to understand the nature of GRBs. They can also help us probe the high-redshift universe (Xu et al. 2021; Hu et al. 2021; Zhao & Xia 2022; Jia et al. 2022; Deng et al. 2023). Our study shows that they are due to the ultra-relativistic effect of highly collimated jets. To be more specific, they mainly result from the variation in the Lorentz factor in different events. It may provide useful insights for a better understanding of these empirical correlations.

Acknowledgments

We would like to thank the anonymous referee for helpful suggestions. This study is supported by the National Natural Science Foundation of China (Grant Nos. 12233002, 12041306, 12147103, U1938201, U2031118, 12273113, 11903019, 11833003), by the National Key R&D Program of China (2021YFA0718500), by National SKA Program of China No. 2020SKA0120300, and by the Youth Innovations and Talents Project of Shandong Provincial Colleges and Universities (Grant No. 201909118).

References

All Tables

All Figures

	Fig. 1. Geometry of our jet model. The upper panel shows the case for an on-axis jet, where θjet ≥ θobs. The lower panel shows an off-axis jet with θjet < θobs.
In the text

	Fig. 2. Time-integrated spectrum for different γ values (marked in the figure). Panel a: Spectrum of on-axis cases. Panel b: Spectrum of off-axis cases. The black star on each curve marks the position of the peak energy (Ep).
In the text

Fig. 3. Numerical results showing the dependence of the prompt emission features on the input parameters. The calculated prompt emission parameters (Ep, Eiso, and Lp) are plotted against the input parameters (γ, B′, ${\gamma }_{\mathrm{inj}}^{\prime }$, |θobs − θjet|, and p) for both the on-axis cases (with cross symbols) and the off-axis cases (with circle symbols). For comparison, the power-law indices derived from our analytical solutions are marked correspondingly in each panel.

In the text

Fig. 4. Monte Carlo simulation results. (a,b) Distribution of simulated bursts on the Ep − Eiso plane and the Ep − Lp plane for the on-axis cases. (c,d) Distribution of simulated bursts on the Ep − Eiso plane and the Ep − Lp plane for off-axis cases. The solid lines show the best-fitting results for the simulated data points, and the dashed lines represent the 3σ confidence level. To conclude, for on-axis bursts, the best-fit Amati relation is $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.495\pm 0.015}$ and the best-fit Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.511\pm 0.016}$. For off-axis bursts, the best-fit Amati relation is $E_{\mathrm{p}}\propto E_{\mathrm{iso}}^{0.282\pm 0.006}$, while the best-fit Yonetoku relation is $E_{\mathrm{p}}\propto L_{\mathrm{p}}^{0.206\pm 0.006}$.

In the text

Fig. 5. Distribution of simulated off-axis bursts as compared with that of on-axis ones. a) All mock bursts on the Ep − Eiso plane. b) All mock bursts on the Ep − Lp plane. The dash-dotted line and the dotted line correspond to the best-fit result for on- and off-axis samples, respectively. The dashed lines represent the corresponding 3σ confidence level. Most of the off-axis bursts are too weak to appear in the two panels here.

In the text

Fig. 6. Distribution of the mock GRBs on the Ep − Eγ plane for the on-axis case. Panel a: on-axis case. Panel b: off-axis case. The solid line shows the best-fit result for the simulated data points, and the dashed lines represent the 3σ confidence level. The on-axis Ghirlanda relation is $E_{\mathrm{p}}\propto E_{\gamma }^{0.459\pm 0.018}$, and the off-axis Ghirlanda relation reads $E_{\mathrm{p}}\propto E_{\gamma }^{0.269\pm 0.008}$.

In the text

	Fig. 7. Distribution of all the observed GRBs. Panel a: Ep vs. Eiso. Panel b: Ep vs. Lp. The dash-dotted lines represent our theoretical results for on-axis bursts, and the thick dotted line shows our theoretical results for off-axis bursts. The 3σ confidence levels are marked correspondingly. The solid lines are the best-fit results for normal GRBs.
In the text