From fb317fbea8ea9435eb76a79d1e29e7399fe55c1f Mon Sep 17 00:00:00 2001
From: Gregory Ashton <gregory.ashton@ligo.org>
Date: Wed, 11 Sep 2019 17:57:20 -0700
Subject: [PATCH] Resolve issues identified in !408
1) Make the constants consistent with those used in lal
2) Deprecate the gps_time_to_gmst function
3) Use the lal.GreenwichMeanSiderealTime(time) method directly
---
bilby/core/utils.py | 14 +++++++++-----
bilby/gw/utils.py | 8 ++++----
2 files changed, 13 insertions(+), 9 deletions(-)
diff --git a/bilby/core/utils.py b/bilby/core/utils.py
index 6f53601f7..4c8305e5a 100644
--- a/bilby/core/utils.py
+++ b/bilby/core/utils.py
@@ -10,6 +10,7 @@ import subprocess
import multiprocessing
from importlib import import_module
import json
+import warnings
import numpy as np
from scipy.interpolate import interp2d
@@ -18,12 +19,11 @@ import pandas as pd
logger = logging.getLogger('bilby')
-
-# Constants: values taken from astropy v3.0.4
+# Constants: values taken from LAL cd65f38ce43cef3a1dec217c060de25caf99bf14
speed_of_light = 299792458.0 # m/s
-parsec = 3.0856775814671916e+16 # m
-solar_mass = 1.9884754153381438e+30 # Kg
-radius_of_earth = 6378100.0 # m
+parsec = 3.085677581491367e+16 # m
+solar_mass = 1.9885469549614615e+30 # Kg
+radius_of_earth = 6378136.6 # m
_TOL = 14
@@ -254,6 +254,10 @@ def gps_time_to_gmst(gps_time):
float: Greenwich mean sidereal time in radians
"""
+ warnings.warn(
+ "Function gps_time_to_gmst deprecated, use "
+ "lal.GreenwichMeanSiderealTime(time) instead",
+ DeprecationWarning)
omega_earth = 2 * np.pi * (1 / 365.2425 + 1) / 86400.
gps_2000 = 630720013.
gmst_2000 = (6 + 39. / 60 + 51.251406103947375 / 3600) * np.pi / 12
diff --git a/bilby/gw/utils.py b/bilby/gw/utils.py
index eb2fce923..dc90e054d 100644
--- a/bilby/gw/utils.py
+++ b/bilby/gw/utils.py
@@ -6,7 +6,7 @@ import numpy as np
from scipy.interpolate import interp1d
import matplotlib.pyplot as plt
-from ..core.utils import (gps_time_to_gmst, ra_dec_to_theta_phi,
+from ..core.utils import (ra_dec_to_theta_phi,
speed_of_light, logger, run_commandline,
check_directory_exists_and_if_not_mkdir,
SamplesSummary)
@@ -87,7 +87,7 @@ def time_delay_geocentric(detector1, detector2, ra, dec, time):
float: Time delay between the two detectors in the geocentric frame
"""
- gmst = gps_time_to_gmst(time)
+ gmst = np.mod(lal.GreenwichMeanSiderealTime(time(time), 2*np.pi)
theta, phi = ra_dec_to_theta_phi(ra, dec, gmst)
omega = np.array([np.sin(theta) * np.cos(phi), np.sin(theta) * np.sin(phi), np.cos(theta)])
delta_d = detector2 - detector1
@@ -120,8 +120,8 @@ def get_polarization_tensor(ra, dec, time, psi, mode):
array_like: A 3x3 representation of the polarization_tensor for the specified mode.
"""
- greenwich_mean_sidereal_time = gps_time_to_gmst(time)
- theta, phi = ra_dec_to_theta_phi(ra, dec, greenwich_mean_sidereal_time)
+ gmst = np.mod(lal.GreenwichMeanSiderealTime(time(time), 2*np.pi)
+ theta, phi = ra_dec_to_theta_phi(ra, dec, gmst)
u = np.array([np.cos(phi) * np.cos(theta), np.cos(theta) * np.sin(phi), -np.sin(theta)])
v = np.array([-np.sin(phi), np.cos(phi), 0])
m = -u * np.sin(psi) - v * np.cos(psi)
--
GitLab