From c2b668fdb36e1af2d26c52494b675d15bf867c01 Mon Sep 17 00:00:00 2001 From: Chad Hanna <chad.hanna@ligo.org> Date: Wed, 17 Oct 2018 19:03:11 -0700 Subject: [PATCH] inspiral_lr: allow the mass model filename to be encoded in XML --- gstlal-inspiral/python/stats/inspiral_lr.py | 30 ++++++++++++++++++--- 1 file changed, 26 insertions(+), 4 deletions(-) diff --git a/gstlal-inspiral/python/stats/inspiral_lr.py b/gstlal-inspiral/python/stats/inspiral_lr.py index a66cb03afc..d85abef0bb 100644 --- a/gstlal-inspiral/python/stats/inspiral_lr.py +++ b/gstlal-inspiral/python/stats/inspiral_lr.py @@ -86,7 +86,7 @@ class LnLRDensity(snglcoinc.LnLRDensity): # SNR, \chi^2 binning definition snr_chi_binning = rate.NDBins((rate.ATanLogarithmicBins(2.6, 26., 300), rate.ATanLogarithmicBins(.001, 0.2, 280))) - def __init__(self, template_ids, instruments, delta_t, min_instruments = 2): + def __init__(self, template_ids, instruments, delta_t, min_instruments = 2, mass_model_filename = None): # # check input # @@ -105,6 +105,7 @@ class LnLRDensity(snglcoinc.LnLRDensity): # self.template_ids = frozenset(template_ids) if template_ids is not None else template_ids + self.mass_model_filename = mass_model_filename self.instruments = frozenset(instruments) self.min_instruments = min_instruments self.delta_t = delta_t @@ -120,6 +121,8 @@ class LnLRDensity(snglcoinc.LnLRDensity): # that one or the other is None to make it possible to # construct generic seed objects providing initialization # data for the ranking statistics. + if self.mass_model_filename is not None and other.mass_model_filename is not None and other.mass_model_filename != self.mass_model_filename: + raise ValueError("incompatible mass model file names") if self.template_ids is not None and other.template_ids is not None and self.template_ids != other.template_ids: raise ValueError("incompatible template IDs") if self.instruments != other.instruments: @@ -131,6 +134,8 @@ class LnLRDensity(snglcoinc.LnLRDensity): if self.template_ids is None and other.template_ids is not None: self.template_ids = frozenset(other.template_ids) + if self.mass_model_filename is None and other.mass_model_filename is not None: + self.mass_model_filename = other.mass_model_filename for key, lnpdf in self.densities.items(): lnpdf += other.densities[key] @@ -171,7 +176,7 @@ class LnLRDensity(snglcoinc.LnLRDensity): self.densities["%s_snr_chi" % event.ifo].count[event.snr, event.chisq / event.snr**2.] += 1.0 def copy(self): - new = type(self)(self.template_ids, self.instruments, self.delta_t, self.min_instruments) + new = type(self)(self.template_ids, self.instruments, self.delta_t, self.min_instruments, self.mass_model_filename) for key, lnpdf in self.densities.items(): new.densities[key] = lnpdf.copy() return new @@ -237,6 +242,9 @@ class LnLRDensity(snglcoinc.LnLRDensity): xml = super(LnLRDensity, self).to_xml(name) # FIXME: switch to .from_pyvalue() when it can accept None xml.appendChild(ligolw_param.Param.build(u"template_ids", u"lstring", ",".join("%d" % template_id for template_id in sorted(self.template_ids)) if self.template_ids is not None else None)) + # FIXME this is not an ideal way to get only one into the file + if self.mass_model_filename is not None: + xml.appendChild(ligolw_param.Param.from_pyvalue(u"mass_model_filename", self.mass_model_filename)) xml.appendChild(ligolw_param.Param.from_pyvalue(u"instruments", lsctables.instrumentsproperty.set(self.instruments))) xml.appendChild(ligolw_param.Param.from_pyvalue(u"min_instruments", self.min_instruments)) xml.appendChild(ligolw_param.Param.from_pyvalue(u"delta_t", self.delta_t)) @@ -251,13 +259,25 @@ class LnLRDensity(snglcoinc.LnLRDensity): def from_xml(cls, xml, name): xml = cls.get_xml_root(xml, name) template_ids = ligolw_param.get_pyvalue(xml, u"template_ids") + # FIXME find a better way to do this: Allow a file to not have + # a mass model filename. This is can happen when e.g., gstlal + # inspiral produces a ranking stat file which just records but + # doesn't evaluate LRs. Then you marginalize it in with a + # create prior diststats file which does have the mass model. + # The iadd method makes sure that the mass model is unique or + # that it doesn't exist. + try: + mass_model_filename = ligolw_param.get_pyvalue(xml, u"mass_model_filename") + except ValueError: + mass_model_filename = None if template_ids is not None: template_ids = frozenset(int(i) for i in template_ids.split(",")) self = cls( template_ids = template_ids, instruments = lsctables.instrumentsproperty.get(ligolw_param.get_pyvalue(xml, u"instruments")), delta_t = ligolw_param.get_pyvalue(xml, u"delta_t"), - min_instruments = ligolw_param.get_pyvalue(xml, u"min_instruments") + min_instruments = ligolw_param.get_pyvalue(xml, u"min_instruments"), + mass_model_filename = mass_model_filename ) for key in self.densities: self.densities[key] = self.densities[key].from_xml(xml, key) @@ -287,8 +307,10 @@ class LnSignalDensity(LnLRDensity): # network self.horizon_history = horizonhistory.HorizonHistories((instrument, horizonhistory.NearestLeafTree()) for instrument in self.instruments) self.population_model_file = population_model_file - self.population_model = inspiral_intrinsics.SourcePopulationModel(self.template_ids, filename = self.population_model_file) + # source population model + # FIXME: introduce a mechanism for selecting the file + self.population_model = inspiral_intrinsics.SourcePopulationModel(self.template_ids, filename = self.mass_model_filename) self.InspiralExtrinsics = inspiral_extrinsics.InspiralExtrinsics(self.min_instruments) def __call__(self, segments, snrs, chi2s_over_snr2s, phase, dt, template_id): -- GitLab