From 638510a5cda8e628f3ba1690b973216ddf833e30 Mon Sep 17 00:00:00 2001
From: ChiWai Chan <chiwai.chan@ligo.org>
Date: Mon, 25 Mar 2019 12:24:59 +0900
Subject: [PATCH] svd_bank_snr.py: change drop snr data -> select snr data

---
 gstlal-inspiral/bin/gstlal_inspiral_calc_snr | 17 ++++++----
 gstlal-inspiral/python/svd_bank_snr.py       | 35 +++++++++++++-------
 2 files changed, 33 insertions(+), 19 deletions(-)

diff --git a/gstlal-inspiral/bin/gstlal_inspiral_calc_snr b/gstlal-inspiral/bin/gstlal_inspiral_calc_snr
index 6d2ea6ebe6..82394b3209 100755
--- a/gstlal-inspiral/bin/gstlal_inspiral_calc_snr
+++ b/gstlal-inspiral/bin/gstlal_inspiral_calc_snr
@@ -136,8 +136,8 @@ def parse_command_line():
 	group.add_option("--outdir", metavar = "directory", type = "str", help = "Output directory for SNR(s) (requires).")
 	group.add_option("--mode", metavar = "method", type = "int", default = 0, help = "The method (0 = LLOID / 1 = FIR) that is used to calculate SNR (default = 0).")
 	group.add_option("--complex", action = "store_true", help = "Choose whether to output the complex snr or not.")
-	group.add_option("--drop-first", metavar = "seconds", type = "int", default = 0, help = "Dropping the initital '--drop-first' seconds of SNR data (default = 0).")
-	group.add_option("--drop-last", metavar = "seconds", type = "int", default = 0, help = "Dropping the last '--drop-last' seconds of SNR data (default = 0).")
+	group.add_option("--start", metavar = "seconds", type = "int", help = "Start SNR time series at GPS time '--start'.")
+	group.add_option("--end", metavar = "seconds", type = "int", help = "End SNR time series at GPS time '--end'.")
 	group.add_option("--output-width", metavar = "bits", type = "int", default = 32, help = "The size of the output data, can only be 32 or 64 bits (default = 32 bits).")
 	group.add_option("--instrument", metavar = "name", help = "The detector from which the --reference-psd and --frame-cache are loaded (require).")
 	parser.add_option_group(group)
@@ -146,6 +146,10 @@ def parse_command_line():
 
 	options, args = parser.parse_args()
 
+	# Check SNR series output
+	if options.start >= options.end:
+		raise ValueError("--start must less than --end.")
+
 	# Setting up GW data
 	gw_data_source_info = datasource.GWDataSourceInfo(options)
 	if options.instrument not in gw_data_source_info.channel_dict.keys():
@@ -247,15 +251,14 @@ if options.mode == 0:
 		verbose = options.verbose
 	)
 
-	for index, snr in enumerate(lloid_snr(COMPLEX = options.complex, row_number = options.row_number, drop_first = options.drop_first, drop_last = options.drop_last)):
+	for index, snr in enumerate(lloid_snr(COMPLEX = options.complex, row_number = options.row_number, start = options.start, end = options.end)):
 		snrdict = {options.instrument: [snr]}
 		svd_bank_snr.write_url(svd_bank_snr.make_xmldoc(snrdict), os.path.join(options.outdir, "%s-SNR_%d-%d-%d.xml.gz" % (options.instrument, index, int(snr.epoch), int(snr.data.length * snr.deltaT))), verbose = options.verbose)
 	#
 	# uncomment to save all snrs in one single XML
 	#
 	#snrdict = {options.instrument : lloid_snr(COMPLEX = options.complex, row_number = options.row_number, drop_first = options.drop_first, drop_last = options.drop_last)}
-	#tempsnr = 
-	#svd_bank_snr.write_url(svd_bank_snr.make_xmldoc(snrdict), os.path.join(options.outdir, "%s-SNR-%d-%d.xml.gz" % (options.instrument, int(snrdict.values()[0].epoch), int(snrdict.values()[0].data.length * snrdict.values()[0].deltaT))), verbose = options.verbose)
+	#svd_bank_snr.write_url(svd_bank_snr.make_xmldoc(snrdict), os.path.join(options.outdir, "%s-SNR-%d-%d.xml.gz" % (options.instrument, int(snrdict.[options.instrument][0].epoch), int(snrdict[options.instrument][0].data.length * snrdict[options.instrument][0].deltaT))), verbose = options.verbose)
 
 elif options.mode == 1:
 	#FIXME: proper handle for latency
@@ -272,5 +275,5 @@ elif options.mode == 1:
 		verbose = options.verbose
 		)
 
-	snrdict = {options.instrument : fir_snr(COMPLEX = options.complex, drop_first = options.drop_first, drop_last = options.drop_last)}
-	svd_bank_snr.write_url(svd_bank_snr.make_xmldoc(snrdict),os.path.join(options.outdir, "%s-SNR-%d-%d.xml.gz" % (options.instrument, int(snrdict.values()[0].epoch), int(snrdict.values()[0].data.length * snrdict.values()[0].deltaT))), verbose = options.verbose)
+	snrdict = {options.instrument : fir_snr(COMPLEX = options.complex, start = options.start, end = options.end)}
+	svd_bank_snr.write_url(svd_bank_snr.make_xmldoc(snrdict),os.path.join(options.outdir, "%s-SNR-%d-%d.xml.gz" % (options.instrument, int(snrdict[options.instrument][0].epoch), int(snrdict[options.instrument][0].data.length * snrdict[options.instrument][0].deltaT))), verbose = options.verbose)
diff --git a/gstlal-inspiral/python/svd_bank_snr.py b/gstlal-inspiral/python/svd_bank_snr.py
index 88135ff271..371928d9f7 100644
--- a/gstlal-inspiral/python/svd_bank_snr.py
+++ b/gstlal-inspiral/python/svd_bank_snr.py
@@ -119,14 +119,25 @@ class SNR_Pipeline(object):
 
 		return tseries
 
-	def get_snr_series(self, COMPLEX = False, row_number = None, drop_first = 0, drop_last = 0):
-		assert drop_first >= 0, "must drop positive number of data"
-		assert drop_last >= 0, "must drop positive number of data"
-		bps = drop_first * int(round(1 / self.snr_info["deltaT"]))
-		bpe = -drop_last * int(round(1 / self.snr_info["deltaT"])) if drop_last != 0 else None
+	def get_snr_series(self, COMPLEX = False, row_number = None, start = None, end = None):
+		gps_start = self.snr_info["epoch"].gpsSeconds + self.snr_info["epoch"].gpsNanoSeconds * 10.**-9
+		gps = gps_start + numpy.arange(len(self.snr_info["data"])) * self.snr_info["deltaT"]
 
-		self.snr_info["epoch"] += drop_first
-		self.snr_info["data"] = self.snr_info["data"][bps:bpe].T
+		if start >= end:
+			raise ValueError("Start time must be less than end time.")
+
+		if start - gps[0] >= 0 and start - gps[-1] <= 0:
+			s = abs(gps - start).argmin()
+		else:
+			raise ValueError("Invalid choice of start time %f." % start)
+
+		if end - gps[0] >= 0 and end - gps[-1] <= 0:
+			e = abs(gps - end).argmin()
+		else:
+			raise ValueError("Invalid choice of end time %f." % end)
+
+		self.snr_info["epoch"] = gps[s]
+		self.snr_info["data"] = self.snr_info["data"][s:e].T
 
 		if row_number is None:
 			temp = []
@@ -136,7 +147,7 @@ class SNR_Pipeline(object):
 				return temp
 			else:
 				for data in self.snr_info["data"]:
-					temp.append(numpy.abs(self.make_series(data)))
+					temp.append(self.make_series(numpy.abs(data)))
 				return temp
 		else:
 			self.snr_info["data"] = self.snr_info["data"][row_number]
@@ -233,8 +244,8 @@ class LLOID_SNR(SNR_Pipeline):
 		self.run(gw_data_source_info.seg)
                 self.snr_info["data"] = numpy.concatenate(numpy.array(self.snr_info["data"]), axis = 0)
 
-	def __call__(self, COMPLEX = False, row_number = 0, drop_first = 0, drop_last = 0):
-		return self.get_snr_series(COMPLEX, row_number, drop_first, drop_last)
+	def __call__(self, COMPLEX = False, row_number = 0, start = None, end = None):
+		return self.get_snr_series(COMPLEX, row_number, start, end)
 
 class FIR_SNR(SNR_Pipeline):
 	def __init__(self, gw_data_source_info, template, instrument, rate, latency, psd = None, psd_fft_length = 32, ht_gate_threshold = float("inf"), veto_segments = None, width = 32, track_psd = False, verbose = False):
@@ -282,8 +293,8 @@ class FIR_SNR(SNR_Pipeline):
 		self.snr_info["data"] = numpy.vectorize(complex)(self.snr_info["data"][:,0], self.snr_info["data"][:,1])
 		self.snr_info["data"].shape = len(self.snr_info["data"]), 1
 
-	def __call__(self, COMPLEX = False, row_number = 0 , drop_first = 0, drop_last = 0):
-		return self.get_snr_series(COMPLEX, row_number, drop_first, drop_last)
+	def __call__(self, COMPLEX = False, row_number = 0 , start = None, end = None):
+		return self.get_snr_series(COMPLEX, row_number, start, end)
 
 #=============================================================================================
 #
-- 
GitLab