20 """This module provides a Python class for generating antenna response
21 functions. By default the class uses its own implementation of the antenna
22 response functions, but it also provides a wrapper to the equivalent
127 from lal
import (LALDetectorIndexLHODIFF, LALDetectorIndexLLODIFF,
128 LALDetectorIndexGEO600DIFF, LALDetectorIndexVIRGODIFF,
129 LALDetectorIndexTAMA300DIFF, LALDetectorIndexKAGRADIFF,
130 LALDetectorIndexLIODIFF, LALDetectorIndexE1DIFF,
131 LALDetectorIndexE2DIFF, LALDetectorIndexE3DIFF,
132 CachedDetectors, LIGOTimeGPS, GreenwichMeanSiderealTime,
133 ComputeDetAMResponse, ComputeDetAMResponseExtraModes,
134 DAYSID_SI, Detector, TranslateHMStoRAD, TranslateDMStoRAD)
136 from .
import git_version
138 __author__ =
"Matthew Pitkin <matthew.pitkin@ligo.org>"
139 __version__ = git_version.verbose_msg
140 __date__ = git_version.date
144 DETMAP = {
'H1': LALDetectorIndexLHODIFF,
145 'H2': LALDetectorIndexLHODIFF,
146 'LHO': LALDetectorIndexLHODIFF,
147 'L1': LALDetectorIndexLLODIFF,
148 'LLO': LALDetectorIndexLLODIFF,
149 'G1': LALDetectorIndexGEO600DIFF,
150 'GEO': LALDetectorIndexGEO600DIFF,
151 'GEO600': LALDetectorIndexGEO600DIFF,
152 'V1': LALDetectorIndexVIRGODIFF,
153 'VIRGO': LALDetectorIndexVIRGODIFF,
154 'T1': LALDetectorIndexTAMA300DIFF,
155 'TAMA': LALDetectorIndexTAMA300DIFF,
156 'TAMA300': LALDetectorIndexTAMA300DIFF,
157 'K1': LALDetectorIndexKAGRADIFF,
158 'KAGRA': LALDetectorIndexKAGRADIFF,
159 'LCGT': LALDetectorIndexKAGRADIFF,
160 'I1': LALDetectorIndexLIODIFF,
161 'LIO': LALDetectorIndexLIODIFF,
162 'E1': LALDetectorIndexE1DIFF,
163 'E2': LALDetectorIndexE2DIFF,
164 'E3': LALDetectorIndexE3DIFF}
168 VAR_RANGES = {
'psi': [0., 2.*np.pi],
169 'time': [0., DAYSID_SI],
170 'ra': [0., 2.*np.pi],
175 def __init__(self, detector, ra, dec, psi=0., times=None, tensor=True,
176 vector=False, scalar=False, use_lal=False, lookup=False,
177 lookuppair=None, bins1=100, bins2=100):
179 Calculate the long-wavelength limit antenna response functions for a given
180 ground-based gravitational wave detector. The response can include tensor,
181 vector, and scalar modes.
183 @param detector: (str) a valid detector name (e.g., 'H1') or
185 @param ra: (array_like) the right ascension of the source in radians,
186 or a string of the format 'hh:mm:ss.s'.
187 @param dec: (array_like) the declination of the source in radians, or a
188 string of the format 'dd:mm:ss.s'.
189 @param psi: (array_like) the polarization angle in radians. Defaults to
191 @param times: (array_like) an array of GPS time values at which to
192 calculate the response function.
193 @param tensor: (bool) set to calculate and store the tensor
194 polarization components (plus and cross). Defaults to True.
195 @param vector: (bool) set to calculate and store the vector
196 polarization components ("x" and "y"). Defaults to False.
197 @param scalar: (bool) set to calculate and store the scalar
198 polarization components (longitudinal and breathing). Defaults to
200 @param use_lal: (bool) set to internally use the
201 XLALComputeDetAMResponse() and XLALComputeDetAMResponseExtraModes()
202 functions. Defaults to False.
203 @param lookup: (bool) set to generate and use a look-up table in a pair
204 of parameters for computing the antenna responses. Defaults to
205 False. If using the look-up table, the arrays of values being
206 "looked-up" must be in ascending order.
207 @param lookuppair: (list, tuple) a list of the two parameters that will
208 be used in the look-up table interpolation. Defaults to
209 <tt>['psi', 'time']</tt> (allowed values are <tt>'ra'</tt>,
210 <tt>'dec'</tt>, <tt>'psi'</tt>, or <tt>'time'</tt>)
211 @param bins1: (int) the number of bins in the grid in the first look-up
212 table parameter. Defaults to 100.
213 @param bins2: (int) the number of bins in the grid in the second
214 look-up table parameter. Defaults to 100.
216 Example usage for tensor polarizations:
218 >>> from lal.antenna import AntennaResponse
219 >>> # compute tensor response for a single time
220 >>> resp = AntennaResponse('H1', ra=1.2, dec=-0.3, psi=2.9,
222 >>> print('Fplus: {}'.format(resp.plus))
224 >>> print('Fcross: {}'.format(resp.cross))
225 Fcross: [-0.79809163]
226 >>> # re-use class to get response at multiple new times
227 >>> resp.compute_response([1010101010., 1234567890.])
228 >>> print('Fplus: {}'.format(resp.plus))
229 Fplus: [ 0.09498567 -0.45495654]
230 >>> print('Fcross: {}'.format(resp.cross))
231 Fcross: [0.1706959 0.21690418]
234 Example usage for tensor, vector and scalar polarizations (at a series
237 >>> import numpy as np
238 >>> times = np.linspace(1000000000.0, 1000086340.0, 1440)
239 >>> resp = AntennaResponse('H1', ra=1.2, dec=-0.3, psi=2.9,
240 ...scalar=True, vector=True, times=times)
242 array([0.32427018, 0.32805983, 0.3318344 , ..., 0.32780195, 0.33157755,
245 array([-0.79809163, -0.79607858, -0.79404097, ..., -0.7962166 ,
246 -0.79418066, -0.79212028])
247 >>> resp.x # vector "x" polarization
248 array([-0.46915186, -0.46773594, -0.46627224, ..., -0.46783399,
249 -0.46637354, -0.46486538])
250 >>> resp.y # vector "y" polarization
251 array([-0.17075718, -0.17475991, -0.17875012, ..., -0.17448742,
252 -0.17847849, -0.18245689])
253 >>> resp.b # scalar "breathing" mode
254 array([0.05365678, 0.05573073, 0.05780282, ..., 0.05558939, 0.05766162,
256 >>> resp.l # scalar "longitudinal mode"
257 array([-0.05365678, -0.05573073, -0.05780282, ..., -0.05558939,
258 -0.05766162, -0.05973181])
290 if lookuppair
is None:
293 self.
set_lookup_pairset_lookup_pair(pair=lookuppair, bins1=bins1, bins2=bins2)
306 Set the detector for which to calculate the antenna response.
308 @param det: (str) a valid detector name.
311 if isinstance(det, Detector):
312 self.
_detector_detector = det.frDetector.prefix
314 elif isinstance(det, str):
315 if det.upper()
not in DETMAP.keys():
316 raise ValueError(
"Detector is not a valid detector name")
323 raise TypeError(
"Detector must be a string or lal.Detector object")
332 Set the lal.Detector.
334 @param det: (str) a valid detector name.
337 if isinstance(det, Detector):
339 elif isinstance(det, str):
341 detector = DETMAP[det.upper()]
343 raise KeyError(
"Key {} is not a valid detector name.".format(det))
345 self.
_laldetector_laldetector = CachedDetectors[detector].response
347 raise TypeError(
"Detector must be a string or lal.Detector object")
356 Set the right ascension.
362 elif isinstance(raval, float)
or isinstance(raval, int):
363 self.
_ra_ra = np.array([raval], dtype=
'float64')
364 elif isinstance(raval, list)
or isinstance(raval, np.ndarray):
365 self.
_ra_ra = np.copy(raval).astype(
'float64')
366 elif isinstance(raval, str):
368 rarad = TranslateHMStoRAD(raval)
369 self.
_ra_ra = np.array([rarad], dtype=
'float64')
371 raise ValueError(
"Could not convert '{}' to a right "
372 "ascension".format(raval))
374 raise TypeError(
"Right ascension must be an array")
389 def dec(self, decval):
399 elif isinstance(decval, float)
or isinstance(decval, int):
400 self.
_dec_dec = np.array([decval], dtype=
'float64')
401 elif isinstance(decval, list)
or isinstance(decval, np.ndarray):
402 self.
_dec_dec = np.copy(decval).astype(
'float64')
403 elif isinstance(decval, str):
405 decrad = TranslateDMStoRAD(decval)
406 self.
_dec_dec = np.array([decrad], dtype=
'float64')
408 raise ValueError(
"Could not convert '{}' to a "
409 "declination".format(decval))
411 raise TypeError(
"Declination must be an array")
429 def psi(self, psival):
431 Set the value of the gravitational wave polarization angle psi.
433 @param psival: (float) the polarization angle (radians)
441 elif isinstance(psival, float)
or isinstance(psival, int):
442 self.
_psi_psi = np.array([psival], dtype=
'float64')
443 elif isinstance(psival, list)
or isinstance(psival, np.ndarray):
444 self.
_psi_psi = np.copy(psival).astype(
'float64')
446 raise TypeError(
"Polarization must be an array")
448 self.
_psi_psi = np.mod(self.
_psi_psi, 2.*np.pi)
460 except AttributeError:
464 def times(self, timearr):
466 Set array of times and GPS times.
474 elif isinstance(timearr, float)
or isinstance(timearr, int):
475 self.
_times_times = np.array([timearr], dtype=
'float64')
476 elif isinstance(timearr, list)
or isinstance(timearr, np.ndarray):
477 self.
_times_times = np.copy(timearr).astype(
'float64')
479 raise TypeError(
"Times must be an array")
483 for i, time
in enumerate(self.
_times_times):
484 gps = LIGOTimeGPS(time)
485 gmstrad = GreenwichMeanSiderealTime(gps)
513 def tensor(self, tensorval):
515 Set whether to include tensor polarizations.
518 if not isinstance(tensorval, bool):
519 raise TypeError(
"Must be boolean value")
521 self.
_tensor_tensor = tensorval
528 def vector(self, vectorval):
530 Set whether to include vector polarizations.
533 if not isinstance(vectorval, bool):
534 raise TypeError(
"Must be boolean value")
536 self.
_vector_vector = vectorval
543 def scalar(self, scalarval):
545 Set whether to include scalar polarizations.
548 if not isinstance(scalarval, bool):
549 raise TypeError(
"Must be boolean value")
551 self.
_scalar_scalar = scalarval
560 Set whether to use LAL antenna response functions.
563 if not isinstance(val, bool):
564 raise TypeError(
"Must be a boolean value")
574 if not isinstance(pair, (list, tuple)):
575 raise TypeError(
'Pair must be a list or tuple')
578 raise ValueError(
'Pair must only contain two values')
581 if val.lower()
not in self.
parametersparameters:
582 raise ValueError(
'Parameter {} in pair not '
583 'recognized'.format(val))
585 self.
_lookup_pair_lookup_pair = [val.lower()
for val
in pair]
599 Set the pair of parameters to use for the look-up table.
605 for val, nbins
in zip(pair, [bins1, bins2]):
606 if not isinstance(nbins, int):
607 raise TypeError(
"Value must be an integer")
610 raise ValueError(
"There must be at least 2 bins")
614 vrange = VAR_RANGES[val]
620 vararray = np.arcsin(np.linspace(vrange[0], vrange[1], nbins))
622 vararray = np.linspace(vrange[0], vrange[1], nbins)
650 Set the 2d look-up table.
655 if not isinstance(val, bool):
656 raise TypeError(
"Value must be a boolean")
663 from scipy.interpolate
import RectBivariateSpline
665 raise ImportError(
"Cannot import scipy")
672 except AttributeError:
673 raise AttributeError(
"A time must be set for the look-up table")
764 return self.
responseresponse[
'plus']
773 self.
responseresponse[
'plus'] = resp
777 return self.
responseresponse[
'cross']
786 self.
responseresponse[
'cross'] = resp
842 Convert one-dimensional arrays into a mesh. The meshes dimensions are
843 ordered as right ascension, declination, psi, and time.
851 ramesh, decmesh, psimesh, timemesh = np.meshgrid(self.
rararara, self.
decdecdecdec,
867 if isinstance(val, np.ndarray):
869 self.
_ra_mesh_ra_mesh = val.squeeze()
887 if isinstance(val, np.ndarray):
907 if isinstance(val, np.ndarray):
922 if isinstance(val, np.ndarray):
933 Compute the detector response.
935 @param times: (array_like) an array of GPS times at which to compute
936 the response function. If not set the times set at initialization
937 of the class, or using the <tt>times</tt> property.
940 if times
is not None:
954 def _compute_response(self):
956 Compute antenna pattern.
964 for i
in range(len(self.
shapeshape)):
965 einsum_indices += eindices[i]
968 einstr1 =
'i{},j{}->ij{}'.format(*3*[einsum_indices])
969 einstr2 =
'ij{},ij->{}'.format(*2*[einsum_indices])
978 mm = np.einsum(einstr1, M, M)
979 mn = np.einsum(einstr1, M, N)
980 nm = np.einsum(einstr1, N, M)
981 nn = np.einsum(einstr1, N, N)
995 mq = np.einsum(einstr1, M, Q)
996 qm = np.einsum(einstr1, Q, M)
997 nq = np.einsum(einstr1, N, Q)
998 qn = np.einsum(einstr1, Q, N)
1004 qq = np.einsum(einstr1, Q, Q)
1009 def _compute_response_lal(self):
1011 Compute antenna pattern using LAL functions.
1016 slen = np.prod(self.
shapeshape)
if len(self.
shapeshape) != 0
else 1
1020 fp = np.zeros(self.
shapeshape)
1021 fc = np.zeros(self.
shapeshape)
1025 antenna_func = ComputeDetAMResponse
1034 for i
in range(slen):
1035 idxs = np.unravel_index(i, self.
shapeshape)
1046 antenna_func = ComputeDetAMResponseExtraModes
1055 fb = np.zeros(self.
shapeshape)
1056 fl = np.zeros(self.
shapeshape)
1057 fx = np.zeros(self.
shapeshape)
1058 fy = np.zeros(self.
shapeshape)
1060 for i
in range(slen):
1061 idxs = np.unravel_index(i, self.
shapeshape)
1088 def _compute_response_lookup(self):
1093 fp = np.zeros(lushape)
1094 fc = np.zeros(lushape)
1097 fb = np.zeros(lushape)
1098 fl = np.zeros(lushape)
1101 fx = np.zeros(lushape)
1102 fy = np.zeros(lushape)
1106 unsorted_idxs =
None
1109 pairs.append(self.
rararara)
1117 unsorted_idxs = np.argsort(np.argsort(mtimes))
1119 pairs.append(mtimes)
1125 sshape = tuple(sshape)
1128 pos = 4*[slice(
None)]
1137 fp[tuple(pos)] = self.
_lookup_func_lookup_func[
'plus'][i, j](*pairs).reshape(sshape)
1138 fc[tuple(pos)] = self.
_lookup_func_lookup_func[
'cross'][i, j](*pairs).reshape(sshape)
1141 fb[tuple(pos)] = self.
_lookup_func_lookup_func[
'b'][i, j](*pairs).reshape(sshape)
1142 fl[tuple(pos)] = self.
_lookup_func_lookup_func[
'l'][i, j](*pairs).reshape(sshape)
1145 fx[tuple(pos)] = self.
_lookup_func_lookup_func[
'x'][i, j](*pairs).reshape(sshape)
1146 fy[tuple(pos)] = self.
_lookup_func_lookup_func[
'y'][i, j](*pairs).reshape(sshape)
1148 if unsorted_idxs
is not None:
1150 self.
plusplusplusplus = fp[:, :, :, unsorted_idxs].squeeze()
1154 self.
bbbb = fb[:, :, :, unsorted_idxs].squeeze()
1155 self.
llll = fl[:, :, :, unsorted_idxs].squeeze()
1158 self.
xxxx = fx[:, :, :, unsorted_idxs].squeeze()
1159 self.
yyyy = fy[:, :, :, unsorted_idxs].squeeze()
1166 self.
bbbb = fb.squeeze()
1167 self.
llll = fl.squeeze()
1170 self.
xxxx = fx.squeeze()
1171 self.
yyyy = fy.squeeze()
1176 def __call__(self, times, ra=None, dec=None, psi=None, detector=None):
1178 Return the antenna response function as a dictionary.
1181 if detector
is not None:
def vector(self, vectorval)
Set whether to include vector polarizations.
def __init__(self, detector, ra, dec, psi=0., times=None, tensor=True, vector=False, scalar=False, use_lal=False, lookup=False, lookuppair=None, bins1=100, bins2=100)
Calculate the long-wavelength limit antenna response functions for a given ground-based gravitational...
def psi(self, psival)
Set the value of the gravitational wave polarization angle psi.
def ra(self, raval)
Set the right ascension.
def compute_response(self, times=None)
Compute the detector response.
def tensor(self, tensorval)
Set whether to include tensor polarizations.
def detector(self, det)
Set the detector for which to calculate the antenna response.
def lookup(self, val)
Set the 2d look-up table.
def _compute_response_lookup(self)
def __call__(self, times, ra=None, dec=None, psi=None, detector=None)
Return the antenna response function as a dictionary.
def _compute_response(self)
def scalar(self, scalarval)
Set whether to include scalar polarizations.
def times(self, timearr)
Set array of times and GPS times.
def dec(self, decval)
Set the declination.
def use_lal(self, val)
Set whether to use LAL antenna response functions.
def laldetector(self, det)
Set the lal.Detector.
def set_lookup_pair(self, pair=['psi', 'time'], bins1=100, bins2=100)
Set the pair of parameters to use for the look-up table.
def lookup_pair(self, pair)
def _compute_response_lal(self)