Joint gravitational wave-short GRB detection of Binary Neutron Star mergers with existing and future facilities

We explore the joint detection prospects of short gamma-ray bursts (sGRBs) and their gravitational wave (GW) counterparts by the current and upcoming high-energy GRB and GW facilities from binary neutron star (BNS) mergers. We consider two GW detector networks: (1) A four-detector network comprising LIGO Hanford, Livingston, Virgo, and Kagra, (IGWN4) and (2) a future five-detector network including the same four detectors and LIGO India (IGWN5). For the sGRB detection, we consider existing satellites Fermi and Swift and the proposed all-sky satellite Daksha. Most of the events for the joint detection will be off-axis, hence, we consider a broad range of sGRB jet models predicting the off-axis emission. Also, to test the effect of the assumed sGRB luminosity function, we consider two different functions for one of the emission models. We find that for the different jet models, the joint sGRB and GW detection rates for Fermi and Swift with IGWN4 (IGWN5) lie within 0.07-0.62$\mathrm{\ yr^{-1}}$ (0.8-4.0$\mathrm{\ yr^{-1}}$) and 0.02-0.14$\mathrm{\ yr^{-1}}$ (0.15-1.0$\mathrm{\ yr^{-1}}$), respectively, when the BNS merger rate is taken to be 320$\mathrm{\ Gpc^{-3}~yr^{-1}}$. With Daksha, the rates increase to 0.2-1.3$\mathrm{\ yr^{-1}}$ (1.3-8.3$\mathrm{\ yr^{-1}}$), which is 2-9 times higher than the existing satellites. We show that such a mission with higher sensitivity will be ideal for detecting a higher number of fainter events observed off-axis or at a larger distance. Thus, Daksha will boost the joint detections of sGRB and GW, especially for the off-axis events. Finally, we find that our detection rates with optimal SNRs are conservative, and noise in GW detectors can increase the rates further.


INTRODUCTION
Binary neutron star (BNS) mergers have long been hypothesized to be associated with short-duration gamma-ray bursts (sGRBs).This is observationally confirmed for the first time by the joint detection of gravitational waves (GW170817, Abbott et al. 2017) and the spatially coincident sGRB (GRB170817A, Savchenko et al. 2017;Goldstein et al. 2017) from a BNS merger.Except for the understanding of the sources of the sGRBs, the extensive multi-wavelength follow-up observations of the afterglow significantly improves our understanding of the sGRB jets and emission processes (Troja et al. 2017(Troja et al. , 2018;;Margutti et al. 2018;Lyman et al. 2018;Alexander et al. 2018;Mooley et al. 2018a,b;Ghirlanda et al. 2019;Troja et al. 2020;Beniamini et al. 2020aBeniamini et al. , 2022)).
An interesting aspect of the joint detection of GW and sGRB is that the exact time of the merger is known.Hence, much information about the jet (such as, jet launching time) and the propagation of the jet through the ejecta can be analyzed.Therefore, joint detection of sGRB and GW from BNS mergers is the ideal probe to understand not only the sources of the sGRBs but also the mechanism of the jet and the emission (Nakar 2007;Berger 2014;Nakar 2020).Unfortunately, after the joint detection of GRB170817A and GW170817, there has been no other confirmed joint sGRB and GW detection, which hinders advancing our understanding of the sGRB sources and the jet mechanism.
Joint GW and sGRB detections depend on the sensitivities and capabilities of both the GW detectors and the satellites.Currently, the fourth observation (O4) run of the LIGO-Virgo-KAGRA collaboration is ongoing (LIGO Scientific Collaboration et al. 2015;Acernese et al. 2015;Aso et al. 2013).In the future, LIGO detectors at Hanford and Livingston will be upgraded to A+ sensitivity (Barsotti et al. 2018), significantly improving their performance.Moreover, the planned GW detector LIGO-India (Iyer et al. 2011;Saleem et al. 2022) will also join the detector network, further improving the prospects of GW detection.
As the sensitivity and range of GW detectors increase significantly in the future, it is likely that the bottleneck for joint detections will soon be the sensitivity of EM detectors like Fermi and Swift.Thus, more sensitive high-energy satellites are required to improve the prospects of joint detections.
One such satellite proposed to study explosive astrophysical sources like GRBs and other high-energy electromagnetic counterparts to GW sources is Daksha (Bhalerao et al. 2022a,b).The Daksha mission will have two high-energy space telescopes with three types of detectors, covering the broad energy range from 1 keV to ∼ 1 MeV.The high sensitivity of Daksha and its near-uniform all-sky sensitivity arise from the Cadmium Zinc Telluride (CZT) detectors covering the medium energy range from 20 − 200 keV, with a median effective area of ∼ 1310 cm2 .Thus, Daksha has an effective area significantly higher than Fermi-GBM and will achieve Swift-BAT like sensitivity over the entire sky.With this high sensitivity, Daksha will clearly be a formidable instrument in the search for and study of high-energy transients.
In this paper, we study the prospects of joint sGRB and GW detection by the existing missions Fermi and Swift and the proposed Daksha mission.We consider multiple sGRB emission models to probe the various proposed models for the prompt emission.There have been discussions about sGRB-like emissions from neutron starblack hole (NSBH) mergers as well, however the theoretical details are still murky, and there is no observational confirmation.Hence, we do not consider NSBH mergers in this study.The plan of the paper is as follows: In §2.1, we provide details of the high-energy satellites and the synthetic BNS population employed here.In §2.2, we describe the various prompt emission models considered in this work, followed by the method used to calculate the flux for the satellites in §2.3.In §3, we compare the properties of different models and calculate the sGRB detection rates for the three satellites.§4 presents the rates and prospects of joint sGRB and GW detection.Finally, in §5, we discuss our main results, compare them with previous works, and present our conclusion.

METHODS
To study the prospects of the joint sGRB and GW detection from BNS mergers, we inject the sources in comoving volume.We then calculate the detection probabilities of the GW and GRB detectors, both separately and jointly.To probe the sensitivity of joint detection rates on the specifics of the GRB jet, we consider a broad range of sGRB jet models.Note that most of the candidate sources for the joint detection will be observed off-axis.Hence, in this work, we limit ourselves only to a set of models that make quantitative predictions for off-axis emission.

The Sources and the Detectors
We inject BNS merger sources up to a luminosity distance of   = 1.6 GPc (comoving distance of ≈1.2 Gpc and  ∼ 0.31 ).We distributed the injected sources uniformly in the source-frame Figure 1.The distribution of the injected BNS merger events over viewing angle, , and luminosity distance   (blue points).The histograms show the corresponding distributions (blue histogram).In the scatter plot, we also show the distribution of the events which are detected over threshold SNR of 8 either by IGWN4 or IGWN5 network (grey points).In the histograms, we show the distribution of these events individually for IGWN4 (dashed grey histogram) and IGWN5 (solid grey histogram).
comoving volume, which translates to a redshift probability distribution of The inclinations of the sources are isotropically distributed.The masses of the neutron stars in the sources are drawn from a normal distribution with a mean of 1.33 M ⊙ and a standard deviation of 0.09 M ⊙ (Özel & Freire 2016).The spins of the neutron stars are drawn uniformly between 0 <  < 0.05, with the upper limit set by the maximum known NS spin (Zhu et al. 2018).The GW waveforms are generated using IMRPhenomPv2 (Schmidt et al. 2012;Hannam et al. 2014;Khan et al. 2016).Note that the total number of sources considered are  inj = 0.3 million, which is over-sampled by a factor of ∼ 130 from the rate of the BNS mergers within the same comoving volume (the median rate ∼ 2316, Abbott et al. 2021c,a).Such over-sampling is done to avoid any effect of small number statistics on our results.Figure 1 shows the distribution of the sources with the observation angle () and luminosity distance (  ).
To determine the joint detection probability of GW and sGRB, we consider the detection prospects of GW sources by two different GW detector networks, comprising of (1) International Gravitational Wave Network-4 or IGWN4 with four detectors: Advanced LIGO Hanford and Livingston (LIGO Scientific Collaboration et al. 2015), Advanced Virgo (Acernese et al. 2015), and design sensitivity Kagra (Aso et al. 2013) and (2) International Gravitational Wave Network-5 or IGWN5 with five detectors: A+ sensitivities of LIGO Hanford and Livingston (Barsotti et al. 2018), A+ sensitivity LIGO-India (Iyer et al. 2011), Advanced Virgo and upgraded Kagra.We use the publicly available 2 power spectral densities (PSDs): aligo_O4high.text,avirgo_O4high_NEW.text, and Table 1.The energy bands, flux thresholds, and effective detection fraction of the high-energy missions. sky denotes the fraction of sky which is not earthocculted for the satellite, and  not_SAA denotes the fraction of the satellite's orbit inside SAA. df , thus, denotes the fraction of merger events that are not missed by the satellite for either being earth-occulted or for SAA-outage.

Mission
Energy band Detection Threshold Sky-coverage SAA-outage fraction Total detection fraction (keV) 10 −8 erg s −1  sky 1 −  not_SAA  df =  sky ×  not_SAA Daksha (two) (a)  20-200 4 1 0.134 0.866 Fermi-GBM (b)  50-300 3 0.7 0.15 0.595 Swift-BAT (c)  15 We calculate the optimal SNR of the sources by simulating waveforms without noise and using PSD models (see §4.1 for a discussion on the possible effect of noise).As the criteria for the detection of GW, we consider the events above a particular threshold value of the network signal-to-noise (SNR net = √︃ det SNR 2 det , where SNR det is the signal-to-noise ratio (SNR) of a single detector) (Pai et al. 2001).In this study, we consider two values of SNR net as ∼ 8 and 6.5 (following Petrov et al. 2022), and we define these SNR net criteria as threshold and sub-threshold detections, respectively.
To determine the prospects of sGRB detection, we consider two already existing high-energy GRB detectors: Fermi-GBM (50 − 300 keV) and Swift-BAT (15 − 150 keV), and the Medium-Energy (ME) detectors of the proposed mission Daksha (20−200 keV).Note that these satellites are low earth-orbit satellites.These satellites suffer from limited sky coverage since the Earth blocks about 30% of the sky for satellites in low earth orbit.For instance, Fermi-GBM can cover roughly about a fraction of 0.7 of the sky.In the case of Daksha, the two satellites on opposite sides of the earth overcome this limitation, making the fraction 1.In addition, all satellites in low earth-orbit are inactive during passage through the South Atlantic Anomaly (SAA).
We account for the effects of limited sky coverage and SAA by assigning an effective detection fraction (  df ) composed of the SAA out-time and sky-coverage fraction.To determine whether a sGRB is detectable by a satellite, we consider the thresholds in terms of the average flux received from the source by each of the detectors.The flux thresholds and the  df values considered for the three satellites are listed in Table 1.

The emission models
The sGRB emission is produced from a relativistic jet launched after the central remnant of the BNS merger collapses into a BH or hypermassive NS.If the emission arises only from the relativistic jet inside the jet opening angle ( 0 ), the emission declines very rapidly beyond the jet axis.However, the jet can have angular structure either from the intrinsic structure of the jet or the interaction of the jet with the surrounding ejecta, forming a cocoon surrounding the jet (Nakar 2020).The angular structure of the jet depends on several factors, such as the jet opening angle, the jet launching time after the merger, the initial energy of the jet, which determines whether the jet can successfully break out of the ejecta.The angular structure of the jet makes the observations beyond the jet core possible.
In this work, we adopt only jets with angular structures beyond the jet cores since we consider the sGRB associated with the GW, which, in most cases, will be observed off-axis.For example, we consider the jets with a Gaussian or power law distribution of the total energy ( or luminosity, ), the energy at the spectral peak ( peak ), and/or the Lorentz factor (Γ) with the angle from the jet axis (hereafter structured jets) and the jets with cocoon structure beyond the jet core (hereafter jet-cocoon).Note that all these models considered here fall under the broad category of structured jets (Nakar 2020).However, in this work, we define the models separately for convenience.

Structured Jets
We adopt two structured jet models.In the first model, the beamingcorrected isotropic energy in the bolometric band,  bol iso , varies as follows:  bol iso is assumed to be constant within the jet core ( bol iso () =  0 at  <  0 ), beyond which it varies as a power law as  bol iso () =  0 (/ 0 ) −  .Furthermore, we also assume that the Lorentz factor distribution has a similar angular profile with Γ 0 at the core.We call this model the "Power-Law Jet" (henceforth PLJ, note that we name the models after the type of angular dependence of the energy for ease of reference).We follow the treatment of Beniamini et al. (2019); Beniamini & Nakar (2019) for the model parameters:  0 = 0.1 rad ≈ 5.7°and  = 4.5.The log 10 ( peak,c ) distribution along the core is drawn from a normal distribution with a median value of 2.7 and a standard deviation of 0.19 (Nava et al. 2011).We then vary the value of  peak with  such that  peak () =  peak,c (Γ()/Γ 0 ), i.e. considering a situation in which the value of the peak energy remains constant in the co-moving frame.We note that when the observed emission is dominated by  ≪  obs , then due to relativistic Doppler beaming, we have  peak,obs =  peak () [( obs − )Γ()] −2 .The event duration was drawn randomly from a distribution of cosmological sGRB (Beniamini et al. 2020b).All the above jet parameters are defined in the observer frame.
We consider another structured jet model with a Gaussian variation of the radiative energy with angle:   () =  0  −  2 /2 2 0 , following the works of Ioka & Nakamura (2019).Unlike the previous model, the parameters here are defined in the comoving frame.It further assumes a power variation law of the Lorentz factor: Γ() = Γ max /(1 + (/ 0 )  ).Following the treatment of Ioka & Nakamura (2019), we use  0 = 0.059 rad ≈ 3.4°.Appropriate relativistic beaming corrections are then applied to obtain the beamingcorrected isotropic emission profile in the observer frame.For the details of the calculation, refer to §3 of their paper.Based on our convention, we call this model the "Gaussian."

Jet-Cocoon model
We adopt the Jet-Cocoon estimates of (Beniamini et al. 2019).The model consists of a top hat jet (with constant bolometric luminosity  0 within the jet core angle  0 ) surrounded by a quassi-spherical cocoon emission with bolometric isotropic equivalent luminosity  ,co .Hence the cocoon model will follow  bol iso () =  0 +  ,co for  obs <  0 ; and  iso () =  ,co for  obs >  0 , where  0 = 0.1 rad ≈ 5.7°.The value of the jet core's  peak is taken to be the same as described above for the PLJ model.The value corresponding to the cocoon is taken as 100 keV.The jet engine durations were taken as mentioned before for the structured jets.For more details, see Beniamini et al. (2019).We call this model the "Jet-Cocoon" model.The interaction between the jet and the ejecta may result in some low-energy jets failing to break out of the ejecta.We refer to such events as "failed jet" events.This results in very low luminosity values, rendering them mostly undetectable.Under fiducial simulation conditions, about 1/3 rd of the jets fail in the estimates of Beniamini et al. (2019).For the details of the jet-ejecta interaction parameters, refer to the paper.We denote the fraction of events having successful jets as  jet , which, in this work, is 2/3.

Flux Calculation
For each of the injected events, we calculate the flux received by the different satellites (e.g., the existing satellites Swift and Fermi and the upcoming satellite Daksha) within their respective energy bands.The flux received is calculated as  =  sat iso /4 2  , where  sat iso is the isotropic luminosity within the energy band of the satellite, and   is the luminosity distance to the source.For all the models except Gaussian, we calculate the energy in the respective energy bands of the satellites ( sat iso ) as: where  1 and  2 are the energy band limits and  is the source redshift.To calculate  (), we assume a comptonized spectral model with a cutoff power law: where   is the cut-off energy defined as:   =  peak /(2 + ) and  is the power-law index.We use  = −0.6,following Abbott et al. (2017) for the analysis of GRB170817A, and the value of  peak at different angles are given by the models (see 2.2).The Gaussian model, instead, uses a Band-function (Band et al. 1993) like spectral energy distribution in the comoving frame: where  0 () = 0.15 keV(1 + (/ 0 ) 0.75 ) with  0 = 0.059 rad.The values of  and  are taken to be 1 and 2.5, respectively (Kaneko et al. 2006).The spectral distribution is used in conjunction with the angular distribution of energy and Lorentz factor (see §2.2.1) to determine the beaming corrected emission energy within a desired energy range, which is given by Equation 14in Ioka & Nakamura (2019).We use this equation to calculate the energy received by the different satellites within their respective energy bands.
Note that the uncertainty in the assumed luminosity function may affect the predicted rates.To study the sensitivity of our results to the choice of the luminosity function, we also use a different luminosity function taking the Gaussian jet as the representative model.For this purpose, we use the isotropic energy distribution of the observed Swift GRBs as reported in Fong et al. (2015, hereafter F15).The distribution takes the form of a Gaussian in the log space with center at log( 0 ) = 51.31erg s −1 and a standard deviation of 0.98.To calculate the corresponding luminosity, we assume a jet activity time of 0.3 seconds to calculate the luminosity for all the events, which corresponds to typical sGRB duration (Kouveliotou et al. 1993).Hereafter, we call this model as the Gaussian-F model.The distribution is biased towards higher energy values since this is designed from the isotropic energy distribution of the observed sGRBs, unlike the WP15 function.
After calculating the flux, we calculate the detection rates of prompt emission by various satellites and also the joint detection rates based on various criteria (see following sections).In the jetcocoon model we have to explicitly account for unsuccessful jets as the cocoon emission might still be detectable in those cases.Therefore we multiply the calculated detection rates with the structured jet models by  jet = 2/3, i.e. the number of successful jet break-outs is assumed to be the same in the different models.This provides a fair comparison between Jet-Cocoon and Structured Jet models.

The model dependence of sGRB luminosity
To calculate the EM detection rates by various satellites, we need to calculate the source fluxes in the appropriate energy bands.We begin by calculating the bolometric luminosities of the injected events using the models described above (Figure 2, left panel).Then we use the source spectrum and source redshift to calculate the flux for each satellite.As an illustration, the right panel of Figure 2 shows fluxes of all simulated sources in the Fermi GBM (50-300 kev) band.
The luminosities of the structured jets (PLJ and Gaussian) are high at smaller angles, reaching the value of ∼ 10 50 erg, which declines with the observer angles.The luminosities inside the jet opening angle are similar for both the models, owing to their intrinsically similar energy structure inside the jet-core.However, beyond the jetcore, the PLJ is brighter than Gaussian Jet as a consequence of the assumed angular energy distributions for the two models.
For the Jet-Cocoon model, the total gamma-ray energy of the jet and the cocoon differ by ∼ 4 orders of magnitudes (Figure 2).This reflects the energy difference between the jet and the cocoon.Moreover, the cocoon shows a nearly constant luminosity over viewing angles since it has nearly isotropic energy distribution.Note that although we consider both the successful and failed jets for this model, Table 2.The sGRB and sGRB and GW joint detection rates for all GW network, mission, and emission model combinations.The above-threshold GW detection rates for BNS events are 15 yr −1 and 148 yr −1 for IGWN4 and IGWN5 respectively, with corresponding sub-threshold rates of 26 yr −1 and 264 yr −1 .The three "sGRB detection rate" columns give the EM rates for the three missions, independent of GW detection.The last six columns give joint detection rates.For instance, the last cell of the table states that for the Gaussian-F model with the F15 luminosity function, Daksha + IGWN5 will detect 8.27 above-threshold events and another 13.22 sub-threshold events per year, and 94% of the above-threshold events will be off-axis events.

Models
Luminosity  the results shown in Figure 2 include only the successful jets.This is due to the fact that if the jet is failed and the cocoon is the source of the emission, the luminosities across all angles are extremely low (⪅ 10 45 erg s −1 ), making them undetectable by the different facilities (e.g., Fermi) for typical source distances.
The median on-axis luminosity of the jet-cocoon model is higher than that of the structured jet models (PLJ and Gaussian, left panel of Figure 2).This is a selection effect.As discussed before, weak jets fail to break out of the ejecta, creating very faint emissions.We discard these as mostly non-detectable.Hence, the models shown in this plot are only ones with successful jets, which in turn means that their on-axis luminosity is higher.For the off-axis emission, PLJ shows the weakest angular dependence among the different models (left panel of Figure 2).
The Gaussian-F model has higher bolometric luminosity than the Gaussian model at smaller viewing angles since the luminosity function is brighter (left panel of Figure 2).Additionally, the events span a larger range of luminosity since the luminosity function is wider.
Note that the Gaussian-F model is consistent with the GRB170817A luminosity since the Gaussian model (Ioka & Nakamura 2019) is designed to explain this particular sGRB, and the luminosity function is drawn from the observation of the sGRBs.
The luminosities in Fermi energy band (right panel in Figure 2) show similar behavior as the bolometric luminosities for most models, except for the fact that the luminosities are scaled down by a factor of ∼ 2 − 5.This is expected since we consider the luminosity in a smaller wavelength range.However, for the PLJ model, the luminosities at larger angles are significantly different from the bolometric luminosities, showing a sharp decline with angle.This is because of the assumed  peak () distribution for this model, which drops sharply with the viewing angle.

EM counterpart detection rates
In this section, we compare the sGRB detection rates for different satellites (existing Fermi and Swift, and upcoming Daksha).We con- sider an event to be detected if the flux from the sGRB meets the detection threshold of the satellite concerned (see Table 1 for the threshold flux values).We calculate the rates considering the median BNS merger event rate of  BNS = 320 yr −1 GPc −3 (Abbott et al. 2021b).We perform a linear scaling from the number of injected events to the number of expected events within the redshifted volume (derived from Equation 1, see Chen et al. 2021 for a detailed discussion on volumes in this context) of the injections.To account for the limited sky coverage of the satellites and SAA outages, we perform a second scaling with a net detection fraction (  df ) for a given satellite (see Table 1).Thus, the final relation for rate calculation becomes: where  det is the number of detections among the injected events,   is the comoving volume, and  lim is the redshift value corresponding to 1.6 Gpc.
The detection rates of the sGRB counterparts (for BNS merger events within   = 1.6 Gpc, the limit of this work) for Fermi, Swift and Daksha lie within the ranges of 3.5 − 13.6, 0.6 − 3.5, and 5.5 − 29.1 yr −1 , respectively, for the different models.In all the cases, the rates for Daksha are factors of ≈ 2 − 4 and ≈ 9 higher than Fermi and Swift, respectively.
The detection rates drop sharply at off-axis viewing angles (right panel of Figure 4).For the structured jet models, the cutoff occurs in the viewing angle range of ∼ 20 − 25°, making larger viewing angle detections rare.Similar trends have also been reported in several other studies (Howell et al. 2019;Saleem 2020;Mohan et al. 2022).For the Jet-Cocoon model, the cutoff occurs at  ∼  0 .The cocoon component, however, allows fairly large off-axis sGRB detections.Such detections, however, are only possible for nearby events (  ≲ 100 MPc).
Figure 3 summarizes the ability of the three satellites to detect the sGRB counterparts of the injected BNS events.Swift has the lowest detection rate among the three satellites due to its low sky coverage (∼ 10%).Between Fermi and Daksha, although the total detection rates are not much different, nevertheless Daksha has better capability for detection of the fainter (off-axis or distant) events due to its higher sensitivity and the full sky-coverage.Thus, Daksha will be able to increase the detection rates of the sGRBs, specially at the fainter ends.
We note that our analysis is limited by the luminosity distance cutoff of 1.6 Gpc ( ∼ 0.3).Hence, the events beyond this distance are missed by our analysis.However, in reality, many bright events with  > 0.3 will be detected by the satellites.For instance, several sGRBs at much higher redshifts are regularly observed by Swift and Fermi (Fong et al. 2015;Poolakkil et al. 2021).Additionally, we do not consider potential sGRBs from NSBH mergers.These make the sGRB rates calculated in this work lower than the true sGRB detection rates by the satellites.This, however, is not a concern for the joint GW and sGRB detection rates, discussed in the next section, as GW detection is limited to much smaller distances with the current as well as near future detector networks.

JOINT GW AND SGRB DETECTION
We now evaluate the rates of joint high-energy and gravitational wave detections of these events using the same criteria as defined above.As discussed before, we consider three missions: Fermi, Swift, and the upcoming twin-satellite Daksha.For both detector networks that we consider (IGWN4 and IGWN5), we follow Petrov et al. (2022) to define an "above threshold" detection criterion to be SNR net = 8.On the other hand, if there is a confident EM detection of any event, the GW threshold can be lowered without compromising the false alarm rate, hence giving a significant joint detection.We refer to such events as "sub-threshold" events.While the computation of an exact network SNR cut-off value for sub-threshold events is beyond the scope of this work, we adopt a fiducial value of SNR net,min = 6.5 for sub-threshold events.Considering the binary neutron star merger rate  BNS = 320 GPc −3 yr −1 as discussed in §3.2, and the 30% downtime for GW detectors ( §2.1), our cutoffs lead to 15 (28) abovethreshold events per year and 148 (264) sub-threshold events per year for IGWN4 (IGWN5), with the corresponding median luminosity distance of ∼ 236 Mpc (∼ 508 Mpc).
For the sensitivity of IGWN4, the overall joint sGRB + GW detection rates are rather low.The total joint detection rates for the existing sGRB satellites Fermi and Swift lie in the range of 0.07 − 0.62 yr −1 and 0.02 − 0.14 yr −1 respectively, for the different models (Table 2, Figure 4).With an increased GW detection horizon with IGWN5, the detection rates increase to 0.8 − 4.0 yr −1 and 0.15 − 1.0 yr −1 for Fermi and Swift, respectively.For both IGWN4 and IGWN5, the fraction of GW events that are jointly detected with prompt emission,  GW , takes value in 0.5 − 4.0% and 0.1 − 1.0% for Fermi and Swift respectively.
As expected, the rate of on-axis / off-axis detections is modeldependent.Off-axis events are geometrically more probable, and indeed most of the joint detections are expected to be off-axis ( obs >  0 ) for the PLJ and Gaussian Jet models.For these models, fraction of jointly detected above-threshold off-axis events (  off ) takes values in the range of ∼ 60 − 95% for Fermi and ∼ 80 − 95% for Swift respectively (Table 2).However, the strong decline in off-axis luminosity for the Jet-Cocoon models makes these hard to detect, with  off ∼ 50% and ∼ 15% for IGWN4 and IGWN5, respectively.The sharp decline in  off from IGWN4 to IGWN5 for Jet-Cocoon model originates from the domination of on-axis sGRB detections from the further away GW-detected events.The higher sensitivity of Swift plays an important role here, as seen in Figure 3 -a significant fraction of these off-axis events fall below the Fermi detection threshold.
To consider a future perspective with a higher sensitivity mission, we examine the detection rates for Daksha, which has a higher volumetric sensitivity for GW170817-like events than other missions (Bhalerao et al. 2022b).As expected, the higher sensitivity and all-sky coverage yields much higher detection rates: 0.2 − 1.3 yr −1 for IGWN4, and 1.3 − 8.3 yr −1 for IGWN5: detecting  GW =∼ 1.0% − 8.4% of the GW events.These rates are factors of ∼ 2 − 9 times higher than Fermi or Swift.The comparable sensitivities of Daksha and Swift lead to similar  off values, but the increased sky coverage gives higher rates.The region spanned by the jointly detected events in the  obs −   space of the merger events increases significantly from IGWN4 to IGWN5 for all the satellites and models (Figure 5).However, the overall increment in the joint detection rates, especially for the fainter off-axis and far away events, is limited by the satellite capability.For example, Swift probes a larger region in  −   space, i.e., can detect fainter ends of the events due to its better sensitivity.Nevertheless, the total number of events detected by Swift is still not high due to its low sky coverage.Fermi, on the other hand, probes much smaller region in this space owing to the higher flux threshold (i.e.lower sensitivity).However, with Daksha, the detection rates will increase mainly within the viewing angle range of 5 − 20°.Such events will serve as a connecting link between cosmological on-axis sGRBs and GRB170817-like events, improving our understanding of the jet structure.
In addition to the satellites, the performance of GW detectors also limits joint detections.To quantify this, we consider the EM-detected sGRBs that are within the GW detection horizons, and calculate the fraction that is also detected in gravitational waves (  sGRB ).We define the GW detection horizon to be the distance at which the cumulative rate attains 99% of the total rate, which is ≈ 523 Mpc and ≈ 1127 Mpc for IGWN4 and IGWN5.For all the missions and emission models,  sGRB takes values between ∼ 40 − 50% and ∼ 55 − 65% for IGWN4 and IGWN5 networks, respectively.These values are lower than 100% due to both the non-zero downtime of the GW detectors and the limitation of detecting signals from the outer margins of the detection volume.The increment in the fraction from IGWN4 to IGWN5 results from improved GW detection efficiency.

Comparing our rates with other works in the literature
In this section, we conduct a comparative analysis of our results with those from previous studies in the literature.We first compare the GW detection rates, followed by the sGRB+GW rates.

GW detection rates
With  BNS = 320 Gpc −3 yr −1 and detection threshold of SNR det = 8, we obtain GW detection rates of 15 yr −1 (123 yr −1 ) with IGWN4 (IGWN5) network.These rates align well with those constrained in other studies employing different population models within the uncertainty limits.For instance, Colombo et al. (2022) obtain a rate of 7.7 yr −1 with a similar  BNS , a stricter threshold detection criterion of SNR net = 12, and 80% detector duty cycle.
One notable difference is that our rates are less than half the value of 34 yr −1 obtained by Petrov et al. (2022).Kiendrebeogo et al. (2023) obtain similar high rates: 36 yr −1 with  BNS ≈ 210 GPc −3 yr −1 and 17 yr −1 with  BNS ≈ 170 GPc −3 yr −1 .While we have slightly different assumptions about the underlying astrophysical merger rate densities, the discrepancy can mainly be attributed to the methodology employed in estimating detectable sources.Both these works emulate the detection methods currently in use in LIGO/Virgo/KAGRA analyses.They employ similar population models as us but inject waveforms into simulated Gaussian noise, with detections based on matched filter SNR.In this scenario, noise can decrease or increase the signal strength, leading to false non-detections or chance detections respectively.Thanks to a uniform distribution of sources in volume, there are a larger number of events just below the detection threshold than above it: hence, this effect increases the number of events detected in their analyses.While these effects accurately mirror the actual data processing in IGWN, they overestimate the number of detections.
On the other hand, we base our estimations on the injected SNR of the sources, computed by simulating waveforms without noise and using power spectral density models.To test if this completely accounts for the discrepancy, we used the Kiendrebeogo et al. ( 2023) data set to calculate the expected SNRs following our methods we found that the rates become consistent with our rates.

Joint sGRB+GW detection rates
Our sGRB+GW joint detection rates calculated in §4 agree well with those estimated in several past works.Recent work by Colombo et al. (2022), provides a joint rate of ∼ 0.17 yr −1 and ∼ 0.03 yr −1 with Fermi and Swift respectively, which are consistent with with IGWN4 rates of this work.A lower fraction of successful jets (  jet ∼ 52%) and stricter threshold detection criteria of SNR net = 12 in their work leads to relatively lower rates.Our rates are also broadly consistent with Howell et al. (2019) after appropriately scaling for the higher assumed  BNS and  jet = 100% used in their work.However, other works by Patricelli et al. (2022) and Saleem (2020) yield significantly higher rates.This may arise from the brighter luminosity function or emission models.Patricelli et al. (2022) included Gaussian noise of GW detectors, which, in light of the previous discussion about GW rates, can also be responsible for their relatively higher joint detection rate estimations.Our prediction of low values of  GW is consistent with most past works (Howell et al. 2019;Colombo et al. 2022;Patricelli et al. 2022).If we use GW detections following the Petrov et al. (2022) or Kiendrebeogo et al. (2023) method, we find that the GW+sGRB detection rates mentioned in Table 2 and elsewhere will approximately double.

DISCUSSION AND CONCLUSION
In this work, we examine the prospects of detecting sGRB counterparts and joint sGRB and GW detection of BNS merger events by three high-energy satellites: the existing missions Fermi and Swift and the upcoming all-sky high sensitivity mission Daksha; for two GW detector networks: current sensitivity 4-detector (HLVK, IGWN4 case) and future A+ sensitivity 5-detector (HLVKI, IGWN5 case) network.Many of the jointly detected events are observed beyond the jet core (off-axis).Hence, we use several jet models predicting off-axis emission to calculate the joint detection rates (power law jet and jet-cocoon models by Beniamini et al. (2019) and Gaussian jet model by Ioka & Nakamura (2019), see §2.2).In addition, to test the effect of the assumed luminosity function on the predicted rates, we use two different luminosity functions, given by Wanderman & Piran (2015) and Fong et al. (2015), for one of the representative models (the Gaussian model, see §2.3).
We use the volumetric BNS merger rate  BNS = 320 GPc −3 yr −1 (Abbott et al. 2021c,a) for calculating joint detector rates.Our predicted joint detection rates for the existing satellites Fermi and Swift with IGWN4 (IGWN5) are 0.07-0.62yr −1 (0.8-4.0 yr −1 ) and 0.02-0.14yr −1 (0.15-1.0 yr −1 ), respectively.These rates increase by a factor of ∼ 2 − 9 for the proposed Daksha mission.The predicted joint detection rates with Daksha with IGWN4 (IGWN5) network and GW detection over threshold SNR lie in the range of 0.2 − 1.3 yr −1 (1.3 − 8.3 yr −1 ).This highlights the need for a more sensitive future mission: The IGWN5 rates imply that Daksha will lead to at least one joint sGRB and GW detection from BNS merger event per year, which is a notable improvement from no event in 6 years (since GRB170817A).Figure 6 summarizes the results of this work (also see Table 2).We note here that our rate estimations with optimal SNRs of the injected BNS merger events are conservative.Comparison with the works of Petrov et al. (2022) and Kiendrebeogo et al. (2023) suggests that Gaussian noise in GW detectors can assist the detection of otherwise GW-faint events, which can potentially increase both the GW and sGRB+GW detection rates by a factor of two (see §4.1).
We show that with future GW detection networks (like IGWN5), joint detection rates at lower viewing angles ( obs ⪅ 2 0 ,  0 being the jet opening angle) increase significantly from IGWN4 to IGWN5.Larger angle detections are, however, limited by satellite performance.Nevertheless, within the detectable viewing angle range, Daksha is expected to perform better than the existing satellites.It will probe the largest region in the space of viewing angle and distance of the merger events with the highest efficiency among the satellites (see Figure 5).The predicted range of off-axis joint detection rate with Daksha and IGWN4 (IGWN5) of 0.09 − 1.2 yr −1 (0.2−7.77 yr −1 ), is higher than that of both Swift and Fermi by factors of 2 − 6.This underscores the need for future missions that combine high sensitivity with all-sky coverage for detecting and characterizing these events, and in turn understanding the physics of compact object mergers and the post-merger radiative processes.

Figure 2 .
Figure 2. The luminosity ranges for all the prompt emission models considered in this work within the bolometric (Left) and Fermi-GBM (Right) energy bands.The line represents the median luminosity, whereas the band represents the interquartile spread.GRB170817A is marked with a red cross.The jet opening angles of the Gaussian model (0.059 rad ≈ 3.4°) and the PLJ and Jet Cocoon models (0.1 rad ≈ 5.7°) are also marked (vertical lines grey and black lines, respectively).

Figure 3 .
Figure 3. Left: The histogram of event rate per flux bin scaled by  df .In other words, the figure shows the distribution of all the events which are not missed by the satellites due to limited sky coverage or SAA outage.Vertical lines mark the detection thresholds of the satellites.Right: The angular distribution of the rate of sGRB detection.Vertical lines mark the jet opening angles ( 0 ).For all the models, Daksha shows the highest detection rates.

Figure 4 .
Figure 4. Cumulative joint detection rate over distance for IGWN4 (solid lines) and IGWN5 (dashed lines) networks.The grey lines show the GW detection rate, and the colored lines show joint sGRB and GW detection rates by the satellites (blue: Daksha, green: Swift and orange: Fermi).The rates show a significant increase from IGWN4 to IGWN5.Daksha performs 2 − 9 times better than the existing missions.

Figure 5 .
Figure5.The distribution of the join sGRB and GW detected events in the viewing angle ( obs ) and luminosity distance (  ) space for IGWN4 (Left) and IGWN5 (Right) cases.Horizontal lines mark the jet opening angles ( 0 ).We span a larger region in  −   space in IGWN5.With Daksha, we span the largest space in  −   with the highest efficiency and detection rates.

Figure 6 .
Figure6.The sGRB and joint sGRB and GW rates for the three missions with different prompt emission models.The rates are computed for  BNS = 320 GPc −3 yr −1 .The sGRB rate is only for events with   < 1.6 Gpc, the distance limit of this work.