Skip to content
Snippets Groups Projects
Commit 507ae6b6 authored by Jacob Peoples's avatar Jacob Peoples Committed by Kipp Cannon
Browse files

gstlal.far: move .add_background_prior(), .add_foreground_prior()

from DistributionStats class to CoincParamsDistributions class
parent a8920b71
No related branches found
No related tags found
No related merge requests found
......@@ -145,11 +145,11 @@ coincparamsdistributions, likelihood_seglists = Far.distribution_stats, Far.live
# calculate injections before writing to disk
if options.synthesize_injections != 0:
coincparamsdistributions.add_foreground_prior(n = options.synthesize_injections, instruments = Far.trials_table.get_sngl_ifos(), verbose = options.verbose)
coincparamsdistributions.raw_distributions.add_foreground_prior(n = options.synthesize_injections, instruments = Far.trials_table.get_sngl_ifos(), verbose = options.verbose)
# add a uniform prior to background, by default 0 is added so it has no effect
if options.background_prior != 0:
coincparamsdistributions.add_background_prior(n = options.background_prior, instruments = Far.trials_table.get_sngl_ifos(), verbose = options.verbose)
coincparamsdistributions.raw_distributions.add_background_prior(n = options.background_prior, instruments = Far.trials_table.get_sngl_ifos(), verbose = options.verbose)
if options.verbose:
print >>sys.stderr, "computing event densities ..."
......
......@@ -100,9 +100,9 @@ options, filenames = parse_command_line()
coincparamsdistributions = far.DistributionsStats()
if options.background_prior != 0:
coincparamsdistributions.add_background_prior(options.background_prior, instruments = options.instrument, verbose = options.verbose)
coincparamsdistributions.raw_distributions.add_background_prior(options.background_prior, instruments = options.instrument, verbose = options.verbose)
coincparamsdistributions.add_foreground_prior(verbose = options.verbose, n = options.synthesize_injection_count, instruments = options.instrument)
coincparamsdistributions.raw_distributions.add_foreground_prior(verbose = options.verbose, n = options.synthesize_injection_count, instruments = options.instrument)
# no segments, this file is independent of time
FAR = far.LocalRankingData(segments.segment(None, None), coincparamsdistributions)
......
......@@ -417,35 +417,8 @@ class ThincaCoincParamsDistributions(snglcoinc.CoincParamsDistributions):
instruments.discard("H2")
return dict(("%s_snr_chi" % event.ifo, (event.snr, event.chisq / event.snr**2)) for event in events if event.ifo in instruments)
#
# Paramter Distributions
#
class DistributionsStats(object):
"""
A class used to populate a ThincaCoincParamsDistribution instance using
event parameter data.
"""
def __init__(self):
self.raw_distributions = ThincaCoincParamsDistributions()
self.smoothed_distributions = ThincaCoincParamsDistributions()
self.likelihood_pdfs = {}
self.target_length = 1000
def __add__(self, other):
out = type(self)()
out.raw_distributions += self.raw_distributions
out.raw_distributions += other.raw_distributions
#FIXME do we also add the smoothed distributions??
return out
def add_single(self, event):
self.raw_distributions.add_background(self.raw_distributions.coinc_params((event,), None))
def add_background_prior(self, n = 1., transition = 10., instruments = None, prefactors_range = (1.0, 10.0), df = 40, verbose = False):
for param, binarr in self.raw_distributions.background_rates.items():
for param, binarr in self.background_rates.items():
new_binarr = rate.BinnedArray(binarr.bins)
# FIXME only works if there is a 1-1 relationship between params and instruments
instrument = param.split("_")[0]
......@@ -469,51 +442,14 @@ class DistributionsStats(object):
# add to raw counts
binarr += new_binarr
# FIXME, an adhoc way of adding glitches, us a signal distribution with bad matches
# FIXME if we keep this, make a function that can be reused for both signal and background
# FIXME: for maintainability, this should be modified to
# use the .add_injection() method of the .raw_distributions
# attribute, but that will slow this down
pfs = numpy.linspace(prefactors_range[0], prefactors_range[1], 10)
for param, binarr in self.raw_distributions.background_rates.items():
new_binarr = rate.BinnedArray(binarr.bins)
# FIXME only works if there is a 1-1 relationship between params and instruments
instrument = param.split("_")[0]
# save some computation if we only requested certain instruments
if instruments is not None and instrument not in instruments:
continue
if verbose:
print >> sys.stderr, "synthesizing background for %s" % param
# Custom handle the first and last over flow bins
snrs = new_binarr.bins[0].centres()
snrs[0] = snrs[1] * .9
snrs[-1] = snrs[-2] * 1.1
chi2_over_snr2s = new_binarr.bins[1].centres()
chi2_over_snr2s[0] = chi2_over_snr2s[1] * .9
chi2_over_snr2s[-1] = chi2_over_snr2s[-2] * 1.1
for i, snr in enumerate(snrs):
for j, chi2_over_snr2 in enumerate(chi2_over_snr2s):
chisq = chi2_over_snr2 * snr**2 * df # We record the reduced chi2
dist = 0
for pf in pfs:
nc = pf * snr**2
v = stats.ncx2.pdf(chisq, df, nc)
if numpy.isfinite(v):
dist += v
dist *= (snr / snrs[0])**-4
if numpy.isfinite(dist):
new_binarr[snr, chi2_over_snr2] += dist
# normalize to the requested count
new_binarr.array *= n / new_binarr.array.sum()
# add to raw counts
binarr += new_binarr
# FIXME, an adhoc way of adding glitches, use a signal distribution with bad matches
self.add_foreground_prior(n = n, prefactors_range = prefactors_range, df = df, instruments = instruments, verbose = verbose)
def add_foreground_prior(self, n = 1., prefactors_range = (0.0, 0.10), df = 40, instruments = None, verbose = False):
# FIXME: for maintainability, this should be modified to
# use the .add_injection() method of the .raw_distributions
# attribute, but that will slow this down
# use the .add_injection() method of self, but that will slow this down
pfs = numpy.linspace(prefactors_range[0], prefactors_range[1], 10)
for param, binarr in self.raw_distributions.injection_rates.items():
for param, binarr in self.injection_rates.items():
new_binarr = rate.BinnedArray(binarr.bins)
# FIXME only works if there is a 1-1 relationship between params and instruments
instrument = param.split("_")[0]
......@@ -521,7 +457,7 @@ class DistributionsStats(object):
if instruments is not None and instrument not in instruments:
continue
if verbose:
print >> sys.stderr, "synthesizing injections for %s" % param
print >> sys.stderr, "synthesizing background/injections for %s" % param
# Custom handle the first and last over flow bins
snrs = new_binarr.bins[0].centres()
snrs[0] = snrs[1] * .9
......@@ -546,6 +482,33 @@ class DistributionsStats(object):
# add to raw counts
binarr += new_binarr
#
# Paramter Distributions
#
class DistributionsStats(object):
"""
A class used to populate a ThincaCoincParamsDistribution instance using
event parameter data.
"""
def __init__(self):
self.raw_distributions = ThincaCoincParamsDistributions()
self.smoothed_distributions = ThincaCoincParamsDistributions()
self.likelihood_pdfs = {}
self.target_length = 1000
def __add__(self, other):
out = type(self)()
out.raw_distributions += self.raw_distributions
out.raw_distributions += other.raw_distributions
#FIXME do we also add the smoothed distributions??
return out
def add_single(self, event):
self.raw_distributions.add_background(self.raw_distributions.coinc_params((event,), None))
def finish(self, verbose = False):
self.smoothed_distributions = self.raw_distributions.copy()
#self.smoothed_distributions.finish(verbose = verbose)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment