-
Notifications
You must be signed in to change notification settings - Fork 348
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Add a module and class for handling sky grids for `pycbc_multi_inspir…
…al` (#4741) * Add basic sky grid module and class * Use SkyGrid in multi_inspiral * Run black * More cleanup * Fix errors * Make detectors and ref time properties of the instance (to be continued…) * Use --instruments and --trigger-time * CLI improvements * Fix error * Add angular units and convention
- Loading branch information
1 parent
1641227
commit fc4c210
Showing
2 changed files
with
173 additions
and
78 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,136 @@ | ||
| """Functionality for handling grids of points in the sky for coherent SNR | ||
| calculation via `pycbc_multi_inspiral`. The main operation to be performed on | ||
| these points is calculating the antenna pattern functions and time delays from | ||
| the Earth center for a network of detectors. | ||
| """ | ||
|
|
||
| import numpy as np | ||
| import h5py | ||
|
|
||
| from pycbc.detector import Detector | ||
|
|
||
|
|
||
| class SkyGrid: | ||
| def __init__(self, ra, dec, detectors, ref_gps_time): | ||
| """Initialize a sky grid from a list of RA/dec coordinates. | ||
| Parameters | ||
| ---------- | ||
| ra: iterable of floats | ||
| Right ascensions for each point in radians, in the interval [0,2π). | ||
| dec: iterable of floats | ||
| Declination for each point in radians, where π/2 is the North pole, | ||
| -π/2 is the South pole, and 0 is the celestial equator. | ||
| detectors: iterable of str | ||
| List of detector names associated with the sky grid, typically the | ||
| detectors for which the grid has been placed or is to be placed. | ||
| The detectors will be used when calculating the antenna pattern | ||
| functions and time delays from Earth center. | ||
| ref_gps_time: float | ||
| Reference GPS time associated with the sky grid. This will be used | ||
| when calculating the antenna pattern functions and time delays from | ||
| Earth center. | ||
| """ | ||
| # We store the points in a 2D array internally, first dimension runs | ||
| # over the list of points, second dimension is RA/dec. | ||
| # Question: should we use Astropy sky positions instead? | ||
| self.positions = np.vstack([ra, dec]).T | ||
| self.detectors = sorted(detectors) | ||
| self.ref_gps_time = ref_gps_time | ||
|
|
||
| def __len__(self): | ||
| """Returns the number of points in the sky grid.""" | ||
| return self.positions.shape[0] | ||
|
|
||
| def __getitem__(self, index): | ||
| """Returns the coordinates of a single point in the grid.""" | ||
| return self.positions[index] | ||
|
|
||
| @property | ||
| def ras(self): | ||
| """Returns all right ascensions in radians, in the interval [0,2π).""" | ||
| return self.positions[:, 0] | ||
|
|
||
| @property | ||
| def decs(self): | ||
| """Returns all declinations in radians, where π/2 is the North pole, | ||
| -π/2 is the South pole, and 0 is the celestial equator.""" | ||
| return self.positions[:, 1] | ||
|
|
||
| @classmethod | ||
| def from_cli(cls, cli_parser, cli_args): | ||
| """Initialize a sky grid from command-line interface, via argparse | ||
| objects. | ||
| """ | ||
| if cli_args.sky_grid is not None: | ||
| if cli_args.ra is not None or cli_args.dec is not None: | ||
| cli_parser.error( | ||
| 'Please provide either a sky grid via --sky-grid or a ' | ||
| 'single sky position via --ra and --dec, not both' | ||
| ) | ||
| return cls.read_from_file(cli_args.sky_grid) | ||
| if cli_args.ra is not None and cli_args.dec is not None: | ||
| return cls( | ||
| [cli_args.ra], | ||
| [cli_args.dec], | ||
| cli_args.instruments, | ||
| cli_args.trigger_time | ||
| ) | ||
| cli_parser.error( | ||
| 'Please specify a sky grid via --sky-grid or a single sky ' | ||
| 'position via --ra and --dec' | ||
| ) | ||
|
|
||
| @classmethod | ||
| def read_from_file(cls, path): | ||
| """Initialize a sky grid from a given HDF5 file.""" | ||
| with h5py.File(path, 'r') as hf: | ||
| ra = hf['ra'][:] | ||
| dec = hf['dec'][:] | ||
| detectors = hf.attrs['detectors'] | ||
| ref_gps_time = hf.attrs['ref_gps_time'] | ||
| return cls(ra, dec, detectors, ref_gps_time) | ||
|
|
||
| def write_to_file(self, path): | ||
| """Writes a sky grid to an HDF5 file.""" | ||
| with h5py.File(path, 'w') as hf: | ||
| hf['ra'] = self.ras | ||
| hf['dec'] = self.decs | ||
| hf.attrs['detectors'] = self.detectors | ||
| hf.attrs['ref_gps_time'] = self.ref_gps_time | ||
|
|
||
| def calculate_antenna_patterns(self): | ||
| """Calculate the antenna pattern functions at each point in the grid | ||
| for the list of GW detectors specified at instantiation. Return a dict, | ||
| keyed by detector name, whose items are 2-dimensional Numpy arrays. | ||
| The first dimension of these arrays runs over the sky grid, and the | ||
| second dimension runs over the plus and cross polarizations. | ||
| """ | ||
| result = {} | ||
| for det_name in self.detectors: | ||
| det = Detector(det_name) | ||
| result[det_name] = np.empty((len(self), 2)) | ||
| for i, (ra, dec) in enumerate(self): | ||
| result[det_name][i] = det.antenna_pattern( | ||
| ra, dec, 0, t_gps=self.ref_gps_time | ||
| ) | ||
| return result | ||
|
|
||
| def calculate_time_delays(self): | ||
| """Calculate the time delays from the Earth center to each GW detector | ||
| specified at instantiation, for each point in the grid. Return a dict, | ||
| keyed by detector name, whose items are 1-dimensional Numpy arrays | ||
| containing the time delays for each sky point. | ||
| """ | ||
| result = {} | ||
| for det_name in self.detectors: | ||
| det = Detector(det_name) | ||
| result[det_name] = np.empty(len(self)) | ||
| for i, (ra, dec) in enumerate(self): | ||
| result[det_name][i] = det.time_delay_from_earth_center( | ||
| ra, dec, t_gps=self.ref_gps_time | ||
| ) | ||
| return result | ||
|
|
||
|
|
||
| __all__ = ['SkyGrid'] |