From 67c262595b2674cddaa8ae140d2fb645c6e416d4 Mon Sep 17 00:00:00 2001 From: Max Isi Date: Wed, 16 Oct 2024 17:27:28 -0400 Subject: [PATCH] detector wip --- src/jimgw/single_event/data.py | 35 ++++++++++++++++-------- src/jimgw/single_event/detector.py | 43 +++++++++++++++++------------- 2 files changed, 48 insertions(+), 30 deletions(-) diff --git a/src/jimgw/single_event/data.py b/src/jimgw/single_event/data.py index 071c86e7..2acd9531 100644 --- a/src/jimgw/single_event/data.py +++ b/src/jimgw/single_event/data.py @@ -1,3 +1,5 @@ +__include__ = ["Data", "PowerSpectrum"] + from abc import ABC, abstractmethod import jax @@ -25,8 +27,6 @@ "V1": "https://dcc.ligo.org/public/0169/P2000251/001/O3-V1_sensitivity_strain_asd.txt", } -_DEF_GWPY_KWARGS = {"cache": True} - class Data(ABC): """ @@ -82,10 +82,10 @@ def has_fd(self) -> bool: """Whether the Fourier domain data has been computed.""" return bool(np.any(self.fd)) - def __init__(self, td: Float[Array, " n_time"], - delta_t: float, + def __init__(self, td: Float[Array, " n_time"] = jnp.array([]), + delta_t: float = 0., epoch: float = 0., - name: Optional[str] = None, + name: str = '', window: Optional[Float[Array, " n_time"]] = None)\ -> None: """Initialize the data class. @@ -113,6 +113,13 @@ def __init__(self, td: Float[Array, " n_time"], self.window = window self.name = name or '' + def __repr__(self): + return f"{self.__class__.__name__}(name='{self.name}', delta_t={self.delta_t}, epoch={self.epoch})" + + def __bool__(self) -> bool: + """Check if the data is empty.""" + return len(self.td) > 0 + def set_tukey_window(self, alpha: float = 0.2) -> None: """Create a Tukey window on the data; the window is stored in the :attr:`window` attribute and only applied when FFTing the data. @@ -201,7 +208,6 @@ def from_gwosc(cls, Whether to cache the data (default: True) **kws: dict, optional Keyword arguments for `gwpy.timeseries.TimeSeries.fetch_open_data` - defaults to {} """ duration = gps_end_time - gps_start_time logging.info(f"Fetching {duration} s of {ifo} data from GWOSC " @@ -211,9 +217,6 @@ def from_gwosc(cls, cache=cache, **kws) return cls(data_td.value, data_td.dt.value, data_td.epoch.value, ifo) - from_gwosc.__doc__ = from_gwosc.__doc__.format(_DEF_GWPY_KWARGS) - - class PowerSpectrum(ABC): name: str values: Float[Array, " n_freq"] @@ -234,13 +237,23 @@ def duration(self) -> Float: """Duration of the data in seconds.""" return 1 / self.delta_f - def __init__(self, values: Float[Array, " n_freq"], - frequencies: Float[Array, " n_freq"], + @property + def sampling_frequency(self) -> Float: + """Sampling frequency of the data in Hz.""" + return self.frequencies[-1] * 2 + + def __init__(self, values: Float[Array, " n_freq"] = jnp.array([]), + frequencies: Float[Array, " n_freq"] = jnp.array([]), name: Optional[str] = None) -> None: self.values = values self.frequencies = frequencies + assert len(self.values) == len(self.frequencies), \ + "Values and frequencies must have the same length" self.name = name or '' + def __repr__(self) -> str: + return f"{self.__class__.__name__}(name='{self.name}', frequencies={self.frequencies})" + def slice(self, f_min: float, f_max: float) -> \ tuple[Float[Array, " n_sample"], Float[Array, " n_sample"]]: """Slice the power spectrum. diff --git a/src/jimgw/single_event/detector.py b/src/jimgw/single_event/detector.py index c4c64de8..fce30d72 100644 --- a/src/jimgw/single_event/detector.py +++ b/src/jimgw/single_event/detector.py @@ -9,6 +9,7 @@ from beartype import beartype as typechecker from scipy.interpolate import interp1d from scipy.signal.windows import tukey +from . import data as jd from jimgw.constants import C_SI, EARTH_SEMI_MAJOR_AXIS, EARTH_SEMI_MINOR_AXIS from jimgw.single_event.wave import Polarization @@ -99,10 +100,8 @@ class GroundBased2G(Detector): Array of noise power spectral density. """ polarization_mode: list[Polarization] - frequencies: Float[Array, " n_sample"] - data: Float[Array, " n_sample"] - psd: Float[Array, " n_sample"] - epoch: Float = 0 + data: jd.Data + psd: jd.PowerSpectrum latitude: Float = 0 longitude: Float = 0 @@ -115,22 +114,26 @@ class GroundBased2G(Detector): def __repr__(self) -> str: return f"{self.__class__.__name__}({self.name})" - def __init__(self, name: str, **kwargs) -> None: + def __init__(self, name: str, latitude: float = 0, longitude: float = 0, + elevation: float = 0, xarm_azimuth: float = 0, + yarm_azimuth: float = 0, xarm_tilt: float = 0, + yarm_tilt: float = 0, modes: str = "pc"): self.name = name - self.latitude = kwargs.get("latitude", 0) - self.longitude = kwargs.get("longitude", 0) - self.elevation = kwargs.get("elevation", 0) - self.xarm_azimuth = kwargs.get("xarm_azimuth", 0) - self.yarm_azimuth = kwargs.get("yarm_azimuth", 0) - self.xarm_tilt = kwargs.get("xarm_tilt", 0) - self.yarm_tilt = kwargs.get("yarm_tilt", 0) - modes = kwargs.get("mode", "pc") + self.latitude = latitude + self.longitude = longitude + self.elevation = elevation + self.xarm_azimuth = xarm_azimuth + self.yarm_azimuth = yarm_azimuth + self.xarm_tilt = xarm_tilt + self.yarm_tilt = yarm_tilt self.polarization_mode = [Polarization(m) for m in modes] - self.frequencies = jnp.array([]) - self.data = jnp.array([]) - self.psd = jnp.array([]) + self.data = jd.Data() + + # self.frequencies = jnp.array([]) + # self.data = jnp.array([]) + # self.psd = jnp.array([]) @staticmethod def _get_arm( @@ -312,6 +315,8 @@ def load_data( self.psd = psd[(freq > f_min) & (freq < f_max)] load_data.__doc__ = load_data.__doc__.format(_DEF_GWPY_KWARGS) + # def load_data(self, data: ) + def fd_response( self, frequency: Float[Array, " n_sample"], @@ -494,7 +499,7 @@ def load_psd( xarm_tilt=-6.195e-4, yarm_tilt=1.25e-5, elevation=142.554, - mode="pc", + modes="pc", ) L1 = GroundBased2G( @@ -506,7 +511,7 @@ def load_psd( xarm_tilt=0, yarm_tilt=0, elevation=-6.574, - mode="pc", + modes="pc", ) V1 = GroundBased2G( @@ -518,7 +523,7 @@ def load_psd( xarm_tilt=0, yarm_tilt=0, elevation=51.884, - mode="pc", + modes="pc", ) detector_preset = {