H.E.S.S. Follow-up Observations of Binary Black Hole Coalescence Events during the Second and Third Gravitational-wave Observing Runs of Advanced LIGO and Advanced Virgo

In 2017, the detection of both gravitational waves (GWs) and electromagnetic radiation from the merger of two neutron stars (NSs), GW170817, revolutionized multimessenger astronomy (Abbott et al. 2017b). The electromagnetic radiation was observed as a short set of flashes of low-energy γ-rays, which is commonly referred to as a short gamma-ray burst (GRB). This single event solved a decades-old mystery of high-energy astrophysics by confirming that some short GRBs are produced by the merger of two compact objects, with at least one being an NS. The other population of GRBs, long GRBs, have long been associated with core-collapse supernovae and are known to arise from the collapses of the most massive stars. In either GRB progenitor scenario, the catastrophic event launches an ultrarelativistic jet. Interactions within the jet produce the highly variable prompt emission, while the jet's later interactions with the surrounding material cause it to decelerate and produce long-lived afterglow emission. This emission is observed as smoothly fading emission when the GRB is viewed near the axis of the jet. In contrast, as was the case in GRB 170817A (Nynka et al. 2018; Troja et al. 2020), if the GRB is observed off-axis, the afterglow emission is instead observed to initially rise before peaking and then fading due to geometric effects.

Object GRB 170817A was well observed over a wide range of photon energies. With a redshift of 0.0098 (Hjorth et al. 2017), corresponding to a luminosity distance of 40 Mpc, its proximity meant that its emission could be observed to much later times; the X-ray afterglow, for instance, was still detectable 1000 days after the initial event (Troja et al. 2020). It was also observed by very high energy (VHE) γ-ray detectors such as the High Energy Stereoscopic System (H.E.S.S.; Abdalla et al. 2017), an array of five imaging atmospheric Cerenkov telescopes that can detect γ-rays with energies greater than tens of GeV. The upper limits on VHE emission from the late-time H.E.S.S. observations provide constraints on the magnetic field of the outflow (Abdalla et al. 2020).

Short GRBs could also be produced by the merger of an NS and a black hole (BH; see, e.g., Foucart 2020 for a recent review). In contrast, a merger of two stellar mass BHs is typically not expected to produce electromagnetic radiation due to the lack of accreting material, which is considered necessary for launching a jet. However, it is essential to verify this assumption with observations. Indeed, the announcement of a γ-ray transient temporally coincident with the GW150914 binary black hole (BBH) merger (Connaughton et al. 2016, 2018) sparked much interest and controversy and inspired explanations such as a circumbinary or remnant disk (Kotera & Silk 2016; Murase et al. 2016; Perna et al. 2016, 2019; Martin et al. 2018) or charged BHs (Liebling & Palenzuela 2016; Zhang 2016; Fraschetti 2018). However, no GRB-like candidate for a BBH merger has been observed since then, although other potential sources of electromagnetic emission have been discussed (Bartos et al. 2017; Stone et al. 2017; Graham et al. 2020).

A single interferometer has a broad antenna pattern and hence poor directional sensitivity; therefore, the localization of transient GW events with only one detector can encompass large parts of the sky. However, the subset of events that are detected by more than two GW instruments are localized to tens of square degrees or less (Pankow et al. 2018). H.E.S.S. has a field of view (FoV) of about 20 deg² and can observe large portions of the GW localization regions using an observational pattern known as tiling (Ashkar et al. 2021), during which different parts of the sky region are observed according to some predefined sequence. During the second and third observing runs (O2 and O3) of the GW detectors Advanced LIGO and Advanced Virgo, H.E.S.S. conducted follow-up observations of six well-localized GW events, chosen so that a sky area corresponding to at least 50% localization probability could be covered in 1 night for BBH events ⁴⁰ and at least 10% could be covered in 1 night for mergers involving an NS. Of these six events, one was the NS–NS merger GW170817, discussed in Abdalla et al. (2017, 2020); one was an NS–BH merger, GW200105 (Abbott et al. 2021b), and the other four were BBH mergers, GW170814 (Abbott et al. 2017a), GW190512_180714, GW190728_064510 (Abbott et al. 2021a), and S200224ca (LIGO Scientific Collaboration & Virgo Collaboration 2020a). Of these last five, the NS–BH binary was only covered by one pointing (or observation run) due to poor weather conditions, amounting to less than 1% of the localization region; since this translates to such a small chance of having observed the correct sky region, we do not discuss the NS–BH event here and instead focus only on the BBH mergers. The same analysis methods can be used for all GW events regardless of the type of progenitor system, so that the four BBH events discussed here are a good test of H.E.S.S.'s sensitivity to VHE γ-ray emission from GW mergers.

The H.E.S.S. observations during O2 and O3 all began with delays of at least a few hours after the original merger events and were conducted by scanning regions of the sky to maximize the amount of the GW localization region probability that could be covered. The rate of GW events is expected to increase in the fourth GW observing run (O4), while the localization uncertainties should greatly decrease (Abbott et al. 2020). In order to prepare for O4, it is important to assess the sensitivity obtained by the current H.E.S.S. observation strategy and determine if modifications are necessary.

In this paper, we present the H.E.S.S. observations of four well-localized BBH mergers. In Section 2, we discuss the details of the observations; then, we describe the analysis of the data and report on the results in Section 3. We search for a pointlike source anywhere in the region covered by the observations and do not find any significant signal, so we calculate integral upper limits on the γ-ray flux. Finally, we discuss these upper limits in Section 4 and conclude in Section 5. Throughout this paper, we focus on GW190728_064510 as an illustrative example and include the equivalent information and figures for the three other GW events in the Appendix. The cosmological parameters used are taken from Ade et al. (2016) with a Hubble constant H₀ = 67.8 ± 0.9 km s⁻¹ Mpc⁻¹ and a matter density parameter Ω_m = 0.308 ± 0.012.

H.E.S.S. is located in the Khomas Highland in Namibia. It consists of four 12 m telescopes with 20 deg² FoVs and one 28 m telescope with a 9 deg² FoV and is capable of detecting VHE γ-rays ranging from a few tens of GeV to 100 TeV in energy.

The H.E.S.S. Transients Alert System (M. Füßling et al. 2021, in preparation) monitors GW alerts via the automated GCN (Gamma-ray burst Coordinate Network (GCN) 2020) notices provided by the LIGO and Virgo Collaborations. These notices contain GW localization maps in HEALPix (Górski et al. 2005; Zonca et al. 2019) format. The pixels of the map contain four layers of information: a layer containing the probability of finding the GW event in the corresponding position in the sky and three layers containing distance information. The BBH events followed by H.E.S.S. are at large distances of several hundreds to thousands of megaparsecs. Current galaxy catalogs, such as GLADE (Dálya et al. 2018), are not complete at these distances and cannot be used to further optimize the follow-up strategy. Therefore, the BBH events are followed by 2D grading schemes that cover the regions in the sky that have the greatest 2D probability integrated inside the H.E.S.S. FoV without taking into account distance information. By default, each region is observed with one run.

The most probable regions to host the GW event are observed first, taking into consideration observation and visibility constraints. Imaging atmospheric Cerenkov telescopes like H.E.S.S. observe during astronomical nights with low levels of light falling into the cameras. To ensure this, the Sun and Moon altitude, Moon phase, and Moon-to-source separation need to be considered. Moreover, due to the absorption of Cerenkov shower light in the atmosphere, the energy threshold of the H.E.S.S. observations depends heavily on the respective zenith angle. In general, low zenith angle observations are favored, and observations are restricted to zenith angles below 60°. This constraint on zenith angle combined with the telescope location on Earth determines the sky visibility for H.E.S.S. at a given time.

The GW follow-ups are adapted to these H.E.S.S. visibility and observation conditions in order to observe the GW events as quickly and efficiently as possible. More details on the pointing patterns of H.E.S.S. GW follow-up observations, the description of the GW follow-up triggering framework, and the GW follow-up tools integrated within the H.E.S.S Transient Follow-up System (M. Füßling et al. 2021, in preparation) are provided in Ashkar et al. (2021).

Event GW170814 (Abbott et al. 2017a) was the first interferometer BBH merger detection that occurred during O2. It was also the first BBH merger whose location permitted H.E.S.S. follow-up observations. Observations were conducted over 3 nights (2017 August 17, 18, and 19) and are presented in Ashkar et al. (2019); as this was the first tiled GW observation, this was the first opportunity to test the developed follow-up algorithms.

During O3, H.E.S.S. performed follow-up observations of three BBH events: GW190512_180714, GW190728_064510 (Abbott et al. 2021a), and S200224ca (LIGO Scientific Collaboration & Virgo Collaboration 2020a). In contrast to O2, the follow-up observations were restricted to a maximum delay of 24 hr and all performed during the same night following the arrival of the GW alert.

Event GW190512_180714 was the first event that was detected by three interferometers, and the first BBH to be well localized in O3 with a favorable zenith angle for H.E.S.S. The observation delay was ∼7.1 hr, and with a total of six observation runs, the total observation duration was ∼2.6 hr. This GW event in particular was used to commission and streamline the H.E.S.S. response to GW alerts during O3.

Event GW190728_064510 was initially classified as a MassGap event (LIGO Scientific Collaboration & Virgo Collaboration 2019a). The H.E.S.S Transient System automatically scheduled follow-up observations starting with a delay of ∼11 hr. A neutrino candidate (IceCube Collaboration 2019) event from IceCube was temporally coincident with the GW event (within 360 s) and consistent with its sky localization. When the H.E.S.S observations started, GW190728_064510 was reclassified as a BBH merger, and the localization map was updated (LIGO Scientific Collaboration & Virgo Collaboration 2019b). A new schedule was automatically computed by the H.E.S.S. Transient System, and 10 runs were taken summing to a total of 5.65 hr of data. Follow-up observations of GW190728_064510 covered parts of the uncertainty region of the neutrino candidate from IceCube.

Event S200224ca (LIGO Scientific Collaboration & Virgo Collaboration 2020a) was one of the well-localized and therefore best followed-up BBH events (according to the Treasure Map ⁴¹ ; Wyatt et al. 2020). H.E.S.S. took three runs covering the larger part of the localization region (LIGO Scientific Collaboration & Virgo Collaboration 2020b). Due to weather conditions, these observations started with ∼3 hr of delay and lasted for 1.4 hr.

A new camera for the 28 m telescope was in the commissioning phase during the S200224ca follow-up. Therefore, to have a consistent reconstruction for all observed GW events, we exclude the 28 m telescope and only use the data from the 12 m telescopes for all GW event follow-up analyses discussed here. We note that the 28 m telescope data were analyzed, and the results are consistent with the results shown in this paper with a slight decrease in energy threshold. We also note that during the GW170814 and GW190512_180714 observations, one of the four 12 m telescopes was not available due to maintenance. The pointing positions for the GW190728_064510 follow up are presented in Table 1 and in Table 4 of the Appendix section for GW170814, GW190512_180714 and S200224ca. Some positions were observed more than once, either due to technical or weather conditions or because there was a hint of a signal in the preliminary real-time analysis.

The data analysis is performed with the model analysis (de Naurois & Rolland 2009), which is a reconstruction method that performs a log-likelihood comparison between the recorded shower images and a semianalytical model of γ-ray showers. As explained in Section 2, the data from the four 12 m telescopes have been analyzed in stereoscopic mode. The background level in the FoV is determined from the data set itself using the ring background technique (Berge et al. 2007). For that, all γ-ray-like events from the runs are accumulated, and a reconstructed map with 0.02° pixel size is created. On each position (pixel) of the map, the ON region is defined as a 0.1° radius circular region. The excess counts are computed by subtracting the background of γ-ray-like events after analysis and selection cuts, determined from a ring region around the ON region and accounting for the radial variation of acceptance. As an example, we show in the left plot of Figure 1 the resulting excess map for the GW190728_064510 GW event. The edges of the FoV are constrained by the statistics of the γ-ray-like background events. We only consider the regions where α N_OFF > 5, where N_OFF is the number of events in the background region and α is the ratio of the effective exposure in the signal region over the effective exposure in the background region. Significance values are then computed following Li & Ma (1983). A significance map is computed as shown in the middle panel of Figure 1. The significance distribution shown in the right panel follows a Gaussian distribution and is consistent with the background-only hypothesis. The same is found for the other three GW events. Therefore, we conclude that the analyses resulted in no significant detection of any γ-ray counterparts to the BBH merger events. These findings have been confirmed with an independent analysis using the Image Pixel-wise fit for Atmospheric Cerenkov Telescopes (Parsons & Hinton 2014) software, yielding consistent results.

Zoom In Zoom Out Reset image size

In order to constrain the VHE γ-ray emission from these BBH mergers at the time of the H.E.S.S. observations, we derived integral upper-limit maps for each GW event for a given energy interval $[{E}_{\min }$, ${E}_{\max }]$ and a given spectrum $\phi {(E)\propto (E/{E}_{0})}^{-{\rm{\Gamma }}}$, with E₀, the pivot energy, set to 1 TeV. We first derive upper-limit maps in the 1–10 TeV energy range, in which H.E.S.S. is most sensitive, and for an E⁻² spectrum following the procedure described in Abdalla et al. (2018). The upper-limit map for GW190728_064510 is presented in Figure 2, and Figure 6 of the Appendix presents the three remaining GW events.

Zoom In Zoom Out Reset image size

We also exploit the widest energy range possible for each observation between ${E}_{\min }$, which is the threshold energy E_th, and a maximum energy ${E}_{\max }$. The energy range is imposed by instrument limitations and varies between runs, as it depends on the zenith angle of the position studied in the sky at the time of observation. Higher zenith angles lead to higher energy thresholds. The energy threshold E_th and maximum energy ${E}_{\max }$ are chosen where the energy reconstruction bias is less than 10%. The E_th is increased if it is below the energy at which the effective area is at 10% of its maximum value. For the latter, a per-pixel effective area constrained energy threshold is considered, and a median over the entire coverage is taken. The resulting energy ranges are presented in Table 2.

Table 2. Spectral Indices (γ) at a GW Event's Corresponding Redshift (Abbott et al. 2019, 2021a) and at E_th Assuming an Intrinsic E⁻² Spectrum

GW Event	Redshift	γ(E = Eth, z = zGW)	Energy Range (TeV)	Coverage	Tstart, Tstop (s)
GW170814	${0.12}_{-0.04}^{+0.03}$	2.73	0.42–34.80	75.4%	2.22 × 105, 4.10 × 105
GW190512_180714	${0.27}_{-0.10}^{+0.09}$	3.50	0.31–38.31	34.5%	1.84 × 104, 2.82 × 104
GW190728_064510	${0.18}_{-0.07}^{+0.05}$	2.98	0.35–26.10	50.8%	4.88 × 104, 7.28 × 104
S200224ca	0.29	3.08	0.24–38.31	62.13%	1.07 × 104, 1.59 × 104

Note. The redshift for S200224ca was estimated from the distance in LIGO Scientific Collaboration & Virgo Collaboration (2020b) using the cosmological parameters from Ade et al. (2016). The energy range used to derive the specific integral upper-limit maps and the corresponding coverage are presented in the fourth and fifth columns, respectively, and the sixth column lists the start and end of the H.E.S.S. observations of the GW event, as calculated from the reported GW merger time.

Download table as:  ASCIITypeset image

For the intrinsic spectrum shape, the same generic E⁻² power-law spectrum is assumed. Extragalactic background light (EBL) in the ultraviolet, optical, and infrared bands interacts with VHE γ-rays via electron–positron pair production absorbing parts of the VHE γ-ray flux. The EBL absorption effects increase with energy. Given the EBL absorption (Franceschini et al. 2008) at the respective redshift of each GW event, we evaluate the equivalent power law at the threshold energy E_th of each GW event and compute the equivalent spectral indices. The E_th is considered due to the fact that low-energy events are expected to be dominant. The attenuation due to EBL causes the observed spectrum to be softer than the intrinsic one. The resulting power-law indices are presented in Table 2, and the effective power-law spectra are used to compute specific upper-limit maps. Varying the intrinsic spectral index between 1.5 and 2.5 leads to a variation of the equivalent spectral index between +16% and −16% for the absorbed emission. This affects the integral upper-limit values between +19.4% and −10.7%. Moreover, the equivalent spectral index shows less than 14% dependence on the EBL model in the lower part of the energy spectrum up to 1 TeV. Propagating this dependency to the upper-limit maps will result in less than 10% variation in the upper-limit values.

The 95% C.L. integral upper limits are computed with the Rolke method (Rolke et al. 2005) at each observed position of the map for all four BBH mergers within the energy ranges defined in Table 2, as explained above using the effective power-law spectra. The maps are presented in Figure 7 in the Appendix for all four GW events. The effective VHE γ-ray coverage of the GW localization region after data selection and analysis is presented in the last column of Table 2.

The H.E.S.S. observations presented in this paper cover large portions of the GW localization regions, as shown in Table 2. The computed upper limits constrain the VHE γ-ray emission from BBH mergers for the first time. No electromagnetic or neutrino counterpart has been confirmed for any of the GW events observed by H.E.S.S.

The models that predict electromagnetic emission from BBH mergers are primarily inspired by the candidate electromagnetic counterpart to GW150914. Because these models focus on the production of X-ray and low-energy γ-ray emission, they do not make many statements about VHE emission. For instance, the scenarios involving a remnant disk make predictions on the electromagnetic emission based on synchrotron radiation from the relativistic jets that are launched after the merger. Assuming that the particle acceleration and radiation of the parent charged particle population occur in the same region, there is a maximum synchrotron photon energy that can be estimated by equating the particle acceleration and emission timescales; for GRB-like conditions, after a few hours, this tends to be ${ \mathcal O }$(1 GeV) (e.g., Abdalla et al. 2021), more than an order of magnitude below the H.E.S.S. energy range and 2 orders of magnitude below the GW event-specific energy ranges discussed here (Table 2). The H.E.S.S. upper limits therefore do not constrain these models, and an additional photon production mechanism—such as inverse Compton emission—would be necessary to generate a potentially detectable VHE signal. Note that, while VHE emission from GRBs has been detected, its origin is not completely understood, and its applicability to BBH mergers is unclear. Given the amount of uncertainty involved, any predictions on the production of VHE emission from BBH mergers hours after the merger event is beyond the scope of this paper.

While the VHE emission from a BBH merger is difficult to predict, we can still compare the energy flux and luminosity upper limits for the BBH events to the measurements of GRBs, given the proposed models that link GRBs and BBH mergers. The energy flux upper limits for the H.E.S.S. GW observations were calculated both with and without correcting for EBL absorption. The EBL-corrected energy flux upper limits F_unabs are used for calculating the luminosity upper limits and are calculated assuming an intrinsic E⁻² spectrum over 1–10 TeV in the source frame. This results in an event-specific energy range of $\tfrac{1}{1+z}$ to $\tfrac{10}{1+z}$ TeV in the observer frame. Note that calculating the intrinsic values in this way in general runs the risk of integrating outside of an instrument's bandpass, as the energy range over which spectral fitting is performed is usually defined by instrumental limitations; however, $\tfrac{1}{1+z}$ to $\tfrac{10}{1+z}$ TeV is safely within the event-specific energy ranges (Table 2) for all of the events discussed here.

We then calculate the isotropic luminosity from the unabsorbed energy flux. In short,

$\begin{eqnarray}&&{L}_{\mathrm{iso}}=4\pi {D}_{{\rm{L}}}^{2}{\int }_{{E}_{1}/(1+z)}^{{E}_{2}/(1+z)}{{EN}}_{0}{\left(\displaystyle \frac{E}{{E}_{0}}\right)}^{-{\rm{\Gamma }}}{dE},\end{eqnarray} \tag{ 1 }$

where E₁ = 1 TeV, E₂ = 10 TeV, Γ = 2, and D_L is the luminosity distance. (Here N₀ encompasses the observation-specific measurements; see Equation (2) of Abdalla et al. 2018.) We use the per-pixel luminosity distance from the latest published GW sky maps as, for a given GW event, the mean luminosity distance can vary by several hundred megaparsecs in the regions observed by H.E.S.S. The EBL-absorbed energy flux upper limits are calculated with the GW event-specific spectra and energy ranges (Table 2). The luminosity upper-limit map for GW190728_064510 is shown in Figure 3, and the maps for the three remaining events are shown in Figure 8 of the Appendix.

Zoom In Zoom Out Reset image size

Among the GRBs that have been detected by VHE telescopes, the two detected by H.E.S.S. (Abdalla et al. 2019, 2021) have measurements hours after the GRB began, which is roughly the same time delay as for the GW events, so we include these two in our comparison sample. In order to enlarge the sample, we also take the set of GRBs with known redshift and temporally extended emission observed by the Large Area Telescope (LAT) on board the Fermi space observatory (Ajello et al. 2019). For the luminosity estimates, we consider the energy flux at late times measured by the LAT in the 100 MeV–100 GeV energy range. Using the spectral index γ_LAT measured by the LAT at these times, we extrapolate the spectrum to calculate the energy flux in the H.E.S.S. energy band, assuming an optimistic scenario in which the index stays the same. We then extrapolate these forward in time using the late-time power-law decay index measured by the LAT and convert the energy fluxes into isotropic luminosities using the redshifts of the GRBs.

We emphasize that these extrapolations should be considered more as guidelines to illustrate a range of potential behavior and should not be considered as strict predictions of VHE emission from BBH mergers. While short GRBs are the ones that are associated with binary mergers, all but one of the LAT GRBs (and both of the H.E.S.S. GRBs) in our sample are long GRBs, which have been observed to have isotropic energy releases that are, on average, 1 or 2 orders of magnitude larger than those of short GRBs (Li et al. 2016). In addition, the extrapolations themselves are rather simplistic and result in overpredictions on the VHE flux. Some of the LAT GRBs have spectral indices of around 2, suggesting that the spectrum must soften somewhere above the LAT energy range—potentially in the H.E.S.S. energy range—to avoid the energy flux diverging, whereas we have assumed that it continues unchanged. In addition, the spectral index should evolve with time as the blast wave decelerates and the densities decrease, whereas we have not included any spectral evolution. To estimate the effect of our assumption and illustrate the magnitude of uncertainty, we also calculate what the energy flux will be if the spectral index softens when going from the LAT to the H.E.S.S. energy range. We find that if the spectrum softens from γ_LAT → γ_LAT + 0.5, the luminosity drops by a factor of around 50; if the spectrum softens to γ_LAT + 1, the luminosity drops by 3 orders of magnitude.

The luminosity comparisons are shown in Figure 4. The upper limits for the four BBH events lie below some of the extrapolations of LAT GRBs—given that the GRBs in general are located at larger distances—and are at a similar level as the VHE-detected GRB 190829A, which is at a similarly low redshift. When comparing the four BBH events, it can be seen that the two closer events (GW170814 and GW190728_064510) have the lowest-luminosity upper limits. The value of the upper limit is not greatly affected by the observation duration; since the observations are tiled, increasing the duration increases the amount of sky coverage rather than the sensitivity to a signal in any particular part of the sky. In comparison, the upper limit for GW170817—taken over the entire observation duration—is 3 orders of magnitude lower, primarily because of this event's proximity but also due to the deeper observations of this well-localized source (Abdalla et al. 2020).

Zoom In Zoom Out Reset image size

By comparing the luminosities, we are comparing the intrinsic properties of the source. In order to compare the sources as they are observed on Earth, we additionally calculate the energy flux upper limits using the EBL-attenuated spectra with the power-law index and energy range given in Table 2. For the VHE-detected GRB 180720B, we calculate the energy flux for the observed spectrum using the simple power-law fit to the data (i.e., the fit that does not correct for the EBL absorption). For GRB 190829A, a power-law fit to the data is not provided for the third night on its own, so we calculate the EBL-attenuated energy flux for each of the 3 nights separately using the constant intrinsic photon index of 2.07 derived by combining the data from all 3 nights. For the LAT GRBs, we take the extrapolated spectra and calculate the effect of the EBL on these spectra if the LAT GRBs were all at the redshift of GW190728_064510, over the specific energy range of GW190728_064510. (Using one of the other BBH events results in a decrease in the energy flux extrapolations by at most 50% for GW190512_180714 and S200224ca and an increase of less than 75% for GW170814.)

The average energy flux upper limits for the GW events are shown in Figure 5, along with the LAT GRB extrapolations and the measured values for the H.E.S.S.-detected GRBs. The third night of GRB 190829A (gray circles) resulted in an energy flux measurement below the GW upper limits, mainly because the GRB measurement was the result of almost 5 hr of observations on a single position. (The difference in the energy ranges does not greatly affect this comparison; assuming the measured GRB 190829A spectrum extends up to 40 TeV increases the absorbed energy flux by 7%, and the increase is negligible for GRB 180720B given its larger distance.) In contrast, for GW observations, no single sky position gets as much exposure as it would in a standard single-position multihour GRB follow-up.

Zoom In Zoom Out Reset image size

The upper limits on the energy flux of the BBH events lie over an order of magnitude above almost all of the LAT GRB extrapolations. In addition, the upper limits for the four events are found at similar levels, although they span a large range in observation delays. In order to test whether BBH mergers produce gamma-ray radiation in the H.E.S.S. energy range at a level similar to these extrapolations of LAT GRBs, the H.E.S.S. upper limits would need to be comparable to these extrapolations. This indicates that either shorter observation delays or deeper, more sensitive observations by H.E.S.S. are required to test this.

The temporal extrapolations of the LAT GRBs in Figures 4 and 5 assume that these are all on-axis events. Indeed, this assumption is valid, given that all of the GRBs in this group had afterglow light curves that decayed with time, whereas off-axis events such as GRB 170817A produce emission that rises at early times (Nynka et al. 2018; Troja et al. 2020). Given the current typical sensitivities of GW detectors, the inclination angles—i.e., the angle between the line of sight and the total angular momentum vector—of most GW events can only be weakly constrained due to the degeneracy between the inclination angle and the source distance (Nissanke et al. 2010; Usman et al. 2019). It can only be constrained by knowledge of the host galaxy, which is difficult to ascertain without an electromagnetic counterpart. Because of this, the inclination angles of the events discussed here are mostly unconstrained and so cannot be used to further predict the behavior of any potential electromagnetic emission.

During O2 and O3, H.E.S.S. observed four BBH mergers: GW170814, GW190512_180714, GW190728_064510, and S200224ca. For these GW events, the delay between merger and H.E.S.S. observation is at least a few hours, due to the necessity of waiting for the part of the sky containing the merger to be at favorable zenith angles above the telescopes. The low rate of sufficiently well-localized GW events in the previous Advanced LIGO and Advanced Virgo runs meant we were unlikely to have observed an event with a favorable time and position. However, this is expected to change in O4 and beyond. By comparing the expectations for O3 and O4 (Abbott et al. 2020), ⁴² we can see that, prior to the start of O3, only around 25% of GW events in a three-interferometer scenario were expected to have 90% credible regions (high-latency off-line analysis), similar to the values for the four events discussed here (90% credible region ⪅100 deg²); from GraceDB (2021), it can be seen that this was an accurate prediction of the low-latency localizations in O3. For an event rate of ${18}_{-12}^{+53}$ yr^–1 (combining binary neutron stars, neutron star-black hole and BBH events for simplicity), this translates to only a handful of events per year with sufficiently small credible regions to pass the requirements for H.E.S.S follow-up. This number increases to 75% of events for the four-interferometer network in O4, and with an estimated event rate of ${90}_{-55}^{+232}$ yr^–1, this means around ${67}_{-41}^{+174}$ mergers yr^–1 would pass the requirements for H.E.S.S. follow-up. Given that around half the sky is observable by H.E.S.S., and taking into account that, on average, H.E.S.S. can observe for ∼6 hr a night (averaging over a year), this results in approximately ${8}_{-5}^{+22}$ events yr^–1 occurring during darktime, i.e., being observable with minimal latency. For well-localized GW events, H.E.S.S. observations can begin within 1 minute. On the GW side, in O3, the preliminary notices arrived with a latency of ≈10 minutes (GraceDB 2021); for O4, this would mean the upper limits in Figures 4 and 5 moving to values of T − T₀ < 1000 s for around eight events yr^–1 (not accounting for bad weather). For a fraction of GW events, this latency will become negligible or even negative (as measured from the merger time) as GW detector sensitivity improves with the implementation and improvement of the early-warning system (Nitz et al. 2020; Sachdev et al. 2020; Magee et al. 2021).

Alternatively, given the more precise localizations expected in O4 and beyond, H.E.S.S. could instead choose to spend time on deeper observations of well-localized GW events. In O4, around 35% of mergers observed by the four-interferometer network are expected to have 90% credible areas of less than 20 deg², which would translate to 50% credible areas of less than a few square degrees. For the expected total merger rate of ${90}_{-55}^{+232}$ yr^–1, given that half of these events will be too far north and roughly half of them will be behind the Sun, this translates to around ${8}_{-5}^{+20}$ events yr^–1 for which the 50% credible region could be covered with only one pointing within 24 hr. This would then allow the energy flux upper limit in Figure 5 to improve by up to a factor of ∼5, depending on the amount of observation time, given that the maximum amount of contiguously available observation time is around 12 hr.

Our observation strategy aims to scan over the GW localization regions, thereby emphasizing sky coverage over sensitivity. Improving the chances of a VHE detection would require reducing the delay of follow-up observations and/or spending more time observing single sky regions. Assuming that the GW events are observed fairly on-axis (i.e., that the temporal decay of any electromagnetic emission fades like GRB afterglows), we find that minimizing the observation delay (i.e., following up GW events that are immediately observable by H.E.S.S.) would have a greater effect on the detectability than reducing the sky coverage and spending more time observing single positions. This is expected to happen naturally in the next observing runs; with the increased rate of GW detections, the number of events with favorable conditions for prompt follow-up observations will also increase. Therefore, we do not expect to fundamentally alter our observing strategy and will continue to prioritize sky coverage.

We provide our upper-limit maps in fits format at https://www.mpi-hd.mpg.de/hfm/HESS/pages/publications/auxiliary/2021_BBH_O2_O3/.

The support of the Namibian authorities and the University of Namibia in facilitating the construction and operation of H.E.S.S. is gratefully acknowledged, as is support by the German Ministry for Education and Research (BMBF), the Max Planck Society, the German Research Foundation (DFG), the Helmholtz Association, the Alexander von Humboldt Foundation, the French Ministry of Higher Education, Research and Innovation, the Centre National de la Recherche Scientifique (CNRS/IN2P3 and CNRS/INSU), the Commissariat a l'energie atomique et aux energies alternatives (CEA), the U.K. Science and Technology Facilities Council (STFC), the Knut and Alice Wallenberg Foundation, the National Science Centre, Poland grant No. 2016/22/M/ST9/00382, the South African Department of Science and Technology and National Research Foundation, the University of Namibia, the National Commission on Research, Science & Technology of Namibia (NCRST), the Austrian Federal Ministry of Education, Science and Research and the Austrian Science Fund (FWF), the Australian Research Council (ARC), the Japan Society for the Promotion of Science, and the University of Amsterdam. We appreciate the excellent work of the technical support staff in Berlin, Zeuthen, Heidelberg, Palaiseau, Paris, Saclay, Tübingen, and in Namibia in the construction and operation of the equipment. This work benefited from services provided by the H.E.S.S. Virtual Organisation, supported by the national resource providers of the EGI Federation.

For GW170814, GW190512_180714 and S200224ca, the pointing positions are presented in Table 4. Figure 6 presents the integral upper limit maps for the same events in the 1–10 TeV range assuming an E⁻² spectrum. Figure 7 presents their integral upper limit maps in the specific energy range and specific spectrum indices presented in Table 2. Figure 8 presents their luminosity integral upper limit maps.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size