diff --git a/gstlal-inspiral/bin/gstlal_inspiral_calc_snr b/gstlal-inspiral/bin/gstlal_inspiral_calc_snr index b449d7ca540eabdd476db2175186f14ea6ef7941..a379a45badbb556f7020cad1303cf7160fe2efc4 100755 --- a/gstlal-inspiral/bin/gstlal_inspiral_calc_snr +++ b/gstlal-inspiral/bin/gstlal_inspiral_calc_snr @@ -43,7 +43,6 @@ Typical Usages: 2. Template options: --table (require) --approximant (require) - --template-duration (require) --sample-rate (default = 2048Hz) --f-low (default = 10) --f-high (optional) @@ -67,7 +66,6 @@ import numpy from optparse import OptionParser, OptionGroup, IndentedHelpFormatter import os -from gstlal import cbc_template_fir from gstlal import datasource from gstlal import pipeparts from gstlal import reference_psd @@ -209,27 +207,16 @@ def parse_command_line(): missing_required_options.append("--outdir") if options.instrument is None: missing_required_options.append("--instrument") - if options.template_duration is None: - missing_required_options.append("--template-duration") if options.template_psd is None: missing_required_options.append("--template-psd") # Raise VauleError is required option(s) is/are missing if missing_required_options: raise ValueError("Missing required option(s) %s" % ", ".join(sorted(missing_required_options))) + xmldoc = ligolw_utils.load_url(options.table, contenthandler = ContentHandler, verbose = options.verbose) template_table = table.get_table(xmldoc, lsctables.SnglInspiralTable.tableName) - row = [row for row in template_table if row.ifo == options.instrument] - if len(row) != 1 : - raise ValueError("Expecting only one template for --instrument=%s or cannot find template for --instrument=%s" %(options.instrument, options.instrument)) - - template_psd = lal.series.read_psd_xmldoc(ligolw_utils.load_url(options.template_psd, contenthandler = lal.series.PSDContentHandler)) - - # work around for building a single whitened template - time_slice = numpy.array([(options.sample_rate, 0, options.template_duration)], dtype = [("rate", "int"),("begin", "float"), ("end", "float")]) - workspace = cbc_template_fir.templates_workspace(template_table, options.approximant, template_psd[options.instrument], options.f_low, time_slice, autocorrelation_length = None, fhigh = options.f_high) - template, autocorrelation, sigma = workspace.make_whitened_template(row[0]) - return options, gw_data_source_info, template, psd + return options, gw_data_source_info, template_table, psd # Unknown mode else: @@ -257,9 +244,16 @@ if options.mode == 0: verbose = options.verbose ) - 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) + if options.row_number is None: + for index, snr in enumerate(lloid_snr(COMPLEX = options.complex, row_number = options.row_number, start = options.start, end = options.end)): + snr.epoch += bank.sngl_inspiral_table[index].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) + else: + lloidsnr = lloid_snr(COMPLEX = options.complex, row_number = options.row_number, start = options.start, end = options.end) + lloidsnr[0].epoch += bank.sngl_inspiral_table[options.row_number].end + snrdict = {options.instrument: lloidsnr} + 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, options.row_number, int(lloidsnr[0].epoch), int(lloidsnr[0].data.length * lloidsnr[0].deltaT))), verbose = options.verbose) # # uncomment to save all snrs in one single XML # @@ -267,6 +261,9 @@ if options.mode == 0: #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: + template_table = template + template, time_offset = svd_bank_snr.FIR_SNR.make_template(template_table, options.template_psd, options.sample_rate, options.approximant, options.instrument, options.f_low, f_high = options.f_high, verbose = options.verbose) + #FIXME: proper handle for latency fir_snr = svd_bank_snr.FIR_SNR( gw_data_source_info, @@ -281,5 +278,7 @@ elif options.mode == 1: verbose = options.verbose ) - snrdict = {options.instrument : fir_snr(COMPLEX = options.complex, start = options.start, end = options.end)} + firsnr = fir_snr(COMPLEX = options.complex, start = options.start, end = options.end) + firsnr[0].epoch += time_offset + snrdict = {options.instrument : firsnr} 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 5dd1cefe99852155c66f52130ba6cf213e8aa37f..85b97217b1285c2889507e9b9657ccd0a26c184c 100644 --- a/gstlal-inspiral/python/svd_bank_snr.py +++ b/gstlal-inspiral/python/svd_bank_snr.py @@ -14,6 +14,7 @@ from gi.repository import GObject, Gst, GstAudio GObject.threads_init() Gst.init(None) +from gstlal import cbc_template_fir from gstlal import datasource from gstlal import lloidparts from gstlal import multirate_datasource @@ -25,11 +26,13 @@ from gstlal import simplehandler import lal from lal import LIGOTimeGPS import lal.series +import lalsimulation as lalsim from ligo.lw import array as ligolw_array from ligo.lw import ligolw from ligo.lw import param as ligolw_param from ligo.lw import utils as ligolw_utils +from ligo.lw import lsctables @ligolw_array.use_in @ligolw_param.use_in @@ -296,6 +299,24 @@ 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 + @staticmethod + def make_template(template_table, template_psd, sample_rate, approximant, instrument, f_low, f_high = None, autocorrelation_length = None, verbose = False): + row = [row for row in template_table if row.ifo == instrument] + if len(row) != 1 : + raise ValueError("Expecting only one template for instrument=%s or cannot find template for instrument=%s" % (instrument, instrument)) + + template_psd = lal.series.read_psd_xmldoc(ligolw_utils.load_url(template_psd, contenthandler = lal.series.PSDContentHandler)) + if instrument not in set(template_psd): + raise ValueError("No such instrument: %s in template psd: (%s)"% (instrument, ", ".join(set(template_psd)))) + + # work around for building a single whitened template + template_duration = lalsim.SimInspiralChirpTimeBound(f_low, row[0].mass1 * lal.MSUN_SI, row[0].mass2 * lal.MSUN_SI, 0., 0.) + time_slice = numpy.array([(sample_rate, 0, template_duration)], dtype = [("rate", "int"),("begin", "float"), ("end", "float")]) + workspace = cbc_template_fir.templates_workspace(template_table, approximant, template_psd[instrument], f_low, time_slice, autocorrelation_length = None, fhigh = f_high) + template, autocorrelation, sigma = workspace.make_whitened_template(row[0]) + + return template, row[0].end + def __call__(self, COMPLEX = False, row_number = 0 , start = None, end = None): return self.get_snr_series(COMPLEX, row_number, start, end)