Central engine of GRB170817A: Neutron star versus Kerr black hole based on multimessenger calorimetry and event timing

The central engine of GRB170817A post-merger to GW170817 is probed by GW-calorimetry and event timing, applied to a post-merger descending chirp which can potentially break the degeneracy in spin-down of a neutron star or black hole remnant by the relatively large energy reservoir in the angular momentum, $E_J$, of the latter according to the Kerr metric. This analysis derives from model-agnostic spectrograms with equal sensitivity to ascending and descending chirps generated by time-symmetric butterfly matched filtering. The sensitivity was calibrated by response curves generated by software injection experiments. The statistical significance for candidate emission from the central engine of GRB170817A is expressed by probabilities of false alarm (PFA; type I errors) derived from an event-timing analysis. PDFs were derived for start-time $t_s$, identified via high-resolution image analyses of the available spectrograms. For merged (H1,L1)-spectrograms of the LIGO detectors, a PFA $p_1$ derives from causality in $t_s$ given GW170817-GRB17081A. A statistically independent confirmation is presented in individual H1 and L1 analyses, in a second PFA $p_2$ of consistency in their respective observations of $t_s$. A combined PFA derives from their product since mean and (respectively) difference in timing are statistically independent. Applied to GW170817-GRB170817A, PFAs of event timing in $t_s$ produce $p_1=8.3\times 10^{-4}$ and $p_2=4.9\times 10^{-5}$ of a post-merger output ${\cal E}_{GW}\simeq 3.5\%M_\odot c^2$ ($p_1p_2=4.1\times 10^{-8}$, equivalent $Z$-score 5.48). ${\cal E}_{GW}$ exceeds $E_J$ of the hyper-massive neutron star in the immediate aftermath of GW170817, yet it is consistent with $E_J$ rejuvenated in delayed gravitational collapse to a Kerr black hole. Similar emission may be expected from energetic core-collapse supernovae producing black holes. (Abbr.)


Introduction
Gamma-ray bursts represent the most relativistic transient events in the sky. Discovered serendipitously (Klebesadel et al. 1973), they are now known to have an astronomical origin in the end point of massive stars, broadly classified as short duration events (SGRBs), some with extended emission (SGRBEEs), and longduration events (LGRBs, reviewed in e.g., van Putten et al. 2014b). The end of life of relatively massive stars is believed to be the origin of neutron stars and stellar mass black holes by gravitational collapse when their nuclear burning phase ceases. Neutron stars can be found in supernova remnants harboring pulsars, as single objects such as Vela and the Crab (Hewish 1970), but also in binaries, notably the Hulse-Taylor system PSR1913+16 (Hulse & Taylor 1975) and the short-period binary PSR J0737-3039 (Burgay et al. 2003).
While pulsar emission gives unambiguous evidence of spinning neutron stars, identifying stellar mass black holes (at similar levels of rigor) appears more challenging, despite their common astronomical origin in core-collapse of massive stars. Follow-ing the seminal discovery of cosmological redshift of gamma-ray bursts in optical follow-up observations (e.g., Bloom et al. 2001) to an X-ray afterglow to GRB970228 by BeppoSAX (Costa et al. 1997), LGRBs have been identified thanks to core-collapse supernovae -energetic events of type Ib/c from relatively more massive stars (Galama et al. 1998;Hjorth et al. 2003;Stanek et al. 2003;Matheson et al. 2003;Modjaz et al. 2006;Guetta & Della Valle 2007;Kelly et al. 2008). On the other hand, GRBs originating in compact binary coalescence (Paczynski 1986) increasingly appear to be of the short-duration variety, as notably demonstrated by GW170817-GRB170817A (Abbott et al. 2017); however, this may be extended to include long-duration events, depending on black hole spin (van Putten & Levinson 2003) or progenitor binary evolution (e.g., Rueda et al. 2021).
The astronomical origin of GRBs identified with the endpoint of stellar evolution leaves the central engine producing the ultra-relativistic outflows responsible for their non-thermal gamma-ray emission to be either a strongly magnetized neutron star (magnetar) or black hole (e.g. Piran 2004;Piran et al. 2019;Piro et al. 2019;Nakar 2020). In contrast to neutron stars, elec-Article number, page 1 of 17 arXiv:2212.03295v1 [astro-ph.HE] 19 Nov 2022 A&A proofs: manuscript no. 42974corrR1 tromagnetic observations alone appear to be insufficient to conclusively distinguish between the two, despite a half-century of multi-wavelength observations covering GRBs by their prompt gamma-ray emission and afterglows in X-rays and optical-radio.
The recent advance of the LIGO-Virgo-KAGRA (LVK) detectors offers unprecedented and independent opportunity to advance this conundrum. Quite generally, an astrophysical black hole central engine powering an extreme transient event will be exposed to surrounding high-density matter. Gravitational radiation thereby produced is expected to be dominant over MeVneutrinos and electromagnetic radiation when powered by it ample angular momentum energy reservoir in its angular momentum (van Putten & Levinson 2003) which, according to the Kerr metric, readily exceeds the same of a neutron star by an order of magnitude. True calorimetry hereby promises the break the degeneracy of neutron stars and black holes as the central engine, upon including total energy output in gravitational radiation (GW-calorimetry).
Observational opportunities may be found in compact binary mergers involving a neutron star (Abbott et al. 2017) and corecollapse supernovae, especially those associated with LGRBs (LSC 2018;van Putten et al. 2019b). At a distance D 40Mpc (Coultier et al. 2017;Cantiello et al. 2018), GW170817 presents a fortuitously nearby event offering a unique prospect to witness the formation of stellar mass black holes in gravitational collapse, directly or delayed, of an intermediate hypermassive neutron star formed in the immediate aftermath of this double neutron star merger (e.g., Murguia-Berthier et al. 2020) -in addition to the new opportunities for probing physical properties of the neutron star progenitors (Abbott et al. 2018;Drago & Pagliara 2018;de Pietri et al. 2019;Bauswein 2019).
In association with GRB170817A, the precursor emission GW170817 makes it feasible to achieve a key objective of multimessenger observations (Acernese et al. 2007). An earlier demonstration of significance in timing across different observational channels is found in SN1987A, with essentially coincident arrival times of an MeV-neutrino burst detected independently by Kamiokande II and IMB, and alongside an optical light curve of SN1987A in the satellite galaxy LMC. For GW170817, independent pickup of the merger signal GW170817 by H1 and L1 is likewise coincident (within the light travel time of 10 ms between H1 and L1), followed by GRB170817A in NGC4993. At the distance of 40 Mpc, however, any further MeV-neutrino emission is a priori undetectable. Quite generally, timing coincidences carry observational significance based on a uniform prior on event timing on scales less than a Hubble time, subject only to instrumental constraints.
Gravitational wave observations of compact binary mergers is realized at observational sensitivity close to the limit defined by strain-noise amplitude of the detectors (detector limit; Abbott et al. 2021). This limit is defined by the theory of ideal matched filtering in the face of finite detector noise exploiting phase-coherence in correlation with model gravitational-wave templates over a large number of wave periods (e.g., Abbott et al. 2020). Ideally, un-modeled searches are pursued at similar sensitivities, even when phase-coherence is limited to intermediate time-scales (between the period of the wave and the total duration of the event). For transient emission featuring ascending or descending frequencies in gravitational-wave emission, broadband spectrograms provide a suitable starting point to the search for transient emission by image analysis. High sensitivity spectrograms can be realized by butterfly filtering over a dense bank of time-symmetric templates of intermediate duration τ with sensitivity exceeding that of Fourier-based spectro- Fig. 1. Schematic of event timing in the first double neutron star merger detected (top): the multimessenger event GW170817 at merger time t m with accompanying GRB170817A at t GRB = t m + t g following a gap t g 1.7 s (bottom). The central engine of GRB170817A formed in this gap by causality, whereby t m < t s < t GRB for the start-time of accompanying gravitational radiation. This gap condition poses a discrete event time problem for the probability p of false alarm (PFA) for a candidate detection to be in the 1.7 s gap between GW1709817-GRB170817A: p = t g /T marked by the (unique) global maximum of an indicator function over an observational time T given the uniform prior on astrophysical event times (Appendix A). When sufficiently large, the energy output E GW marked by t s can break the degeneracy between neutron stars and black holes. grams by over an order of magnitude (van Putten et al. 2014a). It and high-resolution mage analysis can be accelerated on modern graphic processor units (GPUs) with high bandwidth memory (Appendix D). This opens a window to probing the central engine of GRB170817A by GW-calorimetry at sensitivities on par with the pre-cursor GW170817.
The statistical significance of a candidate gravitational-wave event derives from the probability of false alarm (PFA), that is, Type I errors assuming the null-hypothesis of no gravitational waves passing by the LIGO detectors. Multimessenger observations involving the H1 and L1 detectors suggest exploiting both energy and event times t k = t k (t): a Boolean stochastic variable over {0, 1} as a function of real time t indicating an event k at some specific instant t s (e.g., Brown et al. 2004;Donges et al. 2016) (Fig. 1). For a candidate emission feature, PFAs may derive from consistency conditions applied to a probability density (PDF) of t s . Consistency in timing can be contextual or acontextual. Quite generally, this can be approached with and without priors from other observational channels, shown here in merged and individual H1 and L1 analyses, respectively. Event timing of a candidate post-merger emission feature to GW170817 associated with the central engine of GRB170817A is subject to the gap condition ( Fig. 1): given the interval G = t m , t m + t g as an observational prior on its start-time, t s , by causality; C 1 hereby introduces a Booleanvalued statistic, t k (t) = 1 (t G) , and 0 otherwise. Given the uniform prior on astrophysical event times, such a discrete event time carries a PFA equal to t g /T over an observation of duration T as a conventional probability. As elucidated in Appendix A, this is analogous to the odds of winning a bet in roulette provided that the measurement of t s is sufficiently precise. Combining PFAs in one-class classification schemes (e.g. Simonso et al. 2017) is common in image analysis, here spectrograms in §3-4. Quite generally, we differentiate between con-tinuous or discrete distributions (e.g. Block et al. 2006). The product of two p-values of continuous distributions, for instance, signal amplitude, is not a p-value, but can be merged by the celebrated Fisher's combined probability test (Fisher 1932(Fisher , 1948Simonso et al. 2017;Lyone & Wardle 2018). For a pair of statistically independent p-values p i (i = 1, 2), Fisher's method defines an equivalent p-value by equivalent information content I i = − log p i according to χ 2 F = −2 log P, where P = p 1 × p 2 and χ 2 F denotes the χ 2 distribution with F = 4 degrees of freedom. The merged p-value derives from the cumulative density function of χ 2 F . Fisher's method has been widely used and indeed led to a number of potential improvements (Whitlock 2005;Heard & Rubin-Delanchy 2017). When p-values p 1 and p 2 are small, p = P 1 − log P gives a merged p-value that is effectively equivalent to Fisher's (Theiler 2004;Lyone & Wardle 2018). This is not the case in the present discrete event timing analysis of C 1 , where the PFA derives as a conventional probability for (1) that can be merged as such with other PFAs (e.g., Theiler 2004).
In the next section, we zoom in on the complex merger sequence of the double neutron star merger GW170817 with the aim of identifying the central engine -a hyper-massive neutron star or rotating black hole -of the associated GRB170817A by GW-calorimetry and event timing, previously introduced in van Putten & Della Valle (2019); van Putten et al. (2019a). These works are hereafter referred to Paper I and Paper II, respectively.
To this end, we present two statistically independent analysis of H1 and L1 data by model-agnostic spectrograms. Observational results are compared and combined based on highresolution PDF(t s )s realized by modern exascale heterogeneous computing (Appendix D). In particular, a joint PFA from event timing, factored over PFAs derived from PDF(t s )s of each, provides a significantly enhanced level of confidence in the observations.
To start ( §3), we revisit and extend our analysis of merged (H1,L1)-spectrograms in Paper I and II with (a) a one-parameter family of response curves calibrated by signal injection experiments and (b) PDF(t s ) generated over an extended foreground derived from a large number of small time-slides between H1 and L1, satisfying (c) a condition of extremal clustering robust against false positives from data anomalies. Statistical significance is expressed by PFA 1 from event timing by causality in the context of GW170817-GRB170817A (C 1 ) (Appendix A).
Next ( §4), we present a statistically independent analysis of H1 and L1 individually with accompanying PDF(t s ). Statistical significance is expressed by PFA 2 of consistency in timing, t s , by cross-correlating their respective PDF(t s )s (C 2 ). In §5, we discuss our results and the observational consequences for the central engine of GRB170817A with joint PFA factored over PFA 1 and PFA 2 based on statistical independence of mean and difference in event time observations, providing a key advance over previous analysis. We summarize our conclusions ( §6) with an outlook ( §7) on upcoming opportunities to probe merger sequences involving neutron stars and core-collapse supernovae in the Local Universe ( §6).

Rejuvenation in gravitational collapse
The double neutron star (DNS) merger sequence of GW170817 is complex, involving a hyper-massive neutron star (HNS) formed in the immediate aftermath (Abbott et al. 2017;Dai 2019;de Pietri et al. 2019;Murguia-Berthier et al. 2020) with the possibility of (delayed) gravitational collapse to a rotating black hole (BH; Akutsu et al. 2020): with the ellipses referring to additional post-merger emission in MeV-neutrinos and gravitational radiation. In a continuing gravitational collapse -direct or delayed -a black hole will form surrounded by a high-density disk or torus from debris and dynamical mass-ejecta from the merger prior (e.g., Rosswog et al. 1999;Baiotti & Rezzolla 2017;Ciolfi 2020). In Eq. 2 the GW170817-GRB170817A association appears quite secure with independent p-values from continuous random variables: temporal (P temp = 5 × 10 −6 ) and spatial (P spat = 1 × 10 −2 ) consistency in the LIGO H1 and L1-detectors, on the one hand, and GRB170817A detected by Fermi Gamma-ray Burst Monitor (GBM) (Abbott et al. 2017).
In Eq. 2, the potential for breaking aforementioned degeneracy between a neutron star or black hole remnant derives, in part, from a rejuvenation of the energy reservoir in angular momentum, E J , in the process of gravitational collapse. By conservation of specific angular momentum a/M = sin λ in prompt gravitational collapse of the HNS with spin-energy E − J (1/2)JΩ, E + J = 2Mc 2 sin 2 (λ/4) of the black hole (Kerr 1963) produces: Here, the right hand size refers to GW170817, using J = IΩ for the moment of inertia I (2/5)MR 2 of the HNS with gravitational radius R g = GM/c 2 , radius R and angular velocity Ω with Newton's constant G and the velocity of light c. Further enhancement in a/M might derive from hyper-accretion (e.g., Bardeen 1970;Levinson & Globus 2013), although this does not appear to be called for in GW170817. Therefore, in the case of the gravitational collapse of a ∼ 3M hyper-massive neutron star, E J increases by Eq. 3, setting the initial condition for potentially significant post-merger emission from a black hole-torus system.

Calorimetry in gravitational radiation
Indeed, rejuvenation Eq. 3 can safely account for Extended Emission feature in (H1,L1)-spectrogram ( Fig. 3) with (Paper II): Notably, (4) breaks the degeneracy between neutron stars and black holes, since E gw 6.3 × 10 52 erg exceeds the limit of rotational energy of the hypermassive neutron star progenitor, even more so after the initial one second of neutrino cooling and internal dissipation leading to a uniform rotator -a supermassive neutron star (Beniamini & Lu 2021). Further by the observed gravitational-wave frequencies f gw 700 Hz (Fig. 3) would imply a spin frequency well-below the Keplerian frequency limits (Haensel et al. 2009) for quadrupole emission at twice its spinfrequency, limiting its rotational energy to: Thus, Eq. 4 exceeds Eq. 5 by a factor greater than 4, defined by maximal values of mass M and radius R. Yet, Eq. 4 is in quantitative agreement with the output of black hole spin-down against a non-axisymmetric high-density disk or torus (Paper II). Accordingly, GW-calorimetry hints at continuing gravitational collapse to a black hole in Eq. 2. In contrast, calorimetry in the (post-merger) electromagnetic observation of GRB170817A and kilonova AT 2017gfo (Connaughton 2017;Savchenko et al. 2017;Mooley et al. 2018a,b;Ascenzi et al. 2020) shows a combined energy output limited to This output falls woefully short to break the degeneracy between a neutron star of black hole remnant in the chronicle Eq. 2, prompting the need for GW-calorimetry outlined above. Even so, electromagnetic observations do provide us with crucial timing information. Post-merger rejuvenation in Eq. 2 appears delayed by constraints on the lifetime of the hypermassive neutron star (Pooley et al. 2018;Gill et al. 2019;Radice et al. 2018;Lucca & Sagunski 2019), inferred from a "blue" component of the associated kilonova AT 2017gfo (Smartt et al. 2017;Pian et al. 2017 (Granot et al. 2017;Gottlieb et al. 2018;Nakar et al. 2018;Metzger et al. 2018;Xie et al. 2018;Gill et al. 2019;Lazzati et al. 2020;Hamidani et al. 2020).

Sensitivity requirements
Performing GW calorimetry Eq. 4 to post-merger emission requires sensitivity on par with the CBC search for GW170817. Gravitational radiation from non-axisymmetric disks or tori swept up by a rotating black hole lack the phase-coherence over long duration, normally exploited in model-dependent searches applicable to compact binary coalescence (CBC).
Existing un-modeled searches by power-excess methods (Abbott et al. 2017(Abbott et al. , 2019a) fall short by a threshold for detection of E th,gw 6.5M c 2 (Sun & Melatos 2019). This represents a sensitivity to post-merger signals of 0.3% compared to E gw 2.5%M c 2 in GW170817 itself (Abbott et al. 2017), a priori in an excluded zone of the parameter space by the total progenitor mass-energy E DNS 3M c 2 of GW170817.
Here, we apply butterfly matched filtering which aims to capture signals with frequencies slowly wandering in time. Using a dense filter bank of time-symmetric templates, it realizes equal sensitivity to ascending and descending branches relevant to the present study of merger and post-merger emission to Eq. 2, squarely in the allowed zone E GW < E DNS of any post-merger radiation.
As a pass filter for frequencies with finite time rate-ofchange in frequency, constant frequency signals are relatively suppressed (van Putten 2016), while stochastic signals can be detected at sensitivity far beyond the sensitivity of conventional Fourier analysis (van Putten et al. 2014a). A precise characterization for the purpose of the present study follows and is based on detailed signal injection experiments.

Analysis of merged (H1,L1)-spectrograms
In this section, we consider (H1,L1)-spectrograms merged by frequency coincidences, generated by butterfly matched filtering (Paper I and II). Merged (H1,L1)-spectrograms were considered for the following steps (Appendix C): 1) the sensitivity is the same for ascending and descending chirps by time-symmetry of matched filtering templates (this step is un-modeled). 2) parameter estimation is performed by an image scan, quantifying features by a goodness of fit over a family of descending chirps using χ-image analysis (detailed in Papers I and II).
Extending previous work, we analyze merged spectrograms by application of surrogate time-slides of the following type: (S 0 ) foreground in conventional analysis of merged (H1,L1)spectrograms with zero time-slide, ∆t = 0; (S 1 ) extended foreground over small time-slides ∆t less than the intermediate duration τ of our templates; (S 2 ) background from time-slides ∆t greater than τ. We note that τ is considerably greater than the light travel time between H1 and L1.
We extract PDF(t s ) from extended foreground (S 1 ). This is a principle extension to previous work, restricted to foreground (S 0 ), that is, zero time-slide. A high-resolution PDF(t s )s produced by heterogeneous computing (Appendix D) permit accurate parameter estimation, rigorous application of event timing ( §3.2, Appendix A), robustness against false positives by extremal clustering ( §3.3), and enable a direct comparison and combination of results between statistically independent analyses ( §4). Our descending chirp follows a fit to the post-merger feature attributed to radiating E J of the compact remnant, distinct from an ascending chirp in binary coalescence, given by (Paper I and II):

Response curves to descending branches
These curves introduce four parameters: a start frequency, f s , late-time frequency, f 0 , in addition to aforementioned (t s , τ s ) used in the χ-image analysis of spectrograms. Stable estimates derive in goodness of fit over data-segments of ∆t = 7 s, moving over small time-steps of 33 ms (supplementary data, Paper I).
Our response curves are generated by injection experiments cover the two-parameter space of energies and τ s , using a source distance of D = 40 Mpc of GW170817 (Cantiello et al. 2018), f s = 680 Hz, and f 0 = 98 Hz. We identify Eq. 7 in merged (H1,L1)-spectrograms by χimage analysis parameterized by (t s , τ s , f s , f 0 ) at a combined resolution of 0.5 MHz (2 14 steps in (τ s , f s , f 0 ) at 30 Hz in t s ) in a high-frequency scan over 10 9 parameter values over the 2048 s frame of H1L1-data sampled at 4096Hz. Such is applied following a time-slides applied to H1 and L1 (prior to generating spectrograms). Figure 2 shows results produced by 56 injections over E gw in seven steps and τ s in eight steps. The response curves shown are averaged over time-slides covering the light-travel time of 10 ms between H1 and L1 using 41 "tiny" time-slides. (A movie of this injection analysis has been published online. 1 ) In this process, an indicator function χ = χ(t s , τ s , f s , f 0 ) represents a normalized (H1,L1)-correlation strength over strips along family of curves Eq. 7 (Paper I, supplementary data). This is reduced toχ(t s , τ s ) = max ( f s , f 0 ) χ(t s , τ s , f s , f 0 ), namely, upon projecting out ( f s , f 0 ) by maxima over ( f s , f 0 ) for each (t s , τ s ). These projections are necessary by memory limitations, when Eight response curves of the two-level pipeline by χ-image analysis of merged (H1,L1)-spectrograms produced by butterfly matched filtering using a total of 56 signal injection experiments on LIGO O2 data covering GW170817 (top-left panel). These model injections are defined by Eq. 7. The dashed line is the mean ofχ-peaks averaged over (tiny) time-sides ∆t within the light-travel time of 10 ms between H1 and L1 identified with extended emission to GW170817EE in a post-merger descending branch (lower panels). Precision realized is shown by parameter recovery (top-right panel). The parameter "field-of-view" covered for descending chirps covers energies 0 < E gw < 8%M c 2 (7 steps) and timescales of descent τ s (eight steps, in units of seconds) with start frequency f s = 680 Hz for a model source distance D = 40 Mpc given by GW170817. These injection signals are in the shot-noise dominated regime of the LIGO detectors, where noise increases with frequency and hence sensitivity decreases with τ s . Inserted is a plot of observed τ s of post-merger emission to GW170817 (black square), overlaid to the calibration curve (blue curve) to infer uncertainty (along abscissa) from observed scatter (along the ordinate). Sample of a single injection of a combined merger and post-merger signal (1810-1825 s) in the snippet of H1L1-data containing GW170817 (t m = 1842.43 s) (reprinted from Paper II), shown in the lower panel. The present 56 injections Eq. 7 vary over energy E GW and time-scale of descent τ s , keeping merger signal injection fixed at the parameters of GW170817.
scanning over a large number of 2 × 10 9 parameter values per frame of 4096 s at 4096 Hz sampling rate.
Response curves vary with τ s due to shot-noise in the H1 and L1 detectors noise above a few hundred Hz. Small τ s corresponds to fast descent in frequency, putting injections mostly over low frequencies close to f 0 , where detector noise is less and hence sensitivity is high. Conversely, large τ s puts injections close to f s , where shot-noise adversely affects sensitivity.

PFA from event timing based on PDF(t s )
Significance of a candidate event -satisfying extremal clustering ( §3.3, below) -is estimated by causality in formation of a central engine of GRB170917A (Paper I and II): the gap condition (1) t m < t s < t m + t g , t g = 1.7 s, between GW170817 and GRB170817A, schematically indicated in Fig. 1.
Including causality in butterfly matched filtering with template duration τ, we consider the slightly more restrictive gap condition, in the reduced gap G of width t g − τ. Statistical significance is expressed by the PFA: Type I errors in C 1 according to t s satisfying t k (t s ) = 1.
Formally, as per §1, t k (t) : [0, T ] → {0, 1} is a Boolean random variable over our snippet of H1L1-data of duration T = 2048 s. In the present case, one such candidate event (Fig. 4) marked by the (unique) global maximum at event time t = t s in an indicator function χ(t) over [0, T ], satisfying extremal clustering over extended foreground S 1 . Under the null-hypothesis of noise only and no signal, t 1 (t) satisfies a uniform prior on 0 ≤ t ≤ T . The Type I error, t 1 (t) = 1, in finding t s G (or G ) (satisfying C 1 (or C 1 ) by mere chance) is determined by G as a subset of [0, T ]. Upon partitioning of [0, T ] into n = 1204 intervals of size |G| = t g , this Type I error carries a PFA given by the ordinary probability p 1 = 1/n = t g /T . We note that this estimate is conservative compared to the same derived from G , since |G | = t g − τ.
Next, we generated PDF(t s ) by application of χ-image analysis over extended foreground: the start-time t s of a candidate emission feature according to Eq. 7 in our spectrograms. For descending chirps characteristic for spin-down of a compact object (neutron star or black hole) powering GRB170817A, Eq. 7 appears minimal in the "1+3" parameters t s and (nuisance parameters) (τ s , f s , f 0 ). Fig. 3 shows PDF(t s ) to be entirely in the t g = 1.7 s gap between GW170817 and GRB170817A. Under the null-hypothesis -a background of uncorrelated noise in H1 and L1 -over the original snippet of T = 2048 s of H1L1-data, the gap condition C 1 is represented by a Boolean random vari- Fig. 3. PDF of the start time t s (blue) of post-merger gravitational radiation to GW170817 following rejuvenation in gravitational collapse to a rotating black hole, setting the initial data to GRB170817A (top-left panel). The PDF shown (bin size 0.05 s) is superposed to estimates of the lifetime of the progenitor hypermassive neutron star formed in the immediate aftermath of the merger (reviewed in (Murguia-Berthier et al. 2020)). Highlighted is the 1.7 s gap (magenta) between merger time t m = 1842.43 s and the onset of GRB170817A. Post-merger emission to GW170817 appears as a descending chirp in the (H1,L1)-spectrogram merged by frequency coincidences f H1 − f L1 < 10 Hz (lower panels). For the merged spectrogram, parameter estimation ( Fig. 4) derives from time-slides less than the duration τ of the templates (11) against a background generated by larger time-slides. Post-merger gravitational radiation is identified with the formation of the central engine of GRB170817A, whose duration is constrained by the lifetime of black hole spin (bottom panel).
able t k (t) with values that are either "true" or "false," that is: Given G, t s in G according to PDF(t s ) (Appendix A) carries a PFA p 1 equal to Pr(t G) satisfying (Paper I): Causality over a background T hereby carries a PFA with twosided Gaussian-equivalent significance of 3.33 σ and false alarm rate (FAR) of about 1 per month, defined by FAR = PFA/T (e.g., Williams & Huntington 2018).

Extremal clustering in extended foreground
Parameter estimation PDF(t s ) (also PDF(τ s )) for our candidate feature is derived from merged (H1,L1)-spectrograms over an extended foreground (S 1 ) by application of small time-slides in the interval, where τ 0.5 s is the intermediate duration of the templates in butterfly matched filtering. We apply this to the original snippet of T = 2048 s of H1L1-data in 99 steps symmetric about zero, including a neighborhood of zero covered by 41 tiny time-slides within 10 ms. Following Eq. 11, the parameter estimations of candidate features give rise to "clustering" about peaksχ * =χ(t * s , τ * s ) shown in Fig. 4, where the remaining two parameters ( f s , f 0 ) in our χ-image analysis have been projected out by taking maxima. Such results by χ-image analysis of spectrograms for each ∆t. In the extended foreground thus produced, clustering in data points -the output of χ-image analysis of merged (H1,L1)spectrograms -is identified in the box, For the background, the ensembles of peaksχ * are identified over a partitioning of a data in cells of duration w. Each cell contains a total of Q = 161 data points (one for each time-slide). According to Eqs. 11-12, up to 99 may be clustered.
Accordingly, on the extended foreground (S 1 ) defined by domain (11), our two-level pipeline of (H1,L1)-spectrograms followed by χ-image analysis defines a map, whose image may be partially or wholly contained in the box (12). In case of the latter, the normalized cluster size: defined by the ratio of counts in B and counts in I. The limiting case η = 100% will be referred to as extremal, shown in Fig. 4.
(A movie of this time-slide analysis has been published online. 2 ) Absent a signal, t s and τ s are uncorrelated. A uniform probability distribution of the event time t s over the data at hand effectively prevents clustering. In response to a signal, t s and τ s tend to be clustered about some central value associated with a goodness of fit by Eq. 7. Such clustering is representative for correlations between H1 and L1 imparted by a gravitational-wave passing by (Fig. 4). Clustering in spectrograms hereby increases with signal strength (van Putten (2017) and tends to persistent over small time-slides Eq. 11. By this property, the condition of extremal clustering provides a robust guard against false positives.
With essentially no free parameters, the extremal condition over extended foreground (S 1 ) provides a robust alternative to conventional thresholding in amplitude-based criteria over foreground. A common implementation of the latter is a pre-defined threshold of signal-to-noise ratio (S/N) of a continuously valued statistic. Figure 4 with extremal clustering in (t s , τ s ) over Eq. 11 confirms identification based on a global maximumχ * =χ(t * s , τ * s ) at ∆t = 0 reported previously in Paper I, signaling Extended Emission post-merger in Fig. 3. This cluster is extremal following Eq. 13. Importantly, it is unique over the entire H1L1-data snippet of T = 2048 s following a partition over 1024 cells of 2 s cells. We hereby identify: satisfying C 1 -a more precise statement of C 1 . In Eq. 14 and in what follows, results are quoted by central value and STD of the PDF. Estimating PDF(t s ) by Eq. 7 serves to produce a statistic t s in which (τ s , f s , f 0 ) are nuisance parameters. Evaluated over afore-mentioned moving data-segments of ∆T obtains well-defined event timing t s at the maximum ofχ(t s ), similar though not identical to deriving a stable estimate of the duration of light curves of long GRBs in the face of large temporal variability by means of matched filtering (van Putten & Gupta 2009).
Given the uniform prior in timing (in the absence of astrophysical context), Eq. 9 is satisfied with associated probability of false alarm Eq. 10. Tracing the cluster of where the first follows from Eq. 14. Robustness of maximalχ conditional on extremal clustering in Fig. 4 pans out in uniqueness already in a minimally extended foreground of 7 time-slides extending over T 100 T in analyzing 50 frames of H1-L1 data, blindly selected to avoid biases in data-quality. This background analysis increases significance (10) in timing by causality by a factor of T /T , that is, improving (10) by: with a two-sided Gaussian-equivalent significance of 4.46 σ and FAR of about 1 per 782 years. This degree of robustness motivates a further analysis of timing based on individual analysis of H1 and L1 data with no regards to context. Starting point for is generating PDF's in parameter estimation circumventing the use of time-slides, to be applied to each of the H1 and L1 datachannels.

Analysis of individual H1 and L1 spectrograms
Independent of context, statistical significance in timing can be evaluated according to consistency in t s obtained from the mutually independent H1 and L1 observations: within the finite precision of our observation.
To evaluate C 2 , we appeal to PDF 1 (t s ) of H1 and PDF 2 (t 2 ) of L1 and their cross-correlation function. The cross-correlation of the two will have a maximum, whose proximity to zero will be determined by the associated STD. As before, given a total duration of observation, T , this STD normalized to T provides us with a second PFA. This PFA from [t s ] is independent of the PFA(t s ) in the mean of H1 and L1 ( §3). Absent gravitational waves passing by, H1 and L1 output are uncorrelated. Given this null-hypothesis, mean and difference of timing (t s , τ s ) are independent. A PDF t s,H1 , t s,L1 , τ s,H1 , τ s,L1 hereby factorizes as p A (t s , τ s ) p B ([t s ] , [τ s ]). Integration over differences reduces it to p A (t s , τ s ); and to (10) when further integration over τ s . On the other hand, integrating over (t s , τ s ) reduces the same to p B ([t s ] , [τ s ]) independently of the context, that is, with no regards to the multimessenger data GW170817-GRB170817.
The finite dispersion in the cross-correlation of PDF 1 (t s ) and PDF 2 (t s ) is expected to be quite small, albeit greater than the light travel-time between H1 and L1. In butterfly matched filtering over templates of intermediate duration τ, temporal resolution in a single cross-correlation between data and template is governed by τ. However, resolution is enhanced when averaging over a number of trials. To this end, we generated PDF's by nm trials in a search extending over n seeds of our template banks in butterfly matched filtering, while output is branched into m spectrograms defined by a stride m. Here, template seeds refer to the master templates prior to time-slicing over τ and conversion to time-symmetric templates (van Putten et al. 2014a). As in the analysis of the previous section, these calculations are realized by modern heterogeneous computing (Appendix D).
We generate PDFs with a total of nm = 320 elements from n = 10 different template seeds and a stride m = 32. From a single template seed, the output of butterfly filtering comprises L = O 10 8 lines of output of the form: where t i = i/4096 s is the sample time, ρ i denotes the correlation between data and specific template for detector H1 (i = 1) and L1 (i = 2) and f i j refers to the initial ( j = 1) and final ( j = 2) frequency of the template. For a stride m, this output branches into m spectrograms. Image analysis of each hereby produces a PDF of the parameters at hand, thus broadening our search by a factor nm in each of H1 and L1. Inherent to our approach is a bank of templates that is dense rather than orthogonal, and approximate over an intermediate times scale (van Putten et al. 2014a;van Putten 2017). A sample i in (18) is hereby typically covered by a multiplicity of different templates. The average time-stepping: in t i is correspondingly smaller than the sample time 1/4096 s of data; ∆t ik is typically zero or 1/4096 s with STD of about 40 µ s. A stride m = 32 hereby covers about 7 samples. The associated STD in the mean of frequencies is similarly small, i.e., O 10 −2 Hz. Offsets m 0 = 1, 2, · · · , m hereby provide exceptionally highresolution time-stepping in our matched filtering analysis. Figure 5 shows the results of parameter estimation in the (t s , τ s )-plane of H1 and L1 by χ-image analysis applied to the H1-and L1-spectrograms, generated over a large number of trials. Noticeably, our candidate signal appears stronger in H1 than in L1, indicated byχ(t * s ) and the peak in the H1-and L1histograms over t s . Furthermore, the results over  Figure 5 shows PDF 1 (t s ) and PDF 2 (t s ) thus extracted from H1 and L1 with their cross-correlation obtained by time-slide analysis applied to this pair. A key finding is that these PDF's are maximally correlated at a time-slide consistent with zero within the time-slide resolution δt res = 0.025 s. This maximum is global (maximal over all such maxima identified in the 64 segments of 32 s covering our analysis of the H1L1 data of duration T = 2048 s. It thereby carries a PFA: with equivalent two-sided Gaussian-equivalent significance of 4.06 σ and FAR of ∼ 1 per 1.4 year.
In deriving Eq. 20 from |[t s ]| < δt res at the peak of crosscorrelating PDF t s,H1 and PDF t s,L1 , a factor of 2 derives from absolute value of time-differences [t s ] = |t H1 − t L1 | in H1 and L1; another factor of 2 derives from the top-hat distribution in [t s ] based on uniform priors in t H1 and t L1 . Similar to p 1 , consistency in the PDF's of t s is hereby formally cast under the nullhypothesis by the Boolean random variable Y {0, 1} with Y = 1 if |[t s ]| < δt res , Y = 0 otherwise.
The mean, now from the pair of PDFs of t s from the independent H1-and L1-analyses, shows a delay time in gravitational collapse t s − t m = (0.92 ± 0.09) s, consistent with Eq. 15 based on the mean t s in merged (H1,L1)spectrograms.
In astrophysical context, Eq. 21 is consistent with the lifetime 0.98 +0.31 −0.26 s of the progenitor hypermassive neutron star derived from jet propagation times and the mass of blue ejecta (Gill et al. 2019) and its survival time, based on neutrino cooling and internal dissipation of rotational energy (e.g. Beniamini & Lu 2021). For τ s , the PDFs associated with the cluster (21) show τ s = (3.1 ± 0.1) s (H1) and τ s = (2.9 ± 0.2) s (L1) with mean: consistent with Eq. 15 in the previous analysis of merged (H1,L1)-spectrograms. Identified with the lifetime of black hole spin, it appears to constrain the duration T 8−70keV 90 of GRB170817A (Pozanenko et al. 2018). Thus, Eqs. 21-22 both appear to possess natural astrophysical context in the electromagnetic spectrum.

Discussion
The complex merger sequence Eq. 2 is probed by multimessenger calorimetry and event timing (Fig. 1). Following calibration  Table 1. Results on event timing in post-merger emission Eq. 4 to GW170817 identified with rejuvenation Eqs. 3 in Eq. 2. Type I errors derive independently from the gap condition, C 1 (PFA p 1 ), and mutual consistency, C 2 , (PFA p 2 ) in H1L1-data over a duration, T . Central values and uncertainties derive from mean and standard deviation of PDF(t s ) and PDF(τ s ). Uncertainties defined by the width of the respective PDFs reflect amplitude information.  (Fig. 2), a post-merger descending chirp characteristic for spin-down of a compact remnant is analyzed for output E GW and statistical significance by event timing subject to two consistency conditions, C 1,2 , on time-of-onset t s seen in merged and individual H1-and L1-spectrograms ( §3-4). Figure 6 shows a detailed summary.
E GW in Eq. 4 suffices to break the degeneracy between a hyper-massive neutron star and a black hole by exceeding limits in E J of the first. Yet, it is amply accounted for by E J following (delayed) gravitational collapse to a rotating black hole in the af-termath of GW170817 ( §2). Statistical significance is expressed by a PFA factorized over two statistically independent PFAs Eq. 10 in §3 and Eq. 20 in §4, derived from consistency conditions on effectively the mean (C 1 ) and, respectively, a difference (C 2 ) in start-time, t s .

Insensitivity to numerical parameters
Event timing analysis is rather insensitive to numerical parameters given that certain discrete conditions are met in PDF(t s ) when evaluating C 1 (Appendix A) and in considering the difference [t s ] in C 2 as follows. PFA p 1 from C 1 (more conservative than PFA from C 1 ) is determined by G relative to the observational window [0, T ], given that PDF(t s ) is well within G (Fig. 3, Appendix A). p 1 is hereby insensitive to trial factors that might otherwise be encountered when inferred from a statistic in continuous stochastic variables such as the S/N. It likewise does not involve Occam factors, common to likelihood analysis of the same in Bayesian approaches. The PFA p 2 from C 2 is derived from cross-correlating PDF(t s ) based on two individual H1 and L1 analyses. Sampled by nm = 320 trials, they appear to be close to the noise limit of the data. Moreover, any systematic uncertainty in PDF(t s ) cancels out in the difference [t s ] ( §4).
An additional safeguard derives from extremal clustering (100% in Eq, 13). It implies the width of PDF(t s ) to be much smaller than the gap size, t g , of GW170817-GRB170817A (Fig.  3 in A.2). Over extended foreground, it derives as a goodness of fit over (t s , τ s , f s , f 0 ), representing an essentially minimal number of "1+3" parameters for descending chirps. (By count this is the same in a canonical CBC search, that is, taken over two masses and orbital ellipticity for merger times 0 ≤ t m ≤ T .) There is no fine-tuning in PDF(t s ) derived from Eq. 7 with maximaχ =χ(t s ) upon projecting out the nuisance parameters (τ s , f s , f 0 ).
Our focus on Eq. 2 gives a concrete example of multimessenger calorimetry and event timing analysis, summarized in the chronicle Fig. 6. The joint PFA p 1 p 2 = 4.1 × 10 −8 from (discrete) event timing is noticeably on par with the PFA 8.91 × 10 −7 of GW170817-GRB170817A, merging the two p-values of (continuous stochastic variables) in temporal and spatial agreement between H1, L1, and the Fermi GBM reported by LIGO-Virgo (Abbott et al. 2017). A further third PFA derives from consistency in τ s from the individual H1 and L1 observations (Table  1). In this regard, our total PFA by t s is conservative.

Advances over previous results
We present two independent analysis of merged and individual H1 and L1 spectrograms, results of which may be compared for consistency and combined to enhance confidence levels on the basis of high-resolution PDFs ( §3-4). To be precise, in the absence of a gravitational-wave signal, data from H1 and L1 are statistically independent (the working hypothesis in multidetector GW analysis). Under the null-hypothesis, merged and individually, H1 and L1 satisfy uniform priors on event timing. We consider the mean ( §3) and difference ( §4) in event timing expressed in a statistic t s , representing two linearly independent combinations of the H1 and L1 time-series equivalent to a rotation over π/4 in the plane. By unitary of rotations, statistical independence of the H1 and L1 data is preserved in the two statistics t s of §3-4. Table 1 lists the measurement results.
A comparison with previous results derived from foreground (S 0 , time-slide zero) is also opportune. Our PDF(t s ) produces t s 0.92 s (Table 1), revising t s 0.67 s (Paper II).
This revised estimate is illustrative for a limitation of foreground in representing a single sample in PDF(t s ) derived from extended. In the face of finite scatter seen in the plot of t s versus time-slide ∆t (Fig. 4, top right panel), foreground estimates incur a statistical uncertainty on the order of the width of PDF(t s ). Additionally, foreground incurs a systematic error given by the finite difference in signal arrival time between the two detectors.
The PDF(t s )s from merged and individual H1 an L1 analysis permit a direct comparison between the two for their implied PFAs p 1 and p 2 ( Table 1). The inferred PFAs are qualitatively distinct, however, being contextual (by the gap t g = 1.7 s between GW170817-GRB170817) and, respectively, acontextual. Strictly speaking, a direct numerical comparison of the two is not opportune.
Nevertheless, our p 1 > p 2 may appear paradoxical in light of a visibly higher contrast of the descending branch in merged rather than individual spectrograms. Indeed, the width in PDF(t s ) is relatively smaller, approximately by a factor of √ 2 (Table 1), as expected. Increasing the time of observation by ×100, p 1 < p 2 (Eq. (16) with confidence level 4.46σ. In fact, including amplitude information, p 1 decreases below p 2 to 2.5 × 10 −5 (4.2σ, supplementary data to Paper I). The same applied to p 1 provides a combined confidence level of 5.16σ, but this avenue is not pursued further here.
Discrete event timing producing PFAs differently from amplitude-based PFA (Appendix A). Fig. 1 shows the candidate signal to be sufficiently strong for PDF(t s ) to be well within the gap of GW170817-GRB170817A, satisfying extremal clustering (Fig. 4). It hereby satisfies causality -a Boolean valued conclusion. True, in the case at hand, fixes PFA 1 at p 1 = t g /T , where T is the time of observation. A further increase in signal strength, narrowing PDF(t s ), will not change this value.
For uniformity of presentation, PFAs are here given by event timing only, sufficient to derive a significantly improved joint PFA, factored over the independent PFA 1 and PFA 2 , to validate the descending branch.

Black hole spin-down
We identified our estimate of τ s 3 s in Eqs. 15, 22 with the Kelvin-Helmholtz time-scale of black hole spin-down against high-density matter through a torus magnetosphere (van Putten & Levinson 2003;van Putten 2015a): where σ = M T /M in the right hand-side refers to the mass-ratio of torus-to-black hole with radius zR g , R g = GM/c 2 , specialized to the present GRB170817A. This timescale Eq. 23 is derived from the first law of thermodynamics. For a black hole with mass, M, angular momentum, J H , and angular velocity, Ω H , interacting with a torus magnetosphere rotating at the angular velocity Ω T of the inner face of a surrounding torus, τ H in (23) in the limit of small Reynolds stresses in suspended accretion, where L j L T refers to the luminosity L j in a baryon-poor jet (BPJ) along an open magnetic flux tube subtended over a finite polar angle θ H π/2 on the event horizon of the black hole and a luminosity L T −Ω TJ in multimessenger radiation from the surrounding inner disk or torus. Here, the inequality L j L T L GW is seen in discrepant energies in electromagnetic Eq. 6 and, respectively, gravitational radiation Eq. 4.
GRB170817A derived from a BPJ is included in Eq. 6. While GRB170817A is negligible in the total energy budget E + J ( §2), it provides crucial timing information t GRB and T (8−70)keV 90 . It can be attributed to L j derived from a Faraday-induced potential (van Putten 2000;van Putten & Levinson 2003),   (2) showing detailed event timing in Fig. 1. Highlighted is E GW by rejuvenation in gravitational collapse of the hyper-massive neutron star formed in the immediate aftermath of the merger at t s in the 1.7 s gap G between GW170817-GRB170817A, in a Kerr black hole with spin-energy E + J exceeding E − J in its progenitor according to (3). E + J amply accounts for E GW = (3.5 ± 1) %M c 2 , breaking the degeneracy between neutron stars and black holes. of charged particles with angular momentum J p = eA ϕ along surfaces of constant flux A ϕ = Φ/2π subject to frame-dragging by a black hole in its lowest energy state. While L j for θ H π/2 is relatively small, it may have observational consequences for GRBs and UHECRs alike and especially so when intermittent (van Putten & Gupta 2009;Shahmoradi & Nemiroff 2015;van Putten 2015b;Gottlieb et al. 2021). For a discussion on highenergy emission from black holes that are not in the lowest energy state, we refer to Rueda et al. (2022).

Conclusions
From analysis of H1 and L1 observations in merged and, new in this work, individual detector spectrograms producing statistically independent PFAs p 1 and, respectively, p 2 , from their respective high-resolution PDF(t s )s of an extended emission feature, a number of findings emerge: 1. A black hole central engine to GRB17017A based on E GW 3.5%M c 2 in a descending chirp exceeding the limits of the HNS in the immediate aftermath of GW170817 (Eqs. 4, 7, 5), attributed to rejuvenation of E J in gravitational collapse ( §2); 2. Statistical significance with combined PFA of 4 × 10 −8 (equivalent Z-score 5.48) factored over two independent PFAs in Eqs. 10 and 20 derived from event timing subject to C 1,2 (Table 1) (Pozanenko et al. 2018). We refer in particular to Figs. 4,6,and Table 2; 5. The gravitational-wave energy emitted in the combined ascending-descending chirp of GW170817 is about 2% of the total mass-energy of the system.

Outlook
Planned LVK observational runs O4-5 offer potentially important new observational opportunities to probe merger sequences involving neutron stars and central engines of energetic corecollapse supernovae in the Local Universe. EM-GW observations involving GW-calorimery and event timing offer some new tools to identify their nature (Cutler & Thorne 2002), particularly when breaking the degeneracy between (super-or hyper-) massive neutron stars and black holes.
Since black holes have no memory of their progenitor except for total mass and angular momentum, such appears notably opportune for type Ib/c supernovae as the parent population of normal long GRBs (van Putten et al. 2019b) and superluminous supernovae (Dong et al. 2015). While rare, failed GRB-supernovae nevertheless may be luminous in gravitational radiation and more frequent than double neutron star mergers.
For the planned LVK runs O4-5, core-collapse supernovae may be probed in blind or optically triggered searches in the Local Universe over distances comparable to GW170817. If detected, gravitational-wave emission is likely to identify their enigmatic central engine, significantly complementing our understanding of core-collapse events currently limited to SN1987A. In the more distant future, searches for broadband gravitational-wave radiation may put rigorous and model-independent bounds on (non-axisymmetric) mass-motion around supermassive black holes such as SgrA* by the planned Laser Interferometric Space Antenna (LISA; van Putten et al. 2019b).

Appendix A: Discrete event timing versus S/N
GW170817-GRB170817A provides a multimessenger context by the gap G of t g = 1.7 s in between these two events with merger time t m and, respectively, the time-of-onset t GRB (Fig. 1).
For a candidate post-merger emission feature associated with the central engine of GRB170817A, a PFA derives from event timing as an alternative to conventional signal-to-noise (SNR) analysis under the null-hypothesis of stationary detector noise and the uniform prior of astrophysical event timing (H0).
To illustrate this, consider an indicator function χ(t) marking the event time t e of a candidate feature by its global maximum over a finite duration of observation T . Since χ(t e ) is continuous stochastic observable, t e is well-defined and unique. Thus, t e tends to be uniformly distributed over [0, T ], namely: with event time now discrete over n = T/t g bins. The gap condition (1) hereby carries a PFA equal p 1 = t g /T . In contrast, gravitational-wave emission from the putative central engine of GRB170817A (H1) changes our expectation to Pr = p(t e |H1) with a bias to satisfying Eq. 1, provided χ(t) is suitably devised to measure the start-time t e = t s of this emission. As may be expected, p(t e |H1) > p(t e |H0) for t e G depends on the S/N. Figure A.1 shows an elementary Monte Carlo simulation of the probability Pr(t e G) of satisfying the gap condition Eq. 1 as a function of S/N, based on event times, t e , defined by maxima of χ(t) given by the sum of normally distributed noise over n = 37 bins plus a signal at bin 30. 3 For n = 37, this PFA equals the chance of winning a single bet in European roulette conform the probability theory of Blaise Pascal.
A PFA from Eq. 1 rather than based on the S/N has the advantage of being independent of trial factors, though the result is a priori limited by G (relative to T ) and the implementation requires PDF(t e ) to be sufficiently narrow relative to G, that is: Once (A.2) is satisfied, p 1 = Pr 0 is independent of the width of PDF(t s ), which might otherwise depend on trial factors. In this sense, (A.2) represents a discrete timing condition. As an ordinary probability, furthermore, p 1 readily combines with other PFAs.

Appendix B: Whitening H1L1-data
As a pre-processing step applied to LIGO strain data, we applied a band pass filter of 10-1700 Hz followed by whitening to suppress numerous lines, mostly violin modes which appear prominently in the spectrum of the detectors (Fig. B.1) (Abbott et al. 2020). Figure B.1 shows the computation by the standard Welch method. Here, the frequency resolution is slightly less than in Fig. 1 of Paper II by a different partitioning in the time domain. By the time-frequency uncertainty relation, lines hereby vary in width yet with the same total energy by Parseval's Theorem. We note, however, that constant frequencies are suppressed in our butterfly matched filtering (Appendix C).
To this end, strain data are normalized in Fourier domain by amplitude following a partitioning of the spectrum over intervals of B 2048 s snippet of H1L1-data covering GW170817 sampled at a n s = 4096 Hz, we have N = 2 12 [T/s] = 2 23 Fourier coefficients c k (k = 1, 2, · · · N). Its Fourier spectrum is partitioned up to the Nyquist frequency 1 2 n s by intervals B s : (s − 1)M < k ≤ s × M of size B[Hz], each comprising M = (B/n s ) N = BT coefficients; to exemplify this, M = 4096 for B = 2 Hz. Normalization to an essentially flat spectrum is obtained by: over all intervals s = 1, 2, · · · , 1 2 n s /B, where A s = M −1 k B s |c k | is the mean of the absolute values of c k in B s with |B s | = B. Crucially, (B.1) preserves phase in each of the Fourier coefficients. Effectively the same whitening obtains upon setting A s equal to the standard deviation of the c k in B s . This procedure suppresses violin modes in the LIGO detectors, provided B is about 1-10 Hz -B must exceed the maximal width of the violin modes yet be relatively modest to effectively suppress the same. The result is essentially flat spectra (Fig. B.1), while preserving signals of interest such as the relatively longduration merger signal GW170817, evidenced in Fig. 3. Whitening (B.1) is effective also in rendering GW170817 to be directly audible, for instance, using the MatLab function sound.
In the present application to the H1L1-data covering GW170817, we empirically verified in previous work (Paper I) that all results remain essentially unchanged for B = 2 − 16 Hz.

Appendix C: Butterfly matched filtering
Butterfly filtering is an essentially linear filter for signals with frequencies slowly wandering in time parameterized by δ (Fig.  C.1). This is devised by matched filtering over an effectively dense bank of time-symmetric chirp-like templates of duration In the time-frequency diagram, the passing of such signals can be schematically indicated by "butterflies" (Fig. C.1), indicated by a finite slew rate (time rate-of-change) in frequency: for some choice of δ > 0. Constant frequency signals fail to pass through the same, indicating relative suppression. Butterfly matched filtering is hereby distinct from Fourier analysis, that favors passing signals with relatively constant frequency. Butterfly filtering was originally developed to extract broadband spectra of light curves of long GRBs of the BeppoSAX catalog (van Putten et al. 2014a). It identifies Kolmogorov spectrum which extends up to the Nyquist frequency of 1024 Hz of the BeppoSAX sampling rate of 2048 Hz (during the first 8 seconds of long bursts), demonstrating a sensitivity one order of magnitude beyond what is attained by conventional Fourier analysis (Fig. C.1).
The template bank of butterfly matched filtering densely covers a two-parameter range in frequency and time rate-of-change of frequency. These templates are produced by time-slicing of a seed template (Fig. C.2) over aforementioned duration τ. Timesymmetric templates obtain by linear combination with their time-reverse, realizing equal sensitivity to signals whose frequency increases or decreases in time.
Butterfly filtering is ported to LIGO data analysis with no principle change in algorithm. An early demonstration shows the suppression of lines in unwhitened LIGO data, as constant frequency signals fail to correlate with chirp-like templates (Fig.  C.2).
Analysis of large data-sets using a dense bank typically containing 2 18 templates requires heterogeneous computing, offloading matched filtering to graphics processor units (GPUs) with additional procedures to manage limitations of bandwidth between GPU and host CPU given exascale computations involved (van Putten 2017, Appendix C). List of symbols and key words in the two-level pipeline of H1L1 data analysis probing the merger sequence Eq. 2 subject to constraints C 1,2 evaluated by PDF(t s ) produced by butterfly matched filtering and χ-image analysis, augmented by time-slide analysis, implemented by load balanced heterogeneous parallel computing.   Here, the template bank size is expressed in units of 1M templates. A throughput in excess of unity defines faster than realtime analysis. Fig. D.1 illustrates D.2-D.1 for a quad-GPU configuration of an air-cooled node with AMD Radeon VII GPUs featuring η 47% and a cross-correlation rate up to 200 kHz between data and templates over segments of 32 s. . Array sizes exceeding 2 12 cause a drop in performance as matrix transpose exceeds 32 kB of local memory, limited now by 1TB/s bandwidth to Global Memory in HBM2. We note that clButterfly uses a default size N = 2 17 (32 s of data at 4096Hz sampling rate) with batch size over n = 128 per frame of 4096 s. A quad-GPU node hereby features a maximum of about 5 teraFLOPs, circumventing limitations of PCIe bandwidth optimized by pre-and post-callback functions (van Putten 2017). The throughput of this quad-GPU node is shown up to six times real-time over a template bank normalized to size 2 20 (right panel). GPU-CPU communication passes through a post-callback function, limiting data-transfer to tails exceeding a threshold κ = 2 in normalized butterfly matched filtering output, effectively 0.1% of all matched filtering results. Throughput peaks in multi-frame analysis at a cross-correlation rate of about 200kHz.

Appendix D.5: Load balancing by synaptic processing
Embarrassingly parallel computing for tasks comprising clButter f ly (step 3) or clChi (step 5) in §3-4 are specified as a line-items for processing on a multi-LAN heterogeneous computing platform.
We used synaptic parallel processing in which line-items are requested by the GPU-nodes, rather than issued by a server. Before execution, these line-items are updated with the URLs for input and output data on the LAN hosting a node. This process ensures maximal throughput on heterogeneous platforms of nodes with different performance characteristics. Synaptic parallel processing is written in bash for Linux and OSX.