From 130dbda4ddf0e3e9dcbebfb18b12ca5f75115f3d Mon Sep 17 00:00:00 2001
From: Leo Tsukada <leo.tsukada@ligo.org>
Date: Wed, 7 Aug 2019 19:03:34 -0700
Subject: [PATCH] gstlal_inspiral_compute_dtdphideff_cov_matrix : changed to
 pass snr dictioary and psd file name to the output .h5 file.

---
 .../bin/gstlal_inspiral_compute_dtdphideff_cov_matrix | 11 +++++++----
 1 file changed, 7 insertions(+), 4 deletions(-)

diff --git a/gstlal-inspiral/bin/gstlal_inspiral_compute_dtdphideff_cov_matrix b/gstlal-inspiral/bin/gstlal_inspiral_compute_dtdphideff_cov_matrix
index d984e591c5..39e879041a 100755
--- a/gstlal-inspiral/bin/gstlal_inspiral_compute_dtdphideff_cov_matrix
+++ b/gstlal-inspiral/bin/gstlal_inspiral_compute_dtdphideff_cov_matrix
@@ -21,7 +21,6 @@ import numpy, scipy.interpolate
 from lal import series
 from gstlal import templates
 from ligo.lw import utils as ligolw_utils
-import sys
 import h5py
 
 parser = argparse.ArgumentParser(description = 'generate a dt dphi covariance matrix')
@@ -46,11 +45,11 @@ refpsd = args.psd_xml
 rho = {"L1": args.L_snr, "H1": args.H_snr, "V1": args.V_snr, "K1": args.K_snr}
 
 instruments = ("H1", "L1", "V1", "K1")
-psd = series.read_psd_xmldoc(ligolw_utils.load_filename(refpsd, verbose = True, contenthandler = series.PSDContentHandler))
+snr = series.read_psd_xmldoc(ligolw_utils.load_filename(refpsd, verbose = True, contenthandler = series.PSDContentHandler))
 psddict = {}
 for ifo in rho:
-	f = numpy.arange(len(psd[ifo].data.data)) * psd[ifo].deltaF
-	psddict[ifo] = scipy.interpolate.interp1d(f, psd[ifo].data.data)
+	f = numpy.arange(len(snr[ifo].data.data)) * snr[ifo].deltaF
+	psddict[ifo] = scipy.interpolate.interp1d(f, snr[ifo].data.data)
 
 def moment(fmin, fmax, n, m1, m2, s1, s2, ifo, psddict = psddict, delta_f = 0.25):
 	farr = numpy.linspace(0, fmax, fmax / delta_f + delta_f)
@@ -94,6 +93,10 @@ for pair in pairs:
 	transdd[frozenset(pair)] = 1. / numpy.sqrt(sigsqdd[a] + sigsqdd[b])
 
 f = h5py.File(args.output, "w")
+f.create_dataset("psd", data=args.psd_xml)
+h5_snr = f.create_group("SNR")
+for ifo, snr in rho.items():
+    h5_snr.create_dataset(ifo, data=snr)
 h5_transtt = f.create_group("transtt")
 h5_transtp = f.create_group("transtp")
 h5_transpt = f.create_group("transpt")
-- 
GitLab