Comparison of dingo detector projection to research code

[max-dax]

The unit test test_detector_projection_agains_research_code in test_transforms.py compares the detector projection of dingo against that of the research code (non-whitened and non-noisy). We find a relative deviation of 1e-3 between these. I suspect that this is due to (i) using bilby instead of pycbc routines and (ii) using float64 instead of float32 for detector time computation and antenna projection. Please double check that this makes sense.

[stephengreen]

One difference between the two codes is the implementation of Greenwich mean sidereal time, which is called when working out antenna patterns and time delays:

• Bilby uses the lal implementation: https://github.com/lscsoft/bilby/blob/master/bilby/gw/utils.py#L987-L990
• pycbc uses astropy.time.Time: https://pycbc.org/pycbc/latest/html/_modules/pycbc/detector.html

When I tested these two time functions I got differences of order a few parts in 10^6. I'm not sure if this accounts for everything.

[stephengreen]

This difference in the sidereal times is approximately 0.02 ms in the case I looked at.

[stephengreen]

Also why is f_min set to 10 Hz in the test? The data file seems to start at 20 Hz (although it has small values for lower frequencies).

[max-dax]

f_min is currently not used. The only use case of f_min that I see is that the domain object should check that a specified truncation_range is entirely above f_min.

[stephengreen]

This is still somewhat mysterious. The error due to the differences in the Greenwich mean sidereal time does not seem to be big enough to account for the discrepancy.

I also discovered that the bulk of the deviation in the waveforms in L1 can be corrected by a time shift of approximately 5e-6 s. However, the time shift needed to correct the waveform in H1 is much smaller. Is it possible that GNPE is sneaking into here somehow, since this puts L1 and H1 on a different footing?

[stephengreen]

Plotting the data, the time-shift needed to correct the waveforms at low frequencies is different than at high frequencies. The only way I can think of for a frequency-dependent error to be introduced is via the time shift function. Any thoughts?

[mpuerrer]

Hmm, I don't have an explanation. Is it possible to take the polarizations in both cases and by hand project onto the detector? Then you can also swap the gmst's around.

[stephengreen]

I have been doing various experiments of the sort. I don't see how we can get to such a large time-shift based on errors in gmst alone. I assume that the projection will be subdominant to the overall time shift.

[stephengreen]

I carried out your experiment, substituting the pycbc gmst into the bilby routines, and the result was only a little bit changed.

There seems to be a bug somewhere...

[mpuerrer]

Strange, then I guess it must be some further downstream operations, assuming that the loaded polarizations are exactly the same in both codes.

[stephengreen]

Okay, this is weird, but the pycbc and bilby antenna responses differ by a part in 1000:

>>> detector_bilby.antenna_response(parameters_ref['ra'], parameters_ref['dec'], ref_time, parameters_ref['psi'], mode='plus')
0.3053774910228462
>>> detector_bilby.antenna_response(parameters_ref['ra'], parameters_ref['dec'], ref_time, parameters_ref['psi'], mode='cross')
-0.34402264322837994
>>> detector_pycbc.antenna_pattern(parameters_ref['ra'], parameters_ref['dec'], parameters_ref['psi'], ref_time)
(0.30578201025297086, -0.34365917031830795)

This accounts for the bulk of the differences.

[stephengreen]

import numpy as np
import time
import bilby
import lal
import lalsimulation

ifos = bilby.gw.detector.InterferometerList(["L1"])
ifo = ifos[0]
gpstime = 1305303144
trial = 100

fp_lal, fp_bilby, fc_lal, fc_bilby = [], [], [], []
for i in range(trial):
    ra, dec, psi = 2. * np.pi * np.random.uniform(), np.pi * np.random.uniform() - np.pi / 2., np.pi * np.random.uniform()
    fplus, fcross = ifo.antenna_response(ra, dec, gpstime, psi, "plus"), ifo.antenna_response(ra, dec, gpstime, psi, "cross")
    fp_bilby.append(fplus)
    fc_bilby.append(fcross)
    gmst = lal.GreenwichMeanSiderealTime(gpstime)
    fplus, fcross = lal.ComputeDetAMResponse(lalsimulation.DetectorPrefixToLALDetector("L1").response, ra, dec, psi, gmst)
    fp_lal.append(fplus)
    fc_lal.append(fcross)

print(np.std(np.array(fc_bilby) - np.array(fc_lal)))
print(np.std(np.array(fp_bilby) - np.array(fp_lal)))

gives

0.00032833988847648507
0.0002990699021891997

This seems higher than I would expect...

[stephengreen]

And if you do it fo H1 you get

8.917933177280991e-09
8.623888720202266e-09

Do they have the wrong L1 position / orientation???

[stephengreen]

I found the source of the error. In https://github.com/lscsoft/bilby/blob/master/bilby/gw/detector/detectors/L1.interferometer they do not include

xarm_tilt = -0.0003121
yarm_tilt = -0.0006107

These arm tilts are included for H1 (although not for V1, but LAL also has this set to 0). When I add them in, the antenna patterns match those of LAL. This should account for the bulk of the discrepancy with our L1 strain. Probably the difference in gmst accounts for the rest.

Detector references:

• https://www.ligo.org/scientists/GW100916/detectors.txt
• https://dcc.ligo.org/LIGO-T980044/public

[max-dax]

In the research code we use the PyCBC (=Lal) routines, and we found that our results agree better with Lal than with bilby. Could this bug in bilby be the reason? Or is the discrepancy too small?

[stephengreen]

I found that the main problem with the Bilby results was non-convergence of the samples. I think that is probably a different issue, although I am not entirely sure how an error in the detector projection of one part in 1000 would impact final results.

[mpuerrer]

Great job finding this issue in bilby! The error is rather small, but I think it would be good to report this as a bug in bilby so it can be fixed.

[stephengreen]

Okay I created an issue https://git.ligo.org/lscsoft/bilby/-/issues/606

[stephengreen]

I created an MR for Bilby to fix the issue. Until this is incorporated into a Bilby release, let's just set the L1 tilts manually at train time as above.

[stephengreen]

The L1 arm tilts should be correct in Bilby version 1.1.5 or newer.

[stephengreen]

Is this done?