#!/usr/bin/python
#
# Project Librarian: Brandon Piotrzkowski
# Graduate Student
# UW-Milwaukee Department of Physics
# Center for Gravitation & Cosmology
# <brandon.piotrzkowski@ligo.org>
#
# 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 3 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, see <http://www.gnu.org/licenses/>.
#
"""
Script to calculate and upload the joint false alarm rate of the provided candidates
to GraceDB, either the spatiotemporal variant if sky information is provided, or
just temporal if not.
"""
__author__ = "Brandon Piotrzkowski <brandon.piotrzkowski@ligo.org>"
# Global imports.
import numpy as np
from ligo.raven import search
from ligo.gracedb.rest import GraceDb
import argparse
# Command line options.
parser = argparse.ArgumentParser(description='Perform query of GraceDB')
parser.add_argument("-s", "--superevent", metavar="(S,MS)XXXXXXX",
default=None,
help="GraceDB ID of superevent (required)"),
parser.add_argument("-e", "--ext_event", metavar="(E,M)XXXXXXX",
default=None,
help="GraceDB ID of external event (required)"),
parser.add_argument("-w", "--window", nargs=2, metavar="t", type=int,
action="append", default=None,
help="Time window [tl, th] seconds to search around event time (required with th > th)"),
parser.add_argument("-g", "--ext_search", metavar="GRB SubGRB SubGRBTargeted HEN MDC",
default=None,
help="Search of external event."),
parser.add_argument("-S", "--superevent_skymap", metavar="bayestar.fits.gz",
default=None,
help="fits(.gz) filename for superevent sky map"),
parser.add_argument("-E", "--ext_event_skymap", metavar="'glg_healpix_all_bn_v00.fit'",
default=None,
help="fits(.gz) filename for external event sky map"),
parser.add_argument("-u", "--gracedb_url", metavar="https://gracedb.ligo.org/api/",
default='',
help="Gracedb url"),
parser.add_argument("-f", "--far_ext", metavar="1e-6", default=None,
help="False alarm rate of GRB"),
parser.add_argument("-n", "--em_rate", metavar="1e-3", default=None,
help="False alarm rate of GRB"),
parser.add_argument("-t", "--far_gw_thresh", metavar="1.16e-5", default=None,
help="Maximum cutoff for GW FAR considered in the search"),
parser.add_argument("-T", "--far_ext_thresh", metavar="1e-4", default=None,
help="Maximum cutoff for GRB FAR considered in the search"),
parser.add_argument("-m", "--se_moc", dest="se_moc", action="store_true",
help="Assume the superevent sky map is multi-ordered (MOC)."),
parser.add_argument("-M", "--ext_moc", dest="ext_moc", action="store_true",
help="Assume the external sky map is multi-ordered (MOC)."),
parser.add_argument("-r", "--se_ring", dest="se_ring", action="store_true",
help="Assume the superevent map uses RING ordering rather than nested."),
parser.add_argument("-R", "--ext_ring", dest="ext_ring", action="store_true",
help="Assume the external sky map uses RING ordering rather than nested."),
parser.add_argument("-c", "--use_radec", dest="use_radec", action="store_true",
help="Choose to use RA and dec of external sky map if a single pixel sky map."),
parser.add_argument("-p", "--use_preferred_event_skymap", dest="use_preferred_event_skymap", action="store_true",
help="Choose to use the GW sky map in the preferred event rather than the superevent.")
args = parser.parse_args()
print(args)
# Check required options are there
if not args.superevent:
raise AssertionError('superevent graceid not given')
else:
se_id = str(args.superevent)
if not args.ext_event:
raise AssertionError('external event graceid not given')
else:
exttrig_id = str(args.ext_event)
if not args.window:
raise AssertionError('search window not given')
else:
tl, th = int(args.window[0][0]), int(args.window[0][1])
# Load other options
if args.gracedb_url:
gracedb = GraceDb(str(args.gracedb_url))
else:
gracedb = None
ext_search = args.ext_search
se_fitsfile = args.superevent_skymap
ext_fitsfile = args.ext_event_skymap
if args.far_ext:
far_ext = float(args.far_ext)
else:
far_ext = None
if args.em_rate:
em_rate = float(args.em_rate)
else:
em_rate = None
if args.far_gw_thresh:
far_gw_thresh = float(args.far_gw_thresh)
else:
far_gw_thresh = None
if args.far_ext_thresh:
far_ext_thresh = float(args.far_ext_thresh)
else:
far_ext_thresh = None
if se_fitsfile and (ext_fitsfile or args.use_radec):
incl_sky = True
else:
incl_sky = False
results = search.calc_signif_gracedb(
se_id, exttrig_id, tl, th, ext_search=ext_search,
se_fitsfile=se_fitsfile, ext_fitsfile=ext_fitsfile,
incl_sky=incl_sky, gracedb=gracedb,
far_ext=far_ext, em_rate=em_rate,
far_gw_thresh=far_gw_thresh,
far_ext_thresh=far_ext_thresh,
se_moc=args.se_moc, ext_moc=args.ext_moc,
se_nested=not args.se_ring, ext_nested=not args.ext_ring,
use_radec=args.use_radec,
use_preferred_event_skymap=args.use_preferred_event_skymap)
print(results)
