Skip to content
Snippets Groups Projects
Commit 703ced9c authored by Prathamesh Joshi's avatar Prathamesh Joshi
Browse files

gstlal_inspiral_plot_background: this commit removes the functionality of... 

gstlal_inspiral_plot_background: this commit removes the functionality of adding Rankingstat and RankingstatPDF objects for the same template ids. It also prevents loading all dist_stats files in memory at one time. Finally, a bug fix is to identify bin_id from the filename
parent 3cf6424b
No related branches found
No related tags found
1 merge request!361gstlal_inspiral_plot_background_fix
Pipeline #480480 failed
Stage: deps
Stage: level0
Stage: level1
Stage: level2
Stage: docker
Stage: test-gstlal
Stage: test-gstlal-full-build
Stage: test-gstlal-ugly
Stage: test-burst
Stage: test-inspiral
Stage: test-offline
Stage: docs
Stage: nightly-pages
Hide whitespace changes
Inline Side-by-side
gstlal-inspiral/bin/gstlal_inspiral_plot_background
+ 75
− 96
View file @ 703ced9c
......@@ -48,7 +48,6 @@ from optparse import OptionParser
import os
import sqlite3
import sys
from tqdm import tqdm
from ligo.lw import dbtables
......@@ -195,31 +194,13 @@ def parse_command_line():
#
def load_distributions(filenames, verbose = False):
# RankingStat objects go into a dictionary indexed by the template
# bank ID set they are for, with RankinStat objects for the same
# templates being summed. RankingStatPDF objects are summed (and
# therefore are required to be fore the same templates). these
# organizational choices reflect the current use cases for this
# program
rankingstats = {}
rankingstatpdf = None
for n, filename in enumerate(filenames, 1):
if verbose:
print("%d/%d:" % (n, len(filenames)), file=sys.stderr)
this_rankingstat, this_rankingstatpdf = far.parse_likelihood_control_doc(far.ligolw_utils.load_filename(filename, verbose = verbose, contenthandler = far.RankingStat.LIGOLWContentHandler))
if this_rankingstat is not None:
template_ids = this_rankingstat.template_ids
if template_ids in rankingstats:
rankingstats[template_ids] += this_rankingstat
else:
rankingstats[template_ids] = this_rankingstat
if this_rankingstatpdf is not None:
if rankingstatpdf is None:
rankingstatpdf = this_rankingstatpdf
else:
rankingstatpdf += this_rankingstatpdf
for rankingstat in tqdm(rankingstats.values(), desc = "Density estimation", disable = not verbose):
def load_distributions(filename, verbose = False):
# Summing up RankingStat objects and RankingStatPDF objects for
# the same template bank ID set is no longer done. It is assumed
# that the files given as input correspond to different SVD bins,
# otherwise their plots will have the same filenames.
rankingstat, rankingstatpdf = far.parse_likelihood_control_doc(far.ligolw_utils.load_filename(filename, verbose = verbose, contenthandler = far.RankingStat.LIGOLWContentHandler))
if rankingstat is not None:
if options.add_zerolag_to_background:
rankingstat.denominator.lnzerolagdensity = rankingstat.zerolag
rankingstat.finish()
......@@ -228,13 +209,12 @@ def load_distributions(filenames, verbose = False):
seg = (0, 0)
if rankingstatpdf is not None:
seg = rankingstatpdf.segments.extent()
elif rankingstats:
rankingstat = list(rankingstats.values())[0]
elif rankingstat:
try:
seg = rankingstat.segmentlists.extent_all()
except ValueError:
pass
return rankingstats, rankingstatpdf, seg
return rankingstat, rankingstatpdf, seg
def load_search_results(filenames, tmp_path = None, verbose = False):
......@@ -337,14 +317,6 @@ options, filenames = parse_command_line()
os.makedirs(options.output_dir, exist_ok=True)
#
# load input
#
rankingstats, rankingstatpdf, seg = load_distributions(filenames, verbose = options.verbose)
if options.database:
timeslide_ln_lr, zerolag_ln_lr, timeslide_sngls, zerolag_sngls = load_search_results(options.database, tmp_path = options.tmp_space, verbose = options.verbose)
if options.scatter_log_lambdas is not None:
......@@ -360,71 +332,78 @@ else:
#
for bin_index, rankingstat in enumerate(sorted(rankingstats.values(), key = lambda rankingstat: sorted(rankingstat.template_ids))):
# SNR and \chi^2
for instrument in rankingstat.instruments:
for snr_chi_type in ("background_pdf", "injection_pdf", "zero_lag_pdf", "LR"):
if instrument in options.event_snr_dict.keys():
fig = plotfar.plot_snr_chi_pdf(rankingstat, instrument, snr_chi_type, options.max_snr, sngls = sngls, event_snr = options.event_snr_dict[instrument], event_chisq = options.event_chisq_dict[instrument], event_bank_chisq = options.event_bank_chisq_dict[instrument])
else:
fig = plotfar.plot_snr_chi_pdf(rankingstat, instrument, snr_chi_type, options.max_snr, sngls = sngls)
if fig is None:
continue
plotname = datafind.T050017_filename(instrument, "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_%04d_%s_SNRCHI2" % (options.user_tag, bin_index, snr_chi_type.upper()), seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
for num, filename in enumerate(filenames):
if verbose:
print("%d/%d:" % (num, len(filenames)), file=sys.stderr)
rankingstat, rankingstatpdf, seg = load_distributions(filename, verbose = options.verbose)
if rankingstat is not None:
# FIXME bin index is inferred from filename
# instead, store bin_index in rankingstat object and load here
bin_index = int(filename.split('-')[1][:4])
# candidate rates
fig = plotfar.plot_rates(rankingstat)
plotname = datafind.T050017_filename("H1L1V1", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_%04d_RATES" % (options.user_tag, bin_index), seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
# SNR PDFs
if options.plot_snr_snr_pdfs:
# assume PDFs are identical for all template bank bins, and pick
# one set to plot
rankingstat = list(rankingstats.values())[0]
for (instruments, horizon_distances) in sorted(rankingstat.numerator.SNRPDF.snr_joint_pdf_cache.keys(), key = lambda x: sorted(x[1])):
# they're stored as a frozen set of quantized key/value
# pairs, need to unquantize them and get a dictionary back
horizon_distances = rankingstat.numerator.SNRPDF.quantized_horizon_distances(horizon_distances)
fig = plotfar.plot_snr_joint_pdf(rankingstat.numerator.SNRPDF, instruments, horizon_distances, rankingstat.min_instruments, options.max_snr, sngls = sngls)
if fig is not None:
plotname = datafind.T050017_filename(instruments, "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_SNR_PDF_%s" % (options.user_tag, "_".join(["%s_%s" % (k, horizon_distances[k]) for k in sorted(horizon_distances)]) ), seg, options.output_format, path = options.output_dir)
# SNR and \chi^2
for instrument in rankingstat.instruments:
for snr_chi_type in ("background_pdf", "injection_pdf", "zero_lag_pdf", "LR"):
if instrument in options.event_snr_dict.keys():
fig = plotfar.plot_snr_chi_pdf(rankingstat, instrument, snr_chi_type, options.max_snr, sngls = sngls, event_snr = options.event_snr_dict[instrument], event_chisq = options.event_chisq_dict[instrument], event_bank_chisq = options.event_bank_chisq_dict[instrument])
else:
fig = plotfar.plot_snr_chi_pdf(rankingstat, instrument, snr_chi_type, options.max_snr, sngls = sngls)
if fig is None:
continue
plotname = datafind.T050017_filename(instrument, "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_%04d_%s_SNRCHI2" % (options.user_tag, bin_index, snr_chi_type.upper()), seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
# candidate rates
fig = plotfar.plot_rates(rankingstat)
plotname = datafind.T050017_filename("H1L1V1", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_%04d_RATES" % (options.user_tag, bin_index), seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
# SNR PDFs
if options.plot_snr_snr_pdfs:
# assume PDFs are identical for all template bank bins, and plot
# only the first one
if num == 0:
for (instruments, horizon_distances) in sorted(rankingstat.numerator.SNRPDF.snr_joint_pdf_cache.keys(), key = lambda x: sorted(x[1])):
# they're stored as a frozen set of quantized key/value
# pairs, need to unquantize them and get a dictionary back
horizon_distances = rankingstat.numerator.SNRPDF.quantized_horizon_distances(horizon_distances)
fig = plotfar.plot_snr_joint_pdf(rankingstat.numerator.SNRPDF, instruments, horizon_distances, rankingstat.min_instruments, options.max_snr, sngls = sngls)
if fig is not None:
plotname = datafind.T050017_filename(instruments, "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_SNR_PDF_%s" % (options.user_tag, "_".join(["%s_%s" % (k, horizon_distances[k]) for k in sorted(horizon_distances)]) ), seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
# ranking statistic PDFs and CCDFs
if rankingstatpdf is not None:
for Title, which, NAME in (("Noise", "noise", "NOISE"), ("Signal", "signal", "SIGNAL")):
fig = plotfar.plot_likelihood_ratio_pdf(rankingstatpdf, (options.min_log_lambda, options.max_log_lambda), Title, which = which)
plotname = datafind.T050017_filename("COMBINED", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_%s_LIKELIHOOD_RATIO_PDF" % (options.user_tag, NAME), seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
# ranking statistic PDFs and CCDFs
if rankingstatpdf is not None:
for Title, which, NAME in (("Noise", "noise", "NOISE"), ("Signal", "signal", "SIGNAL")):
fig = plotfar.plot_likelihood_ratio_pdf(rankingstatpdf, (options.min_log_lambda, options.max_log_lambda), Title, which = which)
plotname = datafind.T050017_filename("COMBINED", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_%s_LIKELIHOOD_RATIO_PDF" % (options.user_tag, NAME), seg, options.output_format, path = options.output_dir)
if rankingstatpdf is not None and rankingstatpdf.zero_lag_lr_lnpdf.array.any():
fapfar = far.FAPFAR(rankingstatpdf.new_with_extinction())
if zerolag_ln_lr is not None:
xhi = max(zerolag_ln_lr + timeslide_ln_lr)[0]
xhi = 5. * math.ceil(xhi / 5.)
xhi = max(xhi, options.max_log_lambda)
else:
xhi = options.max_log_lambda
fig = plotfar.plot_likelihood_ratio_ccdf(fapfar, (options.min_log_lambda, xhi), observed_ln_likelihood_ratios = zerolag_ln_lr, is_open_box = True)
plotname = datafind.T050017_filename("COMBINED", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_NOISE_LIKELIHOOD_RATIO_CCDF_OPEN_BOX" % options.user_tag, seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
if rankingstatpdf is not None and rankingstatpdf.zero_lag_lr_lnpdf.array.any():
fapfar = far.FAPFAR(rankingstatpdf.new_with_extinction())
if zerolag_ln_lr is not None:
xhi = max(zerolag_ln_lr + timeslide_ln_lr)[0]
xhi = 5. * math.ceil(xhi / 5.)
xhi = max(xhi, options.max_log_lambda)
else:
xhi = options.max_log_lambda
fig = plotfar.plot_likelihood_ratio_ccdf(fapfar, (options.min_log_lambda, xhi), observed_ln_likelihood_ratios = zerolag_ln_lr, is_open_box = True)
plotname = datafind.T050017_filename("COMBINED", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_NOISE_LIKELIHOOD_RATIO_CCDF_OPEN_BOX" % options.user_tag, seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
fig = plotfar.plot_likelihood_ratio_ccdf(fapfar, (options.min_log_lambda, xhi), observed_ln_likelihood_ratios = timeslide_ln_lr, is_open_box = False)
plotname = datafind.T050017_filename("COMBINED", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_NOISE_LIKELIHOOD_RATIO_CCDF_CLOSED_BOX" % options.user_tag, seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
fig = plotfar.plot_likelihood_ratio_ccdf(fapfar, (options.min_log_lambda, xhi), observed_ln_likelihood_ratios = timeslide_ln_lr, is_open_box = False)
plotname = datafind.T050017_filename("COMBINED", "GSTLAL_INSPIRAL_PLOT_BACKGROUND_%s_NOISE_LIKELIHOOD_RATIO_CCDF_CLOSED_BOX" % options.user_tag, seg, options.output_format, path = options.output_dir)
if options.verbose:
print("writing %s" % plotname, file=sys.stderr)
fig.savefig(plotname)
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment