From 8c3ce90b6cc70cfd17b19a05050980659588ecda Mon Sep 17 00:00:00 2001 From: Leo Tsukada <leo.tsukada@ligo.org> Date: Mon, 21 Mar 2022 13:45:22 -0700 Subject: [PATCH 1/2] gstlal_inspiral_bank_splitter : implement 2D-sorting by orthogonalized 0/1PN order phase terms --- .../bin/gstlal_inspiral_bank_splitter | 46 +++++++++++++++++-- 1 file changed, 42 insertions(+), 4 deletions(-) diff --git a/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter b/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter index b32e53b74f..023ac238f5 100755 --- a/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter +++ b/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter @@ -31,6 +31,7 @@ from ligo.lw import utils as ligolw_utils from ligo.lw.utils import process as ligolw_process import lal import numpy, scipy.interpolate +from astropy import constants as const from gstlal import chirptime from gstlal import datafind @@ -71,6 +72,33 @@ class LIGOLWContentHandler(ligolw.LIGOLWContentHandler): pass lsctables.use_in(LIGOLWContentHandler) +def calc_mu(mass1, mass2, spin1z, spin2z, fref=100, mu="mu1"): + """ + Calculate the first orthogonal PN phase coefficient + see https://arxiv.org/abs/2007.09108 + """ + + M = mass1 + mass2 + mchirp = (mass1 * mass2)**0.6/M**0.2 + eta = mass1 * mass2 / M**2 + beta = ((113. * (mass1 / M)**2 + 75. * eta) * spin1z + (113. * (mass2 / M)**2 + 75. * eta) * spin2z) / 12. + + norm = const.G.value * const.M_sun.value / const.c.value**3 + v = numpy.pi * mchirp * fref * norm + psi0 = 3. / 4 / (8 * v)**(5./3) + psi2 = 20. / 9 * (743./ 336 + 11. / 4 * eta) / eta**(0.4) * v**(2./3) * psi0 + psi3 = (4 * beta - 16 * numpy.pi) / eta**0.6 * v * psi0 + + # FIXME : the following linear combinations are taken from the ones in the + # paper above, but this will need to be re-computed with o4 representitive + # psd. + if mu=="mu1": + return 0.974 * psi0 + 0.209 * psi2 + 0.0840 * psi3 + elif mu=="mu2": + return -0.221 * psi0 + 0.823 * psi2 + 0.524 * psi3 + else: + raise ValueError("%s is not implemented, so cannnot be computed." % mu) + def group_templates(templates, n, overlap = 0): """ @@ -191,12 +219,19 @@ if options.output_full_bank_file is not None: outputrows = [] -# Bin by Chi -sngl_inspiral_table.sort(key = lambda row: spawaveform.computechi(row.mass1, row.mass2, row.spin1z, row.spin2z)) +if options.sort_by == "mu": + # First partitioning by mu1 + sngl_inspiral_table.sort(key = lambda row: calc_mu(row.mass1, row.mass2, row.spin1z, row.spin2z, mu="mu1")) +else: + # First partitioning by chi + sngl_inspiral_table.sort(key = lambda row: spawaveform.computechi(row.mass1, row.mass2, row.spin1z, row.spin2z)) for chirows in group_templates(sngl_inspiral_table, len(sngl_inspiral_table) / options.group_by_chi, overlap = 0): def sort_func(row, column = options.sort_by): - return getattr(row, column) + if column=="mu": + return calc_mu(row.mass1, row.mass2, row.spin1z, row.spin2z, mu="mu2") + else: + return getattr(row, column) chirows.sort(key=sort_func, reverse = True if options.sort_by == "template_duration" else False) @@ -211,7 +246,10 @@ for chirows in group_templates(sngl_inspiral_table, len(sngl_inspiral_table) / o # A sort of the groups of templates so that the sub banks are ordered. def sort_func(row_rows_tup, column = options.sort_by): - return getattr(row_rows_tup[0], column) + if column=="mu": + return calc_mu(row_rows_tup[0].mass1, row_rows_tup[0].mass2, row_rows_tup[0].spin1z, row_rows_tup[0].spin2z, mu="mu2") + else: + return getattr(row_rows_tup[0], column) outputrows.sort(key=sort_func, reverse = True if options.sort_by == "template_duration" else False) svd_dirs = [] -- GitLab From 3e15b66a8ad25e9cc3d75f6db694a1c5bbbce055 Mon Sep 17 00:00:00 2001 From: Leo Tsukada <leo.tsukada@ligo.org> Date: Mon, 21 Mar 2022 15:16:29 -0700 Subject: [PATCH 2/2] gstlal_inspiral_bank_splitter : make lal constants and add --group-by-mu for consistentcy --- gstlal-inspiral/bin/gstlal_inspiral_bank_splitter | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter b/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter index 023ac238f5..0f4b4f0637 100755 --- a/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter +++ b/gstlal-inspiral/bin/gstlal_inspiral_bank_splitter @@ -31,7 +31,6 @@ from ligo.lw import utils as ligolw_utils from ligo.lw.utils import process as ligolw_process import lal import numpy, scipy.interpolate -from astropy import constants as const from gstlal import chirptime from gstlal import datafind @@ -83,7 +82,7 @@ def calc_mu(mass1, mass2, spin1z, spin2z, fref=100, mu="mu1"): eta = mass1 * mass2 / M**2 beta = ((113. * (mass1 / M)**2 + 75. * eta) * spin1z + (113. * (mass2 / M)**2 + 75. * eta) * spin2z) / 12. - norm = const.G.value * const.M_sun.value / const.c.value**3 + norm = lal.G_SI * lal.MSUN_SI / lal.C_SI**3 v = numpy.pi * mchirp * fref * norm psi0 = 3. / 4 / (8 * v)**(5./3) psi2 = 20. / 9 * (743./ 336 + 11. / 4 * eta) / eta**(0.4) * v**(2./3) * psi0 @@ -134,6 +133,7 @@ def parse_command_line(): parser.add_option("--approximant", type = "string", action = "append", help = "Must specify an approximant given as mchirp_min:mchirp_max:string") parser.add_option("--f-low", type = "float", metavar = "frequency", help = "Lower frequency cutoff. Required") parser.add_option("--group-by-chi", type = "int", metavar = "N", default = 1, help = "group templates into N groups of chi - helps with SVD. Default 1") + parser.add_option("--group-by-mu", type = "int", metavar = "N", default = 20, help = "group templates into N groups of mu1, the 1PN phase term, to help with SVD. Default 20") parser.add_option("--num-banks", metavar = "str", help = "The number of parallel subbanks per SVD bank. can be given as a list like 1,2,3,4 then it will split up the bank into N groups with M banks each. (required)") # FIXME The second half of this help message is incomprehensible options, filenames = parser.parse_args() @@ -222,10 +222,12 @@ outputrows = [] if options.sort_by == "mu": # First partitioning by mu1 sngl_inspiral_table.sort(key = lambda row: calc_mu(row.mass1, row.mass2, row.spin1z, row.spin2z, mu="mu1")) + num_group = options.group_by_mu else: # First partitioning by chi sngl_inspiral_table.sort(key = lambda row: spawaveform.computechi(row.mass1, row.mass2, row.spin1z, row.spin2z)) -for chirows in group_templates(sngl_inspiral_table, len(sngl_inspiral_table) / options.group_by_chi, overlap = 0): + num_group = options.group_by_chi +for chirows in group_templates(sngl_inspiral_table, len(sngl_inspiral_table) / num_group, overlap = 0): def sort_func(row, column = options.sort_by): if column=="mu": -- GitLab