Skip to content
Snippets Groups Projects
Commit 533ca235 authored by Soichiro Kuwahara's avatar Soichiro Kuwahara
Browse files

gstlal-burst: add gstlal_cherenkov_calc_likelihood

- It calculates the likelihood ratio for zero lag candidates
parent 5a17d2d0
No related branches found
No related tags found
No related merge requests found
dist_bin_SCRIPTS = \
gstlal_cherenkov_burst \
gstlal_cherenkov_calc_likelihood \
gstlal_cherenkov_calc_rank_pdfs \
gstlal_cherenkov_inj \
gstlal_cherenkov_plot_rankingstat \
......
#!/usr/bin/env python3
#
# Copyright (C) 2010--2015 Kipp Cannon, Chad Hanna
# Copyright (C) 2021 Soichiro Kuwawhara
#
# This program is free software; you can redistribute it and/or modify it
# under the terms of the GNU General Public License as published by the
# Free Software Foundation; either version 2 of the License, or (at your
# option) any later version.
#
# This program is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
# Public License for more details.
#
# You should have received a copy of the GNU General Public License along
# with this program; if not, write to the Free Software Foundation, Inc.,
# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
### A program to compute the noise probability distributions of likehood ratios for inspiral triggers
#
# =============================================================================
#
# Preamble
#
# =============================================================================
#
from optparse import OptionParser
import sys
from ligo.lw import ligolw
from ligo.lw import lsctables
from ligo.lw import utils as ligolw_utils
from ligo.lw.utils import process as ligolw_process
from lal.utils import CacheEntry
from gstlal import cherenkov
from gstlal.cherenkov import rankingstat as cherenkov_rankingstat
__author__ = "Soichiro Kuwara <soichiro.kuwahara@ligo.org>"
#
# =============================================================================
#
# Command Line
#
# =============================================================================
#
def parse_command_line():
parser = OptionParser(
description = "Rankngstat calculation program for Cherenkov burst search."
)
parser.add_option("--candidate-cache", metavar = "filename", help = "Also load the candidates from files listed in this LAL cache. See lalapps_path2cache for information on how to produce a LAL cache file.")
parser.add_option("-v", "--verbose", action = "store_true", help = "Be verbose.")
parser.add_option("--ranking-stat", metavar = "filename", help = "Load ranking statistic data from this file.")
options, urls = parser.parse_args()
paramdict = options.__dict__.copy()
if options.candidate_cache is not None:
urls += [CacheEntry(line).url for line in open(options.candidate_cache)]
if not urls:
raise ValueError("must provide some candidate files")
if options.ranking_stat is None:
raise ValueError("must set --ranking-stat")
return options, urls, paramdict
#
# =============================================================================
#
# Main
#
# =============================================================================
#
#
# command line
#
options, urls, paramdict = parse_command_line()
#
# load parameter distribution data
#
rankingstat = cherenkov_rankingstat.marginalize_rankingstat_urls([options.ranking_stat], verbose = options.verbose)
#
# invoke .finish() to apply density estimation kernels and correct the
# normalization.
#
rankingstat.finish()
#
# load zero-lag candidates and histogram their ranking statistics
#
for n, url in enumerate(urls, 1):
if options.verbose:
print("%d/%d: " % (n, len(urls)), end = "", file = sys.stderr)
xmldoc = ligolw_utils.load_url(url, contenthandler = cherenkov_rankingstat.LIGOLWContentHandler, verbose = options.verbose)
process = ligolw_process.register_to_xmldoc(xmldoc, "gstlal_cherenkov_calc_likelihood", paramdict = paramdict)
time_slide_index = lsctables.TimeSlideTable.get_table(xmldoc).as_dict()
sngl_burst_index = dict((row.event_id, row) for row in lsctables.SnglBurstTable.get_table(xmldoc))
coinc_index = {}
for row in lsctables.CoincMapTable.get_table(xmldoc):
if row.table_name == "sngl_burst":
if row.coinc_event_id not in coinc_index:
coinc_index[row.coinc_event_id] = []
coinc_index[row.coinc_event_id].append(sngl_burst_index[row.event_id])
coinc_def_id = lsctables.CoincDefTable.get_table(xmldoc).get_coinc_def_id(search = cherenkov.CherenkovBBCoincDef.search, search_coinc_type = cherenkov.CherenkovBBCoincDef.search_coinc_type, create_new = False)
for coinc in lsctables.CoincTable.get_table(xmldoc):
if coinc.coinc_def_id == coinc_def_id:
coinc.likelihood = rankingstat.ln_lr_from_triggers(coinc_index[coinc.coinc_event_id], time_slide_index[coinc.time_slide_id])
process.set_end_time_now()
ligolw_utils.write_url(xmldoc, url, verbose = options.verbose)
xmldoc.unlink()
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