Skip to content

fix slicing tree data being busted when a subset of HLVK is set

Leo Tsukada requested to merge master-dtdphi-treeslice-fix into master

While investigating retractions in MDC9, I realized the current dtdphi codes return inconsistent results depending on instruments argument given during the initialization of InspiralExtrinsics class object. e.g. for G810140 if instruments=("H1", "L1", "V1", "K1"):

>>> IE = InspiralExtrinsics(instruments=("H1", "L1", "V1", "K1"))
>>> times = {}
>>> phases = {}
>>> snrs = {}
>>> horizons = {}
>>>
>>> horizons['H1'] = 434.1128632208363
>>> horizons['L1'] = 535.4751061681576

...
>>> # G810140
...
>>> # Times
... times['H1'] = 0
>>> times['L1'] =  -0.014195348
...
>>> # Coalescenece Phases
... phases['H1'] = -1.6397246
>>> phases['L1'] = -1.1792845
>>>
>>> # SNRs
... snrs['H1'] = 7.529696
>>> snrs['L1'] = 4.0554938
>>> print( IE.time_phase_snr(times, phases, snrs, horizons) )
[1.1154665e-30]

where as if instruments=("H1", "L1", "V1"), which was used in the production analysis

>>> IE = InspiralExtrinsics(instruments=("H1", "L1", "V1"))
>>> times = {}
>>> phases = {}
>>> snrs = {}
>>> horizons = {}
>>>
>>> horizons['H1'] = 434.1128632208363
>>> horizons['L1'] = 535.4751061681576
...
>>> # G810140
...
>>> # Times
... times['H1'] = 0
>>> times['L1'] =  -0.014195348
...
>>> # Coalescenece Phases
... phases['H1'] = -1.6397246
>>> phases['L1'] = -1.1792845
>>>
>>> # SNRs
... snrs['H1'] = 7.529696
>>> snrs['L1'] = 4.0554938
>>> print( IE.time_phase_snr(times, phases, snrs, horizons) )
[0.0592612]

This boils down to the fact that in TimePhaseSNR class, self.slices is derived from the given instruments set. However, the internal tree data in the dtdphi h5 file is structured for HLVK instrument set. The instrument set (i.e. array shape) that these two attibutes are based on needs to be consistent. Otherwise, that would let the dtdphi codes slice the entire tree data in nonsense way, as follows:
for instruments=("H1", "L1", "V1", "K1")

>>> self.slices
{('H1', 'K1'): [0, 1, 2], ('H1', 'L1'): [3, 4, 5], ('H1', 'V1'): [6, 7, 8], ('K1', 'L1'): [9, 10, 11], ('K1', 'V1'): [12, 13, 14], ('L1', 'V1'): [15, 16, 17]}

for instruments=("H1", "L1", "V1")

>>> self.slices
{('H1', 'L1'): [0, 1, 2], ('H1', 'V1'): [3, 4, 5], ('L1', 'V1'): [6, 7, 8]}

After this fix, these two cases have matched up and the probability for G810140 is now identical.

Edited by Leo Tsukada

Merge request reports