Swift-BAT GUANO follow-up of gravitational-wave triggers in the third LIGO–Virgo–KAGRA observing run

We present results from a search for X-ray/gamma-ray counterparts of gravitational-wave (GW) candidates from the third observing run (O3) of the LIGO–Virgo–KAGRA (LVK) network using the Swift Burst Alert Telescope (Swift-BAT). The search includes 636 GW candidates received in low latency, 86 of which have been confirmed by the offline analysis and included in the third cumulative Gravitational-Wave Transient Catalogs (GWTC-3). Targeted searches were carried out on the entire GW sample using the maximum–likelihood NITRATES pipeline on the BAT data made available via the GUANO infrastructure. We do not detect any significant electromagnetic emission that is temporally and spatially coincident with any of the GW candidates. We report flux upper limits in the 15151515–350350350350 keV band as a function of sky position for all the catalog candidates. For GW candidates where the Swift-BAT false alarm rate is less than 10−3superscript10310^{-3}10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT Hz, we compute the GW–BAT joint false alarm rate. Finally, the derived Swift-BAT upper limits are used to infer constraints on the putative electromagnetic emission associated with binary black hole mergers.

The discovery of gravitational waves (GWs) from coalescing binary black holes (BBH) by the Laser Interferometer Gravitational-Wave Observatory (LIGO) opened a new window to the Universe (Abbott et al., 2016a). In addition to GWs, compact binary mergers with at least one neutron star (NS) component are likely to generate electromagnetic (EM) radiation (e.g., Nakar 2020a; Kyutoku et al. 2021). Coincident detection of EM emission from compact binary mergers provides a complete picture of the merger process and can have huge implications for our understanding of the Universe. Such coincidences play a crucial role in tracing the properties of the source host galaxy (Troja et al., 2017; Alexander et al., 2017), mitigating degeneracies in GW parameter estimation (Abbott et al., 2017a; Hughes & Holz, 2003; Wang & Giannios, 2021), placing constraints on the NS equation of state (Bauswein et al., 2017; Radice et al., 2018), and investigating the expansion rate of the Universe, thereby testing cosmological models (Abbott et al., 2017a, b; Hotokezaka et al., 2019; Nissanke et al., 2013; Schutz, 1986). Additionally, they allow for the measurement of arrival time differences between photons and gravitons, providing limits to the mass of a graviton, exploring potential violations to the equivalence principle and Lorentz invariance (Abbott et al., 2017e).

The joint detection of the first GW event consistent with a binary NS (BNS) coalescence GW170817 (Abbott et al., 2017a), and a coincident short gamma-ray burst GRB 170817A (Goldstein et al., 2017a; Savchenko et al., 2017), accompanied by the optical/infrared kilonova counterpart AT 2017gfo (Arcavi et al., 2017; Coulter et al., 2017; Tanvir et al., 2017; Evans et al., 2017; Pian et al., 2017; Smartt et al., 2017; Drout et al., 2017; Cowperthwaite et al., 2017) and the GRB afterglow (in the X-rays: Troja et al. 2017; Margutti et al. 2017 and radio: Hallinan et al. 2017), together ushered in a new era in the field of multi-messenger astrophysics and forever impacted our comprehension of compact binary coalescences (CBC) involving an EM counterpart. Massive coordinated EM follow-up efforts were dedicated to deeply monitor the error regions derived from the joint sky localization of GW detectors and high-energy satellites, helping to reduce the initial three detector network sky localization from 28 deg2 to within a few arcseconds of the host galaxy NGC 4993 (Abbott et al., 2017d, a). The spectacular spectral and light curve evolution of this transient (Abbott et al., 2017d; Villar et al., 2017) suggested that this explosive event was an active site for r-process nucleosynthesis (Pian et al., 2017; Smartt et al., 2017; Coulter et al., 2017; Drout et al., 2017) (for a detailed review of the multimessenger observations of GW170817, see, e.g., Nakar 2020b; Margutti & Chornock 2021).

The expected EM counterpart emission from BNS or neutron star–black hole (NSBH) mergers can potentially be weak due to various factors such as considerable source distances, an off-axis viewing angle, or limited amount of ejected mass. For the specific case of GW170817, despite a coincident GRB detection, it took nearly half a day to localize the host galaxy and begin observations of the kilonova (Abbott et al., 2017a). Prompt targeted searches around the GW trigger times, leveraging facilities with enhanced localization capabilities, can refine search strategies and assist optical or infrared (IR) facilities in correctly identifying and pursuing transient candidates for subsequent follow-up studies. In addition to prompt searches, Fermi-GBM analysis of triggers from the first and second LIGO–Virgo observing runs showed that targeted offline searches are capable of recovering additional candidate joint events that may be of astrophysical relevance (Hamburg et al., 2020; Pillas et al., 2023). Temporal and spatial coincidence information can be used to derive the joint false alarm rate (FAR). These estimates have the potential to elevate subthreshold triggers in either the GW or GRB domains to the status of an above-threshold candidate detection (Nitz et al., 2019).

Unlike Fermi, Swift has been for a long time incapable of relaying a continuous stream of event mode data to the ground in real time. Such a capability was enabled through GUANO (Gamma-ray Urgent Archiver for Novel Opportunities, described in Section 2) (Tohuvavohu et al., 2020), which recovers event data from the Swift Burst Alert Telescope (BAT, Barthelmy et al. 2005), that then get processed by the Non-Imaging Transient Reconstruction And TEmporal Search (NITRATES, DeLaunay & Tohuvavohu 2022) pipeline (see Section 4) to search for subthreshold transient candidates. In addition to other astronomical transients, such as GRBs, fast radio bursts (FRBs), and high-energy neutrinos, the GUANO-NITRATES infrastructure performs targeted searches on GW events communicated by the LVK Collaboration, to detect possible GRBs associated with CBCs.

The impact and potential of Swift-BAT subthreshold searches are crucial for multi-messenger related goals. Indeed, deeper targeted searches increase the joint detection horizon, thus enhancing the probability of finding weak EM counterparts of CBCs in the hard X-ray domain. Moreover, thanks to the high spatial accuracy enabled by the BAT coded mask, subthreshold searches open the possibility of recovering the position of the candidate EM event at the precision level of a few arcminutes, fundamental to drive the subsequent follow-up with ground and space-based EM facilities.

Currently, the targeted search analysis carried out thanks to GUANO, has enabled the discovery of more than 35 GRBs with arcminute localization. A total of 7 of the detected GRBs have a duration <2absent2<2< 2 s (e.g., DeLaunay et al. 2020; Tohuvavohu et al. 2022a), hence they are potentially associated with CBCs containing at least one NS. GUANO data have also been used for the localization of 29 long GRBs through imaging (e.g., DeLaunay et al. 2021a) and non-imaging analysis techniques (e.g., DeLaunay et al. 2021b; Tohuvavohu et al. 2022b). GRB 220107A, detected during BAT slew and localized with arcminute precision, enabled the first optical redshift measured using GUANO data (DeLaunay et al., 2022a). The arcminute localization of GRB 211106A enabled prompt multiband follow up and led to the discovery of the first afterglow in the millimeter band from a short GRB (Tohuvavohu et al., 2021a). With GUANO, one can additionally recover coarse localization information on GRB-like transients that originate from outside the BAT field of view (FOV; e.g., DeLaunay et al. 2023).

In addition to the application to real-time analysis, the availability of BAT data enables us to perform a systematic, deeper targeted search focused on archival LVK triggers. The goal of this study is to perform such an analysis on all the LVK triggers received during the third LIGO–Virgo observing run, during which the GUANO pipeline started to be fully operational. The run duration was comprised of two segments: O3a, which operated from April 1, 2019, 15:00 UTC to October 1, 2019, 15:00, and O3b which operated from November 1, 2019, 15:00 UTC, to March 27, 2020, 17:00 UTC. The alerts distributed during O3 were reporting the following parameters: FAR, the signal classification (CBC or unmodeled Burst), and the associated astrophysical probabilities. The results of O3 are summarized in Gravitational-Wave Transient Catalog data releases GWTC-2 (Abbott et al., 2021), GWTC-2.1 (Abbott et al., 2024), and GWTC-3 (Abbott et al., 2023).

In this work, we use Swift-BAT observations to carry out offline targeted subthreshold searches for EM counterparts of the GW triggers obtained during O3. The rest of the paper is organized as follows: In Section 2, we describe the Swift-BAT instrument and its new capabilities, and in Section 3, we provide details about the GW trigger sample used for the analysis. In Section 4 we summarize the targeted search method adopted for the analysis. We present the results from our targeted search analysis on the various subcategories of triggers in Section 5, and discuss the scientific interpretation in Section 6.

The Neil Gehrels Swift Observatory (henceforth, Swift) is a GRB-focused mission launched in 2004, with three onboard payloads – the Burst Alert Telescope (BAT), the X-ray Telescope (XRT), and the UltraViolet/Optical Telescope (UVOT) – covering the EM spectrum from hard X-rays and gamma-rays all the way to the optical (Gehrels et al., 2004). The BAT instrument (Barthelmy et al., 2005) is a hard X-ray coded mask imager with a wide FOV that operates in the broad energy band of 15–350 keV. It is the primary instrument that detects GRBs and performs an onboard imaging analysis via a cross-correlation between the spatial pattern of the counts across the detector array and the pattern of the coded mask. The sensitivity of BAT is capable of providing arcminute localizations of GRBs triggered onboard (Gehrels et al., 2004). Due to the lack of continuous downlinking of timing, spatial, and energy information for each detector count (event mode data), carrying out targeted searches offline has not been possible in the past. A new infrastructure, called GUANO, was incorporated into the Swift-BAT operations in 2019. Details of the GUANO operations can be found in Tohuvavohu et al. (2020). From the outset, GUANO has demonstrated that recovering event mode data from astrophysically compelling time windows can enhance the overall transient detection rate and sensitivity of BAT (Tohuvavohu et al., 2020).

This paper focuses on the Swift-BAT subthreshold analysis of a sample of GW triggers with a FAR<<<2 per day, distributed by the LVK Collaboration during O3. The subthreshold GW alerts were received by the EM follow-up groups that were part of a Memorandum of Understanding with the LVK Collaboration. For candidates found with CBC search pipelines (Dal Canton et al., 2021; Sachdev et al., 2019; Messick et al., 2017; Aubin et al., 2021; Nitz et al., 2018; Hooper et al., 2012) and Burst search pipelines (Klimenko et al., 2005, 2016), the alerts contain basic information about the GW FAR, the probability of the candidate being astrophysical (pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT) and trigger time. In the case of CBC candidates, the alerts received in low latency report the pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT split in the four CBC classes: BBH, BNS, NSBH, and Mass Gap. The Mass Gap category includes CBC candidates in which at least one component has a mass in the range [3–5] M⊙subscript𝑀direct-productM_{\odot}italic_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT. Using the GW trigger information received in low latency, we can further perform a search for associations in BAT data.

From the list of 1552 alerts that were communicated via low latency channels, we obtained successful GUANO data dumps for 636 triggers. The GW information of the candidates received in low latency are reported in Table  LABEL:low-latency. The FAR and the pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT classification reported here correspond to the preferred event, namely the one with the highest SNR. Swift-BAT event mode data coincident with the GW trigger time were made available for these triggers in real time via the GUANO data dumps. Post-processing on the data was carried out on this sample from O3 using the NITRATES analysis pipeline (see Section 4).

Out of the 636 low-latency alerts, a total of 86 GW candidates have been confirmed by the offline analysis and included in the GWTC-2.1 and GWTC-3 data releases (Abbott et al., 2024, 2023). Among the 86 confirmed candidates, 14 triggers have pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5 and 72 triggers have pastro<0.5subscript𝑝astro0.5p_{\rm astro}<0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT < 0.5. We indicate the details of the confirmed candidates with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5 and pastro<0.5subscript𝑝astro0.5p_{\rm astro}<0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT < 0.5 in Table 2 and Table LABEL:tab:confirmed-subthresh, respectively. The values of FAR, pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT and CBC Class given in Table 2 and Table LABEL:tab:confirmed-subthresh are derived from offline analyses, as reported in GWTC-2.1 and GWTC-3 data releases, hence are considered more reliable than the values obtained in low latency, reported in Table LABEL:low-latency. The FAR and the pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT classification reported in Tables 2 and LABEL:tab:confirmed-subthresh come from the pipeline that gives the highest value of pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT. For high significance events, if multiple pipelines derive a pastro≃1similar-to-or-equalssubscript𝑝astro1p_{\rm astro}\simeq 1italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT ≃ 1, we select the one with highest SNR. According to these rules, in the case of S200225q, reported in Table 2, the selected pipeline is cWB, but the event is classified as CBC. The CBC class is determined by the highest among pBBHsubscript𝑝BBHp_{\rm BBH}italic_p start_POSTSUBSCRIPT roman_BBH end_POSTSUBSCRIPT, pNSBHsubscript𝑝NSBHp_{\rm NSBH}italic_p start_POSTSUBSCRIPT roman_NSBH end_POSTSUBSCRIPT and pBNSsubscript𝑝BNSp_{\rm BNS}italic_p start_POSTSUBSCRIPT roman_BNS end_POSTSUBSCRIPT. In Tables 2 and LABEL:tab:confirmed-subthresh we report the value of pClasssubscript𝑝Classp_{\rm Class}italic_p start_POSTSUBSCRIPT roman_Class end_POSTSUBSCRIPT defined as max⁡[pBBH,pNSBH,pBNS]subscript𝑝BBHsubscript𝑝NSBHsubscript𝑝BNS\max[p_{\rm BBH},p_{\rm NSBH},p_{\rm BNS}]roman_max [ italic_p start_POSTSUBSCRIPT roman_BBH end_POSTSUBSCRIPT , italic_p start_POSTSUBSCRIPT roman_NSBH end_POSTSUBSCRIPT , italic_p start_POSTSUBSCRIPT roman_BNS end_POSTSUBSCRIPT ]. It is possible that some candidates marked as “Burst” in Table  LABEL:low-latency are then re-classified as “CBC” by the offline analysis. Therefore, for each catalog event the most updated group is the one reported in Tables 2 and LABEL:tab:confirmed-subthresh.

In Figure 1 we show the histograms of the pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT probabilities for all the low-latency CBC candidates processed by GUANO, for a total of 424 candidates, divided into 67 BBH, 130 BNS, 148 NSBH, and 79 Mass Gap events (Fig. 2, left panel). In the offline post-processing of GW candidates, the Mass Gap classification was removed, classifying an object with M>3⁢M⊙𝑀3subscript𝑀direct-productM>3M_{\odot}italic_M > 3 italic_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT as a black hole. This implies that all the CBC events with at least one component mass in the range [3−5]⁢M⊙delimited-[]35subscript𝑀direct-product[3-5]M_{\odot}[ 3 - 5 ] italic_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT are re-distributed into the BBH and NSBH classes. The distribution of the updated pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT classification released by LVK is over-plotted in Figure 1, and the classes are divided in 39 BBH, 22 BNS and 17 NSBH candidates (Fig. 2, right panel).

Figure 1: Distribution of the pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT values for the CBC triggers detected during O3, with available GUANO data dumps. We distinguish with different colors the triggers received in low latency and the ones confirmed by offline analysis.

Figure 2: Left: distribution of the CBC triggers from O3 received in low latency, which had successful GUANO data dumps, divided in the BBH, NSBH, BNS, and Mass Gap classes. Right: analogous distribution for the CBC candidates confirmed by the offline analysis. In the post-processing, the Mass Gap classification has been subsumed into BBH and NSBH.

For the selection of the GW sky localization for each candidate we adopt the following scheme:

Above–threshold GW candidate: The GW candidate is contained in the list of high-probability (pastro\>0.5subscript𝑝astro0.5p\_{\\rm astro}>0.5italic\_p start\_POSTSUBSCRIPT roman\_astro end\_POSTSUBSCRIPT > 0.5) candidates reported in Table 2 of GWTC-2.1 (Abbott et al., [2024](https://arxiv.org/html/2407.12867v1#bib.bib12)) or Table 1 of GWTC-3 (Abbott et al., [2023](https://arxiv.org/html/2407.12867v1#bib.bib11)). The GW sky localizations are downloaded from the parameter estimation data releases of GWTC-2.1 and GWTC-3. We use the results derived from a combination of IMRPhenomXPHM (Pratten et al., [2021](https://arxiv.org/html/2407.12867v1#bib.bib66)) and SEOBNRv4PHM (Ossokine et al., [2020](https://arxiv.org/html/2407.12867v1#bib.bib62)) waveforms (labeled as ‘Mixed’ in the release).

Subthreshold GW candidate: The GW candidate is classified as low-probability (pastro<0.5subscript𝑝astro0.5p\_{\\rm astro}<0.5italic\_p start\_POSTSUBSCRIPT roman\_astro end\_POSTSUBSCRIPT < 0.5) in the offline analysis. The GW sky localization is produced by BAYESTAR (Singer & Price, [2016](https://arxiv.org/html/2407.12867v1#bib.bib72); Singer et al., [2016](https://arxiv.org/html/2407.12867v1#bib.bib73)) and is taken from the GWTC-3 release, which contains both O3b events and updated O3a events. If multiple events are present for a single GW candidate, the pipeline with the highest pastrosubscript𝑝astrop\_{\\rm astro}italic\_p start\_POSTSUBSCRIPT roman\_astro end\_POSTSUBSCRIPT is considered for the selection of the sky localization.

Non-confirmed low-latency GW candidate: The GW candidate has an associated low-latency alert, but the event has not been confirmed by the offline analysis. The GW sky localization is downloaded from GraceDB, selecting the preferred event.

Targeted searches are carried out in real time for all types of transients such as GRBs, FRBs, neutrino events as well as GW triggers, on the event mode data obtained using GUANO. The targeted search pipeline that is currently operational is NITRATES. This is a maximum–likelihood framework that forward models signals through the entire instrument response (DeLaunay & Tohuvavohu, 2022). The BAT responses are created by simulating the photon paths through all detector segments using Geant4, which is a particle-interaction simulator software toolkit (Allison et al., 2016). Unlike the standard BAT responses, the NITRATES responses account for all the detectors on the focal plane, regardless of their coding by the mask. The responses also encode details on the gamma-ray interactions inside the instrument, which carry additional information. This approach enables substantial sensitivity gain compared to the conventional technique of cross-correlation imaging, which translates to a 50% increase in the detection horizon distance for a GRB 170817A-like burst compared to the onboard imaging. Details of the NITRATES response generation and the analysis pipeline are provided in DeLaunay & Tohuvavohu (2022).

A GRB-like transient signal is described using the sky localization and parameters that are specific to the assumed spectral model (the peak energy and spectral slope). This framework then computes the significance of each signal using a test statistic (TS) by maximizing the log-likelihood (LLH) as a function of signal parameters. The likelihood ratio test statistic ΛΛ\Lambdaroman_Λ is used to compare the source signal+background model (described using a set of parameters, ΘsigsuperscriptΘsig\Theta^{\text{sig}}roman_Θ start_POSTSUPERSCRIPT sig end_POSTSUPERSCRIPT that maximizes the LLH) to a background-only model (described by the set of parameters ΘbkgoffsubscriptsuperscriptΘoffbkg\Theta^{\text{off}}_{\text{bkg}}roman_Θ start_POSTSUPERSCRIPT off end_POSTSUPERSCRIPT start_POSTSUBSCRIPT bkg end_POSTSUBSCRIPT, that maximizes the LLH in the off-time window) and is defined as follows (DeLaunay & Tohuvavohu, 2022):

Λ=−2⁢[LLH⁢(Θbkgoff|Non)−LLH⁢(Θsig,Θbkgoff|Non)],Λ2delimited-[]LLHconditionalsubscriptsuperscriptΘoffbkgsuperscriptNonLLHsubscriptΘsigconditionalsubscriptsuperscriptΘoffbkgsuperscriptNon\Lambda=-2[\rm LLH(\Theta^{\text{off}}_{\text{bkg}}|N^{\text{on}})-LLH(\Theta_% {\text{sig}},\Theta^{\text{off}}_{\text{bkg}}|N^{\text{on}})],roman_Λ = - 2 [ roman_LLH ( roman_Θ start_POSTSUPERSCRIPT off end_POSTSUPERSCRIPT start_POSTSUBSCRIPT bkg end_POSTSUBSCRIPT | roman_N start_POSTSUPERSCRIPT on end_POSTSUPERSCRIPT ) - roman_LLH ( roman_Θ start_POSTSUBSCRIPT sig end_POSTSUBSCRIPT , roman_Θ start_POSTSUPERSCRIPT off end_POSTSUPERSCRIPT start_POSTSUBSCRIPT bkg end_POSTSUBSCRIPT | roman_N start_POSTSUPERSCRIPT on end_POSTSUPERSCRIPT ) ] ,

(1)

where Nonsuperscript𝑁onN^{\text{on}}italic_N start_POSTSUPERSCRIPT on end_POSTSUPERSCRIPT is on-time data.

The search pipeline workflow can be summarized as follows:

The event mode data is cleaned and filtered to discard potential glitches and artifacts from cosmic rays, and to flag poorly behaving detectors. Good time intervals (GTIs) are determined, where there is quality data for the analysis.

A time window of 50 s (from the pre- and post-trigger intervals) is identified as the background interval. It is then utilized to model contributions to the background from known bright sources and diffuse sources.

To narrow down the search parameter space, a set of simple analyses are performed to select a list of interesting start times and durations (hereafter referred to as time seeds) as well as portions of the BAT FOV (position seeds).

Finally, the log-likelihoods are computed for all parameters corresponding to the shortlisted time and position seeds.

In essence, the set of signal parameters that maximizes the log-likelihood is the most preferred set of parameters.

The NITRATES likelihood analysis outperforms the onboard mask-weighted imaging analysis by delivering superior sensitivity, given the increased effective area (see Fig. 2 in DeLaunay & Tohuvavohu 2022). At the cost of a significantly increased computational time, this method is capable of delivering arcminute scale localization for events that fall inside the BAT FOV, even when the transient event does not trigger Swift-BAT onboard (Tohuvavohu et al. 2021b; Tohuvavohu 2023; DeLaunay et al. 2022b). The NITRATES pipeline has the ability to distinguish between bursts that come from in and outside the BAT FOV. NITRATES has also accurately localized sufficiently bright bursts outside the FOV (DeLaunay & Tohuvavohu, 2022).

Swift-BAT GUANO was operating during the O3 and was successfully procuring event mode information in response to GW subthreshold triggers (Tohuvavohu et al., 2020). We describe the targeted search analysis that has been carried out using the NITRATES version 0.0.1 which was available in early 2022. The targeted search analysis that was operational in O3 corresponded to a preliminary version of the NITRATES code, that has since undergone several stages of development. The most updated version is publicly available on GitHub.

During O3, for a total of 636 GW triggers, GUANO dumped either 200 s or 90 s of event mode data, for public triggers and for privately communicated triggers, respectively. The choice of the width of the temporal window is made to avoid an overload of downlink data in the process of GUANO data dump. The targeted search pipeline was run in a time window of ±plus-or-minus\pm±20 s centered around the trigger time. The search was carried out on 8 time bins (0.128 s, 0.256 s, 0.512 s, 1.024 s, 2.048 s, 4.096 s, 8.192 s and 16.384 s) and 9 energy bins (between 15–350 keV). The results from the search are reported using the following set of parameters: 1) the maximum TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG describes the statistical significance of a potential detection (see Section 5); 2) Δ⁢LLHoutΔsubscriptLLHout\Delta\rm LLH_{out}roman_Δ roman_LLH start_POSTSUBSCRIPT roman_out end_POSTSUBSCRIPT indicates the preference of the search to a location inside or outside the BAT FOV, and 3) Δ⁢LLHpeakΔsubscriptLLHpeak\Delta\rm LLH_{peak}roman_Δ roman_LLH start_POSTSUBSCRIPT roman_peak end_POSTSUBSCRIPT indicates the confidence of the search in localizing the source to arcminute scales.

The NITRATES search was performed on the ROAR supercomputing cluster on a set of 200 virtual cores for a total of ∼600×800similar-toabsent600800\sim 600\times 800∼ 600 × 800 CPU hours for the entire GW sample.

The targeted search analysis provides a list of top candidates whose spatial, temporal, and spectral parameters maximize the log-likelihood. In order for a candidate to be qualified as a confident detection, we require that the resulting detection significance parameter TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG must exceed the threshold value of 8, corresponding to a FAR ∼4×10−5similar-toabsent4superscript105\sim 4\times 10^{-5}∼ 4 × 10 start_POSTSUPERSCRIPT - 5 end_POSTSUPERSCRIPT Hz. Being a targeted search, the NITRATES analysis can give a false positive with TS>8TS8\sqrt{\text{TS}}>8square-root start_ARG TS end_ARG > 8 with a probability which follows a Poissionian distribution:

P⁢(Ndet≥1)=1−P⁢(Ndet=0)=1−e−FAR×Δ⁢t,𝑃subscript𝑁det11𝑃subscript𝑁det01superscript𝑒FARΔ𝑡P(N_{\text{det}}\geq 1)=1-P(N_{\text{det}}=0)=1-e^{-\text{FAR}\times\Delta t},italic_P ( italic_N start_POSTSUBSCRIPT det end_POSTSUBSCRIPT ≥ 1 ) = 1 - italic_P ( italic_N start_POSTSUBSCRIPT det end_POSTSUBSCRIPT = 0 ) = 1 - italic_e start_POSTSUPERSCRIPT - FAR × roman_Δ italic_t end_POSTSUPERSCRIPT ,

(2)

with Δ⁢t=40Δ𝑡40\Delta t=40roman_Δ italic_t = 40 s being the width of the search window. This leads to a pre-trial p-value of 1.6×10−31.6superscript1031.6\times 10^{-3}1.6 × 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT. Since the NITRATES analysis is performed on all GW triggers with a FAR <2absent2<2< 2 day-1, and considering that there are NGW−search=5subscript𝑁GWsearch5N_{\rm GW-search}=5italic_N start_POSTSUBSCRIPT roman_GW - roman_search end_POSTSUBSCRIPT = 5 independent GW pipelines, the rate of expected false positive candidates with TS>8TS8\sqrt{\text{TS}}>8square-root start_ARG TS end_ARG > 8 is ∼5×2×1.6×10−3similar-toabsent521.6superscript103\sim 5\times 2\times 1.6\times 10^{-3}∼ 5 × 2 × 1.6 × 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT /day ∼similar-to\sim∼ 1/(60 day).

For the entire sample of 636 low-latency triggers processed using NITRATES, we have no candidates that qualify as detection of a signal of astrophysical origin. None of the top candidates within the ±20plus-or-minus20\pm 20± 20 s search window are coincident with the GW triggers. A temporal coincidence with a GW trigger is claimed if the NITRATES search finds a candidate with TS>8TS8\sqrt{\text{TS}}>8square-root start_ARG TS end_ARG > 8 and |t0−tstart|<20subscript𝑡0subscript𝑡start20|t_{0}-t_{\rm start}|<20| italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - italic_t start_POSTSUBSCRIPT roman_start end_POSTSUBSCRIPT | < 20 s, where t0subscript𝑡0t_{0}italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is the GW trigger time and tstartsubscript𝑡startt_{\rm start}italic_t start_POSTSUBSCRIPT roman_start end_POSTSUBSCRIPT is the starting time of the temporal bin with highest ranking statistics. A detailed list of all the NITRATES results for the entire sample analyzed during O3 is provided in Table LABEL:low-latency. We discuss specific false positive candidates in Section 5.1.

If the GW trigger time is included in the time window corresponding to slew mode of BAT, the analysis cannot be performed using NITRATES since the targeted search requires stable attitude information to compute the background. Similarly, some triggers have insufficient exposure time, preventing the NITRATES analysis. In this case neither TS results nor flux upper limits can be computed. As a cut to narrow down the parameter space, the targeted search selects time seeds as described in Section 4. If there are no time seeds that pass the preliminary cuts then there will be no final likelihood computations. Results for these types of triggers are indicated as NFL (No Final Likelihood) in Table LABEL:low-latency. In the case of NFL triggers, though, the flux upper limit can be computed, since a full likelihood analysis is not required.

We did not find any candidate associations from any of the triggers with BAT. However, the targeted search pipeline did result in the detection of six candidates with a significance above the NITRATES detection threshold of TS=8TS8\sqrt{\text{TS}}=8square-root start_ARG TS end_ARG = 8. These candidates were examined to understand our false positive population. S200327j (TS∼22similar-toTS22\sqrt{\text{TS}}\sim 22square-root start_ARG TS end_ARG ∼ 22), S200324ax (TS∼11similar-toTS11\sqrt{\text{TS}}\sim 11square-root start_ARG TS end_ARG ∼ 11), and S200225af (TS∼10.5similar-toTS10.5\sqrt{\text{TS}}\sim 10.5square-root start_ARG TS end_ARG ∼ 10.5) are triggers that occurred during the passage of Swift in the proximity of the South Atlantic Anomaly (SAA). The background characterization becomes unreliable when the spacecraft is either entering or leaving the SAA on account of increased background contamination. In S190919au, a peculiar dip (∼similar-to\sim∼20 s) in the background may have contributed to a false elevation in the signal detection statistic, by causing an under representation of the background rate. For S190919u, we obtain a TS∼8.0similar-toTS8.0\sqrt{\text{TS}}\sim 8.0square-root start_ARG TS end_ARG ∼ 8.0, which corresponds to a FAR of ∼4×10−5similar-toabsent4superscript105\sim 4\times 10^{-5}∼ 4 × 10 start_POSTSUPERSCRIPT - 5 end_POSTSUPERSCRIPT Hz. The presence of such a detection, in a total sample of 636 triggers that were analyzed, is compatible with the expected number of false positives, which corresponds to (4×10−5(4\times 10^{-5}( 4 × 10 start_POSTSUPERSCRIPT - 5 end_POSTSUPERSCRIPT Hz) ×(40\times(40× ( 40 s)×)\times) × 636 ∼1similar-toabsent1\sim 1∼ 1.

S200130ai corresponds to a subthreshold GW trigger at T0=subscript𝑇0absentT_{0}=italic_T start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 2020-01-30T09:59:58 that was identified by the CBC search as a NSBH candidate with a pastro=0.008subscript𝑝astro0.008p_{\rm astro}=0.008italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT = 0.008 and a GW FAR ∼1.8×\sim 1.8\times∼ 1.8 ×10-5 Hz. It was detected using NITRATES at a significance of TS∼16.3similar-toTS16.3\sqrt{\text{TS}}\sim 16.3square-root start_ARG TS end_ARG ∼ 16.3 with a Δ⁢LLHout=−19.68ΔsubscriptLLHout19.68\Delta\rm LLH_{out}=-19.68roman_Δ roman_LLH start_POSTSUBSCRIPT roman_out end_POSTSUBSCRIPT = - 19.68 and Δ⁢LLHpeak=2.14ΔsubscriptLLHpeak2.14\Delta\rm LLH_{peak}=2.14roman_Δ roman_LLH start_POSTSUBSCRIPT roman_peak end_POSTSUBSCRIPT = 2.14, consistent with a sky localization outside the BAT-FOV. The highest log-likelihood candidate, was identified to arise 1.5 s prior to the GW trigger time. Due to the low value of Δ⁢LLHpeakΔsubscriptLLHpeak\rm\Delta LLH_{peak}roman_Δ roman_LLH start_POSTSUBSCRIPT roman_peak end_POSTSUBSCRIPT, we do not have an arcminute-level precision on the sky localization. The candidate was associated with a Fermi trigger 602071201 (GCN 26944, Fermi GBM Team 2020) and was classified as a long GRB. The Fermi localization is RA =137.5absent137.5=137.5= 137.5 deg, Dec =−51.3absent51.3=-51.3= - 51.3 deg, with a statistical uncertainty of 3.5 degrees. The Interplanetary Gamma-Ray Burst Timing Network (IPN) further localized the event in a 3-sigma error box with an area of 1487 arcmin2 and centered at an RA =134.742absent134.742=134.742= 134.742 deg and Dec =−49.627absent49.627=-49.627= - 49.627 deg (Hurley et al., 2020). Although this event presents a temporal coincidence with the GW trigger, on account of the lack of spatial coincidence with the GW location, this event is discarded from being associated with the GW subthreshold trigger. Additionally, this low-latency GW candidate has not been confirmed by offline analyses.

Figure 3: Flux upper limit maps are shown for all the O3 catalog events with a pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5 that were processed successfully using NITRATES. The color bar indicates the upper limit in the 15–350 keV Swift-BAT band as a function of the sky position. The part of the sky in white corresponds to the area covered by the Earth. The solid and dashed contours are the GW 90% and 50% credible levels, respectively.

Figure 4: The 15–350 keV flux upper limit ϕULsubscriptitalic-ϕUL\phi_{\rm UL}italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT derived with NITRATES as a function of εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT, namely the probability that the GW candidate is contained inside the BAT FOV. The plot includes all the GW candidates received during O3, including both Burst and CBC events. With different symbols we distinguish the GW candidates received in low latency and the ones confirmed and included in the O3 catalog, separated in pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5 and pastro<0.5subscript𝑝astro0.5p_{\rm astro}<0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT < 0.5. The dashed and solid red lines are the 50%percent\%% and 90%percent\%% containment regions of the scatter plot, respectively. The one-dimensional histograms of ϕULsubscriptitalic-ϕUL\phi_{\rm UL}italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT and εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT are reported on the sides. The color bar indicates the value of ε⊕subscript𝜀⊕\varepsilon_{\earth}italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT, the probability that the GW candidate is occulted by the Earth.

Since each GW trigger processed in this analysis resulted in a non-detection in Swift-BAT, we estimate the flux upper limits in the following manner. The NITRATES analysis generates rates curves in the 15–350 keV energy band from the GTIs of the filtered event list. The number of active detectors corresponding to each trigger is read out from its respective detector mask file. A linear fit is carried out to the background window of duration 50 s. We then estimate the 5σ𝜎\sigmaitalic_σ count rate and the corresponding uncertainty over the full signal window which has a ±plus-or-minus\pm±20 s duration. This is computed for all the 8 time bins (see Section 4). We further convert the 5σ𝜎\sigmaitalic_σ count rates to flux upper limits, as a function of sky position, in the 15–350 keV band, by convolving different spectral models with the NITRATES responses for each time bin iteration. We select 929 grid points on the sky and interpolate upper limit values for locations in between. We assume the following different spectral templates:

Band function (Band et al., [1993](https://arxiv.org/html/2407.12867v1#bib.bib18)) with a soft template (Epeak\=70subscript𝐸peak70E\_{\\rm peak}=70italic\_E start\_POSTSUBSCRIPT roman\_peak end\_POSTSUBSCRIPT = 70 keV, α\=−1.9𝛼1.9\\alpha=-1.9italic\_α = - 1.9, β\=−3.7𝛽3.7\\beta=-3.7italic\_β = - 3.7)

Band function with a normal template (Epeaksubscript𝐸peakE\_{\\rm peak}italic\_E start\_POSTSUBSCRIPT roman\_peak end\_POSTSUBSCRIPT = 230 keV, α\=−𝛼\\alpha=-italic\_α = -1.0, β\=−𝛽\\beta=-italic\_β = -2.3)

Cutoff power law function with a hard template (Epeak\=subscript𝐸peakabsentE\_{\\rm peak}=italic\_E start\_POSTSUBSCRIPT roman\_peak end\_POSTSUBSCRIPT = 1500 keV, α\=𝛼absent\\alpha=italic\_α = 1.5)

Cutoff power law function that has been used to describe GRB 170817A (Epeaksubscript𝐸peakE\_{\\rm peak}italic\_E start\_POSTSUBSCRIPT roman\_peak end\_POSTSUBSCRIPT = 185 keV, α𝛼\\alphaitalic\_α = 0.62) (Goldstein et al., [2017b](https://arxiv.org/html/2407.12867v1#bib.bib41))

The parameters α𝛼\alphaitalic_α and β𝛽\betaitalic_β correspond to the low-energy and high-energy photon indices of the spectrum, respectively. The first three spectral templates are identical to the ones that are routinely adopted by Fermi-GBM (Goldstein et al., 2016a). In the rest of the paper, all the results are reported assuming a Normal spectral template.

Calling Ω=(RA,Dec)ΩRADec\Omega=(\rm RA,Dec)roman_Ω = ( roman_RA , roman_Dec ) the coordinates variable, for each temporal bin and spectral template we convert the upper limit map ϕUL⁢(Ω)subscriptitalic-ϕULΩ\phi_{\rm UL}(\Omega)italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( roman_Ω ) into a unique marginalized upper limit value:

ϕUL=∫Ω∉Ω⊕ϕUL⁢(Ω)⁢PGW⁢(Ω)⁢dΩ,subscriptitalic-ϕULsubscriptΩsubscriptΩ⊕subscriptitalic-ϕULΩsubscript𝑃GWΩdifferential-dΩ\phi_{\rm UL}=\int_{\Omega\notin\Omega_{\earth}}\phi_{\rm UL}(\Omega)P_{\rm GW% }(\Omega)\rm d\Omega,italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT = ∫ start_POSTSUBSCRIPT roman_Ω ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( roman_Ω ) italic_P start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT ( roman_Ω ) roman_d roman_Ω ,

(3)

where PGW⁢(Ω)subscript𝑃GWΩP_{\rm GW}(\Omega)italic_P start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT ( roman_Ω ) is the posterior probability distribution of the GW sky position. The notation Ω∉Ω⊕ΩsubscriptΩ⊕\Omega\notin\Omega_{\earth}roman_Ω ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT means that the integral is limited to the region of the sky not occulted by the Earth. We report in Table LABEL:low-latency the marginalized flux upper limits for a 1 s time bin, assuming the normal spectral template. In Fig. 3, we provide the sky maps reporting both the flux upper limits as a function of sky position and the GW contours (50%percent\%% and 90%percent\%% credible levels) for the GW candidates with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5.

As additional information, we also report the quantity εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT, which quantifies the probability that the GW source is inside the BAT coded FOV and corresponds to

εin⁢BAT=∫Ω∈ΩinPGW⁢(Ω)⁢dΩ,subscript𝜀inBATsubscriptΩsubscriptΩinsubscript𝑃GWΩdifferential-dΩ\varepsilon_{\rm in\,BAT}=\int_{\Omega\in\Omega_{\rm in}}P_{\rm GW}(\Omega)\rm d\Omega,italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT = ∫ start_POSTSUBSCRIPT roman_Ω ∈ roman_Ω start_POSTSUBSCRIPT roman_in end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT ( roman_Ω ) roman_d roman_Ω ,

(4)

where the integral is limited to the solid angle ΩinsubscriptΩin\Omega_{\rm in}roman_Ω start_POSTSUBSCRIPT roman_in end_POSTSUBSCRIPT, namely the portion of the sky where the BAT partial coding fraction is larger than 0.01. The location of ΩinsubscriptΩin\Omega_{\rm in}roman_Ω start_POSTSUBSCRIPT roman_in end_POSTSUBSCRIPT, i.e., the BAT FOV, is identified by the yellow region in the sky maps of Figure 3. The higher the εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT, the better the BAT covering of the GW error region, and the more constraining the derived upper limit. The flux upper limits as a function of εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT for all GW trigger candidates is shown in Fig. 4. We also indicate with different markers the sample of low-latency triggers, the confirmed list of subthreshold candidates (pastro<0.5subscript𝑝astro0.5p_{\rm astro}<0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT < 0.5), and the above-threshold candidates (pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5). In Table LABEL:low-latency we also report the probability that the GW source is occulted by the Earth, defined as:

ε⊕=∫Ω∈Ω⊕PGW⁢(Ω)⁢dΩ,subscript𝜀⊕subscriptΩsubscriptΩ⊕subscript𝑃GWΩdifferential-dΩ\varepsilon_{\earth}=\int_{\Omega\in\Omega_{\earth}}P_{\rm GW}(\rm\Omega)d\Omega,italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT = ∫ start_POSTSUBSCRIPT roman_Ω ∈ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT ( roman_Ω ) roman_d roman_Ω ,

(5)

where Ω⊕subscriptΩ⊕\Omega_{\earth}roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT is the solid angle subtended by the Earth from the Swift reference system.

We further convert the flux upper limits into luminosity upper limits for all the GW triggers with available information about the distance posterior distribution, namely only triggers identified by CBC searches. The luminosity upper limit in the rest frame band 1 keV–10 MeV is estimated as

LUL=⟨4⁢π⁢DL2⁢k⁢ϕUL⟩,subscript𝐿ULdelimited-⟨⟩4𝜋superscriptsubscript𝐷L2𝑘subscriptitalic-ϕULL_{\rm UL}=\langle 4\pi D_{\rm L}^{2}k\phi_{\rm UL}\rangle,italic_L start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT = ⟨ 4 italic_π italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_k italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ⟩ ,

(6)

where DLsubscript𝐷LD_{\rm L}italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT is extracted from the posterior probability P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ) reported in the GW sky localization files, while k𝑘kitalic_k is the k-correction and corresponds to

k=I⁢[1⁢keV/(1+z),10⁢MeV/(1+z)]I⁢[15⁢keV,350⁢keV],𝑘𝐼1keV1𝑧10MeV1𝑧𝐼15keV350keVk=\frac{I[1\,\mathrm{keV}/(1+z),10\,\mathrm{MeV}/(1+z)]}{I[15\,\mathrm{keV},35% 0\,\mathrm{keV}]},italic_k = divide start_ARG italic_I [ 1 roman_keV / ( 1 + italic_z ) , 10 roman_MeV / ( 1 + italic_z ) ] end_ARG start_ARG italic_I [ 15 roman_keV , 350 roman_keV ] end_ARG ,

(7)

where

I⁢[a,b]=∫abE⁢d⁢Nd⁢E⁢(E)⁢dE,𝐼𝑎𝑏superscriptsubscript𝑎𝑏𝐸d𝑁d𝐸𝐸differential-d𝐸I[a,b]=\int_{a}^{b}E\frac{\mathrm{d}N}{\mathrm{~{}d}E}(E)\mathrm{d}E,italic_I [ italic_a , italic_b ] = ∫ start_POSTSUBSCRIPT italic_a end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_b end_POSTSUPERSCRIPT italic_E divide start_ARG roman_d italic_N end_ARG start_ARG roman_d italic_E end_ARG ( italic_E ) roman_d italic_E ,

(8)

and dN𝑁Nitalic_N/dE𝐸Eitalic_E is the assumed photon spectrum. The band 1 keV–10 MeV is chosen since it is usually adopted to report the bolometric luminosity of GRBs.

The luminosity upper limits as a function of the mean value of the luminosity distance is reported in Figure 5. Similar to what was shown previously, we demarcate the various samples. Candidates with a low latency classification of Mass Gap were later re-distributed to other categories as part of post-processing, which is evident from the Mass Gap panel in Figure 5. As expected and as evident already from Figure 4, we see a clear correlation between the luminosity upper limit and εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT in Figure 5, indicating that the inferred constraints on the EM counterpart are more stringent when the GW probability integrated inside the BAT FOV is higher. Since ϕULsubscriptitalic-ϕUL\phi_{\rm UL}italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT is an upper limit and not a measure coming from a detection, Eq. (6) is an approximated method to convert ϕULsubscriptitalic-ϕUL\phi_{\rm UL}italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT in a luminosity upper limit, averaging over the P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ) distribution provided by the GW analysis. In Appendix A we show a more accurate way to estimate the luminosity upper limit, but we find no relevant differences with respect to the method reported in this section. The Eq. (6) is used only to produce the plots of Fig. 5, but this approximation is not used in Section 6 to perform inference about the EM model parameters. Instead, in Section 6 a reverse process is followed, namely the EM model is used to predict the luminosity, which is then convolved with P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ) to obtain a probability distribution of the flux in the BAT energy band.

Figure 5: Upper limits on the luminosity computed in the rest frame 1 keV–10 MeV energy band, as a function of the mean luminosity distance extracted from the sky localization map of each GW candidate. The color bar indicates the quantity εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT, namely the probability that the GW candidate is located inside the BAT coded FOV. With different symbols, we distinguish the GW candidates received in low latency and the ones confirmed and included in the O3 catalog, separated in pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5 and pastro<0.5subscript𝑝astro0.5p_{\rm astro}<0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT < 0.5.

Figure 6: Distribution in the GW FAR–TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG plane of all the triggers which passed the threshold FAR<GRB10−3{}_{\rm GRB}<10^{-3}start_FLOATSUBSCRIPT roman_GRB end_FLOATSUBSCRIPT < 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT Hz. Triggers in the regions below the black and red dashed lines (marked with a green cross) would have triggered the RAVEN alert system, for CBC and Burst events, respectively. The plot does not include all the triggers that have a TS>8TS8\sqrt{\text{TS}}>8square-root start_ARG TS end_ARG > 8 which turned out to be spurious artifacts.

To calculate the joint GW–BAT FARs for each GW trigger, we elaborate on the methods used to compute the individual FARs and subsequently combine them. To derive the sensitivity of the NITRATES search, time-tagged event data assembled from intervals corresponding to calibration runs and data from before and after known GRB signal times (total exposure time of ∼similar-to\sim∼51 ks) were analyzed. The behavior of the background population and its associated FAR were then identified (see Section 7.6 and Fig. 16 in DeLaunay & Tohuvavohu 2022). We further compute the joint Swift-BAT–GW FAR by combining the BAT FAR (calculated using the method described above) with the GW FAR. A targeted joint FAR threshold routine is constructed as part of the Rapid, on-source VOEvent Coincident Monitor (RAVEN; Abbott et al. 2017c; Urban 2016), which combines the FARs obtained from GWs along with those from BAT and computes the joint temporal as well as the joint spatial FAR. We also specify details of the search pipeline used in the process, Burst or CBC. The joint FAR prescription as reported in the RAVEN documentation, is computed as

FARGRB+GW=ZIΩ⁢[1+ln⁡(ZmaxZ)]subscriptFARGRBGW𝑍subscript𝐼Ωdelimited-[]1subscript𝑍𝑍\text{FAR}_{\rm GRB+GW}=\frac{Z}{I_{\Omega}}\left[1+\ln({\frac{Z_{\max}}{Z}})\right]FAR start_POSTSUBSCRIPT roman_GRB + roman_GW end_POSTSUBSCRIPT = divide start_ARG italic_Z end_ARG start_ARG italic_I start_POSTSUBSCRIPT roman_Ω end_POSTSUBSCRIPT end_ARG [ 1 + roman_ln ( divide start_ARG italic_Z start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT end_ARG start_ARG italic_Z end_ARG ) ]

(9)

where Z𝑍Zitalic_Z is the joint ranking statistic given by,

Z=FARGW⁢FARGRB⁢Δ⁢t,𝑍subscriptFARGWsubscriptFARGRBΔ𝑡Z=\text{FAR}_{\rm GW}\text{FAR}_{\rm GRB}\Delta t,italic_Z = FAR start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT FAR start_POSTSUBSCRIPT roman_GRB end_POSTSUBSCRIPT roman_Δ italic_t ,

(10)

Zmax=FARGW,max⁢FARGRB,max⁢Δ⁢t,subscript𝑍maxsubscriptFARGWmaxsubscriptFARGRBmaxΔ𝑡Z_{\rm max}=\text{FAR}_{\rm GW,max}\text{FAR}_{\rm GRB,max}\Delta t,italic_Z start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT = FAR start_POSTSUBSCRIPT roman_GW , roman_max end_POSTSUBSCRIPT FAR start_POSTSUBSCRIPT roman_GRB , roman_max end_POSTSUBSCRIPT roman_Δ italic_t ,

(11)

and we adopt Δ⁢t=30Δ𝑡30\Delta t=30roman_Δ italic_t = 30 s, FARGW,max=subscriptFARGWmaxabsent\text{FAR}_{\rm GW,max}=FAR start_POSTSUBSCRIPT roman_GW , roman_max end_POSTSUBSCRIPT = 2 day-1 and FARGRB,max=10−3subscriptFARGRBmaxsuperscript103\text{FAR}_{\rm GRB,max}=10^{-3}FAR start_POSTSUBSCRIPT roman_GRB , roman_max end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT Hz. IΩsubscript𝐼ΩI_{\Omega}italic_I start_POSTSUBSCRIPT roman_Ω end_POSTSUBSCRIPT is an integral that quantifies the spatial overlap between the GW localization and the GRB localization (Ashton et al., 2018). Even if the search of subthreshold candidates in NITRATES is done in a temporal window [t0−20s,[t_{0}-20\rm~{}s,~{}[ italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - 20 roman_s ,t0+20s]t_{0}+20\rm~{}s]italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT + 20 roman_s ] around the trigger time t0subscript𝑡0t_{0}italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT, for the RAVEN joint alert the adopted temporal window is [t0−10s,[t_{0}-10\rm~{}s,~{}[ italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - 10 roman_s ,t0+20s]t_{0}+20\rm~{}s]italic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT + 20 roman_s ]. Since none of the BAT candidates analyzed in this work has a confident estimation of the sky localization, we adopt a uniform posterior probability on the full sky for the EM candidate. Hence, by definition, we set IΩ=1subscript𝐼Ω1I_{\Omega}=1italic_I start_POSTSUBSCRIPT roman_Ω end_POSTSUBSCRIPT = 1. The candidate triggers a RAVEN alert when the FARGRB+GW×Nt<FARmaxsubscriptFARGRBGWsubscript𝑁𝑡subscriptFARmax\text{FAR}_{\rm GRB+GW}\times N_{t}<\text{FAR}_{\rm max}FAR start_POSTSUBSCRIPT roman_GRB + roman_GW end_POSTSUBSCRIPT × italic_N start_POSTSUBSCRIPT italic_t end_POSTSUBSCRIPT < FAR start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT, with Ntsubscript𝑁𝑡N_{t}italic_N start_POSTSUBSCRIPT italic_t end_POSTSUBSCRIPT being the trials factor of the joint search and FARmax=subscriptFARmaxabsent\text{FAR}_{\rm max}=FAR start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT = (1/30) day-1 for CBC events and FARmax=subscriptFARmaxabsent\text{FAR}_{\rm max}=FAR start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT = 1 yr-1 for Burst events. The trials factor corresponds to Nt=SGW⁢(SGW−1)subscript𝑁𝑡subscript𝑆GWsubscript𝑆GW1N_{t}=S_{\rm GW}(S_{\rm GW}-1)italic_N start_POSTSUBSCRIPT italic_t end_POSTSUBSCRIPT = italic_S start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT ( italic_S start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT - 1 ), where SGWsubscript𝑆GWS_{\rm GW}italic_S start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT is the number of search GW pipelines, 4 for CBC events and 3 for Burst events. Since the GW pipelines are not fully independent, a realistic value of the trials factor is smaller than the one adopted here, therefore the RAVEN threshold can be considered as a conservative estimate for the significance of a joint detection.

We quote the derived joint FARs along with other trigger-specific details only for those triggers with a FARGRB,max<10−3subscriptFARGRBmaxsuperscript103\text{FAR}_{\rm GRB,max}<10^{-3}FAR start_POSTSUBSCRIPT roman_GRB , roman_max end_POSTSUBSCRIPT < 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT Hz in Table 4. We find that, after rejecting false positives, 2 CBC events pass the joint FAR detection threshold to trigger a RAVEN alert. Specifically, S191110x (TS=7.2TS7.2\sqrt{\text{TS}}=7.2square-root start_ARG TS end_ARG = 7.2) and S200108p (TS=7.4TS7.4\sqrt{\text{TS}}=7.4square-root start_ARG TS end_ARG = 7.4) have a joint FAR of 3.02×10−43.02superscript1043.02\times 10^{-4}3.02 × 10 start_POSTSUPERSCRIPT - 4 end_POSTSUPERSCRIPT yr-1 and 21.321.321.321.3 yr-1, respectively. These values are obtained considering the GW FAR received with the low-latency alert. In the offline analysis of the GW candidates, neither S191110x or S200108p have been confirmed. We therefore conclude that, considering the offline joint analysis of GW and Swift-BAT data, none of the candidates is eligible to claim a significant joint detection.

In Fig. 6 we report the location in the GW FAR - TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG plane of all the candidates that pass the condition FARGW,max<subscriptFARGWmaxabsent\text{FAR}_{\rm GW,max}<FAR start_POSTSUBSCRIPT roman_GW , roman_max end_POSTSUBSCRIPT <2 day-1 and FARGRB,max<10−3subscriptFARGRBmaxsuperscript103\text{FAR}_{\rm GRB,max}<10^{-3}FAR start_POSTSUBSCRIPT roman_GRB , roman_max end_POSTSUBSCRIPT < 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT Hz (i.e., TS≳7greater-than-or-equivalent-toTS7\sqrt{\text{TS}}\gtrsim 7square-root start_ARG TS end_ARG ≳ 7), to be considered for a potential joint alert. The astrophysical origin of all the candidates with TS>8TS8\sqrt{\text{TS}}>8square-root start_ARG TS end_ARG > 8 has been rejected as discussed in Section 5.1, and therefore they are not reported in Fig. 6. The dashed black and red lines mark the separation line for the event to pass the RAVEN alert threshold, for CBC and Burst candidates, respectively. Candidates below those lines would have triggered a RAVEN alert.

In this section, we describe how the upper limits derived from the joint subthreshold search can be used to infer constraints about possible EM emission from the GW candidates. Starting from a model of the EM emission, the luminosity in the BAT band can be estimated, whose value will depend on some internal parameters of the model (λ1,…,λk)subscript𝜆1…subscript𝜆𝑘(\lambda_{1},...,\lambda_{k})( italic_λ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_λ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ). The goal is to explore the model parameter space and test if the estimated flux is in agreement with the upper limit constraints derived in this paper.

For this purpose, a knowledge of the distance of the GW candidate is needed, and the GW sky localization is used to extract the posterior distribution P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ). Since only CBC events have such information, Burst events are not considered in this discussion. For the CBC events, we consider a phenomenological model which describes the probability distribution of the luminosity L𝐿Litalic_L (in the 15151515–350350350350 keV rest-frame)

P⁢(L)=(1−f)⁢δ⁢(L=0)+f⁢Π⁢(L).𝑃𝐿1𝑓𝛿𝐿0𝑓Π𝐿P(L)=(1-f)\delta(L=0)+f\Pi(L).italic_P ( italic_L ) = ( 1 - italic_f ) italic_δ ( italic_L = 0 ) + italic_f roman_Π ( italic_L ) .

(12)

Here, the f𝑓fitalic_f parameter is a proxy for the EM-bright nature of the event, i.e., given a CBC source described by a set of GW parameters θ→GWsubscript→𝜃GW\vec{\theta}_{\rm GW}over→ start_ARG italic_θ end_ARG start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT, f⁢(θ→GW)𝑓subscript→𝜃GWf(\vec{\theta}_{\rm GW})italic_f ( over→ start_ARG italic_θ end_ARG start_POSTSUBSCRIPT roman_GW end_POSTSUBSCRIPT ) corresponds to the probability that the EM luminosity of the source is non-zero. On the other hand, Π⁢(L)Π𝐿\Pi(L)roman_Π ( italic_L ) is the intrinsic luminosity function of the EM transient associated with the specific CBC class. In the case of BNS and NSBH candidates, the assumption on Π⁢(L)Π𝐿\Pi(L)roman_Π ( italic_L ) should be informed by our prior knowledge of the luminosity function of merger-driven GRBs. A detailed study of the impact of this work on our knowledge of merger-driven GRBs luminosity function will be reported in a follow-up paper. In this section, instead, we focus only on the BBH class, for which no strong prior exists for Π⁢(L)Π𝐿\Pi(L)roman_Π ( italic_L ). For simplicity and in order to show the constraining power of our joint subthreshold search, we assume that the EM process associated with BBH, if present, produces a universal, viewing angle-independent luminosity L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT. Therefore, in the scenario specified above, we have (λ1,…,λk)=(f,L0)subscript𝜆1…subscript𝜆𝑘𝑓subscript𝐿0(\lambda_{1},...,\lambda_{k})=(f,L_{0})( italic_λ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_λ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ) = ( italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) and

P⁢(L)=(1−f)⁢δ⁢(L=0)+f⁢δ⁢(L−L0)=P⁢(L;f,L0).𝑃𝐿1𝑓𝛿𝐿0𝑓𝛿𝐿subscript𝐿0𝑃𝐿𝑓subscript𝐿0P(L)=(1-f)\delta(L=0)+f\delta(L-L_{0})=P(L;f,L_{0}).italic_P ( italic_L ) = ( 1 - italic_f ) italic_δ ( italic_L = 0 ) + italic_f italic_δ ( italic_L - italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) = italic_P ( italic_L ; italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) .

(13)

Once the model for the EM emission is specified, the probability distribution of the predicted flux is

P⁢(ϕ)=(1−f)⁢δ⁢(ϕ=0)+f⁢PEM⁢(ϕ),𝑃italic-ϕ1𝑓𝛿italic-ϕ0𝑓subscript𝑃EMitalic-ϕP(\phi)=(1-f)\delta(\phi=0)+fP_{\text{EM}}(\phi),italic_P ( italic_ϕ ) = ( 1 - italic_f ) italic_δ ( italic_ϕ = 0 ) + italic_f italic_P start_POSTSUBSCRIPT EM end_POSTSUBSCRIPT ( italic_ϕ ) ,

(14)

where PEM(ϕ)=P(L0/4πkDL2P_{\text{EM}}(\phi)=P(L_{0}/4\pi kD_{\text{L}}^{2}italic_P start_POSTSUBSCRIPT EM end_POSTSUBSCRIPT ( italic_ϕ ) = italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT / 4 italic_π italic_k italic_D start_POSTSUBSCRIPT L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT) is the flux probability distribution in the assumption that the source is EM bright. Hence, for the i𝑖iitalic_i-th candidate, the probability that the predicted flux is below the estimated upper limit ϕ0,isubscriptitalic-ϕ0𝑖\phi_{0,i}italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT corresponds to

Pi⁢(ϕ<ϕ0,i)=(1−f)+f⁢∫0ϕ0,iPi⁢(ϕ)⁢dϕ,subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0𝑖1𝑓𝑓superscriptsubscript0subscriptitalic-ϕ0𝑖subscript𝑃𝑖italic-ϕdifferential-ditalic-ϕP_{i}(\phi<\phi_{0,i})=(1-f)+f\int_{0}^{\phi_{0,i}}P_{i}(\phi)\mathrm{d}\phi,italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ < italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT ) = ( 1 - italic_f ) + italic_f ∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT end_POSTSUPERSCRIPT italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ ) roman_d italic_ϕ ,

(15)

valid in the limit in which the GW candidate is assumed to be real. Therefore, given a candidate GW with a probability of being astrophysical pastro,i=πisubscript𝑝astro𝑖subscript𝜋𝑖p_{\text{astro},i}=\pi_{i}italic_p start_POSTSUBSCRIPT astro , italic_i end_POSTSUBSCRIPT = italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, there are only three possibilities to have a non-detection in BAT:

The source is not astrophysical, with a probability 1−πi1subscript𝜋𝑖1-\\pi\_{i}1 - italic\_π start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT;

The source is astrophysical, but it is occulted by the Earth, with a probability πi⁢ε⊕subscript𝜋𝑖subscript𝜀⊕\\pi\_{i}\\varepsilon\_{\\earth}italic\_π start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT italic\_ε start\_POSTSUBSCRIPT ⊕ end\_POSTSUBSCRIPT;

The source is astrophysical, it is not occulted by the Earth and the predicted flux by the EM model is below the BAT upper limit, with a probability πi⁢(1−ε⊕)⁢Pi⁢(ϕ<ϕ0,i)subscript𝜋𝑖1subscript𝜀⊕subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0𝑖\\pi\_{i}(1-\\varepsilon\_{\\earth})P\_{i}(\\phi<\\phi\_{0,i})italic\_π start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ( 1 - italic\_ε start\_POSTSUBSCRIPT ⊕ end\_POSTSUBSCRIPT ) italic\_P start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ( italic\_ϕ < italic\_ϕ start\_POSTSUBSCRIPT 0 , italic\_i end\_POSTSUBSCRIPT ).

This allows us to define a non-detection likelihood corresponding to

ℒi=(1−πi)+πi⁢[ε⊕+(1−ε⊕)⁢Pi⁢(ϕ<ϕ0,i)].subscriptℒ𝑖1subscript𝜋𝑖subscript𝜋𝑖delimited-[]subscript𝜀⊕1subscript𝜀⊕subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0𝑖\mathcal{L}_{i}=(1-\pi_{i})+\pi_{i}[\varepsilon_{\earth}+(1-\varepsilon_{% \earth})P_{i}(\phi<\phi_{0,i})].caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = ( 1 - italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) + italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT [ italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT + ( 1 - italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT ) italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ < italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT ) ] .

(16)

In the case of L0→0→subscript𝐿00L_{0}\rightarrow 0italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT → 0, Pi⁢(ϕ<ϕ0,i)→1→subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0𝑖1P_{i}(\phi<\phi_{0,i})\rightarrow 1italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ < italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT ) → 1, so ℒi→1→subscriptℒ𝑖1\mathcal{L}_{i}\rightarrow 1caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT → 1. For very large values of L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT, instead, Pi⁢(ϕ<ϕ0,i)→0→subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0𝑖0P_{i}(\phi<\phi_{0,i})\rightarrow 0italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ < italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT ) → 0 and therefore ℒi→(1−πi)+πi⁢ε⊕→subscriptℒ𝑖1subscript𝜋𝑖subscript𝜋𝑖subscript𝜀⊕\mathcal{L}_{i}\rightarrow(1-\pi_{i})+\pi_{i}\varepsilon_{\earth}caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT → ( 1 - italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) + italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT. This last result shows how, even if the luminosity predicted by the model is exceedingly large, a non-detection can occur if the GW source is not real (1−πi1subscript𝜋𝑖1-\pi_{i}1 - italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT), or if it is real but occulted by the Earth (πi⁢ε⊕subscript𝜋𝑖subscript𝜀⊕\pi_{i}\varepsilon_{\earth}italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT). Since the analysis is focused only on BBH events, we consider only those candidates that have pBBH>pNSBH,pBNSsubscript𝑝BBHsubscript𝑝NSBHsubscript𝑝BNSp_{\rm BBH}>p_{\rm NSBH},p_{\rm BNS}italic_p start_POSTSUBSCRIPT roman_BBH end_POSTSUBSCRIPT > italic_p start_POSTSUBSCRIPT roman_NSBH end_POSTSUBSCRIPT , italic_p start_POSTSUBSCRIPT roman_BNS end_POSTSUBSCRIPT. By definition, pBBH+pNSBH+pBNS=pastrosubscript𝑝BBHsubscript𝑝NSBHsubscript𝑝BNSsubscript𝑝astrop_{\rm BBH}+p_{\rm NSBH}+p_{\rm BNS}=p_{\rm astro}italic_p start_POSTSUBSCRIPT roman_BBH end_POSTSUBSCRIPT + italic_p start_POSTSUBSCRIPT roman_NSBH end_POSTSUBSCRIPT + italic_p start_POSTSUBSCRIPT roman_BNS end_POSTSUBSCRIPT = italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT and typically for the candidates classified as BBH we have that pBBH≫pNSBH,pBNSmuch-greater-thansubscript𝑝BBHsubscript𝑝NSBHsubscript𝑝BNSp_{\rm BBH}\gg p_{\rm NSBH},p_{\rm BNS}italic_p start_POSTSUBSCRIPT roman_BBH end_POSTSUBSCRIPT ≫ italic_p start_POSTSUBSCRIPT roman_NSBH end_POSTSUBSCRIPT , italic_p start_POSTSUBSCRIPT roman_BNS end_POSTSUBSCRIPT. The last condition allows us to consider Eq. (16) still valid if we replace πisubscript𝜋𝑖\pi_{i}italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT with pBBH,isubscript𝑝BBH𝑖p_{\text{BBH},i}italic_p start_POSTSUBSCRIPT BBH , italic_i end_POSTSUBSCRIPT, since the contribution of pNSBH,isubscript𝑝NSBH𝑖p_{\text{NSBH},i}italic_p start_POSTSUBSCRIPT NSBH , italic_i end_POSTSUBSCRIPT and pBNS,isubscript𝑝BNS𝑖p_{\text{BNS},i}italic_p start_POSTSUBSCRIPT BNS , italic_i end_POSTSUBSCRIPT to the non-detection probability is negligible.

Given the definition of Eq. (16), ℒisubscriptℒ𝑖\mathcal{L}_{i}caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT indicates the probability, given a set of (λ1,…,λk)subscript𝜆1…subscript𝜆𝑘(\lambda_{1},...,\lambda_{k})( italic_λ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_λ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ) EM parameters, that the BAT upper limit is not violated, taking into account the possible non-astrophysical origin of the candidate and also the probability that, even if astrophysical, the source is occulted by the Earth and therefore not detectable by Swift. Having a collection of E1,…,ENsubscript𝐸1…subscript𝐸𝑁E_{1},...,E_{N}italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT GW candidates, the posterior distribution of the model parameters can be obtained following the Bayes theorem

P⁢(L0,f|E1,…,EN)𝑃subscript𝐿0conditional𝑓subscript𝐸1…subscript𝐸𝑁\displaystyle P(L_{0},f|E_{1},...,E_{N})italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f | italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT )

(17)

=∏i=1Nℒi⁢π⁢(L0)⁢π⁢(f)/∫∏i=1Nℒi⁢π⁢(L0)⁢π⁢(f)⁢d⁢L0⁢d⁢f,absentsuperscriptsubscriptproduct𝑖1𝑁subscriptℒ𝑖𝜋subscript𝐿0𝜋𝑓superscriptsubscriptproduct𝑖1𝑁subscriptℒ𝑖𝜋subscript𝐿0𝜋𝑓dsubscript𝐿0d𝑓\displaystyle=\prod_{i=1}^{N}\mathcal{L}_{i}\pi(L_{0})\pi(f)\Big{/}\int\prod_{% i=1}^{N}\mathcal{L}_{i}\pi(L_{0})\pi(f)\mathrm{d}L_{0}\mathrm{d}f,= ∏ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_π ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) italic_π ( italic_f ) / ∫ ∏ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_π ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) italic_π ( italic_f ) roman_d italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT roman_d italic_f ,

where π⁢(L0)𝜋subscript𝐿0\pi(L_{0})italic_π ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) and π⁢(f)𝜋𝑓\pi(f)italic_π ( italic_f ) are the prior distributions of L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT and f𝑓fitalic_f. We assume a log-uniform prior for both L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT and f𝑓fitalic_f in the respective intervals 46<log10⁡[L0⁢(erg⁢s−1)]<5346subscript10subscript𝐿0ergsuperscripts15346<\log_{10}[L_{0}(\rm erg~{}s^{-1})]<5346 < roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT [ italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ( roman_erg roman_s start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) ] < 53 and −3<log10⁡(f)<03subscript10𝑓0-3<\log_{10}(f)<0- 3 < roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( italic_f ) < 0. The choice of the prior boundaries are poorly informed by theoretical expectations, which are still affected by large uncertainties. Instead, the priors are chosen on the basis of the typical range of upper limit luminosity derived in this work for BBH events and the total number of candidates considered in this analysis. The constraints reported in the following may strongly depend on the choice of the prior boundaries. Therefore, the final goal of this simulation, more than deriving strong limits on the putative EM model, is to show the predictive power of the present analysis in the context of model inference and how this analysis can improve with the addition of more GW events in the future.

In the specific case of our simulation, we consider all the GW candidates released in GWTC-3 (Abbott et al., 2023), including both the above threshold (pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5) and the subthreshold (pastro<0.5subscript𝑝astro0.5p_{\rm astro}<0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT < 0.5) candidates. For the latter, we emphasize that the classification of the CBC candidate as BBH merger is valid under the condition that the subthreshold GW event is of astrophysical origin. The simulation described in this section is set up in such a way that this assumption is taken into account for the final constraints of the physical parameters. All the low-latency candidates not confirmed by the offline analysis are not included in the simulation. The considered BBH sample with full NITRATES results and available flux upper limits consists of 32323232 events, 12 of which with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5.

Figure 7: Constraints on the two parameters L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT and f𝑓fitalic_f of the model for the putative EM counterpart of BBH mergers. The color map reports the likelihood ℒℒ\mathcal{L}caligraphic_L, for the full analysis including all the O3 catalog events with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5. L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is in units of erg s-1. The thick blue solid and dashed contours indicate the exclusion regions in the [L0,f]subscript𝐿0𝑓[L_{0},f][ italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f ] plane at 90%percent\%% and 50%percent\%% credibility levels, respectively. The magenta solid and dashed lines report the same contours, but for an analysis that includes all O3 catalog events, with no cut in pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT.

Figure 8: Posterior distribution of log10⁡(f)subscript10𝑓\log_{10}(f)roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( italic_f ), including the 50%percent\%% and 90%percent\%% upper limits with dot-dashed and dashed lines, respectively. The function is derived from Fig. 7, marginalizing over L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT.

Figure 9: Posterior distribution of log10⁡(L0)subscript10subscript𝐿0\log_{10}(L_{0})roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ), including the 50%percent\%% and 90%percent\%% upper limits with the dot-dashed and dashed lines, respectively. The function is derived from Fig. 7, marginalizing over f𝑓fitalic_f.

In order to compute numerically the functional behavior of the likelihood, we set up a simulation to evaluate P⁢(L0,f|E1,…,EN)𝑃subscript𝐿0conditional𝑓subscript𝐸1…subscript𝐸𝑁P(L_{0},f|E_{1},...,E_{N})italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f | italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT ) in the full [L0,f]subscript𝐿0𝑓[L_{0},f][ italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f ] plane defined by the prior boundaries. The details of the simulation setup are reported in Appendix B. The results of the simulation for the sample of above threshold BBH candidates are reported in Fig. 7, where the color map indicates the value of ℒ=∏ℒiℒproductsubscriptℒ𝑖\mathcal{L}=\prod\mathcal{L}_{i}caligraphic_L = ∏ caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, normalized by the maximum max⁡(ℒ)ℒ\max{(\mathcal{L})}roman_max ( caligraphic_L ) over the full domain. The contour levels defining the 50%percent\%% and 90%percent\%% exclusion regions are reported as well. For comparison, in Fig. 7 we include also the same contour levels obtained with an analysis that considers all the BBH candidates without imposing any cut on the pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT. The level of constraining power of the analysis can be quantified by defining the fraction of the full parameter space excluded with a credibility level η𝜂\etaitalic_η, corresponding to:

Rη=Iη/Itot,subscript𝑅𝜂subscript𝐼𝜂subscript𝐼totR_{\eta}=I_{\eta}/I_{\rm tot},italic_R start_POSTSUBSCRIPT italic_η end_POSTSUBSCRIPT = italic_I start_POSTSUBSCRIPT italic_η end_POSTSUBSCRIPT / italic_I start_POSTSUBSCRIPT roman_tot end_POSTSUBSCRIPT ,

(18)

where Itot=∫d⁢f⁢d⁢L0subscript𝐼totd𝑓dsubscript𝐿0I_{\text{tot}}=\int\text{d}f\text{d}L_{0}italic_I start_POSTSUBSCRIPT tot end_POSTSUBSCRIPT = ∫ d italic_f d italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT, being the integral extended to the full parameters domain, and with

Iη=∫Sdf⁢dL0,subscript𝐼𝜂subscript𝑆differential-d𝑓differential-dsubscript𝐿0I_{\eta}=\int_{S}\mathrm{d}f\mathrm{d}L_{0},italic_I start_POSTSUBSCRIPT italic_η end_POSTSUBSCRIPT = ∫ start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT roman_d italic_f roman_d italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ,

(19)

corresponding to the dimension of the region S𝑆Sitalic_S of the parameter space excluded with a credibility level η𝜂\etaitalic_η. The analysis performed using only BBH with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5 gives a R90%=22%subscript𝑅percent90percent22R_{90\%}=22\%italic_R start_POSTSUBSCRIPT 90 % end_POSTSUBSCRIPT = 22 %, while the analysis performed with the inclusion of subthreshold BBH candidates gives R90%=10%subscript𝑅percent90percent10R_{90\%}=10\%italic_R start_POSTSUBSCRIPT 90 % end_POSTSUBSCRIPT = 10 %. This result indicates that with the inclusion of GW subthreshold events the analysis allows us to exclude a smaller portion of the parameter space, with respect to an analysis carried out using only events with high pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT. Therefore, even adding further information with the inclusion of more GW events, if these last are likely originated by noise, the full Bayesian analysis is affected by more uncertainty, resulting in an overall less constraining power of the analysis.

Figs. 8 and 9 report the posterior distribution of P⁢(f)𝑃𝑓P(f)italic_P ( italic_f ) and P⁢(L0)𝑃subscript𝐿0P(L_{0})italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ), respectively, obtained as:

P⁢(f)=∫P⁢(L0,f|E1,…,EN)⁢dL0𝑃𝑓𝑃subscript𝐿0conditional𝑓subscript𝐸1…subscript𝐸𝑁differential-dsubscript𝐿0P(f)=\int P(L_{0},f|E_{1},...,E_{N})\mathrm{d}L_{0}italic_P ( italic_f ) = ∫ italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f | italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT ) roman_d italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT

(20)

and

P⁢(L0)=∫P⁢(L0,f|E1,…,EN)⁢df.𝑃subscript𝐿0𝑃subscript𝐿0conditional𝑓subscript𝐸1…subscript𝐸𝑁differential-d𝑓P(L_{0})=\int P(L_{0},f|E_{1},...,E_{N})\mathrm{d}f.italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) = ∫ italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f | italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT ) roman_d italic_f .

(21)

Both P⁢(L0)𝑃subscript𝐿0P(L_{0})italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) and P⁢(f)𝑃𝑓P(f)italic_P ( italic_f ) are normalized such that max⁡[P⁢(L0)]=max⁡[P⁢(f)]=1𝑃subscript𝐿0𝑃𝑓1\max[P(L_{0})]=\max[P(f)]=1roman_max [ italic_P ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) ] = roman_max [ italic_P ( italic_f ) ] = 1. The posterior is reported in magenta and blue for both samples, with and without cut in pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT, respectively. The 50%percent5050\%50 % and 90%percent9090\%90 % upper limits are reported as well. From the shape of the posteriors and the values of the 50%percent5050\%50 % and 90%percent9090\%90 % upper limits, it is evident that both the L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT and f𝑓fitalic_f posteriors only slightly differ from the flat priors, if no cut in pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT is applied. On the other hand, for the above threshold sample, the analysis is more informative and more stringent constraints on the parameters can be obtained, especially for L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT. For the sample with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5, the 50%percent5050\%50 % and 90%percent9090\%90 % upper limits for f𝑓fitalic_f are log10⁡(f50%)=−1.67subscript10subscript𝑓percent501.67\log_{10}(f_{50\%})=-1.67roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( italic_f start_POSTSUBSCRIPT 50 % end_POSTSUBSCRIPT ) = - 1.67 and log10⁡(f90%)=−0.42subscript10subscript𝑓percent900.42\log_{10}(f_{90\%})=-0.42roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( italic_f start_POSTSUBSCRIPT 90 % end_POSTSUBSCRIPT ) = - 0.42, while for L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT are log10⁡[L0,50%⁢(erg s−1)]=47.8subscript10subscript𝐿0percent50superscripterg s147.8\log_{10}[L_{0,50\%}(\text{erg s}^{-1})]=47.8roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT [ italic_L start_POSTSUBSCRIPT 0 , 50 % end_POSTSUBSCRIPT ( erg s start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) ] = 47.8 and log10⁡[L0,90%⁢(erg s−1)]=49.2subscript10subscript𝐿0percent90superscripterg s149.2\log_{10}[L_{0,90\%}(\text{erg s}^{-1})]=49.2roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT [ italic_L start_POSTSUBSCRIPT 0 , 90 % end_POSTSUBSCRIPT ( erg s start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) ] = 49.2, respectively.

In the limit of a collection of triggers which correspond only to non-astrophysical events, i.e., all with πi=0subscript𝜋𝑖0\pi_{i}=0italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 0, the likelihood is constant in the full parameter space, not allowing to infer any constraints on the EM model parameters. On the other hand, if we increase the fraction of confident GW events and we keep fixed the total number N𝑁Nitalic_N, their distance distribution P⁢(DL)𝑃subscript𝐷LP(D_{\text{L}})italic_P ( italic_D start_POSTSUBSCRIPT L end_POSTSUBSCRIPT ) and the derived upper limits, then we obtain that ℒℒ\mathcal{L}caligraphic_L decreases accordingly. This implies that increasing the number of events with πisubscript𝜋𝑖\pi_{i}italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT close to 1, the overall exclusion region in the (λ1,…,λk)subscript𝜆1…subscript𝜆𝑘(\lambda_{1},...,\lambda_{k})( italic_λ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_λ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ) parameter space increases as well. This demonstrates that with the collection of more data, in the limit of a GW detector horizon constant in time, this method allows us to improve incrementally our constraints on the EM models of CBC events. Although, realistically the GW detection horizon will increase with time (Abbott et al., 2020), implying an overall increase of the median values of DLsubscript𝐷LD_{\rm L}italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT of the candidate events. Such an effect increases in turn the values of the luminosity upper limits, increasing as well the values of Pi⁢(ϕ<ϕ0)subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0P_{i}(\phi<\phi_{0})italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ < italic_ϕ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) and hence of ℒℒ\mathcal{L}caligraphic_L. This effect tends to decrease the dimension of the exclusion region of the (λ1,…,λk)subscript𝜆1…subscript𝜆𝑘(\lambda_{1},...,\lambda_{k})( italic_λ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_λ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ) parameter space. Overall, the final outcome of the inclusion of additional GW data, in terms of the constraining power of this analysis, will depend on the simultaneous combined effect of increasing the number of confident events and of increasing the detection horizon.

In order to show how the inclusion of more significant GW candidates can improve the constraining power of the present analysis, we carried out the following simulation. We repeated the same procedure adopted to produce the exclusion regions of Fig. 7, but replacing the real πisubscript𝜋𝑖\pi_{i}italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT with πi=1subscript𝜋𝑖1\pi_{i}=1italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 1 for all the confirmed BBH candidates, hence imposing that they are all significant events. All the values of ϕULsubscriptitalic-ϕUL\phi_{\rm UL}italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT, ε⊕subscript𝜀⊕\varepsilon_{\earth}italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT and P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ) of each candidate are left unchanged. The resulting 50%percent\%% and 90%percent\%% exclusion regions are reported in Fig. 10, with black dashed and solid lines, respectively. The fraction of the 90%percent\%% excluded region increases to a value of R90%=17%subscript𝑅percent90percent17R_{90\%}=17\%italic_R start_POSTSUBSCRIPT 90 % end_POSTSUBSCRIPT = 17 %, clearly demonstrating that, even if the BAT flux upper limit are the same, the increase of confidence about the astrophysical nature of the GW improves our final constraints on the model parameter space. Furthermore, Fig. 10 reports also the 50%percent\%% and 90%percent\%% exclusion regions (with red dashed and solid lines, respectively), obtained as before, but imposing both πi=1subscript𝜋𝑖1\pi_{i}=1italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 1 and ε⊕=0subscript𝜀⊕0\varepsilon_{\earth}=0italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT = 0 for each candidate. This combination corresponds to simulate all real BBH candidates, whose sky localization does not overlap with the sky region covered by the Earth. In this case R90%=35%subscript𝑅percent90percent35R_{90\%}=35\%italic_R start_POSTSUBSCRIPT 90 % end_POSTSUBSCRIPT = 35 %, showing that also the fraction of GW sky posterior occulted by Earth has a significant impact on our final results.

Figure 10: Same as Fig. 7, but simulating all the O3 catalog candidates with an associated pastro=πi=1subscript𝑝astrosubscript𝜋𝑖1p_{\rm astro}=\pi_{i}=1italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT = italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 1. The black dashed and solid lines identify the 50%percent5050\%50 % and 90%percent9090\%90 % exclusion regions, respectively. The red dashed and solid lines have the same meaning, but derived imposing both πi=1subscript𝜋𝑖1\pi_{i}=1italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 1 and ε⊕=0subscript𝜀⊕0\varepsilon_{\earth}=0italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT = 0.

In this work we report the systematic search of signals jointly detected by the LIGO–Virgo interferometers and the Swift-BAT telescope, during the third LVK observing run. Thanks to the prompt availability of BAT data using GUANO and the sensitive targeted search capabilities of the NITRATES pipeline, we conducted deep follow-up searches for EM signals on a sample of 636 GW triggers. The search results did not yield any confident joint detection, allowing us to derive upper limits in the 15–350 keV band. We provide comprehensive details on all analyzed GW triggers along with their NITRATES search statistics. This information can be valuable for calibrating and comparing with other offline targeted search pipelines that are currently operational or may be developed in the future.

In the specific case of the BBH class, the BAT flux upper limits have been used to perform a stacking analysis and to derive constraints on the possible nature of an associated EM emission. As illustrated in Section 6, the presence of several BBH candidates with large values of pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT in our sample, enhances our ability to better constrain the parameter space for EM emission from BH mergers, with minimal assumptions on our prior knowledge of the nature of the emission. The prospect of detecting EM emission from BBH mergers has been debated and discussed in detail in recent times. Particularly, the GBM trigger that accompanied the first BBH merger event, GW150914, has served as a case study to test possible association and potential implications (Abbott et al., 2016b; Connaughton et al., 2016; Goldstein et al., 2016b). Though not likely, there are a number of physical models that have been proposed that could give rise to detectable emission in the gamma-ray band. A summary of the various different models has been discussed in Fletcher et al. (2023) and Veres et al. (2019). The models involve parameters pertaining to potential remnant accretion effects, magnetic field strength, black hole charge and spin, among others (e.g., Loeb 2016; Dai et al. 2017; Woosley 2016; Zhang 2016). The method described in Section 6 can be easily extended to any of these models, provided that the luminosity function Π⁢(L)Π𝐿\Pi(L)roman_Π ( italic_L ) of the putative BBH EM emission is known. Additionally, effects possibly related to the viewing angle dependency of the EM emission can be easily included in this approach. Regarding CBCs containing at least one NS (BNS and NSBH), it was not possible to conduct a similar stacking analysis as the one described for BBH in Section 6, due to the paucity of such events with a large enough value of pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT. Further observations, including the fourth LVK observing run (O4), could lead to the collection of a larger number of BNS and NSBH candidates with moderate values of pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT, giving the possibility to repeat the analysis performed in this paper and to derive informative constraints on the EM emission of these classes and the properties of the associated GRB populations. Data products associated with the present analysis are reported in a separate data release.

O4 commenced on the 24th of May, 2023. The number of significant detections is expected to increase by several times during the entire duration of O4 (Abbott et al., 2018; Petrov et al., 2022). Targeted search results using the GUANO-NITRATES infrastructure are publicly available in real-time. In the case of non-detection of an EM counterpart, the GUANO team reports the 15–350 keV flux upper limit for all the GW triggers classified as significant, via GCN Circulars. Additional enhancements to the likelihood search code have reduced the search latency by a factor of 2, with respect to O3.

Thanks to its sensitivity in the hard X-ray band and the possibility to localize EM transients down to a precision of an arcminute, Swift represents one of the main discovery machines for the detection of EM counterparts of GW transients. This paper shows how the GUANO infrastructure has a deep impact on the multi-messenger science case, in particular for optimally exploiting the sensitivity of the Swift-BAT instrument for the detection of EM counterparts of CBCs detected by the LVK Collaboration. The deep subthreshold search enabled by the NITRATES pipeline sensibly increases the detection horizon of Swift, giving the chance to detect transients also outside the BAT FOV and allowing us to possibly detect faint X-ray/gamma-ray transients associated to relativistic jets observed off-axis, as in the case of GW170817. In the case of a confident joint Swift-GW detection, the GUANO team will promptly disseminate all the information about the EM candidate via GCN Circulars, providing an estimate of the sky localization when available. Moreover, also in the case of non-detection, this paper shows how the upper limits derived from the NITRATES analysis can be combined to have the most sensitive constraints on the EM emission from all the CBC classes. The cumulative collection of non-detection will gradually improve our knowledge of the EM nature of CBCs.

Gayathri Raman, Samuele Ronchini, and Jamie Kennea acknowledge the support of NASA grants 80NSSC19K0408 and 80NSSC22K1498 awarded as part of the NASA Neil Gehrels Swift Observatory Guest Investigator program. Jamie Kennea and James Delaunay acknowledge the support of NASA contract NAS5-0136.

This material is based upon work supported by NSF’s LIGO Laboratory which is a major facility fully funded by the National Science Foundation. The authors also gratefully acknowledge the support of the Science and Technology Facilities Council (STFC) of the United Kingdom, the Max-Planck-Society (MPS), and the State of Niedersachsen/Germany for support of the construction of Advanced LIGO and construction and operation of the GEO 600 detector. Additional support for Advanced LIGO was provided by the Australian Research Council. The authors gratefully acknowledge the Italian Istituto Nazionale di Fisica Nucleare (INFN), the French Centre National de la Recherche Scientifique (CNRS) and the Netherlands Organization for Scientific Research (NWO), for the construction and operation of the Virgo detector and the creation and support of the EGO consortium. The authors also gratefully acknowledge research support from these agencies as well as by the Council of Scientific and Industrial Research of India, the Department of Science and Technology, India, the Science & Engineering Research Board (SERB), India, the Ministry of Human Resource Development, India, the Spanish Agencia Estatal de Investigación (AEI), the Spanish Ministerio de Ciencia, Innovación y Universidades, the European Union NextGenerationEU/PRTR (PRTR-C17.I1), the ICSC - CentroNazionale di Ricerca in High Performance Computing, Big Data and Quantum Computing, funded by the European Union NextGenerationEU, the Comunitat Autonòma de les Illes Balears through the Direcció General de Recerca, Innovació i Transformació Digital with funds from the Tourist Stay Tax Law ITS 2017-006, the Conselleria d’Economia, Hisenda i Innovació, the FEDER Operational Program 2021-2027 of the Balearic Islands, the Conselleria d’Innovació, Universitats, Ciència i Societat Digital de la Generalitat Valenciana and the CERCA Programme Generalitat de Catalunya, Spain, the National Science Centre of Poland and the European Union – European Regional Development Fund; Foundation for Polish Science (FNP), the Polish Ministry of Science and Higher Education, the Swiss National Science Foundation (SNSF), the Russian Science Foundation, the European Commission, the European Social Funds (ESF), the European Regional Development Funds (ERDF), the Royal Society, the Scottish Funding Council, the Scottish Universities Physics Alliance, the Hungarian Scientific Research Fund (OTKA), the French Lyon Institute of Origins (LIO), the Belgian Fonds de la Recherche Scientifique (FRS-FNRS), Actions de Recherche Concertées (ARC) and Fonds Wetenschappelijk Onderzoek – Vlaanderen (FWO), Belgium, the Paris Île-de-France Region, the National Research, Development and Innovation Office Hungary (NKFIH), the National Research Foundation of Korea, the Natural Science and Engineering Research Council Canada, Canadian Foundation for Innovation (CFI), the Brazilian Ministry of Science, Technology, and Innovations, the International Center for Theoretical Physics South American Institute for Fundamental Research (ICTP-SAIFR), the Research Grants Council of Hong Kong, the National Natural Science Foundation of China (NSFC), the Leverhulme Trust, the Research Corporation, the National Science and Technology Council (NSTC), Taiwan, the United States Department of Energy, and the Kavli Foundation. The authors gratefully acknowledge the support of the NSF, STFC, INFN and CNRS for provision of computational resources.

This work was supported by MEXT, JSPS Leading-edge Research Infrastructure Program, JSPS Grant-in-Aid for Specially Promoted Research 26000005, JSPS Grant-in-Aid for Scientific Research on Innovative Areas 2905: JP17H06358, JP17H06361 and JP17H06364, JSPS Core-to-Core Program A. Advanced Research Networks, JSPS Grant-in-Aid for Scientific Research (S) 17H06133 and 20H05639 , JSPS Grant-in-Aid for Transformative Research Areas (A) 20A203: JP20H05854, the joint research program of the Institute for Cosmic Ray Research, University of Tokyo, National Research Foundation (NRF), Computing Infrastructure Project of Global Science experimental Data hub Center (GSDC) at KISTI, Korea Astronomy and Space Science Institute (KASI), and Ministry of Science and ICT (MSIT) in Korea, Academia Sinica (AS), AS Grid Center (ASGC) and the National Science and Technology Council (NSTC) in Taiwan under grants including the Rising Star Program and Science Vanguard Research Program, Advanced Technology Center (ATC) of NAOJ, and Mechanical Engineering Center of KEK.

Additional acknowledgements for support of individual authors may be found in the following document:
https://dcc.ligo.org/LIGO-M2300033/public. For the purpose of open access, the authors have applied a Creative Commons Attribution (CC BY) license to any Author Accepted Manuscript version arising. We request that citations to this article use ’A. G. Abac et al. (LIGO-Virgo-KAGRA Collaboration), …’ or similar phrasing, depending on journal convention.

Matplotlib (Hunter, 2007), SEABORN (Waskom, 2021), NumPy (Harris et al., 2020) and SciPy (Virtanen et al., 2020) were used in the preparation of the manuscript.

• Abbott et al. (2016a) Abbott B. P., et al., 2016a, Physical Review X, 6, 041015
• Abbott et al. (2016b) Abbott B. P., et al., 2016b, Phys. Rev. Lett., 116, 061102
• Abbott et al. (2017a) Abbott B. P., et al., 2017a, Phys. Rev. Lett., 119, 161101
• Abbott et al. (2017b) Abbott B. P., et al., 2017b, Nature, 551, 85
• Abbott et al. (2017c) Abbott B. P., et al., 2017c, ApJ, 841, 89
• Abbott et al. (2017d) Abbott B. P., et al., 2017d, ApJ, 848, L12
• Abbott et al. (2017e) Abbott B. P., et al., 2017e, ApJ, 848, L13
• Abbott et al. (2018) Abbott B. P., et al., 2018, Living Reviews in Relativity, 21, 3
• Abbott et al. (2020) Abbott B. P., et al., 2020, Living Reviews in Relativity, 23, 3
• Abbott et al. (2021) Abbott R., et al., 2021, PhRvX, 11, 021053
• Abbott et al. (2023) Abbott R., et al., 2023, PhRvX, 13, 041039
• Abbott et al. (2024) Abbott R., et al., 2024, PhRvD, 109, 022001
• Alexander et al. (2017) Alexander K. D., et al., 2017, ApJ, 848, L21
• Allison et al. (2016) Allison J., et al., 2016, Nuclear Instruments and Methods in Physics Research A, 835, 186
• Arcavi et al. (2017) Arcavi I., et al., 2017, Nature, 551, 64
• Ashton et al. (2018) Ashton G., et al., 2018, ApJ, 860, 6
• Aubin et al. (2021) Aubin F., et al., 2021, Classical and Quantum Gravity, 38, 095004
• Band et al. (1993) Band D., et al., 1993, ApJ, 413, 281
• Barthelmy et al. (2005) Barthelmy S. D., et al., 2005, Space Sci. Rev., 120, 143
• Bauswein et al. (2017) Bauswein A., Just O., Janka H.-T., Stergioulas N., 2017, ApJ, 850, L34
• Connaughton et al. (2016) Connaughton V., et al., 2016, ApJ, 826, L6
• Coulter et al. (2017) Coulter D. A., et al., 2017, Science, 358, 1556
• Cowperthwaite et al. (2017) Cowperthwaite P. S., et al., 2017, ApJ, 848, L17
• Dai et al. (2017) Dai L., McKinney J. C., Miller M. C., 2017, MNRAS, 470, L92
• Dal Canton et al. (2021) Dal Canton T., Nitz A. H., Gadre B., Cabourn Davies G. S., Villa-Ortega V., Dent T., Harry I., Xiao L., 2021, ApJ, 923, 254
• DeLaunay & Tohuvavohu (2022) DeLaunay J., Tohuvavohu A., 2022, ApJ, 941, 169
• DeLaunay et al. (2020) DeLaunay J., Tohuvavohu A., Kennea J., 2020, General Coordinates Network, 27444, 1
• DeLaunay et al. (2021a) DeLaunay J., Tohuvavohu A., Kennea J. A., Raman G., 2021a, General Coordinates Network, 30130, 1
• DeLaunay et al. (2021b) DeLaunay J., Tohuvavohu A., Kennea J. A., Raman G., 2021b, General Coordinates Network, 30302, 1
• DeLaunay et al. (2022a) DeLaunay J., Tohuvavohu A., Raman G., Kennea J. A., 2022a, General Coordinates Network, 31402, 1
• DeLaunay et al. (2022b) DeLaunay J., Tohuvavohu A., Raman G., Kennea J. A., 2022b, General Coordinates Network, 31402, 1
• DeLaunay et al. (2023) DeLaunay J., Tohuvavohu A., Ronchini S., Raman G., Kennea J. A., Parsotan T., 2023, General Coordinates Network, 34747, 1
• Drout et al. (2017) Drout M. R., et al., 2017, Science, 358, 1570
• Evans et al. (2017) Evans P. A., et al., 2017, Science, 358, 1565
• Fermi GBM Team (2020) Fermi GBM Team 2020, General Coordinates Network, 26944, 1
• Fletcher et al. (2023) Fletcher C., et al., 2023, arXiv e-prints, p. arXiv:2308.13666
• Gehrels et al. (2004) Gehrels N., et al., 2004, ApJ, 611, 1005
• Goldstein et al. (2016a) Goldstein A., Burns E., Hamburg R., Connaughton V., Veres P., Briggs M. S., Hui C. M., The GBM-LIGO Collaboration 2016a, arXiv e-prints, p. arXiv:1612.02395
• Goldstein et al. (2016b) Goldstein A., Burns E., Hamburg R., Connaughton V., Veres P., Briggs M. S., Hui C. M., The GBM-LIGO Collaboration 2016b, arXiv e-prints, p. arXiv:1612.02395
• Goldstein et al. (2017a) Goldstein A., et al., 2017a, ApJ, 848, L14
• Goldstein et al. (2017b) Goldstein A., et al., 2017b, ApJ, 848, L14
• Hallinan et al. (2017) Hallinan G., et al., 2017, Science, 358, 1579
• Hamburg et al. (2020) Hamburg R., et al., 2020, ApJ, 893, 100
• Harris et al. (2020) Harris C. R., et al., 2020, Nature, 585, 357
• Hooper et al. (2012) Hooper S., Chung S. K., Luan J., Blair D., Chen Y., Wen L., 2012, Phys. Rev. D, 86, 024012
• Hotokezaka et al. (2019) Hotokezaka K., Nakar E., Gottlieb O., Nissanke S., Masuda K., Hallinan G., Mooley K. P., Deller A. T., 2019, Nature Astronomy, 3, 940
• Hughes & Holz (2003) Hughes S. A., Holz D. E., 2003, Classical and Quantum Gravity, 20, S65
• Hunter (2007) Hunter J. D., 2007, Computing in Science & Engineering, 9, 90
• Hurley et al. (2020) Hurley K., et al., 2020, General Coordinates Network, 26949, 1
• Klimenko et al. (2005) Klimenko S., Mohanty S., Rakhmanov M., Mitselmakher G., 2005, Phys. Rev. D, 72, 122002
• Klimenko et al. (2016) Klimenko S., et al., 2016, Phys. Rev. D, 93, 042004
• Kyutoku et al. (2021) Kyutoku K., Shibata M., Taniguchi K., 2021, Living Reviews in Relativity, 24, 5
• Loeb (2016) Loeb A., 2016, ApJ, 819, L21
• Margutti & Chornock (2021) Margutti R., Chornock R., 2021, ARA&A, 59, 155
• Margutti et al. (2017) Margutti R., et al., 2017, ApJ, 848, L20
• Messick et al. (2017) Messick C., et al., 2017, Phys. Rev. D, 95, 042001
• Nakar (2020a) Nakar E., 2020a, Phys. Rep., 886, 1
• Nakar (2020b) Nakar E., 2020b, Phys. Rep., 886, 1
• Nissanke et al. (2013) Nissanke S., Holz D. E., Dalal N., Hughes S. A., Sievers J. L., Hirata C. M., 2013, arXiv e-prints, p. arXiv:1307.2638
• Nitz et al. (2018) Nitz A. H., Dal Canton T., Davis D., Reyes S., 2018, Phys. Rev. D, 98, 024050
• Nitz et al. (2019) Nitz A. H., Nielsen A. B., Capano C. D., 2019, ApJ, 876, L4
• Ossokine et al. (2020) Ossokine S., et al., 2020, Phys. Rev. D, 102, 044055
• Petrov et al. (2022) Petrov P., et al., 2022, ApJ, 924, 54
• Pian et al. (2017) Pian E., et al., 2017, Nature, 551, 67
• Pillas et al. (2023) Pillas M., et al., 2023, arXiv e-prints, p. arXiv:2306.04373
• Pratten et al. (2021) Pratten G., et al., 2021, Phys. Rev. D, 103, 104056
• Radice et al. (2018) Radice D., Perego A., Zappa F., Bernuzzi S., 2018, ApJ, 852, L29
• Sachdev et al. (2019) Sachdev S., et al., 2019, arXiv e-prints, p. arXiv:1901.08580
• Salafia et al. (2023) Salafia O. S., Ravasio M. E., Ghirlanda G., Mandel I., 2023, A&A, 680, A45
• Savchenko et al. (2017) Savchenko V., et al., 2017, ApJ, 848, L15
• Schutz (1986) Schutz B. F., 1986, Nature, 323, 310
• Singer & Price (2016) Singer L. P., Price L. R., 2016, Phys. Rev. D, 93, 024013
• Singer et al. (2016) Singer L. P., et al., 2016, ApJ, 829, L15
• Smartt et al. (2017) Smartt S. J., et al., 2017, Nature, 551, 75
• Tanvir et al. (2017) Tanvir N. R., et al., 2017, ApJ, 848, L27
• Tohuvavohu (2023) Tohuvavohu A., 2023, General Coordinates Network, 33132, 1
• Tohuvavohu et al. (2020) Tohuvavohu A., Kennea J. A., DeLaunay J., Palmer D. M., Cenko S. B., Barthelmy S., 2020, ApJ, 900, 35
• Tohuvavohu et al. (2021a) Tohuvavohu A., Raman G., DeLaunay J., Kennea J. A., 2021a, General Coordinates Network, 31049, 1
• Tohuvavohu et al. (2021b) Tohuvavohu A., Raman G., DeLaunay J., Kennea J. A., 2021b, General Coordinates Network, 31049, 1
• Tohuvavohu et al. (2022a) Tohuvavohu A., DeLaunay J., Raman G., Kennea J. A., 2022a, General Coordinates Network, 32167, 1
• Tohuvavohu et al. (2022b) Tohuvavohu A., DeLaunay J., Raman G., Kennea J. A., 2022b, General Coordinates Network, 32375, 1
• Troja et al. (2017) Troja E., et al., 2017, Nature, 551, 71
• Urban (2016) Urban A. L., 2016, PhD thesis, University of Wisconsin, Milwaukee
• Veres et al. (2019) Veres P., Dal Canton T., Burns E., Goldstein A., Littenberg T. B., Christensen N., Preece R. D., 2019, ApJ, 882, 53
• Villar et al. (2017) Villar V. A., et al., 2017, ApJ, 851, L21
• Virtanen et al. (2020) Virtanen P., et al., 2020, Nature Methods, 17, 261
• Wang & Giannios (2021) Wang H., Giannios D., 2021, ApJ, 908, 200
• Waskom (2021) Waskom M. L., 2021, Journal of Open Source Software, 6, 3021
• Woosley (2016) Woosley S. E., 2016, ApJ, 824, L10
• Zhang (2016) Zhang B., 2016, ApJ, 827, L31

A more accurate method to derive the luminosity upper limit should be based on the knowledge of P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ) and P⁢(ϕ)𝑃italic-ϕP(\phi)italic_P ( italic_ϕ ), where ϕitalic-ϕ\phiitalic_ϕ is the flux measured in the BAT energy band. Having only an upper limit, P⁢(ϕ)𝑃italic-ϕP(\phi)italic_P ( italic_ϕ ) can be approximated as

P⁢(ϕ)∝{Π⁢(ϕ),ϕ<ϕUL0,ϕ>ϕUL,proportional-to𝑃italic-ϕcasesΠitalic-ϕitalic-ϕsubscriptitalic-ϕUL0italic-ϕsubscriptitalic-ϕULP(\phi)\propto\begin{cases}\Pi(\phi),&\phi<\phi_{\text{UL}}\\ 0,&\phi>\phi_{\text{UL}}\end{cases},italic_P ( italic_ϕ ) ∝ { start_ROW start_CELL roman_Π ( italic_ϕ ) , end_CELL start_CELL italic_ϕ < italic_ϕ start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT end_CELL end_ROW start_ROW start_CELL 0 , end_CELL start_CELL italic_ϕ > italic_ϕ start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT end_CELL end_ROW ,

(A1)

where Π⁢(ϕ)Πitalic-ϕ\Pi(\phi)roman_Π ( italic_ϕ ) is our prior distribution for the flux. Using the conversion from flux to luminosity L=4⁢π⁢DL2⁢ϕ𝐿4𝜋superscriptsubscript𝐷L2italic-ϕL=4\pi D_{\rm L}^{2}\phiitalic_L = 4 italic_π italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_ϕ, the probability distribution of the luminosity can be computed as

P⁢(L)=P⁢(4⁢π⁢DL2⁢ϕ)∝∫1ϕ⁢Pϕ⁢(ϕ)⁢PDL2⁢(L4⁢π⁢ϕ)⁢dϕ,𝑃𝐿𝑃4𝜋superscriptsubscript𝐷L2italic-ϕproportional-to1italic-ϕsubscript𝑃italic-ϕitalic-ϕsubscript𝑃superscriptsubscript𝐷L2𝐿4𝜋italic-ϕdifferential-ditalic-ϕP(L)=P(4\pi D_{\rm L}^{2}\phi)\propto\int\frac{1}{\phi}P_{\phi}(\phi)P_{D_{\rm L% }^{2}}\left(\dfrac{L}{4\pi\phi}\right)\rm d\phi,italic_P ( italic_L ) = italic_P ( 4 italic_π italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_ϕ ) ∝ ∫ divide start_ARG 1 end_ARG start_ARG italic_ϕ end_ARG italic_P start_POSTSUBSCRIPT italic_ϕ end_POSTSUBSCRIPT ( italic_ϕ ) italic_P start_POSTSUBSCRIPT italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT ( divide start_ARG italic_L end_ARG start_ARG 4 italic_π italic_ϕ end_ARG ) roman_d italic_ϕ ,

(A2)

where Pϕsubscript𝑃italic-ϕP_{\phi}italic_P start_POSTSUBSCRIPT italic_ϕ end_POSTSUBSCRIPT is the flux probability distribution and PDL2subscript𝑃superscriptsubscript𝐷L2P_{D_{\rm L}^{2}}italic_P start_POSTSUBSCRIPT italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT is the probability distribution of DL2superscriptsubscript𝐷L2D_{\rm L}^{2}italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT. In the conversion from flux to luminosity, the k-correction has been omitted, since it introduces a mild dependence on the redshift, which is not relevant for the purposes of this section. The 5σ𝜎\sigmaitalic_σ luminosity upper limit LULsubscript𝐿ULL_{\text{UL}}italic_L start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT can be found imposing that

∫0LULP⁢(L)⁢d⁢L=1−ε5⁢σ,superscriptsubscript0subscript𝐿UL𝑃𝐿d𝐿1subscript𝜀5𝜎\int_{0}^{L_{\text{UL}}}P(L)\text{d}L=1-\varepsilon_{5\sigma},∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_L start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT end_POSTSUPERSCRIPT italic_P ( italic_L ) d italic_L = 1 - italic_ε start_POSTSUBSCRIPT 5 italic_σ end_POSTSUBSCRIPT ,

(A3)

with ε5⁢σ=3×10−7subscript𝜀5𝜎3superscript107\varepsilon_{5\sigma}=3\times 10^{-7}italic_ε start_POSTSUBSCRIPT 5 italic_σ end_POSTSUBSCRIPT = 3 × 10 start_POSTSUPERSCRIPT - 7 end_POSTSUPERSCRIPT. The value of LULsubscript𝐿ULL_{\text{UL}}italic_L start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT has been computed adopting two different assumptions for the flux prior, corresponding to Π⁢(ϕ)∝proportional-toΠitalic-ϕabsent\Pi(\phi)\proptoroman_Π ( italic_ϕ ) ∝ const. and Π⁢(ϕ)∝ϕ−3/2proportional-toΠitalic-ϕsuperscriptitalic-ϕ32\Pi(\phi)\propto\phi^{-3/2}roman_Π ( italic_ϕ ) ∝ italic_ϕ start_POSTSUPERSCRIPT - 3 / 2 end_POSTSUPERSCRIPT, with the latter being inspired by the usual trend followed by GRBs (e.g., Salafia et al. 2023). In both cases, we find that LUL∼3×4⁢π⁢⟨DL2⟩⁢ϕULsimilar-tosubscript𝐿UL34𝜋delimited-⟨⟩superscriptsubscript𝐷L2subscriptitalic-ϕULL_{\text{UL}}\sim 3\times 4\pi\langle D_{\rm L}^{2}\rangle\phi_{\text{UL}}italic_L start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT ∼ 3 × 4 italic_π ⟨ italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT ⟩ italic_ϕ start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT.

In this section we specify the details of the simulation used to compute numerically the ℒ⁢(L0,f|E1,…,EN)ℒsubscript𝐿0conditional𝑓subscript𝐸1…subscript𝐸𝑁\mathcal{L}(L_{0},f|E_{1},...,E_{N})caligraphic_L ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f | italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT ) function. For each simulated GW candidate, the single ℒisubscriptℒ𝑖\mathcal{L}_{i}caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is computed for each pairs of values (L0,n,fm)subscript𝐿0𝑛subscript𝑓𝑚(L_{0,n},f_{m})( italic_L start_POSTSUBSCRIPT 0 , italic_n end_POSTSUBSCRIPT , italic_f start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ). The flux predicted by the EM model is predicted injecting 1000 sources whose luminosity distance is distributed according to P⁢(DL)𝑃subscript𝐷LP(D_{\rm L})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ), derived from the GW localization. The probability Pi⁢(ϕ<ϕ0,i)subscript𝑃𝑖italic-ϕsubscriptitalic-ϕ0𝑖P_{i}(\phi<\phi_{0,i})italic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_ϕ < italic_ϕ start_POSTSUBSCRIPT 0 , italic_i end_POSTSUBSCRIPT ) is derived computing the fraction of cases that have a flux below the sky-averaged BAT upper limit, defined by Eq. (3). The computation of ℒℒ\mathcal{L}caligraphic_L is performed on a 100×100100100100\times 100100 × 100 grid of (L0,n,fm)subscript𝐿0𝑛subscript𝑓𝑚(L_{0,n},f_{m})( italic_L start_POSTSUBSCRIPT 0 , italic_n end_POSTSUBSCRIPT , italic_f start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ). Once the previous steps are performed for all the GW candidates, the final combined likelihood is computed as

ℒ⁢(L0,f|E1,…,EN)=∏iℒi⁢(L0,f).ℒsubscript𝐿0conditional𝑓subscript𝐸1…subscript𝐸𝑁subscriptproduct𝑖subscriptℒ𝑖subscript𝐿0𝑓\mathcal{L}(L_{0},f|E_{1},...,E_{N})=\prod_{i}\mathcal{L}_{i}(L_{0},f).caligraphic_L ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f | italic_E start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_E start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT ) = ∏ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT caligraphic_L start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_f ) .

(B1)

In order to produce the credibility contours in the [L0,n,fm]subscript𝐿0𝑛subscript𝑓𝑚[L_{0,n},f_{m}][ italic_L start_POSTSUBSCRIPT 0 , italic_n end_POSTSUBSCRIPT , italic_f start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ] plane, we adopt the following steps:

ℒℒ\\mathcal{L}caligraphic\_L is normalized such that

∑n,mℒ⁢(L0,n,fm)\=1.subscript𝑛𝑚ℒsubscript𝐿0𝑛subscript𝑓𝑚1\\sum\_{n,m}\\mathcal{L}(L\_{0,n},f\_{m})=1.∑ start\_POSTSUBSCRIPT italic\_n , italic\_m end\_POSTSUBSCRIPT caligraphic\_L ( italic\_L start\_POSTSUBSCRIPT 0 , italic\_n end\_POSTSUBSCRIPT , italic\_f start\_POSTSUBSCRIPT italic\_m end\_POSTSUBSCRIPT ) = 1 .

(B2)

A one-dimensional array ℒ⁢\[xn\]ℒdelimited-\[\]subscript𝑥𝑛\\mathcal{L}\[x\_{n}\]caligraphic\_L \[ italic\_x start\_POSTSUBSCRIPT italic\_n end\_POSTSUBSCRIPT \] is created flattening the two-dimensional grid ℒ⁢(L0,n,fm)ℒsubscript𝐿0𝑛subscript𝑓𝑚\\mathcal{L}(L\_{0,n},f\_{m})caligraphic\_L ( italic\_L start\_POSTSUBSCRIPT 0 , italic\_n end\_POSTSUBSCRIPT , italic\_f start\_POSTSUBSCRIPT italic\_m end\_POSTSUBSCRIPT ), then ℒ⁢\[xn\]ℒdelimited-\[\]subscript𝑥𝑛\\mathcal{L}\[x\_{n}\]caligraphic\_L \[ italic\_x start\_POSTSUBSCRIPT italic\_n end\_POSTSUBSCRIPT \] is sorted in ascending order.

We find the element \[L^0,f^\]\=\[xn∗\]subscript^𝐿0^𝑓delimited-\[\]subscript𝑥superscript𝑛\[\\hat{L}\_{0},\\hat{f}\]=\[x\_{n^{\*}}\]\[ over^ start\_ARG italic\_L end\_ARG start\_POSTSUBSCRIPT 0 end\_POSTSUBSCRIPT , over^ start\_ARG italic\_f end\_ARG \] = \[ italic\_x start\_POSTSUBSCRIPT italic\_n start\_POSTSUPERSCRIPT ∗ end\_POSTSUPERSCRIPT end\_POSTSUBSCRIPT \] such that

∑n\=0n∗ℒ⁢\[xn\]\=λ,superscriptsubscript𝑛0superscript𝑛ℒdelimited-\[\]subscript𝑥𝑛𝜆\\sum\_{n=0}^{n^{\*}}\\mathcal{L}\[x\_{n}\]=\\lambda,∑ start\_POSTSUBSCRIPT italic\_n = 0 end\_POSTSUBSCRIPT start\_POSTSUPERSCRIPT italic\_n start\_POSTSUPERSCRIPT ∗ end\_POSTSUPERSCRIPT end\_POSTSUPERSCRIPT caligraphic\_L \[ italic\_x start\_POSTSUBSCRIPT italic\_n end\_POSTSUBSCRIPT \] = italic\_λ ,

(B3)

where λ𝜆\\lambdaitalic\_λ is the credibility level of the contour.

The contour is drawn imposing ℒ\=ℒ⁢(L^0,f^)ℒℒsubscript^𝐿0^𝑓\\mathcal{L}=\\mathcal{L}(\\hat{L}\_{0},\\hat{f})caligraphic\_L = caligraphic\_L ( over^ start\_ARG italic\_L end\_ARG start\_POSTSUBSCRIPT 0 end\_POSTSUBSCRIPT , over^ start\_ARG italic\_f end\_ARG ).

In this appendix we show an alternative method to compute the non-detection likelihood presented in Section 6. The NITRATES analysis allows us to derive a flux upper limit at a given confidence level for each pixel of the GW sky localization, corresponding to the function ϕUL⁢(RA,Dec)subscriptitalic-ϕULRADec\phi_{\rm UL}(\rm RA,Dec)italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( roman_RA , roman_Dec ) defined in Eq. (3). Then combined probability of being located in the pixel xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and to have a non-detectable EM emission is

P⁢(non-det,xi)∝PGW⁢(xi)⁢P⁢[ϕ<ϕUL⁢(xi)]⁢Δ⁢Ωi,proportional-to𝑃non-detsubscript𝑥𝑖subscript𝑃GWsubscript𝑥𝑖𝑃delimited-[]italic-ϕsubscriptitalic-ϕULsubscript𝑥𝑖ΔsubscriptΩ𝑖P(\text{non-det},x_{i})\propto P_{\text{GW}}(x_{i})P[\phi<\phi_{\rm UL}(x_{i})% ]\Delta\Omega_{i},italic_P ( non-det , italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) ∝ italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) italic_P [ italic_ϕ < italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) ] roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ,

(C1)

where

P⁢[ϕ<ϕUL⁢(xi)]=(1−f)+f⁢∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ,𝑃delimited-[]italic-ϕsubscriptitalic-ϕULsubscript𝑥𝑖1𝑓𝑓superscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕP[\phi<\phi_{\rm UL}(x_{i})]=(1-f)+f\int_{0}^{\phi_{\rm UL}(x_{i})}P(\phi|x_{i% })\text{d}\phi,italic_P [ italic_ϕ < italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) ] = ( 1 - italic_f ) + italic_f ∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) d italic_ϕ ,

(C2)

and Δ⁢ΩiΔsubscriptΩ𝑖\Delta\Omega_{i}roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is the area of the pixel. Here we express P⁢(ϕ|xi)𝑃conditionalitalic-ϕsubscript𝑥𝑖P(\phi|x_{i})italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) as the conditional flux probability distribution, namely the flux probability distribution assuming that the GW source is contained in the pixel xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT. To compute the P⁢(ϕ|xi)𝑃conditionalitalic-ϕsubscript𝑥𝑖P(\phi|x_{i})italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ), for a fixed luminosity L𝐿Litalic_L, the luminosity distance is extracted from the the conditional probability distribution P⁢(DL|xi)𝑃conditionalsubscript𝐷Lsubscript𝑥𝑖P(D_{\rm L}|x_{i})italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ), which is derived from the GW sky localization. Finally the overall non-detection probability is obtained integrating Eq. (C1) over the full sky:

P⁢(non-det|f,L0)=∑xiPGW⁢(xi)⁢P⁢[ϕ<ϕUL⁢(xi)]⁢Δ⁢Ωi=(1−f)+f⁢[ε⊕+∑xi∉Ω⊕PGW⁢(xi)⁢Δ⁢Ωi⁢∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ],𝑃conditionalnon-det𝑓subscript𝐿0subscriptsubscript𝑥𝑖subscript𝑃GWsubscript𝑥𝑖𝑃delimited-[]italic-ϕsubscriptitalic-ϕULsubscript𝑥𝑖ΔsubscriptΩ𝑖1𝑓𝑓delimited-[]subscript𝜀⊕subscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑃GWsubscript𝑥𝑖ΔsubscriptΩ𝑖superscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕP(\text{non-det}|f,L_{0})=\sum_{x_{i}}P_{\text{GW}}(x_{i})P[\phi<\phi_{\rm UL}% (x_{i})]\Delta\Omega_{i}=(1-f)+f\Big{[}\varepsilon_{\earth}+\sum_{x_{i}\notin% \Omega_{\earth}}P_{\text{GW}}(x_{i})\Delta\Omega_{i}\int_{0}^{\phi_{\rm UL}(x_% {i})}P(\phi|x_{i})\text{d}\phi\Big{]},italic_P ( non-det | italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) = ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) italic_P [ italic_ϕ < italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) ] roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = ( 1 - italic_f ) + italic_f [ italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) d italic_ϕ ] ,

(C3)

where we have used that

∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ=1⁢ if ⁢xi∈⊕, and ⁢∑xi∈Ω⊕PGW⁢(xi)⁢Δ⁢Ωi=ε⊕.formulae-sequencesuperscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕ1 if subscript𝑥𝑖⊕ and subscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑃GWsubscript𝑥𝑖ΔsubscriptΩ𝑖subscript𝜀⊕\int_{0}^{\phi_{\rm UL}(x_{i})}P(\phi|x_{i})\text{d}\phi=1\text{ if }x_{i}\in% \earth,\quad\text{ and }\sum_{x_{i}\in\Omega_{\earth}}P_{\text{GW}}(x_{i})% \Delta\Omega_{i}=\varepsilon_{\earth}.∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) d italic_ϕ = 1 if italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ ⊕ , and ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT .

(C4)

The resulting probability of non-detecting any EM emission in correspondence to a GW trigger with a given pastro=πisubscript𝑝astrosubscript𝜋𝑖p_{\rm astro}=\pi_{i}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT = italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is

P⁢(non-det|f,L0,πi)=(1−πi)+πi⁢P⁢(non-det|f,L0).𝑃conditionalnon-det𝑓subscript𝐿0subscript𝜋𝑖1subscript𝜋𝑖subscript𝜋𝑖𝑃conditionalnon-det𝑓subscript𝐿0P(\text{non-det}|f,L_{0},\pi_{i})=(1-\pi_{i})+\pi_{i}P(\text{non-det}|f,L_{0}).italic_P ( non-det | italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT , italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = ( 1 - italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) + italic_π start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_P ( non-det | italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) .

(C5)

Eq. (C3) has to be compared with the method used in Section 6, where instead we used the approximation:

P⁢(non-det|f,L0)=(1−f)+f⁢∫0ϕULP⁢(ϕ)⁢d⁢ϕ,𝑃conditionalnon-det𝑓subscript𝐿01𝑓𝑓superscriptsubscript0subscriptitalic-ϕUL𝑃italic-ϕditalic-ϕP(\text{non-det}|f,L_{0})=(1-f)+f\int_{0}^{\phi_{\rm UL}}P(\phi)\text{d}\phi,italic_P ( non-det | italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ) = ( 1 - italic_f ) + italic_f ∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT end_POSTSUPERSCRIPT italic_P ( italic_ϕ ) d italic_ϕ ,

(C6)

with

ϕUL=∫Ω∉Ω⊕ϕUL⁢(Ω)⁢PGW⁢(Ω)⁢d⁢Ω,subscriptitalic-ϕULsubscriptΩsubscriptΩ⊕subscriptitalic-ϕULΩsubscript𝑃GWΩdΩ\phi_{\rm UL}=\int_{\Omega\notin\Omega_{\earth}}\phi_{\text{UL}}(\Omega)P_{% \text{GW}}(\Omega)\text{d}\Omega,italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT = ∫ start_POSTSUBSCRIPT roman_Ω ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_ϕ start_POSTSUBSCRIPT UL end_POSTSUBSCRIPT ( roman_Ω ) italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( roman_Ω ) d roman_Ω ,

(C7)

and P⁢(ϕ)𝑃italic-ϕP(\phi)italic_P ( italic_ϕ ) is obtained extracting DLsubscript𝐷LD_{\text{L}}italic_D start_POSTSUBSCRIPT L end_POSTSUBSCRIPT from the full sky marginalized distribution P⁢(DL)𝑃subscript𝐷LP(D_{\text{L}})italic_P ( italic_D start_POSTSUBSCRIPT L end_POSTSUBSCRIPT ), corresponding to

P⁢(DL)=∑xiPGW⁢(xi)⁢P⁢(DL|xi)⁢Δ⁢Ωi.𝑃subscript𝐷Lsubscriptsubscript𝑥𝑖subscript𝑃GWsubscript𝑥𝑖𝑃conditionalsubscript𝐷Lsubscript𝑥𝑖ΔsubscriptΩ𝑖P(D_{\rm L})=\sum_{x_{i}}P_{\text{GW}}(x_{i})P(D_{\rm L}|x_{i})\Delta\Omega_{i}.italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT ) = ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) italic_P ( italic_D start_POSTSUBSCRIPT roman_L end_POSTSUBSCRIPT | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT .

(C8)

The two methods give comparable results in the assumption that the following approximation is valid:

∫0ϕULP⁢(ϕ)⁢d⁢ϕ≈ε⊕+∑xi∉Ω⊕PGW⁢(xi)⁢Δ⁢Ωi⁢∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ.superscriptsubscript0subscriptitalic-ϕUL𝑃italic-ϕditalic-ϕsubscript𝜀⊕subscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑃GWsubscript𝑥𝑖ΔsubscriptΩ𝑖superscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕ\int_{0}^{\phi_{\rm UL}}P(\phi)\text{d}\phi\approx\varepsilon_{\earth}+\sum_{x% _{i}\notin\Omega_{\earth}}P_{\text{GW}}(x_{i})\Delta\Omega_{i}\int_{0}^{\phi_{% \rm UL}(x_{i})}P(\phi|x_{i})\text{d}\phi.∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT end_POSTSUPERSCRIPT italic_P ( italic_ϕ ) d italic_ϕ ≈ italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) d italic_ϕ .

(C9)

For completeness, we clarify here the main differences in the two methods.

This is the method used in Section 6 and is based on the following steps:

The marginalized upper limit ϕULsubscriptitalic-ϕUL\\phi\_{\\rm UL}italic\_ϕ start\_POSTSUBSCRIPT roman\_UL end\_POSTSUBSCRIPT is computed over the full sky, weighting by the GW sky localization.

Once L0subscript𝐿0L\_{0}italic\_L start\_POSTSUBSCRIPT 0 end\_POSTSUBSCRIPT is fixed, the flux probability distribution P⁢(ϕ)𝑃italic-ϕP(\\phi)italic\_P ( italic\_ϕ ) is computed extracting randomly DLsubscript𝐷LD\_{\\rm L}italic\_D start\_POSTSUBSCRIPT roman\_L end\_POSTSUBSCRIPT from the P⁢(DL)𝑃subscript𝐷LP(D\_{\\rm L})italic\_P ( italic\_D start\_POSTSUBSCRIPT roman\_L end\_POSTSUBSCRIPT ), the latter corresponding to the posterior distribution of the GW luminosity distance, marginalized over the full sky (excluding the part occulted by the Earth).

The integral ∫0ϕULP⁢(ϕ)⁢d⁢ϕsuperscriptsubscript0subscriptitalic-ϕUL𝑃italic-ϕditalic-ϕ\\int\_{0}^{\\phi\_{\\rm UL}}P(\\phi)\\text{d}\\phi∫ start\_POSTSUBSCRIPT 0 end\_POSTSUBSCRIPT start\_POSTSUPERSCRIPT italic\_ϕ start\_POSTSUBSCRIPT roman\_UL end\_POSTSUBSCRIPT end\_POSTSUPERSCRIPT italic\_P ( italic\_ϕ ) d italic\_ϕ which appears in Eq. ([C6](https://arxiv.org/html/2407.12867v1#A3.E6)) is evaluated counting the fraction of simulated events that have a predicted flux below the sky-averaged upper limit ϕULsubscriptitalic-ϕUL\\phi\_{\\rm UL}italic\_ϕ start\_POSTSUBSCRIPT roman\_UL end\_POSTSUBSCRIPT.

This is the method presented in this appendix and summarized by Eqs. (C3) and (C5), consisting in the following procedure:

A set of sources is injected in space and the distribution follows the volumetric probability distribution of the GW candidate. First the coordinates of the injected source are extracted from the sky localization PGW⁢(RA,Dec)subscript𝑃GWRADecP\_{\\text{GW}}(\\rm RA,Dec)italic\_P start\_POSTSUBSCRIPT GW end\_POSTSUBSCRIPT ( roman\_RA , roman\_Dec ), then for each position the distance is extracted according to the conditional probability P⁢(DL|RA,Dec)𝑃conditionalsubscript𝐷LRADecP(D\_{\\text{L}}|\\rm RA,Dec)italic\_P ( italic\_D start\_POSTSUBSCRIPT L end\_POSTSUBSCRIPT | roman\_RA , roman\_Dec ).

For each injected source, once the luminosity L0subscript𝐿0L\_{0}italic\_L start\_POSTSUBSCRIPT 0 end\_POSTSUBSCRIPT is fixed, the predicted flux is compared with the coordinates-dependent BAT upper limit map ϕUL⁢(RA,Dec)subscriptitalic-ϕULRADec\\phi\_{\\rm UL}(\\rm RA,Dec)italic\_ϕ start\_POSTSUBSCRIPT roman\_UL end\_POSTSUBSCRIPT ( roman\_RA , roman\_Dec ).

We define ρ∉⊕subscript𝜌absent⊕\\rho\_{\\notin\\earth}italic\_ρ start\_POSTSUBSCRIPT ∉ ⊕ end\_POSTSUBSCRIPT the fraction of all the sources injected which are not occulted by the Earth and also have a predicted flux below ϕUL⁢(RA,Dec)subscriptitalic-ϕULRADec\\phi\_{\\rm UL}(\\rm RA,Dec)italic\_ϕ start\_POSTSUBSCRIPT roman\_UL end\_POSTSUBSCRIPT ( roman\_RA , roman\_Dec ). Given this definition, we have:

∑xi∉Ω⊕PGW⁢(xi)⁢Δ⁢Ωi⁢∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ\=(1−ε⊕)⁢ρ∉⊕.subscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑃GWsubscript𝑥𝑖ΔsubscriptΩ𝑖superscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕ1subscript𝜀⊕subscript𝜌absent⊕\\sum\_{x\_{i}\\notin\\Omega\_{\\earth}}P\_{\\text{GW}}(x\_{i})\\Delta\\Omega\_{i}\\int\_{0}^% {\\phi\_{\\rm UL}(x\_{i})}P(\\phi|x\_{i})\\text{d}\\phi=(1-\\varepsilon\_{\\earth})\\rho\_{% \\notin\\earth}.∑ start\_POSTSUBSCRIPT italic\_x start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ∉ roman\_Ω start\_POSTSUBSCRIPT ⊕ end\_POSTSUBSCRIPT end\_POSTSUBSCRIPT italic\_P start\_POSTSUBSCRIPT GW end\_POSTSUBSCRIPT ( italic\_x start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ) roman\_Δ roman\_Ω start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ∫ start\_POSTSUBSCRIPT 0 end\_POSTSUBSCRIPT start\_POSTSUPERSCRIPT italic\_ϕ start\_POSTSUBSCRIPT roman\_UL end\_POSTSUBSCRIPT ( italic\_x start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ) end\_POSTSUPERSCRIPT italic\_P ( italic\_ϕ | italic\_x start\_POSTSUBSCRIPT italic\_i end\_POSTSUBSCRIPT ) d italic\_ϕ = ( 1 - italic\_ε start\_POSTSUBSCRIPT ⊕ end\_POSTSUBSCRIPT ) italic\_ρ start\_POSTSUBSCRIPT ∉ ⊕ end\_POSTSUBSCRIPT .

(C10)

The last equality can be justified considering that, if for each pixel i𝑖iitalic_i we inject Ntot,isubscript𝑁tot𝑖N_{\text{tot},i}italic_N start_POSTSUBSCRIPT tot , italic_i end_POSTSUBSCRIPT sources, we can define ρi=NND,i/Ntot,isubscript𝜌𝑖subscript𝑁ND𝑖subscript𝑁toti\rho_{i}=N_{\text{ND},i}/N_{\rm\text{tot},i}italic_ρ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = italic_N start_POSTSUBSCRIPT ND , italic_i end_POSTSUBSCRIPT / italic_N start_POSTSUBSCRIPT tot , roman_i end_POSTSUBSCRIPT, where NND,isubscript𝑁ND𝑖N_{\text{ND},i}italic_N start_POSTSUBSCRIPT ND , italic_i end_POSTSUBSCRIPT is the fraction of injected sources that are not detected, i.e., with a predicted flux below ϕUL⁢(xi)subscriptitalic-ϕULsubscript𝑥𝑖\phi_{\rm UL}(x_{i})italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ). Therefore:

∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ=ρi.superscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕsubscript𝜌𝑖\int_{0}^{\phi_{\rm UL}(x_{i})}P(\phi|x_{i})\text{d}\phi=\rho_{i}.∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) d italic_ϕ = italic_ρ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT .

(C11)

Let us call Ntotsubscript𝑁totN_{\rm tot}italic_N start_POSTSUBSCRIPT roman_tot end_POSTSUBSCRIPT the total number of sources injected on the full sky. Then we have

Ntot,i=Ntot⁢PGW⁢(xi)⁢Δ⁢Ωi,subscript𝑁totisubscript𝑁totsubscript𝑃GWsubscript𝑥𝑖ΔsubscriptΩ𝑖N_{\rm\text{tot},i}=N_{\rm\text{tot}}P_{\text{GW}}(x_{i})\Delta\Omega_{i},italic_N start_POSTSUBSCRIPT tot , roman_i end_POSTSUBSCRIPT = italic_N start_POSTSUBSCRIPT tot end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ,

(C12)

and therefore

∑xi∉Ω⊕PGW⁢(xi)⁢Δ⁢Ωi⁢∫0ϕUL⁢(xi)P⁢(ϕ|xi)⁢d⁢ϕ=∑xi∉Ω⊕Ntot,iNtot⁢ρi=1Ntot⁢∑xi∉Ω⊕NND,i.subscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑃GWsubscript𝑥𝑖ΔsubscriptΩ𝑖superscriptsubscript0subscriptitalic-ϕULsubscript𝑥𝑖𝑃conditionalitalic-ϕsubscript𝑥𝑖ditalic-ϕsubscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑁tot𝑖subscript𝑁totsubscript𝜌𝑖1subscript𝑁totsubscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑁ND𝑖\sum_{x_{i}\notin\Omega_{\earth}}P_{\text{GW}}(x_{i})\Delta\Omega_{i}\int_{0}^% {\phi_{\rm UL}(x_{i})}P(\phi|x_{i})\text{d}\phi=\sum_{x_{i}\notin\Omega_{% \earth}}\frac{N_{\text{tot},i}}{N_{\rm tot}}\rho_{i}=\dfrac{1}{N_{\rm tot}}% \sum_{x_{i}\notin\Omega_{\earth}}N_{\text{ND},i}.∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT GW end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) roman_Δ roman_Ω start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∫ start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_ϕ start_POSTSUBSCRIPT roman_UL end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT italic_P ( italic_ϕ | italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) d italic_ϕ = ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT divide start_ARG italic_N start_POSTSUBSCRIPT tot , italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_N start_POSTSUBSCRIPT roman_tot end_POSTSUBSCRIPT end_ARG italic_ρ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = divide start_ARG 1 end_ARG start_ARG italic_N start_POSTSUBSCRIPT roman_tot end_POSTSUBSCRIPT end_ARG ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_N start_POSTSUBSCRIPT ND , italic_i end_POSTSUBSCRIPT .

(C13)

Since the total number of injected sources not occulted by Earth are Ntot,∉⊕=(1−ε⊕)⁢Ntotsubscript𝑁totabsent⊕1subscript𝜀⊕subscript𝑁totN_{\rm tot,\notin\earth}=(1-\varepsilon_{\earth})N_{\rm tot}italic_N start_POSTSUBSCRIPT roman_tot , ∉ ⊕ end_POSTSUBSCRIPT = ( 1 - italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT ) italic_N start_POSTSUBSCRIPT roman_tot end_POSTSUBSCRIPT, and using that

ρ∉⊕=1Ntot,∉⊕⁢∑xi∉Ω⊕NND,i,subscript𝜌absent⊕1subscript𝑁totabsent⊕subscriptsubscript𝑥𝑖subscriptΩ⊕subscript𝑁ND𝑖\rho_{\notin\earth}=\frac{1}{N_{\rm tot,\notin\earth}}\sum_{x_{i}\notin\Omega_% {\earth}}N_{\text{ND},i},italic_ρ start_POSTSUBSCRIPT ∉ ⊕ end_POSTSUBSCRIPT = divide start_ARG 1 end_ARG start_ARG italic_N start_POSTSUBSCRIPT roman_tot , ∉ ⊕ end_POSTSUBSCRIPT end_ARG ∑ start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∉ roman_Ω start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_N start_POSTSUBSCRIPT ND , italic_i end_POSTSUBSCRIPT ,

(C14)

we finally recover Eq. (C10).

In order to quantify the difference between the two methods, the following test is performed. Having fixed the two parameters f𝑓fitalic_f and L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT, we compute the likelihood ℒℒ\mathcal{L}caligraphic_L for the two methods and we derive the quantity

εℒ=2⁢abs⁢(ℒ1−ℒ2)ℒ1+ℒ2.subscript𝜀ℒ2abssubscriptℒ1subscriptℒ2subscriptℒ1subscriptℒ2\varepsilon_{\mathcal{L}}=2\frac{\text{abs}(\mathcal{L}_{1}-\mathcal{L}_{2})}{% \mathcal{L}_{1}+\mathcal{L}_{2}}.italic_ε start_POSTSUBSCRIPT caligraphic_L end_POSTSUBSCRIPT = 2 divide start_ARG abs ( caligraphic_L start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT - caligraphic_L start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) end_ARG start_ARG caligraphic_L start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + caligraphic_L start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG .

(C15)

Here we use the subscripts 1 and 2 for the respective methods. Both likelihoods are computed considering only BBH candidates with pastro>0.5subscript𝑝astro0.5p_{\rm astro}>0.5italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT > 0.5. In both cases, the total number of injected sources for each BBH candidate is Ntot=1000subscript𝑁tot1000N_{\rm tot}=1000italic_N start_POSTSUBSCRIPT roman_tot end_POSTSUBSCRIPT = 1000. The distribution of εℒsubscript𝜀ℒ\varepsilon_{\mathcal{L}}italic_ε start_POSTSUBSCRIPT caligraphic_L end_POSTSUBSCRIPT is evaluated sampling randomly f𝑓fitalic_f and L0subscript𝐿0L_{0}italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT, for a total of 100 sampled pairs (f,L0𝑓subscript𝐿0f,L_{0}italic_f , italic_L start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT). We obtain that the median value of εℒsubscript𝜀ℒ\varepsilon_{\mathcal{L}}italic_ε start_POSTSUBSCRIPT caligraphic_L end_POSTSUBSCRIPT is 0.04 and that in ∼80%similar-toabsentpercent80\sim 80\%∼ 80 % of the sampled cases εℒ<0.2subscript𝜀ℒ0.2\varepsilon_{\mathcal{L}}<0.2italic_ε start_POSTSUBSCRIPT caligraphic_L end_POSTSUBSCRIPT < 0.2. Since the difference between the two methods is limited and since the Method 2 is more computationally expensive, all the results are used adopting Method 1.

Table 1: *

List of 636 low latency GW triggers analyzed using NITRATES is shown along with their respective pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT values and 15–350 keV band flux upper limits. The maximum TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG is indicated for all the triggers with successful NITRATES results. Observations corresponding to triggers with insufficient exposure time during the BAT pointing mode do not have valid NITRATES results or flux upper limits. For those triggers that do have NITRATES results but fail to meet the criterion for a full likelihood analysis, the max TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG is indicated as NFL (No Final Likelihood). The GW triggers from the Burst pipeline do not have associated pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT values and are therefore left blank. The fraction of the GW sky posterior distribution inside the BAT coded FOV and the fraction of the GW posterior occulted by the Earth are denoted by εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT and ε⊕subscript𝜀⊕\varepsilon_{\earth}italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT, respectively.

SID

Time

GW FAR

group

pastrosubscript𝑝astrop_{\rm astro}italic_p start_POSTSUBSCRIPT roman_astro end_POSTSUBSCRIPT

Class

TSTS\sqrt{\text{TS}}square-root start_ARG TS end_ARG

Flux UL

εin⁢BATsubscript𝜀inBAT\varepsilon_{\rm in\,BAT}italic_ε start_POSTSUBSCRIPT roman_in roman_BAT end_POSTSUBSCRIPT

ε⊕subscript𝜀⊕\varepsilon_{\earth}italic_ε start_POSTSUBSCRIPT ⊕ end_POSTSUBSCRIPT

(UTC)

(Hz)

(erg cm-2 s-1)

(%)

(%)

S190701ah

2019-07-01T20:33:07

1.92×\times×10-8

CBC

0.934

BBH

6.4

1.58×\times×10-7

99.42

0

S190816i

2019-08-16T13:04:31

1.44×\times×10-8

CBC

0.833

NSBH

-

-

-

-

S190828af

2019-08-28T17:51:02

1.83×\times×10-5

CBC

0.012

BNS

-

-

-

-

S190829p

2019-08-29T13:49:01

2.59×\times×10-6

CBC

0.062

Mass Gap

5.8

2.05×\times×10-7

58.03

0.94

S190830y

2019-08-30T15:07:04

7.70×\times×10-8

CBC

0.382

Mass Gap

5.7

1.45×\times×10-7

78.14

5.94

S190831ai

2019-08-31T18:31:02

8.85×\times×10-6

CBC

0.064

NSBH

5.7

8.38×\times×10-7

1.46

41.91

S190901al

2019-09-01T21:01:03

8.78×\times×10-6

CBC

0.015

NSBH

6.6

1.82×\times×10-7

80.18

1.08

S190901d

2019-09-01T02:56:47

1.01×\times×10-5

CBC

0.018

BNS

6.5

5.33×\times×10-7

13.86

43.34

S190901h

2019-09-01T04:38:54

5.26×\times×10-6

Burst

-

-

6.3

4.14×\times×10-7

24.56

43.03

S190902ao

2019-09-02T20:56:00

7.21×\times×10-6

CBC

0.016

NSBH

6.1

5.41×\times×10-7

9.78

24.09

S190904c

2019-09-04T02:59:52

3.64×\times×10-6

CBC

0.037

Mass Gap

5.2

5.35×\times×10-7

11.35

17.04

S190904p

2019-09-04T12:32:03

3.75×\times×10-6

CBC

0.026

NSBH

-

-

-

-

S190904w

2019-09-04T17:49:01

1.56×\times×10-6

Burst

-

-

5.8

4.53×\times×10-7

27.01

0

S190906ad

2019-09-06T18:33:04

1.41×\times×10-5

CBC

0.010

BNS

5.9

5.43×\times×10-7

15.61

46.69

S190906ag

2019-09-06T19:35:03

3.22×\times×10-6

Burst

-

-

5.5

5.73×\times×10-7

24.48

28.62

S190906ah

2019-09-06T20:05:00

8.90×\times×10-7

CBC

0.101

NSBH

5.7

5.47×\times×10-7

12.62

19.27

S190906s

2019-09-06T15:20:02

4.71×\times×10-6

Burst

-

-

6.4

3.65×\times×10-7

51.54

6.26

S190907n

2019-09-07T14:29:05

2.71×\times×10-6

Burst

-

-

NFL

5.91×\times×10-7

16.71

10.63

S190908az

2019-09-08T21:34:01

4.26×\times×10-7

Burst

-

-

5.7

7.04×\times×10-7

8.59

23.3

S190908e

2019-09-08T02:34:06

4.52×\times×10-6

Burst

-

-

5.5

2.50×\times×10-7

21.39

77.13

S190909ac

2019-09-09T14:13:01

4.54×\times×10-6

Burst

-

-

5.9

3.43×\times×10-7

16.23

21.97

S190909aw

2019-09-09T19:41:05

1.07×\times×10-6

Burst

-

-

5.8

2.29×\times×10-7

70.87

2.4

S190909bd

2019-09-09T21:26:03

8.85×\times×10-6

CBC

0.023

Mass Gap

-

-

-

-

S190909y

2019-09-09T12:49:01

1.66×\times×10-6

CBC

0.134

NSBH

-

-

-

-

S190915ak

2019-09-15T23:57:02

9.74×\times×10-10

CBC

0.990

BBH

5.4

1.33×\times×10-7

87.31

0.17

S190915q

2019-09-15T16:03:01

2.66×\times×10-6

Burst

-

-

5.8

4.81×\times×10-7

11.92

13.01

S190916y

2019-09-16T15:55:01

9.70×\times×10-7

CBC

0.143

BNS

6.9

5.58×\times×10-7

10.2

27.57

S190917ad

2019-09-17T19:14:00

1.47×\times×10-5

CBC

0.013

Mass Gap

6.4

4.47×\times×10-7

4.18

3.5

S190918aa

2019-09-18T19:38:04

6.68×\times×10-7

CBC

0.023

BNS

7.2

1.09×\times×10-7

64.36

16.27

S190919ag

2019-09-19T17:58:02

3.42×\times×10-6

Burst

-

-

7.2

2.77×\times×10-7

14.79

56.16

S190919ak

2019-09-19T18:34:03

1.06×\times×10-5

CBC

0.089

BNS

5.8

3.11×\times×10-7

28.36

35.45

S190919au

2019-09-19T20:39:00

3.24×\times×10-6

Burst

-

-

9.4

6.34×\times×10-7

0.04

42.88

S190919u

2019-09-19T12:13:02

8.18×\times×10-6

Burst

-

-

8.0

3.84×\times×10-7

20.56

23.68

S190920an

2019-09-20T19:09:05

2.77×\times×10-7

Burst

-

-

5.7

3.16×\times×10-7

15.6

36.74

S190920ap

2019-09-20T19:27:04

5.51×\times×10-6

Burst

-

-

5.8

2.95×\times×10-7

11.49

46.04

S190920z

2019-09-20T12:55:04

6.08×\times×10-6

Burst

-

-

6.1

1.02×\times×10-7

52.58

38.59

S190922ag

2019-09-22T15:22:01

2.89×\times×10-6

Burst

-

-

5.4

3.02×\times×10-7

50.45

9.47

S190922aq

2019-09-22T18:08:05

1.63×\times×10-5

CBC

0.013

Mass Gap

6.2

-

-

-

S190923aj

2019-09-23T17:04:04

1.91×\times×10-7

Burst

-

-

-

-

-

-

S190923ak

2019-09-23T17:06:03

2.70×\times×10-6

CBC

0.216

BBH

-

-

-

-

S190923x

2019-09-23T12:19:00

2.10×\times×10-6

CBC

0.036

BNS

-

-

-

-

S190923y

2019-09-23T12:55:59

4.78×\times×10-8

CBC

0.670

NSBH

-

-

-

-

S190926z

2019-09-26T16:47:02

1.27×\times×10-6

Burst

-

-

6.2

3.37×\times×10-7

1.81

72.24

S190927an

2019-09-27T14:58:00

3.60×\times×10-6

CBC

0.038

BNS

6.6

4.59×\times×10-7

0.02

58.85

S190928c

2019-09-28T02:11:45

6.73×\times×10-9

Burst

-

-

6.9

2.84×\times×10-7

0.04

1.53

S190928j

2019-09-28T06:30:16

2.51×\times×10-6

CBC

0.092

BNS

6.4

5.64×\times×10-7

8.9

27.8

S190930s

2019-09-30T13:35:41

3.00×\times×10-9

CBC

0.950

Mass Gap

-

-

-

-

S190930t

2019-09-30T14:34:07

1.54×\times×10-8

CBC

0.74

NSBH

5.9

5.15×\times×10-7

11.69

21.15

S191105d

2019-11-05T13:40:51

6.63×\times×10-6

Burst

-

-

-

-

-

-

S191106r

2019-11-06T18:41:51

3.31×\times×10-6

Burst

-

-

-

-

-

-

S191107o

2019-11-07T16:05:23

5.62×\times×10-6

Burst

-

-

5.6

6.00×\times×10-7

3.02

17.27

S191107t

2019-11-07T18:03:55

4.14×\times×10-6

CBC

0.020

NSBH

6.0

4.35×\times×10-7

0.06

17.98

S191110w

2019-11-10T16:48:32

4.43×\times×10-6

Burst

-

-

6.6

1.55×\times×10-7

27.57

47.27

S191110x

2019-11-10T18:08:42

2.93×\times×10-11

CBC

0.999

Mass Gap

7.2

1.75×\times×10-7

15.22

56.88

S191112n

2019-11-12T04:43:25

1.76×\times×10-5

Burst

-

-

5.3

1.09×\times×10-6

38.68

17.28

S191113aj

2019-11-13T14:28:49

2.31×\times×10-5

CBC

0.005

BNS

6.5

2.18×\times×10-7

16.07

0.01

S191114ad

2019-11-14T12:58:04

1.36×\times×10-5

CBC

0.065

BBH

6.0

2.05×\times×10-7

30.58

22

S191114am

2019-11-14T15:39:15

1.57×\times×10-5

CBC

0.021

BBH

6.3

5.12×\times×10-7

0

1.97

S191114at

2019-11-14T16:16:17

8.13×\times×10-6

CBC

0.008

NSBH

6.7

3.66×\times×10-7

15.85

50.49

S191115be

2019-11-15T23:07:27

1.05×\times×10-5

CBC

0.010

NSBH

6.0

9.94×\times×10-7

3

0.64

S191116ac

2019-11-16T14:21:55

9.04×\times×10-6

CBC

0.015

NSBH

NFL

4.95×\times×10-7

8.53

1.28

S191118n

2019-11-18T07:59:05

5.88×\times×10-6

CBC

0.018

NSBH

6.8

8.70×\times×10-8

85.8

6.04

S191118z

2019-11-18T16:49:55

7.31×\times×10-7

CBC

0.164

BNS

6.2

1.17×\times×10-7

38.57

54.47

S191121bf

2019-11-21T13:13:24

3.70×\times×10-6

Burst

-

-

-

-

-

-

S191121bq

2019-11-21T15:54:12

2.72×\times×10-6

Burst

-

-

5.7

2.73×\times×10-7

38.34

29.35

S191121bt

2019-11-21T16:45:42

2.03×\times×10-5

CBC

0.004

NSBH

5.5

3.06×\times×10-7

6.25

7.14

S191123q

2019-11-23T09:01:14

1.07×\times×10-5

Burst

-

-

5.4

5.84×\times×10-7

2.04

14.18

S191127p

2019-11-27T05:02:27

2.63×\times×10-6

CBC

0.037

Mass Gap

-

-

-

-

S191130q

2019-11-30T07:52:23

8.69×\times×10-6

CBC

0.005

NSBH

5.0

5.55×\times×10-7

1.24

28.9

S191202af

2019-12-02T18:42:26

2.36×\times×10-6

Burst

-

-

-

-

-

-

S191204o

2019-12-04T14:17:13

1.16×\times×10-5

CBC

0.009

NSBH

-

-

-

-

S191204r

2019-12-04T17:15:26

3.06×\times×10-25

CBC

1.000

BBH

NFL

1.21×\times×10-7

86.69

0

S191204t

2019-12-04T18:34:16

1.67×\times×10-6

Burst

-

-

5.1

3.73×\times×10-7

15.06

24.91

S191205ae

2019-12-05T20:56:37

2.83×\times×10-7

Burst

-

-

5.7

6.07×\times×10-7

7.73

3.8

S191205ah

2019-12-05T21:52:08

1.25×\times×10-8

CBC

0.932

NSBH

5.2

3.43×\times×10-7

31.7

7.85

S191206ab

2019-12-06T14:05:21

6.19×\times×10-6

CBC

0.024

Mass Gap

6.2

2.82×\times×10-7

34.03

8.59

S191206an

2019-12-06T17:38:57

1.01×\times×10-6

Burst

-

-

5.7

3.14×\times×10-7

5.28

7.69

S191207o

2019-12-07T10:16:32

1.42×\times×10-5

Burst

-

-

-

-

-

-

S191207u

2019-12-07T12:29:56

1.00×\times×10-5

Burst

-

-

6.1

3.61×\times×10-7

28.24

32.09

S191208b

2019-12-08T02:02:15

9.14×\times×10-6

CBC

0.017

BNS

5.0

4.92×\times×10-7

0.03

41.23

S191209ao

2019-12-09T13:58:21

1.02×\times×10-6

CBC

0.097

Mass Gap

5.6

4.34×\times×10-7

0.37

39.42

S191209ar

2019-12-09T14:32:42

1.08×\times×10-6

Burst

-

-

6.2

2.61×\times×10-7

56.27

8.52

S191212ad

2019-12-12T16:57:39

2.55×\times×10-6

Burst

-

-

7.0

5.43×\times×10-6

0.17

17.49

S191212ap

2019-12-12T19:59:21

1.31×\times×10-6

CBC

0.080

BNS

6.0

4.98×\times×10-7

0.67

4.94

S191212b

2019-12-12T00:31:07

9.49×\times×10-6

Burst

-

-

4.9

4.94×\times×10-7

7.29

39.35

S191212l

2019-12-12T07:57:05

9.31×\times×10-6

CBC

0.021

Mass Gap

7.2

4.24×\times×10-7

11.14

69.95

S191213al

2019-12-13T16:09:04

1.27×\times×10-7

CBC

0.518

NSBH

6.3

3.61×\times×10-7

16.31

35.94

S191213an

2019-12-13T16:58:32

8.60×\times×10-7

Burst

-

-

6.2

2.92×\times×10-7

50.3

30.27

S191213au

2019-12-13T18:44:42

7.84×\times×10-7

Burst

-

-

6.0

1.47×\times×10-7

9.88

87.84

S191213ay

2019-12-13T19:16:25

2.93×\times×10-6

Burst

-

-

6.3

5.88×\times×10-7

0.26

11.37

S191213be

2019-12-13T19:54:22

1.72×\times×10-6

Burst

-

-

5.6

6.80×\times×10-7

1.1

79.01

S191213c

2019-12-13T01:17:45

6.71×\times×10-8

CBC

0.395

NSBH

5.9

4.87×\times×10-8

98.61

0

S191215r

2019-12-15T19:57:29

5.46×\times×10-6

CBC

0.023

BNS

6.3

9.51×\times×10-8

73.13

12.09

S191216ap

2019-12-16T21:33:38

1.13×\times×10-23

CBC

1.000

Mass Gap

6.6

6.61×\times×10-7

0.56

4.99

S191219ak

2019-12-19T17:49:47

1.22×\times×10-6

Burst

-

-

7.2

1.40×\times×10-7

83.49

3.23

S191219an

2019-12-19T18:36:24

2.26×\times×10-6

Burst

-

-

6.0

4.64×\times×10-7

0.01

38.33

S191219ap

2019-12-19T19:52:52

3.14×\times×10-6

CBC

0.067

Mass Gap

6.3

3.05×\times×10-7

21.87

42.93

S191220af

2019-12-20T12:24:14

3.96×\times×10-10

CBC

0.996

BNS

6.4

4.93×\times×10-7

4.63

41.78

S191220al

2019-12-20T14:46:02

1.92×\times×10-6

CBC

0.026

BNS

-

-

-

-

S191220aw

2019-12-20T17:49:42

9.31×\times×10-8

CBC

0.522

NSBH

-

-

-

-

S191221aa

2019-12-21T10:31:37

1.02×\times×10-5

CBC

0.018

BNS

5.6

5.73×\times×10-7

52.7

1.74

S191221al

2019-12-21T14:41:21

2.92×\times×10-6

Burst

-

-

5.8

4.53×\times×10-7

5.41

24.69

S191221ar

2019-12-21T17:12:28

1.22×\times×10-5

CBC

0.011

NSBH

6.7

7.21×\times×10-7

0.72

7.38

S191221v

2019-12-21T08:51:06

1.43×\times×10-6

CBC

0.003

NSBH

6.3

1.33×\times×10-7

97.95

0

S191221w

2019-12-21T09:02:03

1.76×\times×10-5

CBC

0.055

Mass Gap

4.6

5.28×\times×10-7

8.91

6.08

S191222a

2019-12-22T01:34:42

8.95×\times×10-6

CBC

0.020

BNS

6.6

9.11×\times×10-7

0.03

13.45

S191222af

2019-12-22T13:57:46

1.94×\times×10-5

CBC

0.006

NSBH

-

-

-

-

S191222an

2019-12-22T16:30:03

1.79×\times×10-5

CBC

0.037

Mass Gap

-

-

-

-

S191223aj

2019-12-23T15:55:41

2.23×\times×10-5

CBC

0.006

BNS

5.9

2.57×\times×10-7

51.53

7.48

S191223p

2019-12-23T08:22:49

1.54×\times×10-5

CBC

0.008

NSBH

5.7

1.22×\times×10-6

3.89

6.33

S191224p

2019-12-24T05:03:59

8.31×\times×10-6

Burst

-

-

6.1

5.61×\times×10-7

10.96

17.7

S191224x

2019-12-24T11:23:11

1.98×\times×10-5

CBC

0.011

BNS

6.3

5.16×\times×10-7

22.09

8.42

S191225aq

2019-12-25T21:57:15

1.27×\times×10-8

CBC

0.390

Mass Gap

5.7

1.78×\times×10-7

55.64

0.05

S191225e

2019-12-25T02:11:26

9.64×\times×10-6

CBC

0.017

BNS

6.0

2.49×\times×10-7

37.71

36.13

S191225q

2019-12-25T10:30:40

1.69×\times×10-5

CBC

0.006

NSBH

6.0

4.39×\times×10-7

19.53

13.61

S191226ad

2019-12-26T13:31:05

1.19×\times×10-5

CBC

0.011

BNS

5.8

3.58×\times×10-7

9.16

13.79

S191226ae

2019-12-26T14:39:20

1.99×\times×10-6

CBC

0.026

NSBH

5.6

2.49×\times×10-7

21.45

4.19

S191226ai

2019-12-26T17:35:27

2.96×\times×10-7

CBC

0.602

BBH

-

-

-

-

S191226aj

2019-12-26T18:10:18

2.46×\times×10-6

CBC

0.073

Mass Gap

6.6

4.71×\times×10-7

0.26

65.36

S191226ap

2019-12-26T20:33:18

5.04×\times×10-6

Burst

-

-

6.0

4.36×\times×10-7

14.53

51.45

S191226d

2019-12-26T01:40:51

2.24×\times×10-6

Burst

-

-

5.7

4.80×\times×10-7

0.01

43.14

S191226u

2019-12-26T10:24:57

9.17×\times×10-6

Burst

-

-

-

-

-

-

S191227aa

2019-12-27T11:47:25

1.60×\times×10-5

CBC

0.006

NSBH

-

-

-

-

S191227af

2019-12-27T13:06:49

2.79×\times×10-6

Burst

-

-

6.6

9.85×\times×10-8

57.79

35.03

S191227aj

2019-12-27T14:54:31

5.45×\times×10-6

CBC

0.030

Mass Gap

5.9

2.56×\times×10-7

22.12

13.71

S191227al

2019-12-27T15:49:54

1.62×\times×10-5

CBC

0.010

BNS

6.2

2.71×\times×10-7

41.93

20.28

S191227am

2019-12-27T15:51:27

6.04×\times×10-6

Burst

-

-

-

-

-

-

S191227an

2019-12-27T16:10:45

8.60×\times×10-6

CBC

0.014

NSBH

6.3

7.13×\times×10-7

29.95

0.15

S191227ap

2019-12-27T16:35:12

1.68×\times×10-5

Burst

-

-

5.9

2.18×\times×10-7

44.6

27.51

S191227as

2019-12-27T17:29:03

1.32×\times×10-5

Burst

-

-

6.1

2.87×\times×10-7

21.77

35.26

S191227az

2019-12-27T21:55:49

1.21×\times×10-5

CBC

0.008

NSBH

6.9

5.10×\times×10-7

4.48

0.34

S191227bb

2019-12-27T23:04:40

3.86×\times×10-7

Burst

-

-

-

-

-

-

S191227h

2019-12-27T02:58:36

1.94×\times×10-5

CBC

0.034

BBH

6.5

8.67×\times×10-7

9.06

26.08

S191227o

2019-12-27T04:55:18

1.17×\times×10-5

Burst

-

-

6.8

5.77×\times×10-7

14.35

1.91

S191228ac

2019-12-28T13:17:16

1.44×\times×10-5

CBC

0.007

NSBH

7.0

2.43×\times×10-7

16.91

65.16

S191228am

2019-12-28T18:51:01

7.20×\times×10-6

Burst

-

-

-

-

-

-

S191228an

2019-12-28T20:07:11

1.83×\times×10-6

Burst

-

-

6.4

2.49×\times×10-7

48.48

0.41

S191228at

2019-12-28T23:57:39

8.92×\times×10-6

Burst

-

-

4.6

2.99×\times×10-7

36.08

31.05

S191228i

2019-12-28T05:44:50

2.27×\times×10-5

CBC

0.008

Mass Gap

6.5

3.47×\times×10-7

38.65

34.86

S191228q

2019-12-28T08:14:37

5.18×\times×10-6

Burst

-

-

5.8

4.05×\times×10-7

27.39

34.8

S191228u

2019-12-28T09:08:39

9.00×\times×10-6

Burst

-

-

5.4

4.01×\times×10-7

10.1

41.03

S191228w

2019-12-28T09:49:35

1.01×\times×10-5

Burst

-

-

6.4

4.10×\times×10-7

20.77

26.8

S191229ah

2019-12-29T21:50:21

1.23×\times×10-6

Burst

-

-

8.9

7.45×\times×10-7

22.06

23.56

S191229ai

2019-12-29T22:11:21

7.60×\times×10-6

CBC

0.015

NSBH

6.4

1.56×\times×10-7

30.56

60.52

S191229ak

2019-12-29T23:16:09

9.76×\times×10-6

Burst

-

-

5.7

5.81×\times×10-7

7.89

40.58

S191229o

2019-12-29T12:02:34

1.07×\times×10-5

CBC

0.011

NSBH

-

-

-

-

S191230aa

2019-12-30T13:51:30

1.07×\times×10-5

CBC

0.011

NSBH

6.6

4.68×\times×10-7

55.93

0.06

S191230ae

2019-12-30T14:19:12

1.44×\times×10-5

CBC

0.011

Mass Gap

5.8

1.89×\times×10-7

34.17

54.3

S191230at

2019-12-30T21:24:48

1.64×\times×10-5

CBC

0.011

BNS

6.1

5.52×\times×10-7

0.09

48.46

S191230au

2019-12-30T22:04:37

2.12×\times×10-6

Burst

-

-

7.0

2.54×\times×10-7

40.47

20.44

S191230e

2019-12-30T02:40:45

1.09×\times×10-5

Burst

-

-

-

-

-

-

S191230k

2019-12-30T04:10:08

3.17×\times×10-7

Burst

-

-

-

-

-

-

S191230v

2019-12-30T11:19:08

8.75×\times×10-6

CBC

0.013

NSBH

6.8

5.44×\times×10-7

0.36

1.48

S191230y

2019-12-30T13:08:19

1.01×\times×10-5

Burst

-

-

5.7

4.53×\times×10-7

2.13

52.59

S191231ad

2019-12-31T11:45:12

1.33×\times×10-6

Burst

-

-

6.9

1.21×\times×10-7

51.96

34.52

S191231an

2019-12-31T16:59:30

1.67×\times×10-5

CBC

0.009

NSBH

6.2

5.52×\times×10-7

11.42

7.79

S200101o

2020-01-01T14:18:13

1.66×\times×10-5

CBC

0.012

BNS

7.4

5.29×\times×10-7

1.52

61.15

S200102ah

2020-01-02T15:04:48

1.70×\times×10-5

CBC

0.007

NSBH

5.3

5.99×\times×10-7

0

38.48

S200102an

2020-01-02T18:05:23

1.45×\times×10-5

Burst

-

-

6.5

5.94×\times×10-7

1.94

32.62

S200102ar

2020-01-02T19:39:36

2.24×\times×10-5

CBC

0.005

NSBH

6.4

1.77×\times×10-7

73.19

2.04

S200102au

2020-01-02T21:01:56

1.77×\times×10-5

CBC

0.009

Mass Gap

5.6

5.10×\times×10-7

6.1

8.29

S200102k

2020-01-02T06:17:35

5.93×\times×10-6

CBC

0.016

NSBH

-

-

-

-

S200102y

2020-01-02T11:15:25

8.42×\times×10-6

Burst

-

-

-

-

-

-

S200103aa

2020-01-03T12:32:22

6.00×\times×10-6

Burst

-

-

5.7

9.42×\times×10-7

1.26

26.1

S200103am

2020-01-03T16:46:33

1.68×\times×10-5

CBC

0.007

NSBH

-

-

-

-

S200103ao

2020-01-03T18:29:37

1.12×\times×10-5

CBC

0.013

NSBH

7.0

3.23×\times×10-6

1.05

67.67

S200103aw

2020-01-03T22:34:12

2.00×\times×10-5

CBC

0.002

BBH

5.2

5.28×\times×10-7

0.22

25.17

S200103az

2020-01-03T23:31:11

1.32×\times×10-6

CBC

0.078

NSBH

6.6

2.58×\times×10-7

1.42

97.64

S200103r

2020-01-03T09:42:24

1.57×\times×10-6

Burst

-

-

6.4

2.60×\times×10-7

27.35

36.86

S200103t

2020-01-03T10:31:18

4.10×\times×10-6

Burst

-

-

5.9

7.13×\times×10-7

10.63

21.96

S200103v

2020-01-03T10:55:34

1.18×\times×10-5

Burst

-

-

6.4

3.68×\times×10-7

24.33

31.88

S200103z

2020-01-03T11:55:03

1.49×\times×10-5

Burst

-

-

-

-

-

-

S200104aa

2020-01-04T10:04:38

9.99×\times×10-6

Burst

-

-

6.1

4.74×\times×10-7

6.58

43.27

S200104ar

2020-01-04T19:50:42

3.13×\times×10-6

CBC

0.024

NSBH

6.7

5.16×\times×10-7

0.4

26.15

S200104d

2020-01-04T04:13:54

2.99×\times×10-6

Burst

-

-

6.2

3.78×\times×10-7

8.92

44.22

S200104r

2020-01-04T08:05:48

1.70×\times×10-5

CBC

0.010

BNS

6.0

4.83×\times×10-7

0.03

12.34

S200105aj

2020-01-05T18:00:59

1.43×\times×10-5

CBC

0.010

BNS

-

-

-

-

S200105p

2020-01-05T09:03:23

9.27×\times×10-6

CBC

0.019

BNS

6.0

2.87×\times×10-7

18.21

77.63

S200105u

2020-01-05T12:01:59

2.26×\times×10-5

CBC

0.005

NSBH

6.8

8.27×\times×10-7

1.99

41.47

S200105w

2020-01-05T12:48:13

6.15×\times×10-6

CBC

0.012

NSBH

5.6

5.22×\times×10-7

1.4

69.28

S200106ar

2020-01-06T17:48:06

3.03×\times×10-6

Burst

-

-

6.4

1.54×\times×10-7

58.45

21.4

S200106az

2020-01-06T18:50:35

1.39×\times×10-5

Burst

-

-

6.5

3.42×\times×10-7

13.15

69.28

S200106bd

2020-01-06T22:24:59

1.76×\times×10-5

CBC

0.005

BNS

5.8

8.99×\times×10-7

0.02

99.56

S200106f

2020-01-06T01:36:45

1.76×\times×10-6

Burst

-

-

6.0

6.32×\times×10-7

4.58

18.82

S200106i

2020-01-06T03:07:57

1.82×\times×10-5

CBC

0.015

Mass Gap

5.5

4.44×\times×10-7

14.62

28.82

S200106k

2020-01-06T04:37:09

1.86×\times×10-5

CBC

0.009

Mass Gap

7.0

3.30×\times×10-7

17.36

36.7

S200106s

2020-01-06T08:37:45

8.19×\times×10-6

Burst

-

-

-

-

-

-

S200107i

2020-01-07T03:16:26

5.14×\times×10-6

Burst

-

-

5.4

5.52×\times×10-7

0.1

23.85

S200107j

2020-01-07T03:22:04

3.44×\times×10-7

Burst

-

-

6.5

6.19×\times×10-7

0

7.88

S200107m

2020-01-07T04:08:28

4.26×\times×10-6

CBC

0.022

NSBH

5.6

4.71×\times×10-7

13.51

34.62

S200107o

2020-01-07T05:11:52

1.05×\times×10-5

CBC

0.006

NSBH

5.3

6.16×\times×10-7

0.77

24.76

S200108ag

2020-01-08T18:30:05

2.24×\times×10-5

CBC

0.009

BNS

6.6

3.40×\times×10-7

37.19

28.06

S200108ah

2020-01-08T18:42:57

2.16×\times×10-5

CBC

0.008

NSBH

7.1

3.15×\times×10-7

25.82

29.67

S200108an

2020-01-08T23:51:51

1.66×\times×10-5

CBC

0.011

BNS

6.8

4.67×\times×10-7

34.39

12.38

S200108j

2020-01-08T03:42:27

2.24×\times×10-5

CBC

0.007

NSBH

5.7

5.72×\times×10-7

15.67

30.15

S200108l

2020-01-08T04:13:13

2.13×\times×10-6

CBC

0.092

BNS

6.4

3.53×\times×10-7

10.23

41.69

S200108p

2020-01-08T05:20:09

1.93×\times×10-7

CBC

0.470

BNS

7.4

7.83×\times×10-7

2.39

46.3

S200109m

2020-01-09T08:48:21

1.94×\times×10-5

CBC

0.006

NSBH

5.9

4.91×\times×10-7

11.43

3.39

S200109o

2020-01-09T13:51:35

1.44×\times×10-5

Burst

-

-

5.8

1.09×\times×10-7

77.72

9.11

S200109r

2020-01-09T15:30:31

3.14×\times×10-6

Burst

-

-

7.0

5.56×\times×10-7

0

99.4

S200109s

2020-01-09T15:44:38

4.34×\times×10-6

Burst

-

-

7.0

4.75×\times×10-7

22.85

20.23

S200110aa

2020-01-10T11:01:48

2.12×\times×10-5

Burst

-

-

7.1

2.92×\times×10-7

47.55

21.69

S200110d

2020-01-10T01:23:11

1.45×\times×10-5

Burst

-

-

7.3

5.40×\times×10-7

7.42

39.3

S200110e

2020-01-10T02:01:40

6.32×\times×10-6

Burst

-

-

5.8

3.77×\times×10-7

0.77

14.67

S200110m

2020-01-10T05:25:05

1.21×\times×10-6

Burst

-

-

6.7

6.61×\times×10-7

1.48

5.31

S200110q

2020-01-10T05:46:19

4.45×\times×10-6

Burst

-

-

6.4

4.99×\times×10-7

0.76

61.25

S200110s

2020-01-10T06:46:45

1.51×\times×10-5

Burst

-

-

6.9

1.35×\times×10-7

79.96

3.4

S200110t

2020-01-10T07:50:59

1.79×\times×10-5

Burst

-

-

6.5

4.70×\times×10-7

11.22

28.08

S200110v

2020-01-10T08:40:52

1.56×\times×10-5

Burst

-

-

6.5

2.89×\times×10-7

25.33

9.06

S200110z

2020-01-10T10:33:06

6.50×\times×10-6

Burst

-

-

5.8

2.03×\times×10-7

56.99

10.49

S200111ae

2020-01-11T22:02:00

1.61×\times×10-5

CBC

0.009

NSBH

6.2

7.87×\times×10-7

11.26

25.72

S200111j

2020-01-11T06:51:45

2.70×\times×10-6

Burst

-

-

6.9

5.25×\times×10-7

12.19

27.34

S200111s

2020-01-11T15:23:44

1.04×\times×10-5

CBC

0.019

BNS

-

-

-

-

S200111w

2020-01-11T19:00:59

1.36×\times×10-5

CBC

0.011

NSBH

5.7

2.70×\times×10-7

60.18

0.99

S200112ac

2020-01-12T21:29:08

1.60×\times×10-5

CBC

0.004

NSBH

6.0

6.91×\times×10-7

1.12

79.05

S200112e

2020-01-12T09:44:25

1.61×\times×10-5

CBC

0.009

BNS

6.2

6.47×\times×10-7

0.01

58.06

S200113f

2020-01-13T02:14:20

1.79×\times×10-5

CBC

0.010

BNS

5.3

3.12×\times×10-7

18.59

56.9

S200113g

2020-01-13T02:20:40

1.81×\times×10-5

CBC

0.009

BNS

5.6

4.69×\times×10-7

3.35

0.16

S200113n

2020-01-13T09:59:40

1.57×\times×10-5

CBC

0.013

Mass Gap

6.4

6.05×\times×10-7

0.52