Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add ability to cache detector response to marginalize time model and us… #4806

Merged
merged 5 commits into from
Jul 22, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 4 additions & 1 deletion examples/inference/margtime/margtime.ini
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,11 @@
name = marginalized_time
low-frequency-cutoff = 30.0

# This is the sample rate used for the model and marginalization
sample_rate = 4096

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I updated the example to make explicit use of this feature.

marginalize_vector_params = tc, ra, dec, polarization
marginalize_vector_samples = 500
marginalize_vector_samples = 2000

; You shouldn't use phase marginalization if the approximant has
; higher-order modes
Expand Down
3 changes: 2 additions & 1 deletion examples/inference/margtime/run.sh
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
OMP_NUM_THREADS=1 pycbc_inference \
--config-file `dirname "$0"`/margtime.ini \
--nprocesses 2 \
--nprocesses 1 \
--processing-scheme mkl \
--output-file marg_150914.hdf \
--seed 0 \
Expand All @@ -23,4 +23,5 @@ pycbc_inference_plot_posterior \
"primary_mass(mass1, mass2) / (1 + redshift(distance)):srcmass1" \
"secondary_mass(mass1, mass2) / (1 + redshift(distance)):srcmass2" \
ra dec tc inclination coa_phase polarization distance \
--vmin 23.2 \
--z-arg snr
43 changes: 33 additions & 10 deletions pycbc/inference/models/marginalized_gaussian_noise.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
"""

import itertools
import logging
import numpy
from scipy import special

Expand Down Expand Up @@ -207,8 +208,10 @@ class MarginalizedTime(DistMarg, BaseGaussianNoise):
def __init__(self, variable_params,
data, low_frequency_cutoff, psds=None,
high_frequency_cutoff=None, normalize=False,
sample_rate=None,
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What's the advantage of doing this instead of just increasing the sample rate of the data? Also, this might need some documentation, or maybe change the name of the argument to make it clearer what it is. People are going to think this is the data sample rate otherwise.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Increasing the sample rate of the data means that you need to do data conditioning and parts of the generation at a higher sample rate. That's only useful if you actually care about the frequency band at higher frequencies. If you just need better timing resolution that would be extremely slow and memory intensive.

I think the separation is clear as the data sample rate is in the data section. This is only about the model itself. BTW, we already have this distinction for the heterodyne likelihoods when doing sky marginalization. This was just advancing the same interface / feature to the more regular model

**kwargs):

self.sample_rate = float(sample_rate)
self.kwargs = kwargs
variable_params, kwargs = self.setup_marginalization(
variable_params,
Expand Down Expand Up @@ -241,6 +244,14 @@ def __init__(self, variable_params,

self.dets = {}

if sample_rate is not None:
for ifo in self.data:
if self.sample_rate < self.data[ifo].sample_rate:
raise ValueError("Model sample rate was set less than the"
" data. ")
logging.info("Using %s sample rate for marginalization",
sample_rate)

def _nowaveform_loglr(self):
"""Convenience function to set loglr values if no waveform generated.
"""
Expand Down Expand Up @@ -296,8 +307,15 @@ def _loglr(self):
hp[self._kmin[det]:kmax] *= self._weight[det][slc]
hc[self._kmin[det]:kmax] *= self._weight[det][slc]

hp.resize(len(self._whitened_data[det]))
hc.resize(len(self._whitened_data[det]))
# Use a higher sample rate if requested
if self.sample_rate is not None:
tlen = int(round(self.sample_rate *
self.whitened_data[det].duration))
flen = tlen // 2 + 1
hp.resize(flen)
hc.resize(flen)
self._whitened_data[det].resize(flen)

cplx_hpd[det], _, _ = matched_filter_core(
hp,
self._whitened_data[det],
Expand Down Expand Up @@ -325,15 +343,20 @@ def _loglr(self):
for det in wfs:
if det not in self.dets:
self.dets[det] = Detector(det)
fp, fc = self.dets[det].antenna_pattern(
params['ra'],
params['dec'],
params['polarization'],
params['tc'])
dt = self.dets[det].time_delay_from_earth_center(params['ra'],
params['dec'],
params['tc'])

if self.precalc_antenna_factors:
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How is precalc_antenna_factors set to True in a run? Is this an argument for the model? I couldn't find where it would be set.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's already handled by the marginalization parent class. The features were implemented first for other cases, and this just extends the use to the marginalized time model.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok, but I was asking how you set that in the config file and how that's passed to the model initialization. I couldn't figure it out from looking at the code, although I didn't dig too deep.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah, it's automatically populated if you do sky marginalization so it's the same set of settings really. These factors are calculated in any case when the sky grid is set up.

fp, fc, dt = self.get_precalc_antenna_factors(det)
else:
fp, fc = self.dets[det].antenna_pattern(
params['ra'],
params['dec'],
params['polarization'],
params['tc'])
dt = self.dets[det].time_delay_from_earth_center(params['ra'],
params['dec'],
params['tc'])
dtc = params['tc'] + dt

cplx_hd = fp * cplx_hpd[det].at_time(dtc,
interpolate='quadratic')
cplx_hd += fc * cplx_hcd[det].at_time(dtc,
Expand Down
Loading