POSyTIVE -- a GRB population study for the Cherenkov Telescope Array (ICRC-2019)

TL;DR

POSyTIVE is a comprehensive simulation project that predicts gamma-ray burst detection rates with the Cherenkov Telescope Array, aiding in understanding GRB emission mechanisms and optimizing observational strategies.

Contribution

It introduces a realistic GRB population synthesis model calibrated with 40 years of data for CTA detection predictions and spectral analysis planning.

Findings

Predicted GRB detection rates with CTA

Simulated time-resolved spectral observations

Explored theoretical models for GRB prompt and afterglow emission

Abstract

One of the central scientific goals of the next-generation Cherenkov Telescope Array (CTA) is the detection and characterization of gamma-ray bursts (GRBs). CTA will be sensitive to gamma rays with energies from about 20 GeV, up to a few hundred TeV. The energy range below 1 TeV is particularly important for GRBs. CTA will allow exploration of this regime with a ground-based gamma-ray facility with unprecedented sensitivity. As such, it will be able to probe radiation and particle acceleration mechanisms at work in GRBs. In this contribution, we describe POSyTIVE, the POpulation Synthesis Theory Integrated project for very high-energy emission. The purpose of the project is to make realistic predictions for the detection rates of GRBs with CTA, to enable studies of individual simulated GRBs, and to perform preparatory studies for time-resolved spectral analyses. The mock GRB population used by POSyTIVE is calibrated using the entire 40-year dataset of multi-wavelength GRB observations. As part of this project we explore theoretical models for prompt and afterglow emission of long and short GRBs, and predict the expected radiative output. Subsequent analyses are performed in order to simulate the observations with CTA, using the publicly available ctools and Gammapy frameworks. We present preliminary results of the design and implementation of this project.

Full text

1 The POSyTIVE project

Gamma-ray bursts ($\gamma\text{-ray}$; GRBs) are among the most powerful phenomena in the Universe. They are appealing targets for the upcoming Cherenkov Telescope Array (CTA), especially considering the recent detections at very high-energies [1, 2, 3]. Observations with CTA will provide valuable insights into particle acceleration at relativistic shocks and radiative processes at work in GRBs. They will enable the use of GRBs as tools to investigate fundamental physics questions (such as violations of Lorentz invariance), and to probe the extragalactic background light.

To perform preparatory studies in view of CTA and predict the CTA GRB detection rate, we are developing the POSyTIVE222 POSyTIVE, POpulation SYnthesis Theory Integrated model for Very high Emission. project. The project will combine population models for both short and long GRBs with emission models of prompt and afterglow radiation, and simulate the detectability of the predicted emission with the CTA.

This method will provide a library of simulated GRBs, which can then be used to i) test follow up strategies of CTA, taking account of the sensitivity on the appropriate time scales; and ii) study which region of the physical parameter space of the GRB population can be constrained by future CTA observations.

2 The population of short and long GRBs

The synthetic GRB populations are simulated based on a minimal set of intrinsic properties, following [4, 5]. Here, we describe the method for the population of long GRBs. A similar scheme is adopted for the population of short bursts. We make the following assumptions:

• •

broken power law distribution of the rest frame peak energy, $E_{\rm peak}$, of GRBs. The Monte Carlo (MC) sampling of this distribution provides for each simulated GRB an intrinsic $E_{\rm peak,i}$;

• •

GRB redshift distribution prescribing the event rate density (in units of Gpc*-3* yr*-1*). The long GRB formation rate is assumed to be proportional to the cosmic star formation rate [6], with the additional possibility of its evolution with redshift. This provides $z_{\rm i}$;

• •

GRBs follow the empirical correlations between $E_{\rm peak}$ and the isotropic equivalent energy, $E_{\rm iso}$, and luminosity, $L_{\rm iso}$ [7, 8]. These provides $E_{\rm iso,i}$ and $L_{\rm iso,i}$.

In addition, we obtain the distribution of observer frame afterglow onset times from the available measurements and lower limits [9]. This observable allows us to infer $\Gamma$, the bulk Lorentz factor of GRBs during the coasting phase. From the onset time, we obtain $\Gamma_{0,{\rm i}}$ for either the case of a constant density circumburst medium, or for a wind density profile.

In order to calibrate the simulated populations with the observed GRB samples, we assume an intrinsic spectrum. The latter is used to compute the flux and fluence (respectively corresponding to $L_{\rm iso,i}$ and $E_{\rm iso,i}$) in the energy bands corresponding to the instruments, which provide the real GRB sample as constraints. We assume that all GRB spectra are described by the Band function [10]. The peak energy is provided as described above. The low and high energy spectral indices are derived by randomly sampling Gaussian distributions.

For each GRB, whose intrinsic properties are drawn via MC from the above probability density distribution, we derive the observer frame properties. The latter include, e.g., the flux and fluence in a given energy band, $\Delta E$. The free parameters of the simulated population are then constrained by observational data. We explore the parameter space in search of combinations that produce a simulated GRB population, consistent with the following:

The GRB population detected by Fermi-GBM. We further limit this comparison to the bright end of the Fermi population in order to avoid possible selection bases affecting the faint end of the peak flux distribution. We compare and reproduce with the simulated population the distributions of peak flux, fluence, duration, and observer frame peak energy of Fermi-GBM bursts. As an independent check of the simulated population, we also verify the consistency with the peak flux distribution of Swift-BAT detected GRBs. 2. 2.

The complete sample of Swift GRBs (selection criteria described in [11, 12]). These are bright Swift GRBs (with a 15–150 keV peak flux $>2.6$ ph cm*-2* s*-1*), which have measured redshifts and prompt emission properties ($E_{\rm iso}$, $L_{\rm iso}$, $E_{\rm peak}$). We require that the simulated GRB population (with the same peak flux cut as above) matches the distributions of $z$, $E_{\rm iso}$, and $L_{\rm iso}$ of the Swift complete sample. The same Swift sample is also adopted (see §4) for the calibration of the afterglow parameters.

The visibility of each GRB from each CTA site is evaluated using a Python code based on Astroplan [13]333 https://astroplan.readthedocs.io/en/latest/ .. This is an open source Python package that, among other features, determines observability of sets of targets from a user-defined observatory, given an arbitrary set of constraints (i.e., altitude, airmass, moon separation & illumination, etc.). Specifically, for each GRB of the population, the code checks the observability from each CTA site during the 24 hours following the trigger time. It then provides a calculation of the corresponding temporal window of visibility. The constraints to determine the observability that we adopted are the following: i) we require that the Sun is 18$\degree$ below the horizon (astronomical twilight); ii) we constrain the minimum altitude of the target, as 10$\degree$ above the horizon.

3 The GRB calibration sample

As described above, the simulated population of GRBs is compared to the observed properties of real samples of long and short GRBs, the so called BAT6 and SBAT4 samples [11, 12, 14].

Both datasets are compiled by selecting GRBs detected by Swift, having favorable observing conditions for redshift determination from the ground. These are events with low Galactic extinction in the direction of the burst, $A_{V}<0.5$ mag, restricted to those that are bright in the 15-150 keV Swift-BAT energy band. In particular, this last criterion is equivalent to a threshold for the peak photon flux, $P$, measured in the 15–150 keV energy band by the Swift-BAT.

For the BAT6 sample, we select long GRBs with $P\geq 2.6$ ph s*-1* cm*-2*, computed using the Swift-BAT light curves binned with $\delta t=1$ s. For the SBAT4 sample, we select short GRBs with $P\geq 3.5$ ph s*-1* cm*-2*, computed using the Swift-BAT light curves binned with $\delta t=64$ ms. Such a high flux cut ensures that our samples are free from any detector-related threshold bias. With the above criteria, the two samples comprise 99 long GRBs (83% having redshifts [12]) and 16 short GRBs (69% having redshifts [14]).

Being free of selection effects (except flux limits), both samples provide the possibility to compare the rest-frame physical properties of GRB prompt and afterglow emission in an unbiased way [11, 15, 16, 17, 18, 19, 20, 14, 21, 22, 12, 23, 24].

4 The prompt model

To investigate the properties of the prompt emission, we use the numerical code described in [25]. The simulated synthetic population of GRBs provides the following properties for each event: peak energy of the prompt spectrum, isotropic energy, $E_{\rm iso}$, redshift, Lorentz factor, and the duration of the prompt emission. At this stage, we use a large set of simulated spectra, calculated in the comoving frame of the shocked material [25, 26]. These are produced for a range of values of the magnetic field, the dynamical timescale, the electron density, and the minimum electron Lorentz factor of the electron distribution.

Each emitted photon spectrum is computed from the time-dependent evolution of relativistic electrons. We do not include the possible contribution from a population of accelerated relativistic protons. Accelerated electrons in the amplified magnetic field radiate synchrotron photons; these synchrotron photons may be scattered to higher energies by relativistic electrons (inverse Compton). At low energies the photons may be absorbed (synchrotron self-absorption). At high energies, we account for photon-photon annihilation, producing electron-positron pairs.444 The contribution of these pairs to the radiation is not considered in the present version of the code.

The resulting spectra in the comoving frame are transformed to the observer frame using the provided values of bulk Lorentz factor and redshift. The peak energy of the spectrum and the observed flux are then compared with the provided values of the synthetic GRB population. In this way, we determine the parameters of the simulated spectra in the comoving frame, which can reproduce the synthetic GRB population. The corresponding spectral components are derived in the energy band observable by CTA.

5 The afterglow model

To simulate afterglow radiation we consider a standard scenario of synchrotron and synchrotron-self Compton (SSC) radiation from electrons accelerated by the forward shock. To describe the dynamics of the blastwave, we follow the method developed by [27]. Particle acceleration and magnetic field amplification at the shock are described by introducing the parameters $\epsilon_{\rm e}$ and $\epsilon_{\rm B}$; these respectively denote the fraction of shock-dissipated energy that is used to accelerate the electrons and to amplify the magnetic field. Electrons are assumed to be efficiently accelerated into a powr law (PL) energy distribution, $dN(\gamma)/d\gamma\propto\gamma^{-p}$, where $\gamma$ is the electron Lorentz factor.

Synchrotron and SSC radiation spectra and their evolution with time are described following [28, 29]. Possible modifications to both spectral components caused by Klein-Nishina (KN) effects follow the prescriptions in [30].

For each GRB of the synthetic population, the initial Lorentz factor is given by §2. The initial blastwave kinetic energy, $E_{\rm k}$, is inferred from $E_{\rm iso}$, assuming an efficiency, $\eta_{\gamma}=20\%$, for the prompt emission mechanism, $E_{\rm k}=E_{\rm iso}\,(1-\eta_{\gamma})/\eta_{\gamma}$. To predict the afterglow emission for each GRB, we need to specify the values of the remaining model parameters: the density of the external medium and the corresponding radial profile, ($n(r)=n_{0}\,r^{-s}$), $\epsilon_{\rm e}$, $\epsilon_{\rm B}$, and $p$. Given the different progenitors and environments of the two classes of GRBs, we assume $s=0$ for short GRBs, and $s=2$ for long ones. The large uncertainties on the remaining free parameters are such that changing the assumed values has a strong impact and small predictive power on the resulting emission. To overcome this problem, before predicting the emission in the CTA energy range, we calibrate the model parameters to reproduce existing afterglow observations at lower frequencies. We use the calibration sample presented in §3. We find the set of values for the free model parameters which are able to describe optical-to-GeV observations. Preliminary results of the parameter calibration are shown in Fig. §2.

These parameters will be used to predict the radiation in the CTA energy range for the full synthetic GRB population.

6 Preliminary tests of detectability with CTA

To test the detection prospects of GRBs by CTA, two methods are being studied. The first uses a sky map likelihood maximization (as for Fermi-LAT), and has been implemented using the ctools framework [31]. The second method relies on the on-off approach, the analysis traditionally used for ground Cherenkov telescopes; it has been implemented as part of a Gammapy-based code [32].

The tests on the simulation pipelines performed to date take as input time-variable very high-energy spectra, derived from GRB afterglows. The selected GRBs are located at different redshifts. Their spectra are properly modified by interactions with the extragalactic background light, following the model of [33]. The GRB spectra are modelled as simple power laws for a series of fixed time bins in a logarithmic scale. Each time slice is considered as an independent observation, to which the CTA instrument response functions (IRFs) [34] are applied. We consider a point-like source in the center of the field-of-view, a zenith angle of 20$\degree$, and an assumed observation time of 100 s. At each time step observations are stacked, and the detection significance is computed. The time slice in which a 3$\sigma$ or 5$\sigma$ detection with a confidence level of 90% is obtained from $\sim$100 MC trials gives the corresponding detection time.

In the first analysis chain, based on the ctools library, a sky map is produced from the full-enclosure IRF. The latter contains the CTA response over the entire field of view ($\pm 2.5\degree$). This allows the use of a likelihood analysis (the ctlike tool), assuming a constant background in time and a simple PL spectral model. A MC simulation involves generating sky maps with varied pixel counts, assuming Poissonian statistics. The detection significance is computed from the test statistics of the likelihood fit.

The current Gammapy-based pipeline uses the CTA on-axis point-like source IRF, for which the 68% containment radius around a source is intrinsically optimized. For all energy bins and each time slice, the predicted signal and background are thus directly obtained (no sky map is simulated). The detection significance calculation follows [35], with $\alpha=1/5$ (i.e. it is considered that the background is obtained from a region which is 5 times larger than the signal region). The simulation involves fluctuating the signal and background counts, based on Poissonian statistics.

Preliminary tests have been performed on a sub-sample of $\sim 10$ generic GRB light-curves, with the intrinsic GRB flux scaled by various factors (between 1/2 and 1/100), testing the capability of CTA to detect faint GRBs. A preliminary analysis of these simulations shows that the majority of the events ($>50\%$) are detected with $\geq 3\sigma$ significance in both analysis chains. The likelihood analysis gives systematically higher significance values than the on-off method, since it uses the full point spread function information. However, the on-off method is much faster. It might be more suitable for situations in which the IRFs are poorly known or have large systematic uncertainties; for instance, it might be appropriate for online analyses.

The two pipelines will be used on the larger, theoretically-based, GRB population, in order to derive our final estimates for the number of potentially detectable GRBs.

Acknowledgements

This work was conducted in the context of the CTA Transients Working Group. We gratefully acknowledge financial support from the agencies and organizations listed here: http://www.cta-observatory.org/consortium_acknowledgments.