# Intrapulse Spectral Evolution in Photospheric Gamma-Ray Bursts

Filip Alamaa¹,²

¹ Sorbonne Université, CNRS, UMR 7095, Institut d'Astrophysique de Paris (IAP), 98 bis boulevard Arago, 75014 Paris, France; filip.alamaa@iap.fr

² Department of Physics and The Oskar Klein Centre, KTH Royal Institute of Technology, SE-10691 Stockholm, Sweden

Received 2024 April 29; revised 2024 June 28; accepted 2024 July 1; published 2024 September 12

## Abstract

Photons that decouple from a relativistic jet do so over a range of radii, leading to a spreading in arrival times at the observer. Therefore, changes to the comoving photon distribution across the decoupling zone are encoded in the emitted signal. In this paper, we study such spectral evolution occurring across a pulse. We track the radiation from the deep subphotospheric regions all the way to the observed time-resolved signal, accounting for emission at various angles and radii. We assume a simple power-law photon spectrum injection over a range of optical depths and let the photons interact with the local plasma. At high optical depths, we find that the radiation exists in one of three characteristic regimes, two of which exhibit a high-energy power law. Depending on the nature of the injection, this power law can persist to low optical depths and manifest itself during the rise time of the pulse with a spectral index $\beta \approx \alpha - 1$, where $\alpha$ is the low-energy spectral index. The results are given in the context of a gamma-ray burst jet, but are general to optically thick, relativistic outflows.

*Unified Astronomy Thesaurus concepts: Gamma-ray bursts (629); Radiative processes (2055)*

## 1. Introduction

Gamma-ray burst (GRB) jets are initially optically thick. As the jet expands, the density drops and any trapped radiation can start leaking out toward the observer. The distance from the central engine where the ejecta transitions from optically thick to optically thin is called the photospheric radius, $R_{ph}$, and the released radiation is called photospheric emission. Photospheric emission is an inevitable part of the fireball model and many magnetic jet acceleration models (Cavallo & Rees 1978; Goodman 1986; Paczyński 1986; Piran 1999; Spruit et al. 2001), and its role in GRB prompt observations has been studied extensively (Daigne & Mochkovitch 2002; Rees & Mészáros 2005; Ryde 2005; Giannios 2006).

Photospheric emission is inherently probabilistic in nature, since it is related to photons experiencing their last scattering. Assuming a spherically symmetric outflow, Abramowicz et al. (1991) showed that the optical depth toward the observer depends on the angle between the radial direction and the observer line of sight, $\theta$. Pe'er (2008) and Beloborodov (2011) constructed a probability function for a photon's last scattering as a function of $\theta$ and radius from the central engine, $r$. They found that photons decouple across a wide range of optical depths. Lundman et al. (2013) expanded upon these works by relaxing the assumption of spherical symmetry, considering the effect of jet structure and arbitrary viewing angle on the observed spectrum.

Photons that decouple from high optical depths outrun the plasma and are the first to arrive at the observer. They are followed by photons that experienced their last scattering at progressively lower optical depths. The spread in arrival times leads to a short-duration pulse in the observer frame. Time-resolved emission across such a pulse was considered in spherical symmetry by Pe'er & Ryde (2011). Time-resolved emission from the photosphere within a structured jet was

studied by Meng et al. (2019). In addition, Meng et al. (2019) also modeled a long-lasting and varying central engine by letting the bulk Lorentz factor and luminosity evolve with time.

All of the aforementioned works have considered the comoving photon distribution to be in thermodynamic equilibrium with the plasma. In this case, the radiation is affected only by adiabatic cooling as the jet propagates outward. This effect, together with the Doppler factor being angle-dependent, leads to a broadening of the observed time-integrated spectrum, which consists of a superposition of many Doppler-boosted comoving spectra emitted at various radii and angles to the line of sight (Beloborodov 2010).

However, a GRB jet is a highly chaotic system and energy dissipation below the photosphere is expected (Rees & Mészáros 2005). The dissipation may be due to radiation-mediated shocks (Levinson & Bromberg 2008; Bromberg et al. 2011; Levinson 2012; Beloborodov 2017; Samuelsson et al. 2022; Samuelsson & Ryde 2023), magnetic reconnection (Spruit et al. 2001; Drenkhahn & Spruit 2002; Giannios 2006), turbulence (Zrake et al. 2019), shearing flows in a structured jet (Ito et al. 2013; Vyas & Pe'er 2023), or nuclear collisions between protons and neutrons (Beloborodov 2010). Dissipation destroys the thermodynamic equilibrium and distorts the photon distribution, which may appear highly nonthermal. When the dissipation ceases, the photon distribution will start to re-establish a new equilibrium. Accounting for dissipation, the spectral change across the decoupling region can be much more dramatic compared to the scenario with adiabatic cooling only.

In this paper, we study the time-resolved signal, accounting for dissipation and thermalization. In contrast to Meng et al. (2019), we do not model the central engine variation. Instead, we focus on a thin slice of the jet with an observed width $\delta r \sim r/\Gamma^2$, where $\Gamma$ is the jet bulk Lorentz factor. Therefore, the observed spectral evolution presented in this work is solely due to changes of the photon distribution within the jet. A complex GRB light curve would consist of radiation from many such slices, where each slice would create a subpulse with its own evolution.

The dissipation is modeled by simply injecting a spectral power-law distribution of photons across a range of optical

<figure><img src="image_1.png" /></figure>

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
depths. We show that under these conditions, the radiation in the deep subphotospheric region exists in one of three characteristic regimes, which we call the slow, the marginally fast, and the fast Compton regime. These regimes have direct analogs in the slow and fast synchrotron-cooling regimes, hence the chosen names. Each regime has a typical spectral shape and we characterize the relation between the spectral indices in each case. Tracking the photon distribution as a function of optical depth, we generate the time-resolved signal in the observer frame, accounting for contributions from different angles and radii. We show that dissipation at low optical depths can lead to a significant intrapulse spectral evolution in the observer frame.

Throughout this paper, we use the term thermalization to indicate that the radiation tends toward a kinetic equilibrium when interacting with the plasma, with high-energy photons losing energy and low-energy photons gaining energy. Furthermore, the terms downscattering and upscattering are used to describe photons losing or gaining energy in a scattering event, respectively.

The paper is structured as follows. In Section 2, we introduce the general picture. We obtain the radiated energy as a function of optical depth and describe the model used in the paper. In Sections 3 to 5, we follow the evolution of the photon distribution on its journey, from the deep subphotospheric regions (Section 3), across the decoupling zone (Section 4), to the observed signal (Section 5). We discuss our results, with a specific emphasis on assumptions, in Section 6, and we conclude in Section 7.

## 2. General Picture

### 2.1. Radiated Energy as Function of Optical Depth

The radiated energy as a function of optical depth is a combination of the probability of the last scattering at a given angle and optical depth, as well as the comoving radiation energy density at that point. When there is no energy dissipation, only two things are needed obtain the radiated energy as a function of optical depth: (1) a photon-decoupling probability function; and (2) a cooling function that describes the energy losses due to adiabatic expansion.

Throughout the paper, we assume the outflow to be spherically symmetric, i.e., that we are observing the jet on-axis and that the jet properties do not vary with angle within ~1/Γ to the line of sight. In the case of spherical symmetry, the probability function for a photon to decouple as a function of optical depth and angle was derived using radiative transfer in the ultrarelativistic regime by Beloborodov (2011). In this work, we adopt the expression from Appendix C in Samuelsson & Ryde (2023), which is identical to the expression given in Beloborodov (2011), but rewritten as a function of optical depth, τ. It reads

$$f(\tau, \mu') = \frac{1}{4} \left\{ \frac{3}{2} + \frac{1}{\pi} \arctan \left[ \frac{1}{3}(\tau - \tau^{-1}) \right] \right\} \times \exp \left[ -\frac{\tau}{6} \left( 3 + \frac{1-\mu'}{1+\mu'} \right) \right]. \quad (1)$$

Here, $\mu' = \cos(\theta')$, where $\theta'$ is the angle between the radial direction and the line of sight as measured in the frame comoving with the outflow.

The optical depth in Equation (1) is measured for a fluid element along the radial direction toward infinity. If the

electron density decreases with radius as $r^{-2}$, as expected in the GRB coasting phase, then $\tau$ is related to radius as (e.g., Beloborodov 2011)

$$\tau = \frac{R_{\text{ph}}}{r}, \qquad R_{\text{ph}} = \frac{L\sigma_{\text{T}}}{4\pi m_{\text{p}} c^3 \Gamma^3}, \quad (2)$$

where L is the isotropic equivalent luminosity and $\sigma_{\text{T}}$ is the Thomson cross section. Equation (2) assumes a negligible amount of pairs at the photosphere and that all magnetic energy initially present in the jet has been converted into kinetic energy below the photosphere.

Under the assumption of no dissipation, the total energy of the comoving radiation in the subphotospheric region of the jet is only affected by adiabatic cooling. Samuelsson & Ryde (2023) found a “photon-cooling function,” $\phi(\tau)$, which well described the energy loss at high and low optical depths, including the transition across the photosphere. The cooling function is given by

$$\phi(\tau) = \frac{\tau^{2/3} + 0.2}{1.2}, \quad (3)$$

which was found to give very good agreement with simulations. The average photon energy at an arbitrary optical depth, $\bar{\epsilon}(\tau)$, is related to the average photon energy at the photosphere, $\bar{\epsilon}_{\text{ph}}$, as $\bar{\epsilon}(\tau) = \phi(\tau)\bar{\epsilon}_{\text{ph}}$. At high optical depths, $\phi(\tau)$ corresponds to an ideal adiabatic cooling of the photon distribution as $\propto \bar{\tau}^{-2/3}$, expected in the coasting phase of GRBs from the conservation of entropy (e.g., Pe'er 2015). At optical depths $\tau \ll 1$, all photons are free-streaming and $\phi(\tau)$ is constant.

With the photon-cooling function, one can calculate the normalized radiated energy as a function of $\tau$ and $\mu'$ in a spherically symmetric outflow as

$$\frac{E_{\tau, \mu'}}{E} = \frac{f(\tau, \mu') \phi(\tau)}{\iiint f \phi \, d\tau d\mu'}. \quad (4)$$

Above, $E_{\tau, \mu'} = dE/d\tau d\mu'$ and $E = \iint E_{\tau, \mu'} d\tau d\mu'$. This quantity, multiplied by $\tau$, is plotted with the solid green lines in Figure 1 for different observing angles, $\theta$. The angle $\theta$ is related to $\mu'$ via $\mu' = (\mu - \beta)/(1 - \mu\beta)$, where $\mu = \cos(\theta)$ and $\beta = \sqrt{1 - 1/\Gamma^2}$. The multiplication with $\tau$ ensures that the curve height indicates where most energy is radiated.

The angle-integrated normalized radiated energy can be calculated as

$$\frac{E_{\tau}}{E} = \frac{\iint f(\tau, \mu') \phi(\tau) d\mu'}{\iiint f \phi \, d\tau d\mu'}, \quad (5)$$

where $E_{\tau} = dE/d\tau$. This, again multiplied by $\tau$, is shown by a dashed gray line in Figure 1. Since scatterings do not change the energy in the photon bath if the electrons are kept at the Compton temperature, Figure 1 is always accurate as long as dissipation is absent.

It is clear from the figure that photons start free streaming at different depths. Additionally, photons reaching the observer from larger angles to the line of sight are preferentially emitted farther out and contribute less to the total observed emission. Both these results are in agreement with previous works (Pe'er 2008; Beloborodov 2011; Bégué et al. 2013). It is interesting to note that the radiated energy peaks already at $\tau \sim 3$.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Normalized radiated energy as a function of optical depth and angle in a spherically symmetric outflow. The different solid lines show the emission at different angles compared to the line of sight calculated using Equation (4), and the gray dashed line shows the angle-integrated radiated energy calculated using Equation (5). The quantity $E_i$ indicates $dE/d\tau d\mu'$ and $dE/d\tau$ for the solid lines and dashed line, respectively. The shaded region shows where two-thirds of the energy is radiated. One can see that photons are emitted over a wide range of optical depths, and that the radiated energy peaks at $\tau \sim 3$. The figure is generated with $\Gamma = 100$, but is essentially identical as long as $\Gamma > 10$. A similar figure for the photon number can be found in Beloborodov (2011).</figcaption></figure>

## 2.2. Jet Profile

The scenario envisioned in the current paper is given in the schematic in Figure 2. We assume that some dissipation process generates a power-law distribution of photons as a function of energy. The assumption of a power-law spectral distribution is discussed in more detail in Section 6.1.1. These photons are continuously injected into a region that we call the interaction region. We assume a single-zone approximation for the interaction region, implying that the whole region is causally connected. The injection continues throughout the injection zone (red), which spans an optical depth $\Delta\tau_{\text{inj}} = \tau_i - \tau_f$.

All injected photons are accumulated in the interaction region. Furthermore, the interaction region contains only the injected photons, i.e., there are no photons in the interaction region at optical depths $\tau > \tau_i$. This scenario mimics that of a region downstream of a subphotospheric shock occurring at $\tau_i$. In the case of a shock, thermal photons from the upstream are energized as they traverse the shock region, after which they are advected downstream. Thus, all photons in the downstream region have passed through the shock and the downstream itself did not exist before the shock was initiated. The photons in the interaction region interact with the local plasma through thermal Comptonization. Meanwhile, the whole ejecta is moving outward toward the photosphere (fluid elements move from left to right in the cartoon).

Photons have a nonzero probability of decoupling anywhere in the jet, as is evident from Equation (1). However, ~99% of the energy is radiated in the “decoupling zone” (blue), at optical depths $\tau \le \tau_{dc}$, where $\tau_{dc} \approx 10$. The decoupling zone, therefore, spans an optical depth of $\Delta\tau_{dc} \approx 10$. In an ultrarelativistic outflow, the average number of scatterings per photon traveling from r to infinity is roughly equal to the optical depth at r. Thus, on average, photons experience several scatterings while traversing the decoupling zone. Note that we introduce the notion of a decoupling zone to facilitate the

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Schematic showing the scenario envisioned in this paper. The plasma moves from left to right, passing through each zone. In the injection zone between $\tau_i$ and $\tau_f$, some dissipation event continuously injects a power-law distribution of photons. The injected radiation interacts with the local plasma. Photons start leaking out of the plasma in the decoupling zone, at $\tau < \tau_{dc}$. If injection occurs within, or close to, the decoupling zone, the observed spectral shape may change significantly with time.</figcaption></figure>

discussion. However, all results presented are obtained using the probability distribution in Equation (1).

In the following sections, we will follow the evolution of the photon distribution in the interaction region as it traverses the GRB jet.

## 3. Spectral Evolution across the Injection Zone

Many authors have made detailed calculations regarding the evolution of the subphotospheric photon distribution in the presence of various dissipation mechanisms (e.g., Beloborodov 2010; Vurm & Beloborodov 2016; Parsotan et al. 2018; Samuelsson et al. 2022). In this paper, we consider a simple one-zone model, which is continuously injected with photons in a power-law spectral distribution between the energies $\epsilon_{min}$ and $\epsilon_{max}$. Photon energies are measured in units of $m_e c^2$ throughout the paper. Furthermore, we assume that $\epsilon_{max} \le 1$, such that we can neglect relativistic effects in the comoving frame. The slope of the injected spectrum in photon number is $\alpha_{\text{inj}}$, i.e., $N_{\epsilon} \propto \epsilon^{\alpha_{\text{inj}}}$, where $N_{\epsilon}$ is the specific photon number density. We show that under these simplified assumptions, the photon distribution can exist in three characteristic regimes. In this section, we assume that the injection is always active.

### 3.1. Slow, Fast, and Marginally Fast Compton Regimes

In the subphotospheric region of a GRB jet, the photons vastly outnumber the electrons once the initial pairs have recombined (Bromberg et al. 2011).<fnref n="3" /> Under such conditions, the electron-scattering time is much shorter compared to the photon-scattering time. This implies that to a good approximation, the electrons are always kept in a thermal distribution at the local Compton temperature, $\Theta_C$. The Compton temperature is defined as the temperature that satisfies no net energy exchange between the photon and electron populations. Note that due to the difference in scattering length, the photon distribution need not be thermal.

In Compton scattering between an isotropic photon field and thermal electrons, the average relative energy change per scattering is given by (Rybicki & Lightman 1979)

$$ \frac{\Delta\epsilon}{\epsilon} = 4\Theta_C - \epsilon, \qquad (6) $$

<footnote n="3">³ This would not be true in a cold, magnetically accelerated jet.</footnote>
where $\Delta\epsilon$ is the photon energy change in the scattering, $\epsilon$ is the photon energy, and $\Theta_C$ is the electron temperature, all measured in units of electron rest-mass energy $m_e c^2$. Equation (6) is valid for photons with energies $\epsilon \ll 1$. It is clear that high-energy photons lose energy more quickly than low-energy photons gain energy. From the equation above, one can derive that all photons with energies $\epsilon > 1/N_{sc}$ lose at least half of their energy after $N_{sc}$ number of scatterings. As the number of scatterings is roughly equal to the optical depth traversed, all photons with energies $\epsilon > 1/\Delta\tau$ lose a significant amount of energy over an optical depth $\Delta\tau$. The down-scattering process continues until $1/\Delta\tau$ becomes comparable to $4\Theta_C$. This leads to three different possible regimes.

**Slow Compton regime:** $4\Theta_C\Delta\tau < 1$. In this case, only the high-energy photons with energies satisfying $\epsilon > 1/\Delta\tau$ have time to be downscattered. No other photons are significantly affected, regardless of whether their energy is above or below $4\Theta_C$. This leads to a spectrum with a power-law slope $\alpha = \alpha_{inj}$ between $\epsilon_{min}$ and $1/\Delta\tau$. At energies higher than $1/\Delta\tau$, the power law steepens due to the additional Compton losses. Since the relative energy loss is proportional to the photon energy, the power-law slope becomes $\beta = \alpha_{inj} - 1$.<fnref n="4" /> This result is valid only when $\alpha_{inj} < -1$, above which no steady-state solution exists (see the further discussion in Section 6.2). A special case of the slow Compton regime is when $\epsilon_{max} < 1/\Delta\tau$, in which case the spectral shape is simply similar to that of the injected spectrum.

The shape of the photon distribution in this regime is completely analogous to the slow-cooling distribution of charged particles downstream of a collisionless shock, hence the chosen name. In the synchrotron case, charged particles are injected as a power law with index $-p$ above the injection Lorentz factor $\gamma_m$. In the slow-cooling regime, particles below the cooling Lorentz factor $\gamma_c > \gamma_m$ do not lose energy efficiently and, thus, the particle distribution between $\gamma_m$ and $\gamma_c$ has a slope $-p$. Above $\gamma_c$, particles are cooled due to synchrotron losses and the energy-loss timescale is inversely proportional to the particle energy. This steepens the power-law slope above $\gamma_c$ to $-p-1$ (Sari et al. 1998). The same thing occurs here, excepts that it is the photons that are the primary particle population. The break occurs at $\epsilon_c = 1/\Delta\tau$.

**Fast Compton regime:** $4\Theta_C\Delta\tau > 1$. In this case, all high-energy photons have time to downscatter, leading to a pileup around $4\Theta_C$. Furthermore, the Compton y-parameter for the low-energy photons, $y \equiv 4\Theta_C\Delta\tau$, is larger than unity. This implies that all low-energy photons will significantly increase their energy, leading to a hardening of the low-energy part. If the difference between $4\Theta_C$ and $\epsilon_{max}$ is large enough, a power-law segment with $\beta = \alpha_{inj} - 1$ develops in the region $4\Theta_C \ll \epsilon \ll \epsilon_{max}$. The analogy between this scenario and the fast-cooling synchrotron scenario is less clear here, since some photons gain energy in the fast Compton regime, while all charged particles in the synchrotron fast-cooling regime lose energy.

**Marginally fast Compton regime:** $4\Theta_C\Delta\tau \sim 1$. In the intermediate regime, we have $y \sim 1$, meaning that the low-energy photons are marginally affected. Furthermore, a significant pileup of high-energy photons is not yet visible

around $4\Theta_C$. In this case, the spectrum consists of a power law with slope $\alpha \sim \alpha_{inj}$ between $\epsilon_{min}$ and $4\Theta_C \sim 1/\Delta\tau$, a possible minor pileup above $4\Theta_C$, and a power law with slope $\beta = \alpha_{inj} - 1$ between $\ge 4\Theta_C$ and $\epsilon_{max}$.

In Figure 3, on the left-hand side, we show schematics of the three spectra in the different regimes. Below $\epsilon_{min}$, the dotted (dashed) line shows a Rayleigh–Jeans slope (Wien slope). Note that which spectral regime is the relevant one in a given context depends partly on $\Delta\tau$, but also on $\alpha_{inj}$, since $\Theta_C$ depends on the slope of the injected power-law spectrum.

## 3.2. The Influence of Adiabatic Cooling

The discussion in Section 3.1 is valid in a plane-parallel geometry and neglects the influence of adiabatic cooling on the spectrum. The characteristic timescale for adiabatic cooling corresponds to one doubling of the radius, which equals one halving of the optical depth. Thus, if a high-energy photon is injected at an optical depth $\tau_{inj}$, it suffers Compton losses until its energy is $\epsilon_c = 2/\tau_{inj}$ (assuming $2/\tau_{inj} > 4\Theta_C$). Thereafter, Compton losses are negligible with respect to adiabatic cooling. This leads to the power-law spectrum injected at $\tau_{inj}$ having a high-energy cutoff at the photosphere of $\epsilon_{ph}^{\text{cut}} = \epsilon_c \frac{\phi^{(1)}}{\phi(\tau_{inj}/2)}$, which is approximately

$$ \epsilon_{ph}^{\text{cut}} \sim \left( \frac{\tau_{inj}}{2} \right)^{-5/3}, \quad 4\Theta_C < \frac{2}{\tau_{inj}}. \quad (7) $$

Similarly, the Compton y-parameter for a low-energy photon injected at $\tau_{inj}$ should be evaluated as $y \sim 4\Theta_C \frac{\tau_{inj}}{2}$.

## 3.3. Spectral Shape as Function of Optical Depth

The spectral shape as a function of optical depth depends on the nature of the injection. In this paper, we assume the photon number injection to be either constant as a function of radius, $dN_{in}/dr = \text{const.}$, or constant as a function of optical depth, $dN_{in}/d\tau \propto r^2 dN_{in}/dr = \text{const.}$, where $N_{in}$ is the number of injected photons (see Equation (9) below).

If the photon number injection is constant with radius, then at least half of the photons in the interaction region were injected during the last doubling of the radius, i.e., between $2\tau$ and $\tau$. Photons injected at $\tau_{inj} \gg \tau$ are relatively few in number and contribute only with higher-order effects to the spectrum. The current spectral regime can be estimated by evaluating $4\Theta_C(\tau)\Delta\tau$ and comparing it to unity, as in Section 3.1. Note that $\Delta\tau = \tau$ in this case.

If instead the photon number injection is constant with optical depth, at least half of all photons were injected between $\tau_i$ and $\tau_i/2$. Thus, the overall spectral regime can be estimated by the value of $4\Theta_C(\tau_i/2)\Delta\tau$, where $\Delta\tau = \tau_i/2$. If the photon distribution is in the slow or marginally fast Compton regime, i.e., $4\Theta_C(\tau_i/2)\Delta\tau \lesssim 1$, then energy transfer via Compton scatterings has no large effect on the photon distribution at optical depths $\tau < \tau_i/2$. With few photons injected and little thermalization, the photon distribution is dominated by adiabatic cooling below $\tau_i/2$. Since the energy loss for each photon during adiabatic cooling is proportional to its current energy, the spectral shape (in a log-log plot) remains ~constant. In the fast Compton regime, the spectral shape may continue to evolve, depending on the value of the Compton temperature at $\tau = \tau_i/2$.

<footnote n="4">⁴ This result can also be obtained by evaluating the steady-state solution to the Kompaneets Equation (8) in the region $4\Theta_C \ll \epsilon \ll \epsilon_{max}$, ignoring adiabatic cooling.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Left: schematics of the slow, marginally fast, and fast Compton regimes discussed in Section 3.1. Right: simulation results in the three regimes. The red cross and the green star show the values of $4\Theta_c$ and $1/\Delta\tau$ at $\tau_f$, respectively. The dashed blue lines show the theoretical slopes, which are determined solely by $\alpha_{\text{inj}}$ in all instances, except for the low-energy slope in the fast Compton case. In the simulations shown by the solid lines, the injection continues over a doubling of the radius. The black dashed lines have identical input parameters, except that the injection has been active for 10 times longer. The solid lines are normalized to 1, while the dashed lines are normalized to 0.5 for clarity. Parameter values are given in Table 1.</figcaption></figure>

### 3.4. Simulating the Three Regimes

To verify the discussion in the previous sections, we make simulations in each of the three regimes. The spectra are obtained as follows.

The evolution of a photon distribution in interaction with thermal electrons is described by the Kompaneets equation. The Kompaneets equation in an ultrarelativistic outflow and modified to account for adiabatic cooling as the jet expands is
<table><caption>**Table 1**
Parameter Values Used for the Right-hand Panels in Figure 3</caption><thead><tr><th>Parameter</th><th>Slow Comp.</th><th>Marg. Comp.</th><th>Fast Comp.</th></tr></thead><tbody><tr><td>τ<sub>i</sub></td><td>200</td><td>40</td><td>200</td></tr><tr><td>τ<sub>f</sub></td><td>100</td><td>20</td><td>100</td></tr><tr><td>α<sub>inj</sub></td><td>-2.3</td><td>-1.9</td><td>-1.8</td></tr><tr><td>C(π̄)</td><td>π̄<sup>-2</sup></td><td>1</td><td>1</td></tr><tr><td></td><td>τ<sub>f</sub> = 10</td><td>τ<sub>i</sub> = 400</td><td>τ<sub>i</sub> = 2000</td></tr></tbody></table>

**Note.** All runs use $\epsilon_{min} = 10^{-4}$ and $\epsilon_{max} = 1$. The bottom row indicates the parameter value change used for the longer-lasting runs, which are shown by the dashed lines in each panel.

given in Vurm & Beloborodov (2016). Neglecting induced emission and assuming constant $\Gamma$, it is given by (Samuelsson et al. 2022):

$$ \frac{\partial}{\partial \bar{r}}(\bar{r}^2 n_{\epsilon}) = s(\epsilon) + \frac{1}{\epsilon^2} \frac{\partial}{\partial \epsilon} \left[ \frac{\epsilon^4}{\bar{r}^2} \left( \Theta_{\text{C}} \frac{\partial(\bar{r}^2 n_{\epsilon})}{\partial \epsilon} + (\bar{r}^2 n_{\epsilon}) \right) + \frac{2}{3} \frac{\epsilon^3 (\bar{r}^2 n_{\epsilon})}{\bar{r}} \right]. \quad (8) $$

Here, $\bar{r}$ is the normalized radius, $\bar{r} = r/R_{ph}$, $n_{\epsilon}$ is the photon occupation number, related to the specific photon number density as $N_{\epsilon} \propto \epsilon^2 n_{\epsilon}$, and $s(\epsilon)$ is a source term. The first term in brackets accounts for the upscattering of low-energy photons, while the second term accounts for the downscattering of high-energy photons. The last term in brackets accounts for an ideal adiabatic cooling of the photon distribution as $\propto \bar{r}^{-2/3}$. The Kompaneets equation is solved using the solver from Chang & Cooper (1970).

The source term injects a power-law distribution of photons as

$$ s(\epsilon) \propto \begin{cases} C(\bar{r}) \epsilon^{p_{\text{inj}}}, & \epsilon_{\text{min}} < \epsilon < \epsilon_{\text{max}}, \\ 0, & \text{otherwise.} \end{cases} \quad (9) $$

The source term is only nonzero at optical depths between $\tau_i$ and $\tau_f$. Note that $p_{inj}$ is related to $\alpha_{inj}$ as $\alpha_{inj} = p_{inj} + 2$. The term $C(\bar{r})$ determines the dependency on $\bar{r}$. If $C(\bar{r}) = 1$, the injection is constant as a function of $r$. If, instead, $C(\bar{r}) = \bar{r}^{-2}$, the injection is constant as a function of $\tau$.

The Compton temperature is given by

$$ 4\Theta_{\text{C}} = \frac{\int \epsilon^4 n_{\epsilon} d\epsilon}{\int \epsilon^3 n_{\epsilon} d\epsilon}, \quad (10) $$

obtained by solving Equation (8) in steady state without injection and adiabatic cooling. Since the Compton temperature is a function of the shape of the photon spectrum, it must be evaluated continuously.

In Figure 3, on the right-hand side, we show the results of six simulation runs, whose parameter values are given in Table 1. For all six simulations, the spectrum is evaluated just at the end of the injection zone, at $\tau_f$.

In each of the three panels, two simulation runs are shown. The solid lines show simulations that have run for one doubling of the radius. The red cross and the green star indicate the values of $4\Theta_{\text{C}}$ and $1/\Delta\tau$, respectively, at the end of the simulation. It is clear that the simulated spectrum in each regime closely resembles the corresponding idealized

schematic shown on the left. Adiabatic cooling is accounted for in the simulation runs, but since the runs only last for one doubling of the radius, its effect is small. The dashed blue lines, with their respective slopes, are shown for clarity. In these cases, the relevant $\Delta\tau$ is $\tau_i - \tau_f$.

To study the effect of longer-lasting injection, we make three additional runs with identical parameter sets, but where the injection has continued for 10 times longer. These runs are shown by the black dashed lines in each panel. In each case, the two simulations have the same global characteristics, verifying the discussion in Section 3.3. However, the detailed spectral shapes are different. Specifically, the relations between the power-law slopes are modified. The relevant $\Delta\tau$ for the simulations with longer-lasting injection are $\Delta\tau = \tau_i/2$ for the slow Compton simulation (injection is proportional to $\tau$) and $\Delta\tau = \tau_f$ for the other two simulations (injection is proportional to $r$).

## 4. Spectral Evolution across the Decoupling Zone

The decoupling zone is the region of optical depth wherein escaping photons contribute significantly to the observed fluence. It is therefore possible to detect spectral change across the decoupling zone in the observed time-resolved signal.<fnref n="5" /> The photon distribution can experience significant thermalization or additional dissipation across the decoupling zone. This leads to three different scenarios.

(1) *No dissipation and no high-energy photons.* Only photons with energies above $2/\Delta\tau_{dc}$ (~100 keV) in the comoving frame have time to downscatter while traversing the decoupling zone (see Sections 3.1 and 3.2). If no such photons exist, and there is no additional dissipation, the comoving radiation is only affected by adiabatic cooling. Note that the Compton y-parameter for the low-energy photons is necessarily less than unity, since $4\Theta_{\text{C}}$ cannot be larger than the maximum photon energy. The observer will see a spectral shape that is constant in time, but that decreases in energy continuously across the pulse. How large the total decrease is depends on the signal-to-background level. If we imagine the observer can distinguish the signal coming from optical depths of $\tau \sim 10$ until $\tau \sim 0.3$, then the observed spectrum would decrease in energy by a factor of $\sim \phi(10)/\phi(0.3) \sim 7$ across the pulse.

(2) *No dissipation but high-energy photons.* In this case, the highest-energy photons have time to downscatter. This leads to a spectral softening with time of the high-energy part of the observed spectrum, shifting the upper cutoff to lower energies. The downscattering is efficient until an optical depth $\tau = \tau_{dc}/2 \sim 5$, after which adiabatic cooling will dominate. The radiated energy as a function of $\tau$ peaks at $\tau_{peak} \sim 3$ (see Figure 1). Thus, we expect a rapid softening of the high-energy part during the rise time of the pulse, which stabilizes and remains ~constant during the decay. The observed upper cutoff energy at the light-curve peak can be estimated as

$$ E_{\text{cut}}^{\text{obs}} \sim \frac{2}{\Delta\tau_{\text{dc}}} \left(\frac{\tau_{\text{dc}}/2}{\tau_{\text{peak}}}\right)^{-2/3} \frac{\Gamma m_{\epsilon} c^2}{1+z} \sim \frac{7 \text{ MeV}}{1+z} \Gamma_2, \quad (11) $$

where $\Gamma_2 \equiv \Gamma/10^2$ and $z$ is the cosmological redshift.

<footnote n="5">⁵ To observe such change requires that the pulse duration is longer than the time resolution of the observing telescope.</footnote>
<table><caption>**Table 2**
Parameter Values Used in Figure 4</caption><thead><tr><th>Parameter</th><th>Run A</th><th>Run B</th></tr></thead><tbody><tr><td>τ<sub>i</sub></td><td>50</td><td>50</td></tr><tr><td>τ<sub>f</sub></td><td>10</td><td>3</td></tr><tr><td>α<sub>inj</sub></td><td>-1.9</td><td>-1.6</td></tr><tr><td>C(¯r)</td><td>¯r<sup>-2</sup></td><td>1</td></tr><tr><td>z</td><td>2</td><td>0</td></tr><tr><td>t<sub>E</sub>/t<sub>dyn</sub></td><td>1</td><td>1</td></tr></tbody></table>

Note. Both runs use $\epsilon_{min} = 10^{-4}$ and $\epsilon_{max} = 1$.

If the Compton temperature is exceptionally high, a slight hardening of the low-energy spectrum could be observed. In addition, we expect the previous behavior of the hard-to-soft evolution of the peak energy due to the adiabatic cooling of the spectrum as a whole. This scenario is demonstrated by Run A in Section 5.

(3) Additional dissipation in the decoupling zone. In this case, the evolution of the observed spectrum will depend on the nature of the dissipation. If the injection is constant with radius, i.e., $C(\bar{r}) = 1$ in the current framework, then the photon distribution will increasingly start to resemble the injected spectrum at smaller optical depths. This is because the total photon number in the interaction region doubles with each doubling of the radius, while the number of scatterings in between each doubling halves. At some point, the radiation is likely to enter the slow-cooling Compton regime even if the injected power-law slope is hard, and when $\tau < 2/\epsilon_{max}$, not even the highest-energy photons have time to downscatter. This scenario is demonstrated by Run B in Section 5 (although $\epsilon_{max} = 1$ in Run B, so there are always some photons being downscattered).

If the dissipation is proportional to $\tau$, the magnitude of the spectral change in the observer frame depends on when the injection started. If $\tau_i \gg \tau_{dc}$, only a small fraction of the total photon number is injected in the decoupling zone and the situation will resemble (1) or (2). Note that in the case of additional dissipation at low optical depths, the shape of the radiated energy as a function of $\tau$ given in Figure 1 will be modified.

## 5. Evolution of the Observed Spectrum

In this section, we study the evolution of the observed signal in two simulated pulses, whose parameters are given in Table 2. The modeling of the comoving photon distribution is done in the same way as in Section 3.4. In Run A, the photon injection stops just at the beginning of the decoupling zone. Since $\epsilon_{max}$ is high, Run A illustrates the case of (2) no dissipation but high-energy photons, as discussed in Section 4. In Run B, photon injection continues until $\tau = 3$. Thus, this run illustrates one possible scenario of (3) additional dissipation in the decoupling zone.

How we obtain the observed signal is outlined in the Appendix. In short, the observed signal is obtained using the code Raylease. We assume that the central engine emits a relativistic outflow during an active period between $t = 0$ s and $t = t_E$. A photon that experiences its last scattering at an optical

depth $\tau$ and an angle $\mu'$ reaches the observer at

$$t_{\text{obs}} = (1 + z) \left[ t_e + \frac{R_{\text{ph}}}{c\beta\Gamma^2\tau(1 + \beta\mu')} \right], \quad (12)$$

where $t_e$ is the emission time from the central engine of the shell from which the photon decoupled. The observer time $t_{\text{obs}}$ is given in relation to an imaginary observed photon, emitted at the line of sight, at $t_e = 0$ and $r = 0$.

The observed spectral flux at the observer time $t_{\text{obs}}$, $F_e(t_{\text{obs}})$, is given by Equation (A5). Since we consider a thin part of the jet with a width $\delta r \sim r/\Gamma^2$, it implies that $t_E$ is short. Therefore, it is safe to assume that $\Gamma$, the central engine luminosity, $L$, and the spectral shape are all independent of $t_e$. In this case, Equation (A5) simplifies to Equation (A6), which is what we use in this paper. The photon distribution is tracked from $\tau_i$ until well above the photosphere, and 1000 snapshots of the comoving photon distribution at linear intervals of $\tau$ have been used to generate the observed time-resolved signal.

### 5.1. Intrapulse Spectral Evolution

In the left panels in Figure 4, four snapshot spectra are shown at different intervals of $t_{\text{obs}}/t_{\text{var}}$ for the two runs. Here, $t_{\text{var}} = t_{\text{dyn}}(1 + z)$, and the dynamical time $t_{\text{dyn}}$ is given by

$$t_{\text{dyn}} = \frac{R_{\text{ph}}}{2c\Gamma^2} = 0.2 \text{ s } L_{53} \Gamma_{2}^{-5}, \quad (13)$$

where $L_{53} \equiv L/10^{53} \text{ erg s}^{-1}$. In both runs, we assume that radial spreading of the plasma causes a perceived jet active period of $t_E = t_{\text{dyn}}$ (this assumption is discussed briefly in Section 6.1). This implies a constant width of the considered slice as $\delta r = R_{\text{ph}}/2\Gamma^2$. The black dashed line in both panels shows the time-integrated spectrum, which is generated by first obtaining the spectral fluence via time integration of 100 snapshots of the observed spectrum, evenly log-spaced in time between $t_{\text{obs}}/t_{\text{var}} = 0.03$ and $t_{\text{obs}}/t_{\text{var}} = 50$, and then dividing the integrated fluence with $t_{\text{var}}$. For Run A (top), the spectra are shown in the observer frame with an assumed redshift $z = 2$. In the case of Run B (bottom), the spectra are instead shown in the central engine frame with $L_e = F_e \times 4\pi d_l^2/(1 + z)$, where $d_l$ is the luminosity distance.<fnref n="6" /> The purple shading shows the energy sensitivity of the Fermi Gamma-ray Burst Monitor (GBM), with a darker color corresponding to higher sensitivity (Meegan et al. 2009).

The first photons to reach the observer originate from the leading parts of the slice ($t_e = 0$) along the line of sight ($\mu' = 1$). By comparison with Equations (12) and (13), one obtains when $\Gamma \gg 1$ that photons observed at $t_{\text{obs}}/t_{\text{var}} = 0.1$ decoupled at $\tau = 10$, i.e., just at the beginning of the decoupling zone. For Run A, the injection is proportional to the optical depth and the spectral regime can be evaluated by estimating $4\Theta_C(\tau_i/2)\Delta\tau$, with $\Delta\tau = \tau_i/2$ (see Section 3.3). When the injection is proportional to the radius, as in Run B, the spectral regime at $\tau = 10$ is instead governed by the value of $4\Theta_C(\tau)\Delta\tau$ with $\Delta\tau = 10$. From the simulations, we find $4\Theta_C(\tau_i/2)\Delta\tau = 0.7$ for Run A and $4\Theta_C(\tau)\Delta\tau = 0.9$ for Run B. Thus, in both runs, the radiation is in the marginally fast

<footnote n="6">⁶ The factor $(1+z)^{-1}$ accounts for the fact that we use the luminosity distance together with the spectral flux (e.g., Hogg 1999).</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Observed spectral evolution for Run A (top) and Run B (bottom). Left: four snapshot spectra from early time (yellow) to late time (red). The gray dashed line shows the time-integrated spectrum. The purple shading shows the Fermi-GBM energy window, with darker shading implying higher sensitivity. Right: evolution of the peak energy and the spectral indices as calculated theoretically (dashed lines) and as obtained by fitting a Band function to 30 mock data sets (points with error bars). Each mock data set has been forward-folded through the Fermi-GBM response matrix, as explained in the main text. The flux in each bin is shown by the purple shading. In Run A, the fitted high-energy power law softens with time as comoving photons are downscattered. The peak energy decreases steadily due to adiabatic cooling. In Run B, the peak energy initially increases due to the injection of high-energy photons at low optical depths. The parameter values are given in Table 2.</figcaption></figure>

Compton regime and we expect the relation $\beta \approx \alpha - 1$ at early times.

At later times, the observed signal is a mix of the leading and trailing parts of the slice, emitting at lower and higher optical depths, respectively, as well as emission from various angles. For Run A, the downscattering of the high-energy radiation in the comoving frame leads to a softening of the high-energy spectrum with time. At the peak of the light curve, the spectrum exhibits a hardening at ~2 MeV. This corresponds to the cutoff energy estimated in Equation (11). However, the contribution from trailing parts at higher optical depths transforms the cutoff in the comoving frame to a hardening in the observer frame.

For Run B, the additional dissipation in the decoupling zone changes the spectral evolution. As more and more photons are injected at lower optical depths, the comoving photon distribution starts to increasingly resemble the injected distribution (see Section 4). In the case of Run B, this leads to a rising peak energy, which only starts to decrease due to adiabatic cooling once injection ceases at $\tau = 3$. The high-energy power law also hardens and the spectral peak broadens during the light-curve rise.

Apart from the parameter values given in Table 2, the simulations assume $L = 10^{53}$ erg s$^{-1}$, $\Gamma = 100$, $R_0 = 10^{10}$ cm. and $\Gamma_0 = 4$, where $R_0$ and $\Gamma_0$ are the radius and Lorentz factor at the base of the jet, respectively. Note that these four parameters only affect the absolute value of the luminosity and the variability time, none of which change the observed spectral shape as a function of $t_{obs}/t_{var}$.

## 5.2. Peak Energy and Spectral Indices

To investigate the changing spectrum more quantitatively, in this section we study the evolution of the peak energy, $E_p$, the low-energy spectral index, $\alpha$, and the high-energy spectral index, $\beta$. We will calculate these quantities theoretically, as well as estimate what would be obtained in actual data.

The theoretical peak energy is calculated as the maximum of the $\epsilon F_{\epsilon}$ spectrum at each observer time. The spectral indices are subsequently obtained as

$$ \alpha, \beta = \frac{\log[N_{\epsilon}(\epsilon_{\alpha,\beta} + \delta\epsilon)] - \log[N_{\epsilon}(\epsilon_{\alpha,\beta})]}{\log(\epsilon_{\alpha,\beta} + \delta\epsilon) - \log(\epsilon_{\alpha,\beta})}. \quad (14) $$
Here, $\epsilon_\alpha$ and $\epsilon_\beta$ are the energies at which the spectral indices $\alpha$ and $\beta$ are estimated, respectively, $\delta\epsilon$ is a small increment in energy, and $N_\epsilon = F_\epsilon/\epsilon$. The evolutions of $\alpha, \beta$, and $E_p$ thus calculated are shown in the right-hand side of Figure 4 with the dashed green, red, and yellow lines, respectively. For the figure, we use $\epsilon_\beta = 5E_p$ in Run A, $\epsilon_\beta = 3E_p$ in Run B, and $\epsilon_\alpha = E_p/10$ in both runs.

To obtain how the spectral evolution would be observed by a real telescope, we do as follows. First, we generate 30 spectra as a function of time for both Run A and Run B. To account for detector and background effects, we forward-fold each spectrum through the GBM response matrix. This implies that the fitting procedure can only account for the signal within the energy range of the GBM (8 keV to 40 MeV; Meegan et al. 2009). For each spectrum, a mock data set is generated. The strength of the signal over the background in the mock data is determined by requiring a specific signal-to-noise ratio (S/N) at the peak of the light curve, $S$. In other words, the mock data set in the brightest time bin has $S/N = S$. The S/Ns of all other time bins are given by the relative luminosity at that specific time. Each mock data set is subsequently fitted with a Band function (Band et al. 1993). The Band function describes a smoothly broken power-law photon spectrum using four parameters: $\alpha, \beta, E_p$, and an overall normalization. The analysis is performed in the Multi-Mission Maximum Likelihood framework (Vianello et al. 2015).

The best-fit values obtained for $\alpha, \beta$, and $E_p$ from the Band fits in Run A and Run B are shown by points with error bars in Figure 4. For the two figures, we used $S = 1000$. The very high S/N was chosen so that the spectral evolution is more clearly visible. The error bars indicate $1\sigma$.

The evolution of the parameters confirms the discussion in the subsections above.

**High-energy spectral index, $\beta$.** In both cases, the early spectra have $\beta \approx \alpha - 1$. For Run A, the downscattering of photons in the comoving frame leads to a shortening of the high-energy power-law segment. Part of the power law persists until late times, as is visible from the theoretically calculated value of $\beta$, which stays more or less constant. However, the fitted values of $\beta$ soften quickly during the pulse rise. The high S/N of the mock data set means that the cutoff is “seen” in the fit, making the fitted values of $\beta$ very low. At late times, the Band function is not a good fit for the input spectrum, which is better described at high energies by a high-energy power law plus an exponential cutoff.

In Run B, $\beta$ hardens initially due to the continued injection. The fitted value of $\beta$ remains high, since the high-energy cutoff is outside of the GBM sensitivity window, and, thus, unaccounted for in the fit. After the peak of the pulse, most of the comoving radiation spectra have lost their high-energy tails. This leads to a rapid decrease of the flux above the peak energy in the observed spectrum. When the cutoff enters the GBM energy window at $t_{obs}/t_{var} \sim 1.2$, the fitted value of $\beta$ drops dramatically.

**Peak energy, $E_p$.** Due to the smooth curvature of the Band function, together with the low value of $\alpha$, the fitted value of $E_p$ is always a factor of $\sim 2$ lower than that of the input spectrum in Run A. However, the general evolution is the same: with no injection of new photons in the decoupling zone, $E_p$ decreases continuously across the pulse. The rate at which the peak energy decreases stalls somewhat between $t_{obs}/t_{var} = 0.7$ and $t_{obs}/t_{var} = 1.2$. The duration of the stalling is related to $t_E$, with

higher values of $t_E$ corresponding to longer-lasting stalling (this is true when the properties of the central engine are constant).

For Run B, the continued injection of photons in the decoupling zone means that $E_p$ increases from its initial value. This spectral change is rapid enough to overcome the decrease of the peak energy due to adiabatic cooling. The peak energy starts decreasing at $t_{obs}/t_{var} \sim 1.33$, which corresponds to line-of-sight emission from the trailing parts of the slice at $\tau_f = 3$. When injection stops, adiabatic cooling starts to dominate and the observed peak energy drops.

**Low-energy spectral index, $\alpha$.** For both runs, the theoretical and fitted values of $\alpha$ stay roughly constant throughout the pulse. This implies that the Compton temperature in the decoupling zone is not high enough for any visible upscattering of the low-energy photons.

It is clear from the figure that the best-fit Band parameters follow the theoretical estimates quite well, with the notable exception of $\beta$ in Run A. However, it is important to remember that the fitted values obtained depend on the S/N, redshift, energy window of the current detector, and spectral model used in the fitting. Another key factor is that the duration of the pulse in the observer frame must be longer than the chosen temporal bin size. Otherwise, a time-resolved spectral analysis is impossible.

### 5.3. Light Curve

In Figure 5, the observed light curves for the two pulses are shown. The gray dashed line shows the bolometric light curve, while the blue and red lines show the light curves in the Fermi NaI (8–1000 keV) and the Fermi BGO (0.4–40 MeV) energy ranges, respectively. The pulse profile is asymmetric. The duration of the rise time is related to the value of $t_E$, with smaller $t_E$ producing shorter rise times. Due to adiabatic cooling of the comoving radiation, the BGO light curve decreases quicker than the NaI light curve after the peak.

Apart from the light curve, we also plot the hardness as a function of time. The hardness is calculated as the ratio of the flux between 100 and 300 keV to the flux between 50 and 100 keV. The hardness evolution depends mostly on the position of the peak energy. In Run A, the hardness decreases continuously once the peak energy drops below 300 keV at $t_{obs}/t_{var} \sim 0.25$. In Run B, $E_p$ remains much larger than 300 keV throughout the pulse. Thus, the hardness remains constant after the first early evolution.

The light curve in Run A at very late times ($t_{obs}/t_{var} \gg 1$) is shown in Figure 6. In Figure 6, the light curve in the Swift-XRT energy range (0.3–10 keV; Gehrels et al. 2004) is also shown. When $t_{obs}/t_{var} \gg 1$, the bolometric flux decreases as $F_{bol} \propto t^{-2}$, while the decrease of the peak energy stalls. These results are in line with those presented in Pe’er & Ryde (2011). The very-late-time light curve for Run B is qualitatively similar.

The very-late-time behavior of the pulse is interesting with regard to the early steep decay phase observed in X-rays in many GRBs (Nousek et al. 2006). The evolutions of the luminosity and peak energy in Figure 6 are quite similar to those observed at the end of GRB pulses (Ronchini et al. 2021; Tak et al. 2023). However, given the small values of $t_{var}$ for canonical GRB parameters, the early steep decay likely requires a different explanation in photospheric models of GRBs (Pe’er et al. 2006; Hascoët et al. 2012; Alamaa et al. 2024), unless $\Gamma$ is small.
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Observed light curves for Run A (top) and Run B (bottom). The dashed gray line shows the bolometric light curve, while the blue and red lines show the light curves in the NaI and BGO energy bands, respectively. The dotted-dashed line in the top (bottom) panel shows a scaled version of the light curve in the BGO (NaI) energy band, to more clearly show the difference between the two light curves. Due to adiabatic cooling, the BGO light curve decreases more rapidly after the peak. The hardness as a function of time is shown on the right y-axis.</figcaption></figure>

# 6. Discussion

## 6.1. Assumptions

There are many simplifying assumptions in the current work. In this section, we list a few of them and discuss how they may affect our results.

### 6.1.1. Power-law-injected Spectrum

Simulations of GRB jets find that a lot of energy is dissipated below the photosphere (Lazzati et al. 2009; Gottlieb et al. 2019). Several dissipation mechanisms are plausible and many of them generate a power-law distribution of photons, e.g., bulk Comptonization in a radiation-mediated shock (Ito et al. 2018; Lundman et al. 2018; Samuelsson et al. 2022), shear interactions in a structured jet (Ito et al. 2013; Vyas & Pe’er 2023), pair cascades as a result of nuclear collisions (Beloborodov 2010), and synchrotron emission from charged particles accelerated at a collisionless subshock (Lundman & Beloborodov 2019; Levinson 2020; Rudolph et al. 2024).

<figure><img src="image_7.png" /><figcaption>**Figure 6.** Very-late-time light curve for Run A. Once $t_{obs}/t_{var} \gg 1$, the bolometric flux decreases with time as $t^{-2}$. After the light-curve peak, the peak energy first decreases roughly as $t^{-1}$. However, the decrease stalls around $t_{obs}/t_{var} \sim 100$. Note that due to the different scales on the y-axes, the gray dashed line indicating a $t^{-2}$ decay is valid only for the flux.</figcaption></figure>

Thus, the choice of a power-law injection spectrum is motivated.

However, the assumption of a power law that remains constant in time as the jet evolves is a strong one. If high-energy particles are present, as in the case of pair cascades and synchrotron emission, inverse Compton scattering on the charged particles can generate complicated radiation profiles (Pe’er et al. 2005; Vurm et al. 2013). Highly relativistic radiation-mediated shocks generate pairs and are subject to Klein–Nishina effects (Ito et al. 2018; Lundman et al. 2018). Furthermore, the injection is likely to evolve, specifically at low optical depths (Levinson 2012).

In this paper, we have neglected a lot of physical complexity in favor of a very simple model that more clearly highlights the main points, which are the three spectral regimes at high optical depths and the evolving observed signal. According to the discussion in Section 3.3, the spectral regime of the radiation is roughly determined during one single doubling of the radius, across which the injection may be ~constant. The two examples in Figure 4 clearly show that the observed signal evolves due to changes of the comoving radiation field. That the observed signal is affected by the changing local radiation in the jet is a general result. However, the details regarding how the observed signal is going to evolve will depend on the considered dissipation scenario. Thus, precise predictions for each case require dedicated studies using physical models.

### 6.1.2. Other Assumptions

The photon distribution is modeled using the Kompaneets Equation (8). However, at small optical depths, the radiation field becomes increasingly anisotropic (Beloborodov 2011). Equation (8) then overestimates the energy transfer between photons and electrons (Vurm & Beloborodov 2016), which leads to a quicker thermalization. If one accounts for this, the high-energy power law in Run A would survive for a longer time.

The interaction region was assumed to be devoid of photons at optical depths $\tau > \tau_i$. As mentioned in Section 2.2, this choice was inspired by the region downstream of a
subphotospheric shock. However, in other dissipation scenarios, the energized photons may exist in the same local region as an initial thermal photon population at low temperature. The thermal population would decrease the Compton temperature in the interaction region, which in turn would increase the parameter space for the slow and the marginally fast Compton regimes. Such a situation may be interesting to investigate in a future work.

We ignored Klein–Nishina suppression of the cross section. However, this does not significantly affect the result. The angle-integrated Klein–Nishina cross section is only lower by a factor of ~2 compared to the Thomson cross section at $\epsilon = 1$. Including the Klein–Nishina corrections, the highest-energy photons would arrive earlier.

We assumed a constant $\Gamma$ throughout the jet evolution. In reality, the high internal energy could lead to a re-acceleration of the plasma and a second fireball. If the outflow accelerates in the decoupling zone, the radiation would suffer from greater adiabatic losses. However, the increased Doppler boost would cancel out this effect in the observer frame (Beloborodov 2013) and the evolution of the peak energy would be the same. The observed light curve would be slightly different compared to that shown in Figure 5.

The spectra shown in Figure 3 are very “clean” compared to the complicated radiation profiles found in several other works (e.g., Pe’er et al. 2005; Vurm et al. 2013; Ahlgren et al. 2015). The main reason is the lack of relativistic electrons in the present work. Indeed, since the electrons are always kept at the Compton temperature, and we assumed $\epsilon_{max} \le 1$, the electrons are always nonrelativistic. That $\epsilon_{max} \le 1$ also ensures that there is no pair production. With relativistic electrons, the radiation profile becomes more complicated due to inverse Compton scattering creating bumps in the spectrum. However, relativistic charged particles are difficult to produce below the photosphere in GRBs (Levinson & Nakar 2020) and simpler spectral shapes similar to the ones presented in Figure 3 are, therefore, expected.

We assumed $t_E = t_{dyn}$. This is based on radial spreading of the plasma, which causes even an originally infinitesimally thin shell to have a width of $\delta r \sim R_{ph}/\Gamma^{-2}$ at the photosphere (e.g., Pe’er 2015). With an assumption of constant $\Gamma$, this causes a perceived jet activity period of $t_E \sim t_{dyn}$. However, $t_E$ can be both smaller and larger. If dissipation is due to a subphotospheric shock, compression of the plasma would lead to $t_E \sim t_{dyn}/3$ (Samuelsson & Ryde 2023). If the jet instead varies on timescales that are much longer than $t_{dyn}$, the pulse duration would increase.

Finally, we assumed injection of high-energy photons up to $\epsilon_{max} = 1$. If $\epsilon_{max} \le 0.1$, then no high-energy photons have time to downscatter in the decoupling zone. In this case, observed spectral evolution, apart from adiabatic cooling, requires additional dissipation in the decoupling zone.

## 6.2. $\alpha-\beta$ Relation

The relation $\beta = \alpha - 1$ is a theoretical prediction that applies to radiation in the marginally fast Compton regime. Due to adiabatic cooling in the comoving frame and the observed signal being a superposition of many comoving spectra, this relation is going to be approximate in the observer frame, $\beta \approx \alpha - 1$. The best chance of detecting it is at early times, during the rise time of the pulse.

As mentioned in Section 3.1, the $\alpha-\beta$ relation is valid only when $\alpha_{inj} < -1$. When $\alpha_{inj} \gtrsim -1$, the Compton temperature is so high that the photon distribution in the interaction region quickly forms a pronounced Wien peak around $\Theta_C$.$^7$ This assumes that $\epsilon_{max}$ is high. (indeed, $4\Theta_C \approx \epsilon_{max}/2$ for the injected photon distribution when $\alpha_{inj} = -1$; see Equation (10)).

Thus, a clear prediction is that the $\alpha-\beta$ relation should only be observed when the detected low-energy slope is quite soft, $\alpha \lesssim -1$. Bursts with harder $\alpha$ must be in the fast Compton regime, in which case we expect a cutoff above the peak.⁷ There could still be some signature at higher energies above the cutoff, as in the bottom panels of Figure 3. If such a feature is detected, it would be a clear sign of Comptonization, as the spectral shape in the bottom panels of Figure 3 is quite complex and unique.

On the other hand, when the injection is very soft, $\alpha_{inj} < -2$, the $\nu F_{\nu}$ peak corresponds to $\epsilon_{min}$ (the top two panels in Figure 3). In this scenario, one has an observed spectrum with a very hard $\alpha$, a low $E_{p}$, and a broken power law at high energies, whose indices are related as $\beta_2 = \beta_1 - 1$.

## 6.3. Early-time Spectra

During the pulse decay, emissions at different optical depths and angles reach the observer simultaneously. In contrast, the first observed photons are all emitted within the same local region of the jet. The early-time spectrum is, therefore, the spectrum that most resembles the comoving one, Doppler-boosted into the observer frame. As such, it should be similar to one of the typical spectra in the three Compton regimes in Figure 3, given that dissipation has occurred at low optical depths.

This part of the pulse can be very difficult to detect due to the low count rate. However, if detected, it can give information regarding the physical properties of the comoving radiation, the dynamics of the jet, and the details of subphotospheric dissipation. Furthermore, it allows for observation of the jet at a radius of up to an order of magnitude smaller the photospheric radius.

## 6.4. Pulse Duration

The observed duration is highly sensitive to the value of $\Gamma$ (see Equation (13)). The duration of the light curve shown for Run A in Figure 5 is ~1 s for the chosen parameters. If $\Gamma$ is increased to 300, the observed duration shrinks to only 5 ms, and for $\Gamma = 50$ the duration is ~40 s. Since the Lorentz factor also scales the Doppler boost, pulses with shorter duration should have higher peak energies on average, scaling as $\delta t_{obs} \propto E_p^{-5}$.

Derebi-Bégué et al. (2022) suggested that GRBs that exhibit an X-ray plateau may have low Lorentz factors, in the order of a few tens. Therefore, this sample may be interesting to investigate as the lower Lorentz factor can facilitate a time-resolved analysis.

## 6.5. Dissipation Radius and Central Engine Activity

Imagine a scenario where most of the subphotospheric energy dissipation takes place over a small distance, equal to

<footnote>$^7$ This is true in the current framework, when $\epsilon_{max} \lesssim 1$ and with the current model assumptions. It may not always be the case—see, e.g., Beloborodov (2010) and Vyas & Pe’er (2023).</footnote>
one doubling of the radius or less. In this case, $\Delta\tau = \tau_i - \tau_f \le \tau_f$. As long as $\alpha_{inj} > -2$ and $\epsilon_{max} > 1/\Delta\tau$, the comoving peak energy equals $\epsilon_p \approx \max[4\Theta_C, 1/\Delta\tau]$ at the end of the dissipation zone (adiabatic cooling is small over one doubling of the radius). Since $1/\tau_f < 1/\Delta\tau$, all injected photons have time to be downscattered before the ejecta reaches the photosphere and the observed spectrum is well described by a cutoff power-law function.

As long as $4\Theta_C$ is not much larger than $1/\Delta\tau$ at $\tau_f$, then $\tau_f$ can be estimated from the peak energy of the observed cutoff power-law spectrum as

$$ \tau_f \sim 2 \left[ \frac{5\Gamma}{3(1+z)} \frac{m_e c^2}{E_p} \right]^{3/5} \approx 110 \Gamma_2^{3/5} E_{p,2}^{-3/5} (1+z)^{-3/5}, \quad (15) $$

where we used Equation (7), a Doppler factor of $5\Gamma/3$ (Samuelsson & Ryde 2023), and $E_{p,2} = E_p/(100\,\text{keV})$. With a measured value of $t_{var}$ from the light curve, the final dissipation radius can be estimated as $R_f = 2c\Gamma^2 t_{var}/(1+z)\tau_f$, where we used $R_f = R_{ph}/\tau_f$ and Equation (13).

If the dissipation is due to an internal collision, the final dissipation radius is related to the central engine duration $\Delta t_{eng}$ as $R_f \approx 2c\Gamma^2\Delta t_{eng}/\psi$, where $\psi$, equal to the ratio of the Lorentz factors of the fast- and slow-moving materials, is of the order a few.<fnref n="8" /> Then the central engine duration can be estimated from the observables as

$$ \Delta t_{eng} \sim \frac{t_{var}}{(1+z)\tau_f}, \quad (16) $$

where $\tau_f$ is given by Equation (15).

### 6.6. Comparison to Shock Breakout Emission

Shock breakout is a likely explanation for the emission in low-luminosity GRBs and for the prompt emission in GRB 170817 (e.g., Campana et al. 2006; Waxman et al. 2007; Nakar & Sari 2012; Bromberg et al. 2018; Gottlieb et al. 2018; Beloborodov et al. 2020). Lundman & Beloborodov (2021) simulated the time-resolved observations from the shock breakout in GRB 170817 using a hydrodynamic code coupled to Monte Carlo photons (radshock; Lundman et al. 2018). They also find a time-resolved spectrum whose shape evolves with time. However, the evolution is quite different compared to the two runs studied in this work.

There are at least two reasons for the different evolution. First, their injection is very hard, leading to $4\Theta_C \lesssim \epsilon_{max}$, in which case the high-energy signature is indistinguishable from the cutoff. Second, the dissipation is dynamic. The shock travels in a speed gradient within the dynamical neutron star ejecta and speeds up once it gets close to breakout, which modifies the dissipation.

## 7. Conclusion

In this paper, we have studied how the evolution of the comoving radiation in the region of last scattering in an optically thick relativistic jet affects the time-resolved observed spectrum. As photons leave the plasma over a range of optical

depths (Figure 1), the average photon experiences multiple interactions within the decoupling zone. The photons that decouple at high optical depths are the first to reach the observer. We may therefore expect an evolving observed signal if dissipation or significant thermalization occurs in the decoupling zone.

To model the energy dissipation, we assumed a simple injection of a power-law photon spectrum across a range of optical depths (the injection zone; see Figure 2). All photons were accumulated in an interaction region, where they were allowed to interact with the local plasma. No other photons apart from those injected existed in the interaction region. We assumed that the photons vastly outnumbered the electrons, such that the electrons were always kept in a thermal distribution at the Compton temperature, $\Theta_C$. In this case, we found that the radiation at high optical depths exists in one of three different characteristic regimes, depending on the value of $4\Theta_C\Delta\tau$, where $\Delta\tau$ is the duration of the injection.

In the slow Compton regime ($4\Theta_C\Delta\tau < 1$) and the marginally fast Compton regime ($4\Theta_C\Delta\tau \sim 1$), the spectrum consists of two power laws. Low-energy photons, $\epsilon < 1/\Delta\tau$, do not have time to be significantly affected by Compton scattering. Thus, the slope of the low-energy power law coincides with the slope of the injected spectrum, $\alpha_{inj}$. In contrast, the high-energy photons lose energy efficiently to the electrons, generating a high-energy power law with index $\alpha_{inj} - 1$. In the fast Compton regime, photons pile up around $4\Theta_C$. The different spectral shapes in the three regimes are given in Figure 3.

The photon distribution in the interaction region was tracked from the onset of injection to far beyond the photosphere, accounting for thermalization and adiabatic cooling. The observed signal as a function of time was obtained via Equation (A6), which correctly accounts for the probability of emission at various optical depths and angles to the line of sight. We focused on a thin slice of the jet, whose emission created an observed pulse of typical duration $\delta t_{obs} \sim t_{var} = 0.2\,\text{s}(1+z)L_{53}\Gamma_2^{-5}$. A complex GRB light curve would consist of many such subpulses. We assumed spherical symmetry and a bulk Lorentz factor $\Gamma \gg 1$.

The evolution of the observed signal across the pulse depends on the comoving photon distribution in the decoupling zone, and we classified three different scenarios. In the case of no dissipation and no high-energy photons (comoving energies $\gtrsim 100\,\text{keV}$), the observed spectral shape is constant, but the whole spectrum decreases in energy by a factor of $\sim 7$ across the pulse due to adiabatic cooling. In the case of no dissipation but high-energy photons, the high-energy photons have time to downscatter within the decoupling zone. This leads to a spectral softening of the high-energy index, which is most pronounced during the rise time of the pulse. This scenario is shown in the two top panels in Figure 4. In the case of additional dissipation within the decoupling zone, the observed signal depends on the nature of the dissipation. One such example is shown in the two bottom panels in Figure 4.

Finally, we discussed the importance of early-time spectra. The early spectra consist only of photons originating from the same local region within the jet and they are emitted at a fairly high optical depth, $\tau \sim 10$. If the comoving radiation is in the marginally fast Compton regime at this time, one obtains a relation between the high- and low-energy spectral indices as $\beta \approx \alpha - 1$. These early spectra require bright GRBs to be seen. However, they would be highly illuminating if detected, since

<footnote n="8">⁸ This estimate assumes that the comoving densities in the fast and slow moving material are roughly equal (Samuelsson et al. 2022).</footnote>
they can give us valuable insights into the jet behavior before it becomes optically thin.

Acknowledgments

I thank the anonymous reviewer for many valuable suggestions. F.A. is supported by the Swedish Research Council (Vetenskapsrådet, 2022-00347). I thank Frédéric Daigne, Robert Mochkovitch, Felix Ryde, and Christoffer Lundman for fruitful discussions. This research has made use of the High Energy Astrophysics Science Archive Research Center (HEASARC) Online Service at the NASA/Goddard Space Flight Center (GSFC). In particular, I thank the GBM team for providing the tools and data.

Appendix

Obtaining the Observed Signal

To obtain the time-resolved observed spectrum, two things are needed. The first is the probability distribution for any given photon to escape the outflow as a function of optical depth and angle. In a spherically symmetric, ultrarelativistic outflow, this is given by Equation (1). The second is the comoving specific photon number density at each point in space and time, $N_{\epsilon}$. Under the assumption of spherical symmetry, the photon density is assumed only to vary with optical depth and time, i.e., $N_{\epsilon} = N_{\epsilon}(\tau, t)$.

A.1. Comoving Photon Distribution

If there is no dissipation and no high-energy photons in the decoupling zone, the photon distribution is only affected by adiabatic cooling. In this case, it is enough to know the spectral shape at a single optical depth and how the spectrum is affected by adiabatic cooling as a function of $\tau$ to obtain $N_{\epsilon}$. This is the special case implicitly assumed in Samuelsson & Ryde (2023). In Samuelsson & Ryde (2023), the photon distribution as a function of optical depth was obtained using the spectrum at the photosphere, $N_{\epsilon ph}$, together with the photon-cooling function given in Equation (3).

To obtain the photon distribution as a function of optical depth in the general case, including the effects of scatterings and potential dissipation, it is necessary to track $N_{\epsilon}$ during the jet evolution. The modeling of the photon distribution is done as in Section 3.4. The distribution is followed from the onset of the injection zone until well above the photosphere and snapshot spectra are saved at even intervals.

Equation (8) assumes an ideal adiabatic cooling $\propto r^{-2/3}$. However, the cooling changes close to the photosphere once photons start to decouple. To account for this, we stop the idealized cooling at $\tau_{\text{cool}} = 1.2^{3/2}$ (see Equation (3) and the discussion in Samuelsson & Ryde 2023). This prescription gives the correct average energy at the photosphere. However, it slightly overestimates the amount of cooling experienced when $\tau \gtrsim 1$ and slightly underestimates it when $\tau \lesssim 1$. To correct for this, the energy grid in each saved snapshot is shifted such that the mean photon energy $\bar{\epsilon}(\tau)$ is given by

$$ \bar{\epsilon}(\tau) = \begin{cases} \bar{\epsilon}^*(\tau)\phi(\tau) & \text{if } \tau \le \tau_{\text{cool}}, \\ \bar{\epsilon}^*(\tau)\frac{\phi(\tau)}{(\tau/\tau_{\text{cool}})^{2/3}} & \text{if } \tau > \tau_{\text{cool}}, \end{cases} \quad (A1) $$

where $\bar{\epsilon}^*(\tau)$ is the mean photon energy in the snapshot before the shift. This gives a very small correction to the observed signal.

A.2. Observed Flux

Consider a relativistic outflow consisting of many thin consecutive shells. We assume that each of the embedded shells acts independently and focus on one specific shell with a total number of photons, $N_{\gamma}$. The observed photon number originating from a volume element $dV$ relates to the emitted photon number as (Pe'er & Ryde 2011; Lundman et al. 2013)

$$ dN_{\text{obs}} = \frac{(1+z)^2}{d_l^2} \frac{dN_{\text{em}}}{d\Omega_v} = \frac{D^2(1+z)^2}{d_l^2} N_{\gamma} f(V, \Omega_v') dV. \quad (A2) $$

Here, $d_l$ is the luminosity distance to the source, $dN_{\text{em}}/d\Omega_v$ is the emitted photon number into a solid angle $d\Omega_v$, $D = [\Gamma(1 - \beta\mu)]^{-1} = \Gamma(1 + \beta\mu')$ is the Doppler factor, and $f(V, \Omega_v')$ is the probability function for an arbitrary volume element and viewing angle. The term $(1+z)^2$ is the correction for redshift, which is necessary to include, since we are using the luminosity distance but the equation is for the photon number. The term $D^2$ in the second equality comes from the solid viewing angle transformation $d\Omega_v = d\Omega_v'/D^2$. For a spherically symmetric outflow, this simplifies to

$$ dN_{\text{obs}} = \frac{N_{\gamma}(1+z)^2}{4\pi d_l^2} f(\tau, \mu') d\tau d\mu'. \quad (A3) $$

A photon that makes its last scattering at a time $t$, radius $r$, and angle $\mu$ compared to the radial direction reaches the observer at a time $\tilde{t}_{\text{obs}}/(1+z) = t + (D - r\mu)/c$. Here, $D$ is the distance to the source, which is slightly different compared to the luminosity distance $d_l$. Since the arrival times for all photons are delayed by a factor $(1+z)D/c$, we define a new observer time as $t_{\text{obs}} = \tilde{t}_{\text{obs}} - (1+z)D/c$. With this definition, a photon emitted at a time $t=0$ s directly from the central engine ($r=0$ cm) would reach the observer at $t_{\text{obs}}=0$ s (such a photon is called the trigger photon in Pe'er 2008).

If we assume that all photons in the shell that decouple between an optical depth $\tau$ and $\tau+d\tau$ do so instantaneously when the shell reaches the optical depth $\tau$, then the time dependence can be described by a Dirac $\delta$-function. A shell that was launched from the central engine at a time $t_e$ and that travels with a velocity $v$ reaches an optical depth $\tau$ at time $t = t_e + R_{\text{ph}}/\tau v$, where we used $\tau = R_{\text{ph}}/r$. Using this, and the fact that $(1 - \beta\mu) = [\Gamma^2(1 + \beta\mu')]^{-1}$, one obtains

$$ dN_{\text{obs}}(t_{\text{obs}}) = \frac{N_{\gamma}(1+z)}{4\pi d_l^2} f(\tau, \mu') d\tau d\mu' \\ \times \delta\left(t_e + \frac{R_{\text{ph}}}{c\beta\Gamma^2\tau(1+\beta\mu')} - \frac{t_{\text{obs}}}{1+z}\right). \quad (A4) $$

The total spectral flux in the observer frame from the shell, $F_{\epsilon}(t_{\text{obs}})$, is obtained by multiplying the above formula with $\epsilon'^{N_{\epsilon}'}$ and integrating over the entire emitting volume. Here, $\epsilon' = (1+z)\epsilon/D$ is the photon energy as measured in the comoving frame. Note that the term $N_{\gamma}f(r, \mu')d\tau d\mu$ already ensures a correct number of decoupling photons. Therefore, $N_{\epsilon}$
should be a probability distribution, normalized as
$$ \int \mathcal{N}_{e'} d\epsilon' = 1. $$

Putting it all together, and inserting $N_\gamma$ from Section A.5, we obtain the observed flux from a series of multiple consecutive shells emitted between $t_e = 0$ s and $t = t_E$, assuming spherical symmetry as

$$ F_{e}(t_{\text{obs}}) = \frac{1+z}{4\pi d_{l}^{2}} \int_{0}^{t_{E}} dt_{e} \int_{-1}^{1} d\mu' \int_{0}^{\infty} d\tau \frac{L}{4kT'_{0}\Gamma_{0}} \\ \times f(\tau, \mu') \epsilon' \mathcal{N}_{e'} \delta \left( t_{e} + \frac{R_{\text{ph}}}{c\beta\Gamma^{2}\tau(1+\beta\mu')} - \frac{t_{\text{obs}}}{1+z} \right), \quad (\text{A5}) $$

where $k$ is the Boltzmann constant and $T'_{0}$ and $\Gamma_{0}$ are the temperature (given in Equation (A12)) and Lorentz factor at the base of the jet, respectively.

### A.3. Steady Source

A detailed study regarding the time-resolved observed signal for different jet emission profiles—i.e., different functional shapes for, e.g., $L(t_e)$ and $\Gamma(t_e)$—was performed by Meng et al. (2019). Here, we make the simplifying assumption that the jet properties do not vary in time. We further assume that the steady jet properties lead to a comoving spectrum, $N_\epsilon$, that is steady in time, such that $N_\epsilon = N_\epsilon(\tau)$.

When the central object is a steady source, the coefficient $L/4kT'_{0}\Gamma_{0}$ in Equation (A5) can be evaluated outside the integral. Note that even in this case, one still needs to perform the integral over $dt_e$ because of its dependence in the $\delta$-function. Using the $\delta$-function to remove the integral over $dt_e$, one obtains

$$ F_{e}(t_{\text{obs}}) = \frac{1+z}{4\pi d_{l}^{2}} \frac{L}{4kT'_{0}\Gamma_{0}} \int_{-1}^{1} d\mu' \int_{0}^{\infty} d\tau \\ \times f(\tau, \mu') \epsilon' \mathcal{N}_{e'} H(t_{e,\delta}) H(t_{E} - t_{e,\delta}), \quad (\text{A6}) $$

where $t_{e,\delta}$ is the value of $t_e$ that satisfies the $\delta$-function in Equation (A5):

$$ t_{e,\delta} = \frac{t_{\text{obs}}}{1+z} - \frac{R_{\text{ph}}}{c\beta\Gamma^2\tau(1+\beta\mu')}. \quad (\text{A7}) $$

In Equation (A6), $H$ is the Heaviside step function, which implies that

$$ H(t_{e,\delta})H(t_{E} - t_{e,\delta}) = \begin{cases} 1 & \text{if } 0 \le t_{e,\delta} \le t_{E}, \\ 0 & \text{otherwise.} \end{cases} \quad (\text{A8}) $$

The Heaviside functions are included so that only regions emitted during the jet activity time contribute to the observed flux. Note that Equation (A6) is only valid when $L$ and $\Gamma$ are constant in time, otherwise there will be extra terms when the $\delta$-function is used to remove the integral over $t_e$.

### A.4. Raylease

Raylease is the name of the numerical code used to obtain the observed time-resolved signal. Raylease takes as input $L(t_e)$, $\Gamma(t_e)$, and $N_\epsilon(\tau, t_e)$ and gives as output $F_e(t_{\text{obs}})$, by solving Equation (A5). One also needs to specify $R_0$ and $\Gamma_0$, which in principle could be functions of $t_e$ as well. In the simplified case of constant jet properties considered in this paper, Raylease solves Equation (A6) for specified values of

$L$, $\Gamma$, $R_0$, $\Gamma_0$, and $N_\epsilon(\tau)$. As the photon distribution is unspecified before injection begins in the current framework, the upper limit in the integral over the optical depth in Equation (A6) is set to $\tau_i$.

Currently, Raylease assumes a spherically symmetric jet, but it can quite easily be extended to include jet structure and an arbitrary viewing angle in the future.

### A.5. Photon Number in an Infinitesimal Shell

This section derives the number of photons in each unit solid angle. The derivation assumes that the net diffusion between different solid angles is negligible. The subscript “0” used throughout this section implies that the quantity is measured at the base of the jet, where the internal pressure starts efficient bulk acceleration. For the derivation in this section, the Lorentz factor $\Gamma$ is allowed to vary with radius.

The number of photons per unit solid angle in an infinitesimally thin shell is given by

$$ \frac{dN_{\gamma}}{d\Omega} = \frac{d\dot{N}_{\gamma}}{d\Omega} dt_{e}. \quad (\text{A9}) $$

Assuming a negligible photon production once bulk acceleration of the outflow begins (see below), the number of photons follows the conservation of particle number as $dN_\gamma = n'_\gamma r^2 \Gamma d\Omega dr$, which gives $d\dot{N}_\gamma/d\Omega = n'_\gamma r^2 \Gamma \beta c$, where $\beta c = dr/dt$. Similarly, we obtain $d\dot{M}/d\Omega = m_p n'_p r^2 \Gamma \beta c$, where $n'_p$ is the baryon number density, which is assumed to be dominated by protons.

Define the specific entropy per baryon, $\eta(\Omega)$, such that $dL/d\Omega = \eta(\Omega)c^2 d\dot{M}/d\Omega$, where $dL/d\Omega$ is the outflow power per unit solid angle. The specific enthalpy, $h'$, and the Lorentz factor, $\Gamma$, both of which vary with radius, have a product that is independent of radius, as $h'\Gamma = \eta$ (e.g., Levinson & Nakar 2020). At the base of the jet, the Lorentz factor is small and the entropy is dominated by the pressure from the radiation and the pairs $h'_{0} \approx 4p_{0}/n'_{p,0}m_{p}c^{2}$, where $p_{0}$ is the total pressure. Assuming a thermal equilibrium between the pairs and the radiation, such that they can be described by a common temperature $T'_{0}$, the pressure is given by $p_{0} \approx n'_{q,0}kT'_{0}$, where $n'_{q,0} = n'_{\gamma,0} + 2n'_{\pm,0}$ is the photon plus pair density. As the temperature drops, the pairs all recombine until $n'_{q} = n'_{\gamma}$. Assuming that recombination is the only source of photon production after the early hot phase, the ratio $n'_{q}/n'_{p}$ is conserved along streamlines and equal to $n'_{q,0}/n'_{p,0}$. Hence, we get

$$ \begin{aligned} n'_{\gamma} &= n'_{p} \frac{\eta m_{p} c^{2}}{4kT'_{0}\Gamma_{0}} = n'_{p} \frac{(dL/d\Omega)m_{p}}{(d\dot{M}/d\Omega)4kT'_{0}\Gamma_{0}} \\ &= \frac{(dL/d\Omega)}{4kT'_{0}\Gamma_{0} r^{2}\Gamma\beta c}, \end{aligned} \quad (\text{A10}) $$

which gives

$$ \frac{dN_{\gamma}}{d\Omega} = \frac{dL}{d\Omega} \frac{1}{4kT'_{0}\Gamma_{0}} dt_{e}. \quad (\text{A11}) $$

Given that pairs and radiation are in thermodynamic equilibrium at the base of the jet, $T'_{0}(\Omega)$ is given by
(Levinson & Nakar 2020)

$$T'_{0}(\Omega) = \left[ \frac{dL}{d\Omega} \frac{3}{11} \frac{1}{R_{0}^{2} ac\Gamma_{0}^{2}\beta_{0}} \right]^{1/4} \quad (A12)$$

ORCID iDs

Filip Alamaa https://orcid.org/0000-0001-7414-5884

References

Abramowicz, M. A., Novikov, I. D., & Paczynski, B. 1991, ApJ, 369, 175

Ahlgren, B., Larsson, J., Nymark, T., Ryde, F., & Pe'er, A. 2015, MNRAS, 454, L31

Alamaa, F., Daigne, F., & Mochkovitch, R. 2024, A&A, 683, A30

Band, D., Matteson, J., Ford, L., et al. 1993, ApJ, 413, 281

Bégué, D., Siutsou, I. A., & Vereshchagin, G. V. 2013, ApJ, 767, 139

Beloborodov, A. M. 2010, MNRAS, 407, 1033

Beloborodov, A. M. 2011, ApJ, 737, 68

Beloborodov, A. M. 2013, ApJ, 764, 157

Beloborodov, A. M. 2017, ApJ, 838, 125

Beloborodov, A. M., Lundman, C., & Levin, Y. 2020, ApJ, 897, 141

Bromberg, O., Mikolitzky, Z., & Levinson, A. 2011, ApJ, 733, 85

Bromberg, O., Tchekhovskoy, A., Gottlieb, O., Nakar, E., & Piran, T. 2018, MNRAS, 475, 2971

Campana, S., Mangano, V., Blustin, A. J., et al. 2006, Natur, 442, 1008

Cavallo, G., & Rees, M. J. 1978, MNRAS, 183, 359

Chang, J., & Cooper, G. 1970, JCoPh, 6, 1

Daigne, F., & Mochkovitch, R. 2002, MNRAS, 336, 1271

Dereli-Bégué, H., Pe'er, A., & Ryde, F. 2022, NatCo, 13, 5611

Drenkhahn, G., & Spruit, H. C. 2002, A&A, 391, 1141

Gehrels, N., Chincarini, G., Giommi, P., et al. 2004, ApJ, 611, 1005

Giannios, D. 2006, A&A, 457, 763

Goodman, J. 1986, ApJL, 308, L47

Gottlieb, O., Levinson, A., & Nakar, E. 2019, MNRAS, 488, 1416

Gottlieb, O., Nakar, E., Piran, T., & Hotokezaka, K. 2018, MNRAS, 479, 588

Hascoët, R., Daigne, F., Mochkovitch, R., & Vennin, V. 2012, MNRAS, 421, 525

Hogg, D. W. 1999, arXiv:astro-ph/9905116

Ito, H., Levinson, A., Stern, B. E., & Nagataki, S. 2018, MNRAS, 474, 2828

Ito, H., Nagataki, S., Ono, M., et al. 2013, ApJ, 777, 62

Lazzati, D., Morsony, B. J., & Begelman, M. C. 2009, ApJL, 700, L47

Levinson, A. 2012, ApJ, 756, 174

Levinson, A. 2020, PhRvE, 102, 063210

Levinson, A., & Bromberg, O. 2008, PhRvL, 100, 131101

Levinson, A., & Nakar, E. 2020, PhR, 866, 1

Lundman, C., & Beloborodov, A. M. 2019, ApJ, 879, 83

Lundman, C., & Beloborodov, A. M. 2021, ApJL, 907, L13

Lundman, C., Beloborodov, A. M., & Vurm, I. 2018, ApJ, 858, 7

Lundman, C., Pe'er, A., & Ryde, F. 2013, MNRAS, 428, 2430

Meegan, C., Lichti, G., Bhat, P. N., et al. 2009, ApJ, 702, 791

Meng, Y.-Z., Liu, L.-D., Wei, J.-J., Wu, X.-F., & Zhang, B.-B. 2019, ApJ, 882, 26

Nakar, E., & Sari, R. 2012, ApJ, 747, 88

Nousek, J. A., Kouveliotou, C., Grupe, D., et al. 2006, ApJ, 642, 389

Paczyński, B. 1986, ApJL, 308, L43

Parsotan, T., López-Cámara, D., & Lazzati, D. 2018, ApJ, 869, 103

Pe'er, A. 2008, ApJ, 682, 463

Pe'er, A. 2015, AdAst, 2015, 907321

Pe'er, A., Mészáros, P., & Rees, M. J. 2005, ApJ, 635, 476

Pe'er, A., Mészáros, P., & Rees, M. J. 2006, ApJ, 652, 482

Pe'er, A., & Ryde, F. 2011, ApJ, 732, 49

Piran, T. 1999, PhR, 314, 575

Rees, M. J., & Mészáros, P. 2005, ApJ, 628, 847

Ronchini, S., Oganesyan, G., Branchesi, M., et al. 2021, NatCo, 12, 4040

Rudolph, A., Tamborra, I., & Gottlieb, O. 2024, ApJL, 961, L7

Rybicki, G. B., & Lightman, A. P. 1979, Radiative Processes in Astrophysics (New York: Wiley)

Ryde, F. 2005, ApJL, 625, L95

Samuelsson, F., Lundman, C., & Ryde, F. 2022, ApJ, 925, 65

Samuelsson, F., & Ryde, F. 2023, ApJ, 956, 42

Sari, R., Piran, T., & Narayan, R. 1998, ApJL, 497, L17

Spruit, H. C., Daigne, F., & Drenkhahn, G. 2001, A&A, 369, 694

Tak, D., Uhm, Z. L., Racusin, J., et al. 2023, ApJ, 949, 110

Vianello, G., Lauer, R. J., Younk, P., et al. 2015, arXiv:1507.08343

Vurm, I., & Beloborodov, A. M. 2016, ApJ, 831, 175

Vurm, I., Lyubarsky, Y., & Piran, T. 2013, ApJ, 764, 143

Vyas, M. K., & Pe'er, A. 2023, ApJL, 943, L3

Waxman, E., Mészáros, P., & Campana, S. 2007, ApJ, 667, 351

Zrake, J., Beloborodov, A. M., & Lundman, C. 2019, ApJ, 885, 30