From a5f9f2d569e3463ff049109cc183a047ec69c490 Mon Sep 17 00:00:00 2001
From: Chad Hanna <chad.hanna@comp-hd-002.gwave.ics.psu.edu>
Date: Wed, 29 May 2019 19:59:25 -0400
Subject: [PATCH] add metric overlap code
---
gstlal-ugly/bin/Makefile.am | 3 ++
.../bin/gstlal_inspiral_add_metric_overlaps | 21 ++++++++
.../bin/gstlal_inspiral_metric_overlap | 52 +++++++++++++++++++
.../bin/gstlal_inspiral_metric_overlap_dag | 44 ++++++++++++++++
4 files changed, 120 insertions(+)
create mode 100755 gstlal-ugly/bin/gstlal_inspiral_add_metric_overlaps
create mode 100755 gstlal-ugly/bin/gstlal_inspiral_metric_overlap
create mode 100755 gstlal-ugly/bin/gstlal_inspiral_metric_overlap_dag
diff --git a/gstlal-ugly/bin/Makefile.am b/gstlal-ugly/bin/Makefile.am
index 27d32d0889..eeea5cffeb 100644
--- a/gstlal-ugly/bin/Makefile.am
+++ b/gstlal-ugly/bin/Makefile.am
@@ -11,8 +11,11 @@ dist_bin_SCRIPTS = \
gstlal_harmonic_mean_psd \
gstlal_mock_data_server \
gstlal_ilwdify \
+ gstlal_inspiral_add_metric_overlaps \
gstlal_inspiral_best_coinc_file \
gstlal_inspiral_check_livetimes \
+ gstlal_inspiral_metric_overlap \
+ gstlal_inspiral_metric_overlap_dag \
gstlal_inspiral_treebank \
gstlal_inspiral_treebank_dag \
gstlal_inspiral_bankviz \
diff --git a/gstlal-ugly/bin/gstlal_inspiral_add_metric_overlaps b/gstlal-ugly/bin/gstlal_inspiral_add_metric_overlaps
new file mode 100755
index 0000000000..bfb88fe396
--- /dev/null
+++ b/gstlal-ugly/bin/gstlal_inspiral_add_metric_overlaps
@@ -0,0 +1,21 @@
+#!/usr/bin/python
+import sys
+import numpy
+import argparse
+import h5py
+
+
+parser = argparse.ArgumentParser()
+parser.add_argument("--out-h5-file", required = True, help = "provide the output hdf5 file name")
+parser.add_argument("infile", nargs="+", help = "provide the output hdf5 file name")
+
+args = parser.parse_args()
+
+output = []
+for f in args.infile:
+ h5f = h5py.File(f, 'r')
+ output += list(numpy.array(h5f["overlaps"]))
+output = numpy.array(output)
+h5f = h5py.File(args.out_h5_file, 'w')
+h5f.create_dataset('overlaps', data = output)
+h5f.close()
diff --git a/gstlal-ugly/bin/gstlal_inspiral_metric_overlap b/gstlal-ugly/bin/gstlal_inspiral_metric_overlap
new file mode 100755
index 0000000000..14a590ca8a
--- /dev/null
+++ b/gstlal-ugly/bin/gstlal_inspiral_metric_overlap
@@ -0,0 +1,52 @@
+#!/usr/bin/python
+import sys
+from gstlal import metric as metric_module
+from ligo.lw import ligolw
+from ligo.lw import utils as ligolw_utils
+from ligo.lw import lsctables
+import numpy
+import argparse
+import h5py
+
+class LIGOLWContentHandler(ligolw.LIGOLWContentHandler):
+ pass
+lsctables.use_in(LIGOLWContentHandler)
+
+parser = argparse.ArgumentParser()
+parser.add_argument("--psd-xml-file", help = "provide a psd xml file")
+parser.add_argument("--bank-file", help = "provide the bank file for which overlaps will be calculated")
+parser.add_argument("--out-h5-file", required = True, help = "provide the output hdf5 file name")
+parser.add_argument("--start-row", type=int, help = "The starting row to calculate overlaps")
+parser.add_argument("--num-rows", type=int, help = "The number of rows to calculate overlaps")
+
+args = parser.parse_args()
+
+g_ij = metric_module.Metric(
+ args.psd_xml_file,
+ coord_func = metric_module.m1_m2_s1z_s2z_func,
+ duration = 1.0, # FIXME!!!!!
+ flow = 10,
+ fhigh = 1024,
+ approximant = "IMRPhenomD"
+)
+
+xmldoc = ligolw_utils.load_filename(args.bank_file, verbose = True, contenthandler = LIGOLWContentHandler)
+sngl_inspiral_table = lsctables.SnglInspiralTable.get_table(xmldoc)
+
+vec1s = numpy.array([[row.template_id, row.mass1, row.mass2, row.spin1z, row.spin2z] for row in sngl_inspiral_table[args.start_row:args.start_row+args.num_rows]])
+vec2s = numpy.array([[row.template_id, row.mass1, row.mass2, row.spin1z, row.spin2z] for row in sngl_inspiral_table])
+
+output = []
+
+for n, vec1 in enumerate(vec1s):
+ g, det = g_ij(vec1[1:])
+ def match(vec2, vec1 = vec1, g = g):
+ return (vec1[0], vec2[0], g_ij.metric_match(g, vec1[1:], vec2[1:]))
+ thisoutput = [row for row in map(match, vec2s) if row[2] > 0.9]
+ print n, len(thisoutput)
+ output += thisoutput
+
+output = numpy.array(output)
+h5f = h5py.File(args.out_h5_file, 'w')
+h5f.create_dataset('overlaps', data = output)
+h5f.close()
diff --git a/gstlal-ugly/bin/gstlal_inspiral_metric_overlap_dag b/gstlal-ugly/bin/gstlal_inspiral_metric_overlap_dag
new file mode 100755
index 0000000000..71c5252115
--- /dev/null
+++ b/gstlal-ugly/bin/gstlal_inspiral_metric_overlap_dag
@@ -0,0 +1,44 @@
+#!/usr/bin/python
+import sys
+from ligo.lw import ligolw
+from ligo.lw import utils as ligolw_utils
+from ligo.lw import lsctables
+import numpy
+import argparse
+import h5py
+from gstlal import pipeline
+from gstlal import dagparts
+
+class LIGOLWContentHandler(ligolw.LIGOLWContentHandler):
+ pass
+lsctables.use_in(LIGOLWContentHandler)
+
+parser = argparse.ArgumentParser()
+parser.add_argument("--psd-xml-file", help = "provide a psd xml file")
+parser.add_argument("--bank-file", help = "provide the bank file for which overlaps will be calculated")
+args = parser.parse_args()
+
+xmldoc = ligolw_utils.load_filename(args.bank_file, verbose = True, contenthandler = LIGOLWContentHandler)
+sngl_inspiral_table = lsctables.SnglInspiralTable.get_table(xmldoc)
+
+try:
+ os.mkdir("logs")
+except:
+ pass
+dag = dagparts.DAG("metric_overlap")
+
+overlapJob = dagparts.DAGJob("gstlal_inspiral_metric_overlap", condor_commands = {"want_graceful_removal":"True", "kill_sig":"15", "accounting_group":"ligo.prod.o3.cbc.uber.gstlaloffline"})
+addJob = dagparts.DAGJob("gstlal_inspiral_add_metric_overlaps", condor_commands = {"want_graceful_removal":"True", "kill_sig":"15", "accounting_group":"ligo.prod.o3.cbc.uber.gstlaloffline"})
+
+num = 1000
+overlapnodes = []
+# FIXME dont hardcode 3345408, it comes from number of tiles in TimePhaseSNR
+for start in range(0, len(sngl_inspiral_table), num):
+ stop = start + num
+ overlapnodes.append(dagparts.DAGNode(overlapJob, dag, parent_nodes = [], opts = {"start-row":str(start), "num-rows":str(num)}, output_files = {"out-h5-file":"%s/metric_overlaps_%d_%d.h5" % (overlapJob.output_path, start, num)}, input_files = {"psd-xml-file": args.psd_xml_file, "bank-file": args.bank_file}))
+
+addnode = dagparts.DAGNode(addJob, dag, parent_nodes = overlapnodes, input_files = {"": [n.output_files["out-h5-file"] for n in overlapnodes]})
+
+dag.write_sub_files()
+dag.write_dag()
+dag.write_script()
--
GitLab