From e24180e1dc93d0e8d1218d9d7b3488d78aca29e0 Mon Sep 17 00:00:00 2001 From: Yifan Wang Date: Wed, 7 Feb 2024 16:37:33 +0100 Subject: [PATCH 01/29] fix for injection file created from pycbc_create_injections (#4623) * fix for injection file created from pycbc_create_injections * Update bin/pycbc_convertinjfiletohdf Co-authored-by: Gareth S Cabourn Davies --------- Co-authored-by: Gareth S Cabourn Davies --- bin/pycbc_convertinjfiletohdf | 18 ++++++++---------- 1 file changed, 8 insertions(+), 10 deletions(-) diff --git a/bin/pycbc_convertinjfiletohdf b/bin/pycbc_convertinjfiletohdf index 9c703d51f98..f0557fa03b2 100755 --- a/bin/pycbc_convertinjfiletohdf +++ b/bin/pycbc_convertinjfiletohdf @@ -23,6 +23,7 @@ or the LIGO HDF format) into a PyCBC hdf injection format import argparse import numpy import h5py +import shutil import pycbc from pycbc.inject import InjectionSet, CBCHDFInjectionSet @@ -179,14 +180,11 @@ else: injclass = LVKNewStyleInjectionSet(args.injection_file) inj_file.close() -if injclass is None: - # Assume this is a PyCBC format - injclass = InjectionSet(args.injection_file) - data = {} - for key in xinj.table[0].__slots__: - data[str(key)] = numpy.array([getattr(t, key) for t in xinj.table]) -else: +if injclass is not None: + # The injection file is of known type data = injclass.pack_data_into_pycbc_format_input() - -samples = FieldArray.from_kwargs(**data) -CBCHDFInjectionSet.write(args.output_file, samples) + samples = FieldArray.from_kwargs(**data) + CBCHDFInjectionSet.write(args.output_file, samples) +else: + # Assume the injection is a HDF file (PyCBC format), only make a copy + shutil.copy(args.injection_file, args.output_file) \ No newline at end of file From e1984b6b6a8e302e543ac7365143958d8a1df5a4 Mon Sep 17 00:00:00 2001 From: Gareth S Cabourn Davies Date: Fri, 9 Feb 2024 16:37:18 +0000 Subject: [PATCH 02/29] Slightly incomplete example in frame documentation (#4624) --- docs/frame.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/frame.rst b/docs/frame.rst index 974f3c81261..0966c1abafc 100644 --- a/docs/frame.rst +++ b/docs/frame.rst @@ -33,7 +33,7 @@ Reading a frame file The ``pycbc.frame`` module provides methods for reading these files into ``TimeSeries`` objects as follows:: >>> from pycbc import frame - >>> data = frame.read_frame('G-G1_RDS_C01_L3-1049587200-60.gwf', 'G1:DER_DATA_H') + >>> data = frame.read_frame('G-G1_RDS_C01_L3-1049587200-60.gwf', 'G1:DER_DATA_H', 1049587200, 1049587200 + 60) Here the first argument is the path to the frame file of interest, while the second lists the `data channel` of interest whose data exist within the file. From 02f77b6d2b0b8fd2a62903f02356954d81412252 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Cl=C3=A9ment=20Robert?= Date: Fri, 9 Feb 2024 17:49:33 +0100 Subject: [PATCH 03/29] MNT: upgrade upload/download-artifact GitHub actions (#4595) --- .github/workflows/basic-tests.yml | 7 ++++--- .github/workflows/distribution.yml | 10 ++++++---- .github/workflows/inference-workflow.yml | 4 ++-- .github/workflows/search-workflow.yml | 4 ++-- .github/workflows/tmpltbank-workflow.yml | 2 +- .github/workflows/workflow-tests.yml | 2 +- 6 files changed, 16 insertions(+), 13 deletions(-) diff --git a/.github/workflows/basic-tests.yml b/.github/workflows/basic-tests.yml index d962026e09a..29739899814 100644 --- a/.github/workflows/basic-tests.yml +++ b/.github/workflows/basic-tests.yml @@ -51,7 +51,7 @@ jobs: tox -e py-inference - name: store documentation page if: matrix.test-type == 'docs' && matrix.python-version == '3.8' - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: documentation-page path: _gh-pages @@ -61,9 +61,10 @@ jobs: if: github.ref == 'refs/heads/master' && github.event_name == 'push' steps: - name: retrieve built documentation - uses: actions/download-artifact@v2 + uses: actions/download-artifact@v4 with: - name: documentation-page + pattern: documentation-page-* + merge-multiple: true - name: debug run: | mkdir _gh-pages diff --git a/.github/workflows/distribution.yml b/.github/workflows/distribution.yml index 321b81ff4d2..e59fb75bc87 100644 --- a/.github/workflows/distribution.yml +++ b/.github/workflows/distribution.yml @@ -29,8 +29,9 @@ jobs: CIBW_BUILD: cp38-* cp39-* cp310-* cp311-* CIBW_SKIP: "*musllinux*" CIBW_ARCHS_MACOS: x86_64 arm64 - - uses: actions/upload-artifact@v2 + - uses: actions/upload-artifact@v4 with: + name: wheels-${{ matrix.os }} path: ./wheelhouse/*.whl deploy_pypi: name: Build and publish Python 🐍 distributions 📦 to PyPI @@ -44,13 +45,14 @@ jobs: uses: actions/setup-python@v4 with: python-version: 3.8 - - uses: actions/download-artifact@v2 + - uses: actions/download-artifact@v4 with: - path: ./ + pattern: wheels-* + merge-multiple: true - name: build pycbc for pypi run: | python setup.py sdist - mv artifact/* dist/ + mv *.whl dist/ - name: Publish distribution 📦 to PyPI if: github.event_name == 'push' && startsWith(github.ref, 'refs/tags') uses: pypa/gh-action-pypi-publish@master diff --git a/.github/workflows/inference-workflow.yml b/.github/workflows/inference-workflow.yml index 16e5b9cc65d..68644bd663d 100644 --- a/.github/workflows/inference-workflow.yml +++ b/.github/workflows/inference-workflow.yml @@ -48,12 +48,12 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: logs path: gw_output/submitdir/work - name: store result page - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: results path: html diff --git a/.github/workflows/search-workflow.yml b/.github/workflows/search-workflow.yml index 18219411078..c8a5bb9f7bb 100644 --- a/.github/workflows/search-workflow.yml +++ b/.github/workflows/search-workflow.yml @@ -55,12 +55,12 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: logs path: output/submitdir/work - name: store result page - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: results path: html diff --git a/.github/workflows/tmpltbank-workflow.yml b/.github/workflows/tmpltbank-workflow.yml index bf29c9051d0..21917cebf7f 100644 --- a/.github/workflows/tmpltbank-workflow.yml +++ b/.github/workflows/tmpltbank-workflow.yml @@ -51,7 +51,7 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: logs path: output/submitdir/work diff --git a/.github/workflows/workflow-tests.yml b/.github/workflows/workflow-tests.yml index 754e25f8f11..3f39ef96ea0 100644 --- a/.github/workflows/workflow-tests.yml +++ b/.github/workflows/workflow-tests.yml @@ -50,7 +50,7 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: logs-${{matrix.test-type}} path: examples/workflow/generic/${{matrix.test-type}}/submitdir/work From ddaa94dec654b40428ea2116f34531c8b887a8fc Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Fri, 9 Feb 2024 16:51:26 +0000 Subject: [PATCH 04/29] Standardize logging in bin/inference (#4621) * inference: standardize verbose argument Only change scripts that already had a verbose argument * revert logging change in pycbc_inference_plot_gelman_rubin Partially reverts 9b7fa3d68e45a65a544f55d8a221ac44949ad6a5 --- bin/inference/pycbc_inference | 3 +-- bin/inference/pycbc_inference_extract_samples | 4 +--- bin/inference/pycbc_inference_model_stats | 3 +-- bin/inference/pycbc_inference_monitor | 2 +- bin/inference/pycbc_inference_plot_acceptance_rate | 4 +--- bin/inference/pycbc_inference_plot_acf | 4 +--- bin/inference/pycbc_inference_plot_acl | 3 +-- bin/inference/pycbc_inference_plot_dynesty_run | 4 +--- bin/inference/pycbc_inference_plot_dynesty_traceplot | 3 +-- bin/inference/pycbc_inference_plot_geweke | 6 ++---- bin/inference/pycbc_inference_plot_inj_recovery | 3 +-- bin/inference/pycbc_inference_plot_movie | 2 +- bin/inference/pycbc_inference_plot_posterior | 3 +-- bin/inference/pycbc_inference_plot_pp | 3 +-- bin/inference/pycbc_inference_plot_prior | 4 +--- bin/inference/pycbc_inference_plot_samples | 3 +-- bin/inference/pycbc_inference_pp_table_summary | 3 +-- bin/inference/pycbc_inference_table_summary | 3 +-- 18 files changed, 19 insertions(+), 41 deletions(-) diff --git a/bin/inference/pycbc_inference b/bin/inference/pycbc_inference index a4c8dd50649..d0cfd7b953a 100644 --- a/bin/inference/pycbc_inference +++ b/bin/inference/pycbc_inference @@ -40,10 +40,9 @@ from pycbc.workflow import configuration # command line usage parser = argparse.ArgumentParser(usage=__file__ + " [--options]", description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=__version__, help="Prints version information.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging messages.") # output options parser.add_argument("--output-file", type=str, required=True, help="Output file path.") diff --git a/bin/inference/pycbc_inference_extract_samples b/bin/inference/pycbc_inference_extract_samples index a404c6b5846..ffbf89ee996 100644 --- a/bin/inference/pycbc_inference_extract_samples +++ b/bin/inference/pycbc_inference_extract_samples @@ -71,7 +71,7 @@ def isthesame(current_val, val): parser = ResultsArgumentParser(defaultparams='all', autoparamlabels=False, description=__doc__) - +pycbc.add_common_pycbc_options(parser) parser.add_argument("--output-file", type=str, required=True, help="Output file to create.") parser.add_argument("--force", action="store_true", default=False, @@ -85,8 +85,6 @@ parser.add_argument("--skip-groups", default=None, nargs="+", "to write all groups if only one file is provided, " "and all groups from the first file except " "sampler_info if multiple files are provided.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Be verbose") opts = parser.parse_args() diff --git a/bin/inference/pycbc_inference_model_stats b/bin/inference/pycbc_inference_model_stats index 8a8512cb4b8..263b6868cd7 100644 --- a/bin/inference/pycbc_inference_model_stats +++ b/bin/inference/pycbc_inference_model_stats @@ -37,6 +37,7 @@ from pycbc.workflow import WorkflowConfigParser parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--input-file", type=str, required=True, help="Input HDF file to read. Can be either an inference " "file or a posterior file.") @@ -51,8 +52,6 @@ parser.add_argument("--nprocesses", type=int, default=1, parser.add_argument("--config-file", nargs="+", type=str, default=None, help="Override the config file stored in the input file " "with the given file(s).") -parser.add_argument("--verbose", action="count", default=0, - help="Print logging messages.") parser.add_argument('--reconstruct-parameters', action="store_true", help="Reconstruct marginalized parameters") diff --git a/bin/inference/pycbc_inference_monitor b/bin/inference/pycbc_inference_monitor index b4311762ca2..090bdf6a9ff 100644 --- a/bin/inference/pycbc_inference_monitor +++ b/bin/inference/pycbc_inference_monitor @@ -6,13 +6,13 @@ import os, sys, time, glob, shutil, logging, argparse, pycbc import pycbc.workflow as wf parser = argparse.ArgumentParser() +pycbc.add_common_pycbc_options(parser) parser.add_argument('--inference-file', help="The name of the inference file") parser.add_argument('--check-interval', default=10, type=int, help="Polling interval to check for file changes (s)") parser.add_argument('--output-dir', help="Output plots / results directory") parser.add_argument('--allow-failure', action='store_true', help="Allow for a failure in plot generation") -parser.add_argument('--verbose', action='store_true') wf.add_workflow_command_line_group(parser) args = parser.parse_args() diff --git a/bin/inference/pycbc_inference_plot_acceptance_rate b/bin/inference/pycbc_inference_plot_acceptance_rate index b3f1e9aafac..fe0b81555b4 100644 --- a/bin/inference/pycbc_inference_plot_acceptance_rate +++ b/bin/inference/pycbc_inference_plot_acceptance_rate @@ -34,6 +34,7 @@ parser = argparse.ArgumentParser( usage="pycbc_inference_plot_acceptance_rate [--options]", description="Plots histogram of the fractions of steps " "accepted by walkers.") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--input-file", type=str, required=True, help="Path to input HDF file.") parser.add_argument('--version', action='version', version=__version__, @@ -56,9 +57,6 @@ parser.add_argument("--temps", type=int, nargs="+", default=None, # add number of bins for histogram parser.add_argument("--bins", type=int, default=10, help="Specify number of bins for the histogram plot.") -# verbose option -parser.add_argument("--verbose", action="store_true", default=False, - help="") # parse the command line opts = parser.parse_args() diff --git a/bin/inference/pycbc_inference_plot_acf b/bin/inference/pycbc_inference_plot_acf index d38645be49c..41ab706395a 100644 --- a/bin/inference/pycbc_inference_plot_acf +++ b/bin/inference/pycbc_inference_plot_acf @@ -38,9 +38,7 @@ from pycbc.inference.sampler import samplers parser = io.ResultsArgumentParser(skip_args='thin-interval', description="Plots autocorrelation function " "from inference samples.") -# verbose option -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging info.") +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=__version__, help='show version number and exit') # output plot options diff --git a/bin/inference/pycbc_inference_plot_acl b/bin/inference/pycbc_inference_plot_acl index 348dbb10763..860a9810713 100644 --- a/bin/inference/pycbc_inference_plot_acl +++ b/bin/inference/pycbc_inference_plot_acl @@ -39,8 +39,7 @@ parser = io.ResultsArgumentParser(skip_args=['thin-interval', 'temps'], description="Histograms autocorrelation " "length per walker from an MCMC " "sampler.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging info.") +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=__version__, help='show version number and exit') # output plot options diff --git a/bin/inference/pycbc_inference_plot_dynesty_run b/bin/inference/pycbc_inference_plot_dynesty_run index 04a1aaf0371..6e2c8d03afb 100644 --- a/bin/inference/pycbc_inference_plot_dynesty_run +++ b/bin/inference/pycbc_inference_plot_dynesty_run @@ -34,6 +34,7 @@ parser = argparse.ArgumentParser( usage="pycbc_inference_plot_dynesty_run [--options]", description="Plots various figures showing evolution of " "dynesty run.") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--input-file", type=str, required=True, help="Path to input HDF file.") parser.add_argument('--version', action='version', version=__version__, @@ -42,9 +43,6 @@ parser.add_argument('--version', action='version', version=__version__, parser.add_argument("--output-file", type=str, required=True, help="Path to output plot.") -parser.add_argument("--verbose", action="store_true", default=False, - help="") - # parse the command line opts = parser.parse_args() diff --git a/bin/inference/pycbc_inference_plot_dynesty_traceplot b/bin/inference/pycbc_inference_plot_dynesty_traceplot index c328313ec94..c7e352f91b7 100644 --- a/bin/inference/pycbc_inference_plot_dynesty_traceplot +++ b/bin/inference/pycbc_inference_plot_dynesty_traceplot @@ -36,6 +36,7 @@ parser = argparse.ArgumentParser( description="Plots various figures showing evolution of " "dynesty run.") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--input-file", type=str, required=True, help="Path to input HDF file.") parser.add_argument('--version', action='version', version=__version__, @@ -44,8 +45,6 @@ parser.add_argument('--version', action='version', version=__version__, parser.add_argument("--output-file", type=str, required=True, help="Path to output plot.") -parser.add_argument("--verbose", action="store_true", default=False, - help="") # parse the command line opts = parser.parse_args() diff --git a/bin/inference/pycbc_inference_plot_geweke b/bin/inference/pycbc_inference_plot_geweke index 16c04463932..755b2f3268a 100644 --- a/bin/inference/pycbc_inference_plot_geweke +++ b/bin/inference/pycbc_inference_plot_geweke @@ -32,11 +32,9 @@ from pycbc.inference import (io, geweke, option_utils) # add options to command line parser = io.ResultsArgumentParser(skip_args=['walkers']) -# program-specific +pycbc.add_common_pycbc_options(parser) -# verbose option -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging info.") +# program-specific parser.add_argument('--version', action='version', version=__version__, help='show version number and exit') diff --git a/bin/inference/pycbc_inference_plot_inj_recovery b/bin/inference/pycbc_inference_plot_inj_recovery index b65dd77c4ac..fc51655bdd8 100644 --- a/bin/inference/pycbc_inference_plot_inj_recovery +++ b/bin/inference/pycbc_inference_plot_inj_recovery @@ -21,12 +21,11 @@ from pycbc.results import save_fig_with_metadata # parse command line parser = io.ResultsArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=__version__, help="Prints version information.") parser.add_argument("--output-file", required=True, type=str, help="Path to save output plot.") -parser.add_argument("--verbose", action="store_true", - help="Allows print statements.") parser.add_argument("--percentiles", nargs=2, type=float, default=[5, 95], help="Percentiles to use as limits.") option_utils.add_scatter_option_group(parser) diff --git a/bin/inference/pycbc_inference_plot_movie b/bin/inference/pycbc_inference_plot_movie index edc04be8c00..1f3a4b83c88 100644 --- a/bin/inference/pycbc_inference_plot_movie +++ b/bin/inference/pycbc_inference_plot_movie @@ -110,6 +110,7 @@ def integer_logspace(start, end, num): # the frame number/step options skip_args = ['thin-start', 'thin-interval', 'thin-end', 'iteration'] parser = io.ResultsArgumentParser(description=__doc__, skip_args=skip_args) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=__version__, help="show version number and exit") # make frame number and frame step mutually exclusive @@ -135,7 +136,6 @@ parser.add_argument("--log-steps", action="store_true", default=False, parser.add_argument("--output-prefix", type=str, required=True, help="Output path and prefix for the frame files " "(without extension).") -parser.add_argument('--verbose', action='store_true') parser.add_argument('--dpi', type=int, default=200, help='Set the dpi for each frame; default is 200') parser.add_argument("--nprocesses", type=int, default=None, diff --git a/bin/inference/pycbc_inference_plot_posterior b/bin/inference/pycbc_inference_plot_posterior index e9df4e79048..7d507a7ab92 100644 --- a/bin/inference/pycbc_inference_plot_posterior +++ b/bin/inference/pycbc_inference_plot_posterior @@ -52,13 +52,12 @@ use('agg') # add options to command line parser = io.ResultsArgumentParser() +pycbc.add_common_pycbc_options(parser) # program-specific parser.add_argument("--version", action="version", version=__version__, help="Prints version information.") parser.add_argument("--output-file", type=str, required=True, help="Output plot path.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Be verbose") parser.add_argument("--plot-prior", nargs="+", type=str, help="Plot the prior on the 1D marginal plots using the " "given config file(s).") diff --git a/bin/inference/pycbc_inference_plot_pp b/bin/inference/pycbc_inference_plot_pp index f976a1f5317..5e6490499c8 100644 --- a/bin/inference/pycbc_inference_plot_pp +++ b/bin/inference/pycbc_inference_plot_pp @@ -42,12 +42,11 @@ from pycbc import __version__ # parse command line parser = io.ResultsArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=__version__, help='show version number and exit') parser.add_argument("--output-file", required=True, type=str, help="Path to save output plot.") -parser.add_argument("--verbose", action="store_true", - help="Allows print statements.") parser.add_argument("--injection-hdf-group", default="injections", help="HDF group that contains injection values. " "Default is 'injections'.") diff --git a/bin/inference/pycbc_inference_plot_prior b/bin/inference/pycbc_inference_plot_prior index 73d39b2d481..436d0c6f0e6 100644 --- a/bin/inference/pycbc_inference_plot_prior +++ b/bin/inference/pycbc_inference_plot_prior @@ -46,6 +46,7 @@ def cartesian(arrays): parser = argparse.ArgumentParser( usage="pycbc_inference_plot_prior [--options]", description="Plots prior distributions.") +pycbc.add_common_pycbc_options(parser) # add input options parser.add_argument("--config-files", type=str, nargs="+", required=True, help="A file parsable by " @@ -75,9 +76,6 @@ parser.add_argument("--nsamples", type=int, default=10000, "plotting. Default is 10000.") parser.add_argument("--output-file", type=str, required=True, help="Path to output plot.") -# verbose option -parser.add_argument("--verbose", action="store_true", default=False, - help="") parser.add_argument("--version", action="version", version=__version__, help="show version number and exit") diff --git a/bin/inference/pycbc_inference_plot_samples b/bin/inference/pycbc_inference_plot_samples index 4ece1a64ac0..651c79e5178 100644 --- a/bin/inference/pycbc_inference_plot_samples +++ b/bin/inference/pycbc_inference_plot_samples @@ -34,8 +34,7 @@ import sys # command line usage parser = argparse.parser = io.ResultsArgumentParser( skip_args=['chains', 'iteration']) -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging info.") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=__version__, help="show version number and exit") parser.add_argument("--chains", nargs='+', default=None, diff --git a/bin/inference/pycbc_inference_pp_table_summary b/bin/inference/pycbc_inference_pp_table_summary index 25d41539d79..23195eb8495 100644 --- a/bin/inference/pycbc_inference_pp_table_summary +++ b/bin/inference/pycbc_inference_pp_table_summary @@ -31,8 +31,7 @@ from pycbc.inference.io import (ResultsArgumentParser, results_from_cli, injections_from_cli) parser = ResultsArgumentParser(description=__doc__) -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging info.") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--output-file", type=str, required=True, help="Path to output plot.") add_injsamples_map_opt(parser) diff --git a/bin/inference/pycbc_inference_table_summary b/bin/inference/pycbc_inference_table_summary index 261715e35da..6e29fb636ff 100644 --- a/bin/inference/pycbc_inference_table_summary +++ b/bin/inference/pycbc_inference_table_summary @@ -29,8 +29,7 @@ from pycbc.inference.io import ResultsArgumentParser, results_from_cli parser = ResultsArgumentParser( description="Makes a table of posterior results.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging info.") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--output-file", type=str, required=True, help="Path to output plot.") parser.add_argument("--print-metadata", nargs="+", metavar="PARAM[:LABEL]", From e2c35eff8367fa4700588a60f80737949d8fd254 Mon Sep 17 00:00:00 2001 From: Jacob Buchanan <83317255+jakeb245@users.noreply.github.com> Date: Sat, 10 Feb 2024 13:46:53 -0500 Subject: [PATCH 05/29] Sky error fix (#4502) --- bin/pygrb/pycbc_pygrb_plot_injs_results | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/bin/pygrb/pycbc_pygrb_plot_injs_results b/bin/pygrb/pycbc_pygrb_plot_injs_results index 22915b945cd..e79494f23e3 100644 --- a/bin/pygrb/pycbc_pygrb_plot_injs_results +++ b/bin/pygrb/pycbc_pygrb_plot_injs_results @@ -100,7 +100,7 @@ def load_mass_data(input_file_handle, key, tag): return data -def load_sky_error_data(input_file_handle, key, tag, trig_data): +def load_sky_error_data(input_file_handle, key, tag): """Extract data related to sky_error from raw injection and trigger data""" if tag == 'missed': # Missed injections are assigned null values @@ -110,8 +110,8 @@ def load_sky_error_data(input_file_handle, key, tag, trig_data): inj['ra'] = input_file_handle[tag+'/ra'][:] inj['dec'] = input_file_handle[tag+'/dec'][:] trig = {} - trig['ra'] = np.rad2deg(trig_data['network/longitude'][:]) - trig['dec'] = np.rad2deg(trig_data['network/latitude'][:]) + trig['ra'] = input_file_handle['network/ra'][:] + trig['dec'] = input_file_handle['network/dec'][:] data = np.arccos(np.cos(inj['dec'] - trig['dec']) - np.cos(inj['dec']) * np.cos(trig['dec']) * (1 - np.cos(inj['ra'] - trig['ra']))) @@ -128,7 +128,7 @@ easy_keys = ['distance', 'mass1', 'mass2', 'polarization', 'spin2_azimuthal', 'spin2_polar', 'dec', 'ra', 'phi_ref'] # Function to extract desired data from an injection file -def load_data(input_file_handle, keys, tag, trig_data): +def load_data(input_file_handle, keys, tag): """Create a dictionary containing the data specified by the list of keys extracted from an injection file""" @@ -148,7 +148,7 @@ def load_data(input_file_handle, keys, tag, trig_data): data_dict[key] = load_incl_data(input_file_handle, key, tag) elif key == 'sky_error': data_dict[key] = load_sky_error_data(input_file_handle, key, - tag, trig_data) + tag) except KeyError: #raise NotImplemented(key+' not allowed: returning empty entry') logging.info(key+' not allowed yet') @@ -247,10 +247,10 @@ max_reweighted_snr = max(trig_data['network/reweighted_snr'][:]) # Post-process injections # ======================= # Extract the necessary data from the missed injections for the plot -missed_inj = load_data(f, [x_qty, y_qty], 'missed', trig_data) +missed_inj = load_data(f, [x_qty, y_qty], 'missed') # Extract the necessary data from the found injections for the plot -found_inj = load_data(f, [x_qty, y_qty], 'found', trig_data) +found_inj = load_data(f, [x_qty, y_qty], 'found') # Injections found surviving vetoes found_after_vetoes = found_inj if 'found_after_vetoes' not in f.keys() else f['found_after_vetoes'] From 0f43c7b1218e8dbfebfd6e5a6e27d061108bb112 Mon Sep 17 00:00:00 2001 From: Tito Dal Canton Date: Tue, 13 Feb 2024 12:48:09 +0100 Subject: [PATCH 06/29] Fix broken --version option (#4630) * Fix broken --version option, add CI test, cleanup and simplify related code * Fix bug preventing .so files from being found * Fix missing import, typo, more codeclimate * pycbc_coinc_time needs dqsegdb; deduplicate --help commands --- pycbc/_version.py | 79 ++++++++++++++++++++++++++------------- pycbc/_version_helper.py | 70 +++++++++++++++++----------------- tools/pycbc_test_suite.sh | 14 ++++++- 3 files changed, 102 insertions(+), 61 deletions(-) diff --git a/pycbc/_version.py b/pycbc/_version.py index 31ebf7f4fea..617cf4a523e 100644 --- a/pycbc/_version.py +++ b/pycbc/_version.py @@ -13,27 +13,40 @@ # You should have received a copy of the GNU General Public License along # with this program; if not, write to the Free Software Foundation, Inc., # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + """ This modules contains a function to provide an argparse action that reports extremely verbose version information for PyCBC, lal, and lalsimulation. """ -import os, sys +import os +import sys +import glob import argparse import inspect import subprocess + def print_link(library): - err_msg = "Could not execute runtime linker to determine\n" + \ - "shared library paths for library:\n " + library + "\n" - FNULL = open(os.devnull, 'w') + err_msg = ( + "Could not execute runtime linker to determine\n" + f"shared library paths for library:\n {library}\n" + ) try: - link = subprocess.check_output(['ldd', library], - stderr=FNULL) + # Linux + link = subprocess.check_output( + ['ldd', library], + stderr=subprocess.DEVNULL, + text=True + ) except OSError: try: - link = subprocess.check_output(['otool', '-L', library], - stderr=FNULL) + # macOS + link = subprocess.check_output( + ['otool', '-L', library], + stderr=subprocess.DEVNULL, + text=True + ) except: link = err_msg except: @@ -41,42 +54,58 @@ def print_link(library): return link +def get_lal_info(module, lib_glob): + """Return a string reporting the version and runtime library information + for a LAL Python import. + """ + module_path = inspect.getfile(module) + version_str = ( + module.git_version.verbose_msg + + "\n\nImported from: " + module_path + + "\n\nRuntime libraries:\n" + ) + possible_lib_paths = glob.glob( + os.path.join(os.path.dirname(module_path), lib_glob) + ) + for lib_path in possible_lib_paths: + version_str += print_link(lib_path) + return version_str + + class Version(argparse.Action): - """ print the pycbc, lal and lalsimulation versions """ + """Subclass of argparse.Action that prints version information for PyCBC, + LAL and LALSimulation. + """ def __init__(self, nargs=0, **kw): super(Version, self).__init__(nargs=nargs, **kw) - def __call__(self, parser, namespace, values, option_string=None): - import pycbc - version_str="--- PyCBC Version --------------------------\n" + \ - pycbc.version.git_verbose_msg + \ + + version_str = ( + "--- PyCBC Version --------------------------\n" + + pycbc.version.git_verbose_msg + "\n\nImported from: " + inspect.getfile(pycbc) + ) version_str += "\n\n--- LAL Version ----------------------------\n" try: import lal.git_version - lal_module = inspect.getfile(lal) - lal_library = os.path.join( os.path.dirname(lal_module), - '_lal.so') - version_str += lal.git_version.verbose_msg + \ - "\n\nImported from: " + lal_module + \ - "\n\nRuntime libraries:\n" + print_link(lal_library) except ImportError: version_str += "\nLAL not installed in environment\n" + else: + version_str += get_lal_info(lal, '_lal*.so') version_str += "\n\n--- LALSimulation Version-------------------\n" try: import lalsimulation.git_version - lalsimulation_module = inspect.getfile(lalsimulation) - lalsimulation_library = os.path.join( os.path.dirname(lalsimulation_module), - '_lalsimulation.so') - version_str += lalsimulation.git_version.verbose_msg + \ - "\n\nImported from: " + lalsimulation_module + \ - "\n\nRuntime libraries:\n" + print_link(lalsimulation_library) except ImportError: version_str += "\nLALSimulation not installed in environment\n" + else: + version_str += get_lal_info(lalsimulation, '_lalsimulation*.so') print(version_str) sys.exit(0) + + +__all__ = ['Version'] diff --git a/pycbc/_version_helper.py b/pycbc/_version_helper.py index 694895a4b93..4dde5936e5c 100644 --- a/pycbc/_version_helper.py +++ b/pycbc/_version_helper.py @@ -64,30 +64,29 @@ def call(command, stdout=subprocess.PIPE, stderr=subprocess.PIPE, if p.returncode != 0 and on_error == 'raise': raise GitInvocationError('Failed to run "%s"' % " ".join(command)) - out = out.decode('utf-8') + out = out.decode('utf-8').strip() if returncode: - return out.strip(), p.returncode - else: - return out.strip() + return out, p.returncode + return out def get_build_name(git_path='git'): """Find the username of the current builder """ - name,retcode = call(('git', 'config', 'user.name'), returncode=True) + name, retcode = call(('git', 'config', 'user.name'), returncode=True) if retcode: name = "Unknown User" - email,retcode = call(('git', 'config', 'user.email'), returncode=True) + email, retcode = call(('git', 'config', 'user.email'), returncode=True) if retcode: email = "" - return "%s <%s>" % (name, email) + return f"{name} <{email}>" def get_build_date(): """Returns the current datetime as the git build date """ - return time.strftime('%Y-%m-%d %H:%M:%S +0000', time.gmtime()) + return time.strftime(r'%Y-%m-%d %H:%M:%S +0000', time.gmtime()) def get_last_commit(git_path='git'): @@ -96,12 +95,15 @@ def get_last_commit(git_path='git'): Returns a tuple (hash, date, author name, author e-mail, committer name, committer e-mail). """ - hash_, udate, aname, amail, cname, cmail = ( - call((git_path, 'log', '-1', - '--pretty=format:%H,%ct,%an,%ae,%cn,%ce')).split(",")) - date = time.strftime('%Y-%m-%d %H:%M:%S +0000', time.gmtime(float(udate))) - author = '%s <%s>' % (aname, amail) - committer = '%s <%s>' % (cname, cmail) + hash_, udate, aname, amail, cname, cmail = call(( + git_path, + 'log', + '-1', + r'--pretty=format:%H,%ct,%an,%ae,%cn,%ce' + )).split(",") + date = time.strftime(r'%Y-%m-%d %H:%M:%S +0000', time.gmtime(float(udate))) + author = f'{aname} <{amail}>' + committer = f'{cname} <{cmail}>' return hash_, date, author, committer @@ -111,8 +113,7 @@ def get_git_branch(git_path='git'): branch_match = call((git_path, 'rev-parse', '--symbolic-full-name', 'HEAD')) if branch_match == "HEAD": return None - else: - return os.path.basename(branch_match) + return os.path.basename(branch_match) def get_git_tag(hash_, git_path='git'): @@ -122,26 +123,26 @@ def get_git_tag(hash_, git_path='git'): '--tags', hash_), returncode=True) if status == 0: return tag - else: - return None + return None + def get_num_commits(): return call(('git', 'rev-list', '--count', 'HEAD')) + def get_git_status(git_path='git'): """Returns the state of the git working copy """ status_output = subprocess.call((git_path, 'diff-files', '--quiet')) if status_output != 0: return 'UNCLEAN: Modified working tree' - else: - # check index for changes - status_output = subprocess.call((git_path, 'diff-index', '--cached', - '--quiet', 'HEAD')) - if status_output != 0: - return 'UNCLEAN: Modified index' - else: - return 'CLEAN: All modifications committed' + # check index for changes + status_output = subprocess.call((git_path, 'diff-index', '--cached', + '--quiet', 'HEAD')) + if status_output != 0: + return 'UNCLEAN: Modified index' + return 'CLEAN: All modifications committed' + def determine_latest_release_version(): """Query the git repository for the last released version of the code. @@ -152,22 +153,21 @@ def determine_latest_release_version(): tag_list = call((git_path, 'tag')).split('\n') # Reduce to only versions - tag_list = [t[1:] for t in tag_list if t.startswith('v')] + re_magic = re.compile(r"v\d+\.\d+\.\d+$") + tag_list = [t[1:] for t in tag_list if re_magic.match(t)] - # Determine if indeed a tag and store largest + # find latest version latest_version = None latest_version_string = None - re_magic = re.compile("\d+\.\d+\.\d+$") for tag in tag_list: - # Is this a version string - if re_magic.match(tag): - curr_version = distutils.version.StrictVersion(tag) - if latest_version is None or curr_version > latest_version: - latest_version = curr_version - latest_version_string = tag + curr_version = distutils.version.StrictVersion(tag) + if latest_version is None or curr_version > latest_version: + latest_version = curr_version + latest_version_string = tag return latest_version_string + def generate_git_version_info(): """Query the git repository information to generate a version module. """ diff --git a/tools/pycbc_test_suite.sh b/tools/pycbc_test_suite.sh index 08ef955dc49..6e1f745194e 100755 --- a/tools/pycbc_test_suite.sh +++ b/tools/pycbc_test_suite.sh @@ -31,7 +31,7 @@ fi if [ "$PYCBC_TEST_TYPE" = "help" ] || [ -z ${PYCBC_TEST_TYPE+x} ]; then # check that all executables that do not require # special environments can return a help message - for prog in `find ${PATH//:/ } -maxdepth 1 -name 'pycbc*' -print 2>/dev/null | egrep -v '(pycbc_live_nagios_monitor|pycbc_make_offline_grb_workflow|pycbc_mvsc_get_features|pycbc_upload_xml_to_gracedb|pycbc_coinc_time)'` + for prog in `find ${PATH//:/ } -maxdepth 1 -name 'pycbc*' -print 2>/dev/null | egrep -v '(pycbc_live_nagios_monitor|pycbc_make_offline_grb_workflow|pycbc_mvsc_get_features|pycbc_upload_xml_to_gracedb|pycbc_coinc_time)' | sort | uniq` do echo -e ">> [`date`] running $prog --help" $prog --help &> $LOG_FILE @@ -45,6 +45,18 @@ if [ "$PYCBC_TEST_TYPE" = "help" ] || [ -z ${PYCBC_TEST_TYPE+x} ]; then echo -e " Pass." fi done + # also check that --version works for one executable + echo -e ">> [`date`] running pycbc_inspiral --version" + pycbc_inspiral --version &> $LOG_FILE + if test $? -ne 0 ; then + RESULT=1 + echo -e " FAILED!" + echo -e "---------------------------------------------------------" + cat $LOG_FILE + echo -e "---------------------------------------------------------" + else + echo -e " Pass." + fi fi if [ "$PYCBC_TEST_TYPE" = "search" ] || [ -z ${PYCBC_TEST_TYPE+x} ]; then From c272f833f03ed5f04a99f6579b6f00edb6d759ac Mon Sep 17 00:00:00 2001 From: Arthur Tolley <32394213+ArthurTolley@users.noreply.github.com> Date: Wed, 14 Feb 2024 18:14:34 +0000 Subject: [PATCH 07/29] [pycbc_live] Simplify and fix how CLI options are passed to the SNR optimizer (#4628) * adding non-optimizer specific options to args_to_string * correct formatting * simplifying if snr opt seed * Adding extra-opts arg * updating options in live example run.sh * restoring deleted space * removing redundant default * moving all snr optimizer options to snr_opt_extra_opts * updating argument help descriptions * removing snr_opt options from pycbc live * removing seed option from example * removing args_to_string * Actually, even simpler --------- Co-authored-by: Tito Dal Canton --- bin/pycbc_live | 12 ++++++------ examples/live/run.sh | 22 +++++++++++++--------- pycbc/live/snr_optimizer.py | 31 +++++-------------------------- 3 files changed, 24 insertions(+), 41 deletions(-) diff --git a/bin/pycbc_live b/bin/pycbc_live index 646d4bebf40..8a808ca5a4f 100755 --- a/bin/pycbc_live +++ b/bin/pycbc_live @@ -106,6 +106,7 @@ class LiveEventManager(object): self.enable_gracedb_upload = args.enable_gracedb_upload self.run_snr_optimization = args.run_snr_optimization self.snr_opt_label = args.snr_opt_label + self.snr_opt_options = args.snr_opt_extra_opts self.gracedb = None # Keep track of which events have been uploaded @@ -121,7 +122,6 @@ class LiveEventManager(object): if platform.system() != 'Darwin': available_cores = len(os.sched_getaffinity(0)) self.fu_cores = available_cores - analysis_cores - self.optimizer = args.snr_opt_method if self.fu_cores <= 0: logging.warning( 'Insufficient number of CPU cores (%d) to ' @@ -133,10 +133,6 @@ class LiveEventManager(object): else: # To enable mac testing, this is just set to 1 self.fu_cores = 1 - # Convert SNR optimizer options into a string - self.snr_opt_options = snr_optimizer.args_to_string(args) - else: - self.snr_opt_options = None if args.enable_embright_has_massgap: if args.embright_massgap_max < self.mc_area_args['mass_bdary']['ns_max']: @@ -987,6 +983,11 @@ parser.add_argument('--snr-opt-timeout', type=int, default=400, metavar='SECONDS help='Maximum allowed duration of followup process to maximize SNR') parser.add_argument('--snr-opt-label', default='SNR_OPTIMIZED', help='Label to apply to snr-optimized GraceDB uploads') +parser.add_argument('--snr-opt-extra-opts', + help='Extra options to pass to the optimizer subprocess. Example: ' + '--snr-opt-extra-opts "--snr-opt-method differential_evolution ' + '--snr-opt-di-maxiter 50 --snr-opt-di-popsize 100 ' + '--snr-opt-seed 42 --snr-opt-include-candidate "') parser.add_argument('--enable-embright-has-massgap', action='store_true', default=False, help='Estimate HasMassGap probability for EMBright info. Lower limit ' @@ -1010,7 +1011,6 @@ Coincer.insert_args(parser) SingleDetSGChisq.insert_option_group(parser) mchirp_area.insert_args(parser) livepau.insert_live_pastro_option_group(parser) -snr_optimizer.insert_snr_optimizer_options(parser) args = parser.parse_args() diff --git a/examples/live/run.sh b/examples/live/run.sh index 82b79db5a24..8857c5d996e 100755 --- a/examples/live/run.sh +++ b/examples/live/run.sh @@ -195,10 +195,11 @@ python -m mpi4py `which pycbc_live` \ --src-class-eff-to-lum-distance 0.74899 \ --src-class-lum-distance-to-delta -0.51557 -0.32195 \ --run-snr-optimization \ ---snr-opt-di-maxiter 50 \ ---snr-opt-di-popsize 100 \ ---snr-opt-include-candidate \ ---snr-opt-seed 42 \ +--snr-opt-extra-opts \ + "--snr-opt-method differential_evolution \ + --snr-opt-di-maxiter 50 \ + --snr-opt-di-popsize 100 \ + --snr-opt-include-candidate " \ --sngl-ifar-est-dist conservative \ --single-newsnr-threshold 9 \ --single-duration-threshold 7 \ @@ -210,11 +211,14 @@ python -m mpi4py `which pycbc_live` \ # If you would like to use the pso optimizer, change --optimizer to pso # and include these arguments while removing other optimizer args. # You will need to install the pyswarms package into your environment. -# --snr-opt-pso-iters 5 \ -# --snr-opt-pso-particles 250 \ -# --snr-opt-pso-c1 0.5 \ -# --snr-opt-pso-c2 2.0 \ -# --snr-opt-pso-w 0.01 \ +# --snr-opt-extra-opts \ +# "--snr-opt-method pso \ +# --snr-opt-pso-iters 5 \ +# --snr-opt-pso-particles 250 \ +# --snr-opt-pso-c1 0.5 \ +# --snr-opt-pso-c2 2.0 \ +# --snr-opt-pso-w 0.01 \ +# --snr-opt-include-candidate " \ # note that, at this point, some SNR optimization processes may still be # running, so the checks below may ignore their results diff --git a/pycbc/live/snr_optimizer.py b/pycbc/live/snr_optimizer.py index 9e6b80481a3..34ab4cd2de2 100644 --- a/pycbc/live/snr_optimizer.py +++ b/pycbc/live/snr_optimizer.py @@ -328,13 +328,13 @@ def insert_snr_optimizer_options(parser): opt_opt_group.add_argument('--snr-opt-include-candidate', action='store_true', help='Include parameters of the candidate event in the initialized ' - 'array for the optimizer. Only relevant for --optimizer pso or ' - 'differential_evolution') + 'array for the optimizer. Only relevant for --snr-opt-method pso ' + 'or differential_evolution') opt_opt_group.add_argument('--snr-opt-seed', default='42', help='Seed to supply to the random generation of initial array to ' - 'pass to the optimizer. Only relevant for --optimizer pso or ' - 'differential_evolution. Set to ''random'' for a random seed') + 'pass to the optimizer. Only relevant for --snr-opt-method pso ' + 'or differential_evolution. Set to ''random'' for a random seed') # For each optimizer, add the possible options for optimizer, option_subdict in option_dict.items(): @@ -345,7 +345,7 @@ def insert_snr_optimizer_options(parser): option_name = f"--snr-opt-{optimizer_name}-{opt_name}" opt_opt_group.add_argument(option_name, type=float, - help=f'Only relevant for --optimizer {optimizer}: ' + + help=f'Only relevant for --snr-opt-method {optimizer}: ' + opt_help_default[0] + f' Default = {opt_help_default[1]}') @@ -381,24 +381,3 @@ def check_snr_optimizer_options(args, parser): key_name = f'snr_opt_{optimizer_name}_{key}' if not getattr(args, key_name): setattr(args, key_name, value[1]) - - -def args_to_string(args): - """ - Convert the supplied arguments for SNR optimization config into - a string - this is to be used when running subprocesses - """ - argstr = f'--snr-opt-method {args.snr_opt_method} ' - optimizer_name = args.snr_opt_method.replace('_', '-') - if optimizer_name == 'differential-evolution': - optimizer_name = 'di' - for opt in option_dict[args.snr_opt_method]: - key_name = f'snr_opt_{optimizer_name}_{opt}' - option_value = getattr(args, key_name) - # If the option is not given, don't pass it and use default - if option_value is None: - continue - option_fullname = f'--snr-opt-{optimizer_name}-{opt}' - argstr += f'{option_fullname} {option_value} ' - - return argstr From 0b11dfc32b99e3e24f3d3a0c64e0b4ea45a24e9e Mon Sep 17 00:00:00 2001 From: maxtrevor <65971534+maxtrevor@users.noreply.github.com> Date: Wed, 14 Feb 2024 22:23:49 -0500 Subject: [PATCH 08/29] Added tables to dq results pages (#4633) * Added tables to dq results pages * Address Tom's comments and fix bin_trigger_rates_dq * Fix decimal formatting * Final change to formatting --- bin/all_sky_search/pycbc_bin_trigger_rates_dq | 16 ++++ bin/plotting/pycbc_page_dq_table | 54 +++++++++++++ bin/plotting/pycbc_page_template_bin_table | 75 +++++++++++++++++++ bin/plotting/pycbc_plot_dq_flag_likelihood | 5 +- .../pycbc_make_offline_search_workflow | 16 ++-- pycbc/workflow/plotting.py | 29 ++++++- 6 files changed, 186 insertions(+), 9 deletions(-) create mode 100644 bin/plotting/pycbc_page_dq_table create mode 100644 bin/plotting/pycbc_page_template_bin_table diff --git a/bin/all_sky_search/pycbc_bin_trigger_rates_dq b/bin/all_sky_search/pycbc_bin_trigger_rates_dq index 25234a687f0..ef069d9a6b1 100644 --- a/bin/all_sky_search/pycbc_bin_trigger_rates_dq +++ b/bin/all_sky_search/pycbc_bin_trigger_rates_dq @@ -46,6 +46,18 @@ logging.info('Start') ifo, flag_name = args.flag_name.split(':') +if args.gating_windows: + gate_times = [] + with h5.File(args.trig_file, 'r') as trig_file: + logging.info('Getting gated times') + try: + gating_types = trig_file[f'{ifo}/gating'].keys() + for gt in gating_types: + gate_times += list(trig_file[f'{ifo}/gating/{gt}/time'][:]) + gate_times = np.unique(gate_times) + except KeyError: + logging.warning('No gating found in trigger file') + trigs = SingleDetTriggers( args.trig_file, ifo, @@ -134,6 +146,7 @@ with h5.File(args.output_file, 'w') as f: frac_dt = abs(segs) / livetime dq_rates[state] = frac_eff / frac_dt bin_grp['dq_rates'] = dq_rates + bin_grp['num_triggers'] = len(trig_times_bin) # save dq state segments for dq_state, segs in dq_state_segs_dict.items(): @@ -142,7 +155,10 @@ with h5.File(args.output_file, 'w') as f: starts, ends = segments_to_start_end(segs) dq_grp['segment_starts'] = starts dq_grp['segment_ends'] = ends + dq_grp['livetime'] = abs(segs) f.attrs['stat'] = f'{ifo}-dq_stat_info' + f.attrs['sngl_ranking'] = args.sngl_ranking + f.attrs['sngl_ranking_threshold'] = args.stat_threshold logging.info('Done!') diff --git a/bin/plotting/pycbc_page_dq_table b/bin/plotting/pycbc_page_dq_table new file mode 100644 index 00000000000..83b39b01569 --- /dev/null +++ b/bin/plotting/pycbc_page_dq_table @@ -0,0 +1,54 @@ +#!/usr/bin/env python +""" Make a table of dq state information +""" +import sys +import argparse +import h5py as h5 +import numpy as np + +import pycbc +import pycbc.results +from pycbc.version import git_verbose_msg as version + +parser = argparse.ArgumentParser() +parser.add_argument("--version", action="version", version=version) +parser.add_argument('--ifo', required=True) +parser.add_argument('--dq-file', required=True) +parser.add_argument('--verbose', action='count') +parser.add_argument('--output-file') +args = parser.parse_args() + +dq_states = { + 'dq_state_0': 'Clean', + 'dq_state_1': 'DQ Flag', + 'dq_state_2': 'Autogating', +} + +f = h5.File(args.dq_file, 'r') +grp = f[args.ifo]['dq_segments'] + +livetimes = [] +total_livetime = 0 +for dq_state in dq_states: + livetime = grp[dq_state]['livetime'][()] + livetimes.append(livetime) + total_livetime += livetime +livetimes.append(total_livetime) + +frac_livetimes = [lt / total_livetime for lt in livetimes] +state_names = list(dq_states.values()) + ['Total'] +columns = [state_names, livetimes, frac_livetimes] +columns = [np.array(c) for c in columns] +col_names = ['DQ State', 'Livetime', '% of Livetime'] + +format_strings = [None, '0.0', '0.00%'] + +html_table = pycbc.results.html_table(columns, col_names, + page_size=len(state_names), + format_strings=format_strings) +title = f'{args.ifo} DQ State Livetimes' +caption = 'Table of DQ state livetimes' + +pycbc.results.save_fig_with_metadata( + str(html_table), args.output_file, title=title, + caption=caption, cmd=' '.join(sys.argv)) diff --git a/bin/plotting/pycbc_page_template_bin_table b/bin/plotting/pycbc_page_template_bin_table new file mode 100644 index 00000000000..ef2f7ffb2aa --- /dev/null +++ b/bin/plotting/pycbc_page_template_bin_table @@ -0,0 +1,75 @@ +#!/usr/bin/env python +""" Make a table of template bin information +""" +import sys +import argparse +import h5py as h5 +import numpy as np + +import pycbc +import pycbc.results +from pycbc.version import git_verbose_msg as version + +parser = argparse.ArgumentParser() +parser.add_argument("--version", action="version", version=version) +parser.add_argument('--ifo', required=True) +parser.add_argument('--dq-file', required=True) +parser.add_argument('--verbose', action='count') +parser.add_argument('--output-file') +args = parser.parse_args() + +f = h5.File(args.dq_file, 'r') +grp = f[args.ifo]['bins'] +bin_names = list(grp.keys()) + +sngl_ranking = f.attrs['sngl_ranking'] +sngl_thresh = f.attrs['sngl_ranking_threshold'] + +livetime = 0 +seg_grp = f[args.ifo]['dq_segments'] +for k in seg_grp.keys(): + livetime += seg_grp[k]['livetime'][()] + +num_templates = [] +num_triggers = [] +total_templates = 0 +total_triggers = 0 +for bin_name in bin_names: + bin_grp = grp[bin_name] + + n_tmp = len(bin_grp['tids'][:]) + num_templates.append(n_tmp) + total_templates += n_tmp + + n_trig = bin_grp['num_triggers'][()] + num_triggers.append(n_trig) + total_triggers += n_trig + +bin_names.append('Total') +num_triggers.append(total_triggers) +num_templates.append(total_templates) +frac_triggers = [n / total_triggers for n in num_triggers] +frac_templates = [n / total_templates for n in num_templates] +trigger_rate = [n / livetime for n in num_triggers] + +col_names = ['Template Bin', 'Number of Templates', '% of Templates', + 'Number of Loud Triggers', '% of Loud Triggers', + 'Loud Trigger Rate (Hz)'] +columns = [bin_names, num_templates, frac_templates, + num_triggers, frac_triggers, trigger_rate] +columns = [np.array(c) for c in columns] + +format_strings = [None, '#', '0.000%', '#', '0.0%', '0.00E0'] + +html_table = pycbc.results.html_table(columns, col_names, + page_size=len(bin_names), + format_strings=format_strings) +title = f'{args.ifo} DQ Template Bin Information' +caption = 'Table of information about template bins ' \ + + 'used for DQ trigger rate calculations. ' \ + + 'Loud triggers are defined as those with ' \ + + f'{sngl_ranking} > {sngl_thresh}.' + +pycbc.results.save_fig_with_metadata( + str(html_table), args.output_file, title=title, + caption=caption, cmd=' '.join(sys.argv)) diff --git a/bin/plotting/pycbc_plot_dq_flag_likelihood b/bin/plotting/pycbc_plot_dq_flag_likelihood index a1d7beed5c1..ed4a6a1231f 100644 --- a/bin/plotting/pycbc_plot_dq_flag_likelihood +++ b/bin/plotting/pycbc_plot_dq_flag_likelihood @@ -77,7 +77,10 @@ ax2.set_ylabel('DQ Log Likelihood Penalty') ax2.set_ylim(numpy.log(ymin), numpy.log(ymax)) # add meta data and save figure -plot_title = f'{ifo} DQ Trigger Rates' +plot_title = f'{ifo}:{args.dq_label} DQ Trigger Rates' + +ax.set_title(plot_title) + plot_caption = 'The log likelihood correction \ during during each dq state for each template bin.' pycbc.results.save_fig_with_metadata( diff --git a/bin/workflows/pycbc_make_offline_search_workflow b/bin/workflows/pycbc_make_offline_search_workflow index a7fc7daf5c8..ac96e5a9815 100755 --- a/bin/workflows/pycbc_make_offline_search_workflow +++ b/bin/workflows/pycbc_make_offline_search_workflow @@ -561,19 +561,25 @@ for key in final_bg_files: # DQ log likelihood plots for dqf, dql in zip(dqfiles, dqfile_labels): + ifo = dqf.ifo + outdir = rdir[f'data_quality/{dqf.ifo}_DQ_results'] + # plot rates when flag was on - outdir = rdir[f'data_quality/{dqf.ifo}_{dql}_trigger_rates'] wf.make_dq_flag_trigger_rate_plot(workflow, dqf, dql, outdir, tags=[dql]) -# DQ background binning plot -ifo_bins = workflow.cp.get_subsections('bin_templates') -for ifo in ifo_bins: + # make table of dq segment info + wf.make_dq_segment_table(workflow, dqf, outdir, tags=[dql]) + + # plot background bins background_bins = workflow.cp.get_opt_tags( 'bin_templates', 'background-bins', tags=[ifo]) - wf.make_template_plot(workflow, hdfbank, rdir['data_quality'], + wf.make_template_plot(workflow, hdfbank, outdir, bins=background_bins, tags=[ifo, 'dq_bins'] + bank_tags) + # make table of template bin info + wf.make_template_bin_table(workflow, dqf, outdir, tags=[ifo]) + ############################## Setup the injection runs ####################### splitbank_files_inj = wf.setup_splittable_workflow(workflow, [hdfbank], diff --git a/pycbc/workflow/plotting.py b/pycbc/workflow/plotting.py index 165786c5715..cdf71b4119a 100644 --- a/pycbc/workflow/plotting.py +++ b/pycbc/workflow/plotting.py @@ -534,7 +534,6 @@ def make_singles_plot(workflow, trig_files, bank_file, veto_file, veto_name, def make_dq_flag_trigger_rate_plot(workflow, dq_file, dq_label, out_dir, tags=None): tags = [] if tags is None else tags makedir(out_dir) - files = FileList([]) node = PlotExecutable(workflow.cp, 'plot_dq_flag_likelihood', ifos=dq_file.ifo, out_dir=out_dir, tags=tags).create_node() @@ -543,5 +542,29 @@ def make_dq_flag_trigger_rate_plot(workflow, dq_file, dq_label, out_dir, tags=No node.add_opt('--ifo', dq_file.ifo) node.new_output_file_opt(dq_file.segment, '.png', '--output-file') workflow += node - files += node.output_files - return files + return node.output_files[0] + + +def make_dq_segment_table(workflow, dq_file, out_dir, tags=None): + tags = [] if tags is None else tags + makedir(out_dir) + node = PlotExecutable(workflow.cp, 'page_dq_table', ifos=dq_file.ifo, + out_dir=out_dir, tags=tags).create_node() + node.add_input_opt('--dq-file', dq_file) + node.add_opt('--ifo', dq_file.ifo) + node.new_output_file_opt(dq_file.segment, '.html', '--output-file') + workflow += node + return node.output_files[0] + + +def make_template_bin_table(workflow, dq_file, out_dir, tags=None): + tags = [] if tags is None else tags + makedir(out_dir) + node = PlotExecutable(workflow.cp, 'page_template_bin_table', + ifos=dq_file.ifo, out_dir=out_dir, + tags=tags).create_node() + node.add_input_opt('--dq-file', dq_file) + node.add_opt('--ifo', dq_file.ifo) + node.new_output_file_opt(dq_file.segment, '.html', '--output-file') + workflow += node + return node.output_files[0] From df3ba441115ea4f8ec2946673c04e800f26dd421 Mon Sep 17 00:00:00 2001 From: Jacob Buchanan <83317255+jakeb245@users.noreply.github.com> Date: Thu, 15 Feb 2024 05:11:47 -0500 Subject: [PATCH 09/29] PyGRB: Update postprocessing functions (#4550) * Add mapping function * Implement mapping in trig combiner * Fix "template_id" dataset being written too early * Force template_id to be integer * Try adding bank files to workflow * Remove fixme * Only use one bank_file * Add template mapping to inj finder * Small change * mapping function works with trig file object * Small change * Remove typo * Add bank file opt to inj_finder * Add template masses to multi_inspiral output * sort_trigs updates * Extract trig properties * Add old imports back for CodeClimate * Remove unused bestnr opts * Update pycbc/results/pygrb_postprocessing_utils.py Co-authored-by: Francesco Pannarale * Codeclimate * Remove masses from multi inspiral output * Correct segment loading names * Add NoneType handling to _slide_vetoes function * Indentation fix * Add 's' back in * Fix docstring(?) * Codeclimate * Codeclimate * Update pycbc/results/pygrb_postprocessing_utils.py Co-authored-by: Tito Dal Canton * Uses event_ids in sort_trigs to avoid confusion * Add possibility of multiple banks (and NotImplementedError) * Remove enumerate and fix indexing issue * Check for single bank earlier * Simplify column name check * Use zip() * Update pycbc/results/pygrb_postprocessing_utils.py Co-authored-by: Tito Dal Canton * Update pycbc/results/pygrb_postprocessing_utils.py Co-authored-by: Tito Dal Canton --------- Co-authored-by: Francesco Pannarale Co-authored-by: Tito Dal Canton --- bin/pygrb/pycbc_grb_inj_finder | 21 ++++- bin/pygrb/pycbc_grb_trig_combiner | 19 +++++ bin/pygrb/pycbc_make_offline_grb_workflow | 9 +- pycbc/results/pygrb_postprocessing_utils.py | 95 +++++++++++---------- pycbc/workflow/grb_utils.py | 14 +-- 5 files changed, 103 insertions(+), 55 deletions(-) diff --git a/bin/pygrb/pycbc_grb_inj_finder b/bin/pygrb/pycbc_grb_inj_finder index 623e45c231e..230e09d3bc9 100644 --- a/bin/pygrb/pycbc_grb_inj_finder +++ b/bin/pygrb/pycbc_grb_inj_finder @@ -42,6 +42,7 @@ from ligo.segments.utils import fromsegwizard from pycbc import __version__ from pycbc.inject import InjectionSet from pycbc.io import FieldArray +from pycbc.results.pygrb_postprocessing_utils import template_hash_to_id __author__ = "Duncan Macleod " @@ -96,6 +97,8 @@ def read_hdf5_triggers(inputfiles, verbose=False): elif name.startswith("network") and \ NETWORK_IFO_EVENT_ID_REGEX.match(name): return + elif "template_id" in name: + return else: shape = obj.shape dtype = obj.dtype @@ -126,7 +129,7 @@ def read_hdf5_triggers(inputfiles, verbose=False): # copy dataset contents for filename in inputfiles: with h5py.File(filename, 'r') as h5in: - # Skipping emtpy files + # Skipping empty files if len(h5in['network']['end_time_gc'][:]): for dset in datasets: data = h5in[dset][:] @@ -174,6 +177,11 @@ parser.add_argument( nargs="+", help="path(s) to input injection file(s)", ) +parser.add_argument( + "--bank-file", + required=True, + help="full template bank file for mapping hashes to ids", +) parser.add_argument( "-o", "--output-dir", @@ -329,12 +337,21 @@ with h5py.File(outfilename, "w") as h5out: for x in out_trigs: xgroup = h5out.create_group(x) for col in out_trigs[x]: + if "template_id" == col: + continue xgroup.create_dataset( col, data=out_trigs[x][col], compression="gzip", compression_opts=9, ) - + # map and write template ids + ids = template_hash_to_id(h5out, args.bank_file) + h5out.create_dataset( + "network/template_id", + data=ids, + compression="gzip", + compression_opts=9, + ) if args.verbose: print("Found/missed written to {}".format(outfilename)) diff --git a/bin/pygrb/pycbc_grb_trig_combiner b/bin/pygrb/pycbc_grb_trig_combiner index 84cdf0b4ca5..e1687a3c237 100644 --- a/bin/pygrb/pycbc_grb_trig_combiner +++ b/bin/pygrb/pycbc_grb_trig_combiner @@ -36,6 +36,7 @@ from ligo import segments from ligo.segments.utils import fromsegwizard from pycbc import __version__ +from pycbc.results.pygrb_postprocessing_utils import template_hash_to_id __author__ = "Duncan Macleod " @@ -55,6 +56,8 @@ def _init_combined_file(h5file, attributes, datasets, **compression_kw): h5file.attrs.update(attributes) # create datasets for dset, (shape, dtype) in datasets.items(): + if "template_id" in dset: + continue h5file.create_dataset(dset, shape=shape, dtype=dtype, **compression_kw) @@ -207,6 +210,11 @@ def merge_hdf5_files(inputfiles, outputfile, verbose=False, **compression_kw): else: ifo_event_id_pos = None + # Template ID becomes unhelpful in the merged file. + # This will be handled separately. + if "template_id" in dset: + continue + # merge dataset position[dset] += _merge_dataset( h5in[dset], h5out[dset], position[dset], @@ -218,6 +226,12 @@ def merge_hdf5_files(inputfiles, outputfile, verbose=False, **compression_kw): while search_datasets: dataset_names.remove(search_datasets.pop()) + template_ids = template_hash_to_id(trigger_file=h5out, bank_path=args.bank_file) + h5out.create_dataset("network/template_id", + data=template_ids, + compression="gzip", + compression_opts=9) + if verbose: print("Merged triggers written to {}".format(outputfile)) return outputfile @@ -423,6 +437,11 @@ parser.add_argument( metavar="TRIGGER FILE", help="read in listed trigger files", ) +parser.add_argument( + "--bank-file", + required=True, + help="full template bank file for mapping hashes to ids", +) parser.add_argument( "-o", "--output-dir", diff --git a/bin/pygrb/pycbc_make_offline_grb_workflow b/bin/pygrb/pycbc_make_offline_grb_workflow index 3dd1ee9d831..59160f4f79c 100644 --- a/bin/pygrb/pycbc_make_offline_grb_workflow +++ b/bin/pygrb/pycbc_make_offline_grb_workflow @@ -273,9 +273,14 @@ all_files.extend(datafind_veto_files) # TEMPLATE BANK AND SPLIT BANK bank_files = _workflow.setup_tmpltbank_workflow(wflow, sciSegs, datafind_files, df_dir) +# Check there are not multiple banks +if len(bank_files) > 1: + raise NotImplementedError("Multiple banks not supported") +full_bank_file = bank_files[0] +# Note: setup_splittable_workflow requires a FileList as input splitbank_files = _workflow.setup_splittable_workflow(wflow, bank_files, df_dir, tags=["inspiral"]) -all_files.extend(bank_files) +all_files.append(full_bank_file) all_files.extend(splitbank_files) # INJECTIONS @@ -411,7 +416,7 @@ if post_proc_method == "PYGRB_OFFLINE": # The content and structure of pp_files are described in the definition # of _workflow.setup_pygrb_pp_workflow pp_files = _workflow.setup_pygrb_pp_workflow(wflow, pp_dir, seg_dir, - sciSegs[ifo][0], inspiral_files, + sciSegs[ifo][0], full_bank_file, inspiral_files, injs, inj_insp_files, inj_tags) sec_name = 'workflow-pygrb_pp_workflow' if not wflow.cp.has_section(sec_name): diff --git a/pycbc/results/pygrb_postprocessing_utils.py b/pycbc/results/pygrb_postprocessing_utils.py index 8760a025224..531f7124a40 100644 --- a/pycbc/results/pygrb_postprocessing_utils.py +++ b/pycbc/results/pygrb_postprocessing_utils.py @@ -310,16 +310,17 @@ def _get_id_numbers(ligolw_table, column): # ============================================================================= # Function to build a dictionary (indexed by ifo) of time-slid vetoes # ============================================================================= -def _slide_vetoes(vetoes, slide_dict_or_list, slide_id): +def _slide_vetoes(vetoes, slide_dict_or_list, slide_id, ifos): """Build a dictionary (indexed by ifo) of time-slid vetoes""" # Copy vetoes - slid_vetoes = copy.deepcopy(vetoes) - - # Slide them - ifos = vetoes.keys() - for ifo in ifos: - slid_vetoes[ifo].shift(-slide_dict_or_list[slide_id][ifo]) + if vetoes is not None: + slid_vetoes = copy.deepcopy(vetoes) + # Slide them + for ifo in ifos: + slid_vetoes[ifo].shift(-slide_dict_or_list[slide_id][ifo]) + else: + slid_vetoes = {ifo: segments.segmentlist() for ifo in ifos} return slid_vetoes @@ -461,36 +462,39 @@ def sort_trigs(trial_dict, trigs, slide_dict, seg_dict): sorted_trigs = {} # Begin by sorting the triggers into each slide - # New seems pretty slow, so run it once and then use deepcopy - tmp_table = glsctables.New(glsctables.MultiInspiralTable) for slide_id in slide_dict: - sorted_trigs[slide_id] = copy.deepcopy(tmp_table) - for trig in trigs: - sorted_trigs[int(trig.time_slide_id)].append(trig) + sorted_trigs[slide_id] = [] + for slide_id, event_id in zip(trigs['network/slide_id'], + trigs['network/event_id']): + sorted_trigs[slide_id].append(event_id) for slide_id in slide_dict: # These can only *reduce* the analysis time curr_seg_list = seg_dict[slide_id] # Check the triggers are all in the analysed segment lists - for trig in sorted_trigs[slide_id]: - if trig.end_time not in curr_seg_list: + for event_id in sorted_trigs[slide_id]: + index = numpy.flatnonzero(trigs['network/event_id'] == event_id)[0] + end_time = trigs['network/end_time_gc'][index] + if end_time not in curr_seg_list: # This can be raised if the trigger is on the segment boundary, # so check if the trigger is within 1/100 of a second within # the list - if trig.get_end() + 0.01 in curr_seg_list: + if end_time + 0.01 in curr_seg_list: continue - if trig.get_end() - 0.01 in curr_seg_list: + if end_time - 0.01 in curr_seg_list: continue err_msg = "Triggers found in input files not in the list of " err_msg += "analysed segments. This should not happen." raise RuntimeError(err_msg) # END OF CHECK # - # The below line works like the inverse of .veto and only returns trigs - # that are within the segment specified by trial_dict[slide_id] - sorted_trigs[slide_id] = \ - sorted_trigs[slide_id].vetoed(trial_dict[slide_id]) + # Keep triggers that are in trial_dict + sorted_trigs[slide_id] = [event_id for event_id in + sorted_trigs[slide_id] + if trigs['network/end_time_gc'][ + trigs['network/event_id'] == event_id][0] + in trial_dict[slide_id]] return sorted_trigs @@ -498,8 +502,7 @@ def sort_trigs(trial_dict, trigs, slide_dict, seg_dict): # ============================================================================= # Extract basic trigger properties and store them as dictionaries # ============================================================================= -def extract_basic_trig_properties(trial_dict, trigs, slide_dict, seg_dict, - opts): +def extract_basic_trig_properties(trial_dict, trigs, slide_dict, seg_dict): """Extract and store as dictionaries time, SNR, and BestNR of time-slid triggers""" @@ -507,16 +510,6 @@ def extract_basic_trig_properties(trial_dict, trigs, slide_dict, seg_dict, sorted_trigs = sort_trigs(trial_dict, trigs, slide_dict, seg_dict) logging.info("Triggers sorted.") - # Local copies of variables entering the BestNR definition - chisq_index = opts.chisq_index - chisq_nhigh = opts.chisq_nhigh - null_thresh = list(map(float, opts.null_snr_threshold.split(','))) - snr_thresh = opts.snr_threshold - sngl_snr_thresh = opts.sngl_snr_threshold - new_snr_thresh = opts.newsnr_threshold - null_grad_thresh = opts.null_grad_thresh - null_grad_val = opts.null_grad_val - # Build the 3 dictionaries trig_time = {} trig_snr = {} @@ -524,21 +517,12 @@ def extract_basic_trig_properties(trial_dict, trigs, slide_dict, seg_dict, for slide_id in slide_dict: slide_trigs = sorted_trigs[slide_id] if slide_trigs: - trig_time[slide_id] = numpy.asarray(slide_trigs.get_end()).\ - astype(float) - trig_snr[slide_id] = numpy.asarray(slide_trigs.get_column('snr')) + trig_time[slide_id] = trigs['network/end_time_gc'][slide_trigs] + trig_snr[slide_id] = trigs['network/coherent_snr'][slide_trigs] else: trig_time[slide_id] = numpy.asarray([]) trig_snr[slide_id] = numpy.asarray([]) - trig_bestnr[slide_id] = get_bestnrs(slide_trigs, - q=chisq_index, - n=chisq_nhigh, - null_thresh=null_thresh, - snr_threshold=snr_thresh, - sngl_snr_threshold=sngl_snr_thresh, - chisq_threshold=new_snr_thresh, - null_grad_thresh=null_grad_thresh, - null_grad_val=null_grad_val) + trig_bestnr[slide_id] = trigs['network/reweighted_snr'][slide_trigs] logging.info("Time, SNR, and BestNR of triggers extracted.") return trig_time, trig_snr, trig_bestnr @@ -699,7 +683,7 @@ def construct_trials(seg_files, seg_dict, ifos, slide_dict, vetoes): seg_buffer.coalesce() # Construct the ifo-indexed dictionary of slid veteoes - slid_vetoes = _slide_vetoes(vetoes, slide_dict, slide_id) + slid_vetoes = _slide_vetoes(vetoes, slide_dict, slide_id, ifos) # Construct trial list and check against buffer trial_dict[slide_id] = segments.segmentlist() @@ -804,3 +788,24 @@ def get_coinc_snr(trigs_or_injs, ifos): coinc_snr = numpy.sqrt(snr_sum_square) return coinc_snr + + +def template_hash_to_id(trigger_file, bank_path): + """ + This function converts the template hashes from a trigger file + into 'template_id's that represent indices of the + templates within the bank. + Parameters + ---------- + trigger_file: h5py File object for trigger file + bank_file: filepath for template bank + """ + with h5py.File(bank_path, "r") as bank: + hashes = bank['template_hash'][:] + ifos = [k for k in trigger_file.keys() if k != 'network'] + trig_hashes = trigger_file[f'{ifos[0]}/template_hash'][:] + trig_ids = numpy.zeros(trig_hashes.shape[0], dtype=int) + for idx, t_hash in enumerate(hashes): + matches = numpy.where(trig_hashes == t_hash) + trig_ids[matches] = idx + return trig_ids diff --git a/pycbc/workflow/grb_utils.py b/pycbc/workflow/grb_utils.py index 8f96f793d39..2d2113516b3 100644 --- a/pycbc/workflow/grb_utils.py +++ b/pycbc/workflow/grb_utils.py @@ -233,8 +233,8 @@ def get_sky_grid_scale( return out -def setup_pygrb_pp_workflow(wf, pp_dir, seg_dir, segment, insp_files, - inj_files, inj_insp_files, inj_tags): +def setup_pygrb_pp_workflow(wf, pp_dir, seg_dir, segment, bank_file, + insp_files, inj_files, inj_insp_files, inj_tags): """ Generate post-processing section of PyGRB offline workflow """ @@ -258,7 +258,7 @@ def setup_pygrb_pp_workflow(wf, pp_dir, seg_dir, segment, insp_files, job_instance = exe_class(wf.cp, "trig_combiner") # Create node for coherent no injections jobs node, trig_files = job_instance.create_node(wf.ifos, seg_dir, segment, - insp_files, pp_dir) + insp_files, pp_dir, bank_file) wf.add_node(node) pp_outs.append(trig_files) @@ -285,7 +285,7 @@ def setup_pygrb_pp_workflow(wf, pp_dir, seg_dir, segment, insp_files, if inj_tag in f.tags[1]]) node, inj_find_file = job_instance.create_node( tag_inj_files, tag_insp_files, - pp_dir) + bank_file, pp_dir) wf.add_node(node) inj_find_files.append(inj_find_file) pp_outs.append(inj_find_files) @@ -321,7 +321,7 @@ def __init__(self, cp, name): self.num_trials = int(cp.get('trig_combiner', 'num-trials')) def create_node(self, ifo_tag, seg_dir, segment, insp_files, - out_dir, tags=None): + out_dir, bank_file, tags=None): node = Node(self) node.add_opt('--verbose') node.add_opt("--ifo-tag", ifo_tag) @@ -331,6 +331,7 @@ def create_node(self, ifo_tag, seg_dir, segment, insp_files, node.add_input_list_opt("--input-files", insp_files) node.add_opt("--user-tag", "PYGRB") node.add_opt("--num-trials", self.num_trials) + node.add_input_opt("--bank-file", bank_file) # Prepare output file tag user_tag = f"PYGRB_GRB{self.trigger_name}" if tags: @@ -385,13 +386,14 @@ class PycbcGrbInjFinderExecutable(Executable): def __init__(self, cp, exe_name): super().__init__(cp=cp, name=exe_name) - def create_node(self, inj_files, inj_insp_files, + def create_node(self, inj_files, inj_insp_files, bank_file, out_dir, tags=None): if tags is None: tags = [] node = Node(self) node.add_input_list_opt('--input-files', inj_insp_files) node.add_input_list_opt('--inj-files', inj_files) + node.add_input_opt('--bank-file', bank_file) ifo_tag, desc, segment = filename_metadata(inj_files[0].name) desc = '_'.join(desc.split('_')[:-1]) out_name = "{}-{}_FOUNDMISSED-{}-{}.h5".format( From 3fc2ad2bc53ce5a4d890d4d2acf0abecf448bbba Mon Sep 17 00:00:00 2001 From: Alex Nitz Date: Thu, 15 Feb 2024 07:20:43 -0500 Subject: [PATCH 10/29] Revert "MNT: upgrade upload/download-artifact GitHub actions (#4595)" (#4629) This reverts commit 02f77b6d2b0b8fd2a62903f02356954d81412252. --- .github/workflows/basic-tests.yml | 7 +++---- .github/workflows/distribution.yml | 10 ++++------ .github/workflows/inference-workflow.yml | 4 ++-- .github/workflows/search-workflow.yml | 4 ++-- .github/workflows/tmpltbank-workflow.yml | 2 +- .github/workflows/workflow-tests.yml | 2 +- 6 files changed, 13 insertions(+), 16 deletions(-) diff --git a/.github/workflows/basic-tests.yml b/.github/workflows/basic-tests.yml index 29739899814..d962026e09a 100644 --- a/.github/workflows/basic-tests.yml +++ b/.github/workflows/basic-tests.yml @@ -51,7 +51,7 @@ jobs: tox -e py-inference - name: store documentation page if: matrix.test-type == 'docs' && matrix.python-version == '3.8' - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: documentation-page path: _gh-pages @@ -61,10 +61,9 @@ jobs: if: github.ref == 'refs/heads/master' && github.event_name == 'push' steps: - name: retrieve built documentation - uses: actions/download-artifact@v4 + uses: actions/download-artifact@v2 with: - pattern: documentation-page-* - merge-multiple: true + name: documentation-page - name: debug run: | mkdir _gh-pages diff --git a/.github/workflows/distribution.yml b/.github/workflows/distribution.yml index e59fb75bc87..321b81ff4d2 100644 --- a/.github/workflows/distribution.yml +++ b/.github/workflows/distribution.yml @@ -29,9 +29,8 @@ jobs: CIBW_BUILD: cp38-* cp39-* cp310-* cp311-* CIBW_SKIP: "*musllinux*" CIBW_ARCHS_MACOS: x86_64 arm64 - - uses: actions/upload-artifact@v4 + - uses: actions/upload-artifact@v2 with: - name: wheels-${{ matrix.os }} path: ./wheelhouse/*.whl deploy_pypi: name: Build and publish Python 🐍 distributions 📦 to PyPI @@ -45,14 +44,13 @@ jobs: uses: actions/setup-python@v4 with: python-version: 3.8 - - uses: actions/download-artifact@v4 + - uses: actions/download-artifact@v2 with: - pattern: wheels-* - merge-multiple: true + path: ./ - name: build pycbc for pypi run: | python setup.py sdist - mv *.whl dist/ + mv artifact/* dist/ - name: Publish distribution 📦 to PyPI if: github.event_name == 'push' && startsWith(github.ref, 'refs/tags') uses: pypa/gh-action-pypi-publish@master diff --git a/.github/workflows/inference-workflow.yml b/.github/workflows/inference-workflow.yml index 68644bd663d..16e5b9cc65d 100644 --- a/.github/workflows/inference-workflow.yml +++ b/.github/workflows/inference-workflow.yml @@ -48,12 +48,12 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: logs path: gw_output/submitdir/work - name: store result page - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: results path: html diff --git a/.github/workflows/search-workflow.yml b/.github/workflows/search-workflow.yml index c8a5bb9f7bb..18219411078 100644 --- a/.github/workflows/search-workflow.yml +++ b/.github/workflows/search-workflow.yml @@ -55,12 +55,12 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: logs path: output/submitdir/work - name: store result page - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: results path: html diff --git a/.github/workflows/tmpltbank-workflow.yml b/.github/workflows/tmpltbank-workflow.yml index 21917cebf7f..bf29c9051d0 100644 --- a/.github/workflows/tmpltbank-workflow.yml +++ b/.github/workflows/tmpltbank-workflow.yml @@ -51,7 +51,7 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: logs path: output/submitdir/work diff --git a/.github/workflows/workflow-tests.yml b/.github/workflows/workflow-tests.yml index 3f39ef96ea0..754e25f8f11 100644 --- a/.github/workflows/workflow-tests.yml +++ b/.github/workflows/workflow-tests.yml @@ -50,7 +50,7 @@ jobs: find submitdir/work/ -type f -name '*.tar.gz' -delete - name: store log files if: always() - uses: actions/upload-artifact@v4 + uses: actions/upload-artifact@v2 with: name: logs-${{matrix.test-type}} path: examples/workflow/generic/${{matrix.test-type}}/submitdir/work From 08bb91c51cdd545171b675fd913d8e4b55849c22 Mon Sep 17 00:00:00 2001 From: Jacob Buchanan <83317255+jakeb245@users.noreply.github.com> Date: Thu, 15 Feb 2024 15:46:37 -0500 Subject: [PATCH 11/29] Update PyGRB efficiency script for HDF files (#4562) * efficiency changes * Delete old comment to retrigger checks * Update bin/pygrb/pycbc_pygrb_efficiency Co-authored-by: Francesco Pannarale * Update bin/pygrb/pycbc_pygrb_efficiency * Don't extract file contents (part 1) * Don't extract file contents (part 2) * Don't extract file contents (part 3) * Don't extract file contents (end?) * Remove missed_inj * Add FIXME * Pull template masses from template bank * Use trigger masses instead of injection masses for trigger mchirp * Avoid h5py error * dataset name correction * Avoid numpy error * Try fixing loudest onsource bestnr finding * Fix issues * Precalculate template chirp masses * Clean up --------- Co-authored-by: Francesco Pannarale --- bin/pygrb/pycbc_pygrb_efficiency | 177 ++++++++++++++----------------- 1 file changed, 81 insertions(+), 96 deletions(-) diff --git a/bin/pygrb/pycbc_pygrb_efficiency b/bin/pygrb/pycbc_pygrb_efficiency index 9c5b48e4660..b3755c228c0 100644 --- a/bin/pygrb/pycbc_pygrb_efficiency +++ b/bin/pygrb/pycbc_pygrb_efficiency @@ -24,6 +24,7 @@ import sys import os import logging +import h5py import matplotlib.pyplot as plt from matplotlib import rc import numpy as np @@ -34,10 +35,7 @@ from pycbc.detector import Detector from pycbc import init_logging from pycbc.results import save_fig_with_metadata from pycbc.results import pygrb_postprocessing_utils as ppu -try: - from glue.ligolw import lsctables -except ImportError: - pass +from pycbc.conversions import mchirp_from_mass1_mass2 plt.switch_backend('Agg') rc("image") @@ -98,18 +96,26 @@ parser.add_argument("-g", "--glitch-check-factor", action="store", "exclusion efficiencies this value is multiplied " + "to the offsource around the injection trigger to " + "determine if it is just a loud glitch.") +parser.add_argument("--found-missed-file", action="store", type=str, + default=None, help="Location of found/missed " + "injections and trigger file") parser.add_argument("--injection-set-name", action="store", type=str, default="", help="Name of the injection set to be " + "used in the plot title.") parser.add_argument("-C", "--cluster-window", action="store", type=float, default=0.1, help="The cluster window used " + "to cluster triggers in time.") +parser.add_argument("--bank-file", action="store", type=str, + help="Location of the full template bank used.") ppu.pygrb_add_injmc_opts(parser) ppu.pygrb_add_bestnr_opts(parser) opts = parser.parse_args() init_logging(opts.verbose, format="%(asctime)s: %(levelname)s: %(message)s") +# Load bank file +bank_file = h5py.File(opts.bank_file, 'r') + # Check options do_injections = opts.found_missed_file is not None @@ -165,8 +171,8 @@ ifos, vetoes = ppu.extract_ifos_and_vetoes(trig_file, opts.veto_files, # Load triggers, time-slides, and segment dictionary logging.info("Loading triggers.") -trigs = ppu.load_xml_table(trig_file, lsctables.MultiInspiralTable.tableName) -logging.info("%d triggers loaded.", len(trigs)) +trigs = ppu.load_triggers(trig_file, vetoes) +logging.info("%d triggers loaded.", trigs['network/event_id'].len()) logging.info("Loading timeslides.") slide_dict = ppu.load_time_slides(trig_file) logging.info("Loading segments.") @@ -181,7 +187,7 @@ logging.info("%d trials generated.", total_trials) # Extract basic trigger properties and store as dictionaries trig_time, trig_snr, trig_bestnr = \ - ppu.extract_basic_trig_properties(trial_dict, trigs, slide_dict, segment_dict, opts) + ppu.extract_basic_trig_properties(trial_dict, trigs, slide_dict, segment_dict) # Calculate BestNR values and maximum time_veto_max_bestnr = {} @@ -219,6 +225,17 @@ offsource_trigs.reverse() max_bestnr, _, full_time_veto_max_bestnr =\ ppu.max_median_stat(slide_dict, time_veto_max_bestnr, trig_bestnr, total_trials) +# ========================== +# Calculate template chirp masses from bank +# ========================== +logging.info('Reading template chirp masses') +with h5py.File(opts.bank_file, 'r') as bank_file: + template_mchirps = mchirp_from_mass1_mass2( + bank_file['mass1'][:], + bank_file['mass2'][:] + ) + + # ======================= # Load on source triggers @@ -226,47 +243,38 @@ max_bestnr, _, full_time_veto_max_bestnr =\ if onsource_file: # Get trigs - on_trigs = ppu.load_xml_table(onsource_file, "multi_inspiral") - logging.info("%d onsource triggers loaded.", len(on_trigs)) + on_trigs = ppu.load_triggers(onsource_file, vetoes) + logging.info("%d onsource triggers loaded.", on_trigs['network/event_id'].len()) - # Separate off chirp mass column - on_mchirp = on_trigs.get_column('mchirp') + # Calculate chirp mass values + on_mchirp = template_mchirps[on_trigs['network/template_id']] # Set loudest event arrays - #loud_on_bestnr_trigs = None loud_on_bestnr = 0 - #loud_on_fap = 1 - # Calculate BestNRs and record loudest trig by BestNR + # Retrieve BestNRs and record loudest trig by BestNR. + # Get indices of loudest triggers and pick the loudest. if on_trigs: - on_trigs_bestnrs = ppu.get_bestnrs(on_trigs, - q=chisq_index, n=chisq_nhigh, - null_thresh=null_thresh, - snr_threshold=snr_thresh, - sngl_snr_threshold=sngl_snr_thresh, - chisq_threshold=new_snr_thresh, - null_grad_thresh=null_grad_thresh, - null_grad_val=null_grad_val) - loud_on_bestnr_trigs, loud_on_bestnr = \ - np.asarray([[x, y] for y, x in sorted(zip(on_trigs_bestnrs, on_trigs), - key=lambda elem: elem[0], - reverse=True)])[0] + loud_on_bestnr_idx = np.argmax(on_trigs['network/reweighted_snr']) + loud_on_bestnr = np.max(on_trigs['network/reweighted_snr']) + + # Convert to float for output + loud_on_bestnr = float(loud_on_bestnr) # If the loudest event has bestnr = 0, there is no event at all! if loud_on_bestnr == 0: - loud_on_bestnr_trigs = None + loud_on_bestnr_idx = None loud_on_fap = 1 logging.info("Onsource analysed.") - if loud_on_bestnr_trigs: - trig = loud_on_bestnr_trigs + if loud_on_bestnr_idx is not None: num_trials_louder = 0 tot_off_snr = np.array([]) for slide_id in slide_dict: - num_trials_louder += sum(time_veto_max_bestnr[slide_id] > \ + num_trials_louder += sum(time_veto_max_bestnr[slide_id] > loud_on_bestnr) - tot_off_snr = np.concatenate([tot_off_snr,\ + tot_off_snr = np.concatenate([tot_off_snr, time_veto_max_bestnr[slide_id]]) #fap_test = sum(tot_off_snr > loud_on_bestnr)/total_trials loud_on_fap = num_trials_louder/total_trials @@ -274,7 +282,7 @@ if onsource_file: else: tot_off_snr = np.array([]) for slide_id in slide_dict: - tot_off_snr = np.concatenate([tot_off_snr,\ + tot_off_snr = np.concatenate([tot_off_snr, time_veto_max_bestnr[slide_id]]) med_snr = np.median(tot_off_snr) #loud_on_fap = sum(tot_off_snr > med_snr)/total_trials @@ -287,55 +295,38 @@ if do_injections: sites = [ifo[0] for ifo in ifos] # Triggers and injections recovered in some form: discard ones at vetoed times - found_trigs = ppu.load_injections(found_file, vetoes) - found_injs = ppu.load_injections(found_file, vetoes, sim_table=True) + # This injs object contains the information about found/missed injections AND triggers + injs = ppu.load_triggers(found_missed_file, vetoes) logging.info("Missed/found injections/triggers loaded.") - # Extract columns of found injections - found_inj = {} - found_inj['ra'] = np.asarray(found_injs.get_column('longitude')) - found_inj['dec'] = np.asarray(found_injs.get_column('latitude')) - found_inj['time'] = np.asarray(found_injs.get_column('geocent_end_time')) +\ - np.asarray(found_injs.get_column('geocent_end_time_ns')*\ - 10**-9) - found_inj['dist'] = np.asarray(found_injs.get_column('distance')) - - # Extract columns of triggers - found_trig = {} - found_trig['mchirp'] = np.asarray(found_trigs.get_column('mchirp')) - found_trig['bestnr'] = ppu.get_bestnrs(found_trigs, q=chisq_index, - n=chisq_nhigh, - null_thresh=null_thresh, - snr_threshold=snr_thresh, - sngl_snr_threshold=sngl_snr_thresh, - chisq_threshold=new_snr_thresh, - null_grad_thresh=null_grad_thresh, - null_grad_val=null_grad_val) + # Calculate quantities not included in trigger files, such as chirp mass + found_trig_mchirp = template_mchirps[injs['network/template_id']] + # Construct conditions for injection: # 1) found louder than background, - zero_fap = np.zeros(len(found_injs)).astype(bool) - zero_fap_cut = found_trig['bestnr'] > max_bestnr + zero_fap = np.zeros(len(injs['network/end_time_gc'])).astype(bool) + zero_fap_cut = injs['network/reweighted_snr'][:] > max_bestnr zero_fap = zero_fap | (zero_fap_cut) # 2) found (bestnr > 0) but not louder than background (non-zero FAP) - nonzero_fap = ~zero_fap & (found_trig['bestnr'] != 0) + nonzero_fap = ~zero_fap & (injs['network/reweighted_snr'] != 0) # 3) missed after being recovered (i.e., vetoed) are not used here # missed = (~zero_fap) & (~nonzero_fap) # Non-zero FAP triggers (g_ifar) g_ifar = {} - g_ifar['bestnr'] = found_trig['bestnr'][nonzero_fap] + g_ifar['bestnr'] = injs['network/reweighted_snr'][nonzero_fap] g_ifar['stat'] = np.zeros([len(g_ifar['bestnr'])]) for ix, (mc, bestnr) in \ - enumerate(zip(found_trig['mchirp'][nonzero_fap], g_ifar['bestnr'])): + enumerate(zip(found_trig_mchirp[nonzero_fap], g_ifar['bestnr'])): g_ifar['stat'][ix] = (full_time_veto_max_bestnr > bestnr).sum() g_ifar['stat'] = g_ifar['stat'] / total_trials # Set the sigma values - inj_sigma = found_trigs.get_sigmasqs() + inj_sigma = {ifo: injs[f'{ifo}/sigmasq'][:] for ifo in ifos} # If the sigmasqs are not populated, we can still do calibration errors, # but only in the 1-detector case for ifo in ifos: @@ -354,10 +345,10 @@ if do_injections: for ifo in ifos: antenna = Detector(ifo) f_resp[ifo] = ppu.get_antenna_responses(antenna, - found_inj['ra'], found_inj['dec'], - found_inj['time']) + injs['found/ra'][:], injs['found/dec'][:], + injs['found/tc'][:]) - inj_sigma_mult = (np.asarray(list(inj_sigma.values())) *\ + inj_sigma_mult = (np.asarray(list(inj_sigma.values())) * np.asarray(list(f_resp.values()))) inj_sigma_tot = inj_sigma_mult[0, :] @@ -368,19 +359,14 @@ if do_injections: for ifo in ifos: inj_sigma_mean[ifo] = ((inj_sigma[ifo]*f_resp[ifo])/inj_sigma_tot).mean() - logging.info("%d found injections analysed.", len(found_injs)) - - # Missed injections (ones not recovered at all) - missed_injs = ppu.load_injections(found_missed_file, vetoes, sim_table=True) + logging.info("%d found injections analysed.", len(injs['found/tc'])) - # Process missed injections 'missed_inj' - missed_inj = {} - missed_inj['dist'] = np.asarray(missed_injs.get_column('distance')) + # Process missed injections (injs['missed']) - logging.info("%d missed injections analysed.", len(missed_injs)) + logging.info("%d missed injections analysed.", len(injs['missed/tc'])) # Create new set of injections for efficiency calculations - total_injs = len(found_injs) + len(missed_injs) + total_injs = len(injs['found/distance']) + len(injs['missed/distance']) long_inj = {} long_inj['dist'] = stats.uniform.rvs(size=total_injs) * (upper_dist-lower_dist) +\ upper_dist @@ -388,10 +374,10 @@ if do_injections: logging.info("%d long distance injections created.", total_injs) # Set distance bins and data arrays - dist_bins = zip(np.arange(lower_dist, upper_dist + (upper_dist-lower_dist),\ - (upper_dist-lower_dist)/num_bins),\ - np.arange(lower_dist, upper_dist + (upper_dist-lower_dist),\ - (upper_dist-lower_dist)/num_bins) +\ + dist_bins = zip(np.arange(lower_dist, upper_dist + (upper_dist-lower_dist), + (upper_dist-lower_dist)/num_bins), + np.arange(lower_dist, upper_dist + (upper_dist-lower_dist), + (upper_dist-lower_dist)/num_bins) + (upper_dist-lower_dist)/num_bins) dist_bins = list(dist_bins) @@ -405,7 +391,7 @@ if do_injections: found_on_bestnr[key] = np.zeros(num_dist_bins_plus_one) # Construct FAP list for all found injections - inj_fap = np.zeros(len(found_injs)) + inj_fap = np.zeros(len(injs['found/distance'])) inj_fap[nonzero_fap] = g_ifar['stat'] # Calculate the amplitude error @@ -428,18 +414,17 @@ if do_injections: # NOTE: the loop on num_mc_injs would fill up the *_inj['dist_mc']'s at the # same time, so filling them up sequentially will vary the numbers a little # (this is an MC, order of operations matters!) - found_inj['dist_mc'] = ppu.mc_cal_wf_errs(num_mc_injs, found_inj['dist'], + found_inj_dist_mc = ppu.mc_cal_wf_errs(num_mc_injs, injs['found/distance'], cal_error, wav_err, max_dc_cal_error) - missed_inj['dist_mc'] = ppu.mc_cal_wf_errs(num_mc_injs, missed_inj['dist'], + missed_inj_dist_mc = ppu.mc_cal_wf_errs(num_mc_injs, injs['missed/distance'], cal_error, wav_err, max_dc_cal_error) long_inj['dist_mc'] = ppu.mc_cal_wf_errs(num_mc_injs, long_inj['dist'], cal_error, wav_err, max_dc_cal_error) - logging.info("MC injection set distributed with %d iterations.",\ + logging.info("MC injection set distributed with %d iterations.", num_mc_injs) # Check injections against on source - more_sig_than_onsource = np.ndarray(len(found_injs)) if onsource_file: more_sig_than_onsource = (inj_fap <= loud_on_fap) else: @@ -447,29 +432,29 @@ if do_injections: distance_count = np.zeros(len(dist_bins)) - found_trig['max_bestnr'] = np.empty(len(found_trig['mchirp'])) - found_trig['max_bestnr'].fill(max_bestnr) + found_trig_max_bestnr = np.empty(len(injs['network/event_id'])) + found_trig_max_bestnr.fill(max_bestnr) - max_bestnr_cut = (found_trig['bestnr'] > found_trig['max_bestnr']) + max_bestnr_cut = (injs['network/reweighted_snr'] > found_trig_max_bestnr) # Check louder than on source - found_trig['loud_on_bestnr'] = np.empty(len(found_trig['mchirp'])) + found_trig_loud_on_bestnr = np.empty(len(injs['network/event_id'])) if onsource_file: - found_trig['loud_on_bestnr'].fill(loud_on_bestnr) + found_trig_loud_on_bestnr.fill(loud_on_bestnr) else: - found_trig['loud_on_bestnr'].fill(med_snr) - on_bestnr_cut = found_trig['bestnr'] > found_trig['loud_on_bestnr'] + found_trig_loud_on_bestnr.fill(med_snr) + on_bestnr_cut = injs['network/reweighted_snr'] > found_trig_loud_on_bestnr # Check whether injection is found for the purposes of exclusion # distance calculation. # Found: if louder than all on source # Missed: if not louder than loudest on source found_excl = on_bestnr_cut & (more_sig_than_onsource) & \ - (found_trig['bestnr'] != 0) + (injs['network/reweighted_snr'] != 0) # If not missed, double check bestnr against nearby triggers near_test = np.zeros((found_excl).sum()).astype(bool) - for j, (t, bestnr) in enumerate(zip(found_inj['time'][found_excl],\ - found_trig['bestnr'][found_excl])): + for j, (t, bestnr) in enumerate(zip(injs['found/tc'][found_excl], + injs['network/reweighted_snr'][found_excl])): # 0 is the zero-lag timeslide near_bestnr = trig_bestnr[0]\ [np.abs(trig_time[0]-t) < cluster_window] @@ -486,10 +471,10 @@ if do_injections: # Loop over the distance bins for j, dist_bin in enumerate(dist_bins): # Construct distance cut - found_dist_cut = (dist_bin[0] <= found_inj['dist_mc'][k, :]) &\ - (found_inj['dist_mc'][k, :] < dist_bin[1]) - missed_dist_cut = (dist_bin[0] <= missed_inj['dist_mc'][k, :]) &\ - (missed_inj['dist_mc'][k, :] < dist_bin[1]) + found_dist_cut = (dist_bin[0] <= found_inj_dist_mc[k, :]) &\ + (found_inj_dist_mc[k, :] < dist_bin[1]) + missed_dist_cut = (dist_bin[0] <= missed_inj_dist_mc[k, :]) &\ + (missed_inj_dist_mc[k, :] < dist_bin[1]) long_dist_cut = (dist_bin[0] <= long_inj['dist_mc'][k, :]) &\ (long_inj['dist_mc'][k, :] < dist_bin[1]) @@ -530,7 +515,7 @@ if do_injections: # Calculate error bars for efficiency/distance plots and datafiles # using max BestNR of background - yerr_low_mc, yerr_high_mc, fraction_mc = efficiency_with_errs(\ + yerr_low_mc, yerr_high_mc, fraction_mc = efficiency_with_errs( found_max_bestnr['mc'], num_injections['mc'], num_mc_injs=num_mc_injs) yerr_low_no_mc, yerr_high_no_mc, fraction_no_mc = efficiency_with_errs(\ found_max_bestnr['no_mc'], num_injections['no_mc']) From 9d965fa9aac645d0c83f8adc7aa4af3fc24c5f77 Mon Sep 17 00:00:00 2001 From: Gareth S Cabourn Davies Date: Fri, 16 Feb 2024 11:38:35 +0000 Subject: [PATCH 12/29] Bug: Live combine fits, different bin edge sizes causes an error (#4636) * Update pycbc_live_combine_single_fits * logging warn -> warning deprecation --- bin/live/pycbc_live_combine_single_fits | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/bin/live/pycbc_live_combine_single_fits b/bin/live/pycbc_live_combine_single_fits index a0638ecb036..c2847a1e774 100644 --- a/bin/live/pycbc_live_combine_single_fits +++ b/bin/live/pycbc_live_combine_single_fits @@ -84,10 +84,11 @@ for f in args.trfits_files: same_conf = (fits_f.attrs['sngl_ranking'] == sngl_rank and fits_f.attrs['fit_threshold'] == fit_thresh and fits_f.attrs['fit_function'] == fit_func + and fits_f['bins_lower'].size == bl.size and all(fits_f['bins_lower'][:] == bl) and all(fits_f['bins_upper'][:] == bu)) if not same_conf: - logging.warn( + logging.warning( "Found a change in the fit configuration, skipping %s", f ) @@ -119,8 +120,8 @@ for f in args.trfits_files: counts_all[ifo].append(ffi['counts'][:]) alphas_all[ifo].append(ffi['fit_coeff'][:]) if any(np.isnan(ffi['fit_coeff'][:])): - logging.warn("nan in %s, %s", f, ifo) - logging.warn(ffi['fit_coeff'][:]) + logging.warning("nan in %s, %s", f, ifo) + logging.warning(ffi['fit_coeff'][:]) # Set up the date array, this is stored as an offset from the first trigger time of # the first file to the last trigger of the file From 320ac7f610ae9cffdf686b5dc6176adfd5eb67f4 Mon Sep 17 00:00:00 2001 From: Tito Dal Canton Date: Mon, 19 Feb 2024 17:18:02 +0100 Subject: [PATCH 13/29] Improve how pycbc_live handles localization-only detectors (#4635) * Improve how pycbc_live handles localization-only detectors LOD = localization-only detectors * Do not initialize single-detector FAR estimators for LOD * Only use triggering detectors for PSD var statistic * Correct the GraceDB annotation to explain which detectors contributed to the FAR and list the LOD * Only list triggering detectors in background dump filename * Code reformatting * Fixes * Fixes * Use attributes of event manager instead * Improve GraceDB annotation * Fix * Ignore loc-only detectors for rank > 0 * Further clarify which detectors are used only for localization * Main loop should go through the loc-only detectors only in rank 0 --- bin/pycbc_live | 119 ++++++++++++++++++++++++++++++------------------- 1 file changed, 74 insertions(+), 45 deletions(-) diff --git a/bin/pycbc_live b/bin/pycbc_live index 8a808ca5a4f..606f081cf9a 100755 --- a/bin/pycbc_live +++ b/bin/pycbc_live @@ -82,7 +82,13 @@ class LiveEventManager(object): def __init__(self, args, bank): self.low_frequency_cutoff = args.low_frequency_cutoff self.bank = bank - self.skymap_only_ifos = [] if args.skymap_only_ifos is None else list(set(args.skymap_only_ifos)) + # all interferometers involved in the analysis, whether for generating + # candidates or doing localization only + self.ifos = set(args.channel_name.keys()) + # interferometers used for localization only + self.skymap_only_ifos = set(args.skymap_only_ifos or []) + # subset of the interferometers allowed to produce candidates + self.trigg_ifos = self.ifos - self.skymap_only_ifos # Figure out what we are supposed to process within the pool of MPI processes self.comm = mpi.COMM_WORLD @@ -117,7 +123,7 @@ class LiveEventManager(object): if self.run_snr_optimization: # preestimate the number of CPU cores that we can afford giving # to followup processes without slowing down the main search - bg_cores = len(tuple(itertools.combinations(ifos, 2))) + bg_cores = len(tuple(itertools.combinations(self.trigg_ifos, 2))) analysis_cores = 1 + bg_cores if platform.system() != 'Darwin': available_cores = len(os.sched_getaffinity(0)) @@ -502,7 +508,7 @@ class LiveEventManager(object): logging.info('computing followup data for coinc') coinc_ifos = coinc_results['foreground/type'].split('-') followup_ifos = set(ifos) - set(coinc_ifos) - followup_ifos = list(followup_ifos | set(self.skymap_only_ifos)) + followup_ifos = list(followup_ifos | self.skymap_only_ifos) double_ifar = coinc_results['foreground/ifar'] if double_ifar < args.ifar_double_followup_threshold: @@ -537,13 +543,20 @@ class LiveEventManager(object): logging.info('Coincident candidate! Saving as %s', fname) # Verbally explain some details not obvious from the other info - comment = ('Trigger produced as a {} coincidence. ' - 'FAR is based on all listed detectors.
' - 'Two-detector ranking statistic: {}
' - 'Followup detectors: {}') - comment = comment.format(ppdets(coinc_ifos), - args.ranking_statistic, - ppdets(followup_ifos)) + comment = ( + 'Trigger produced as a {} coincidence.
' + 'Two-detector ranking statistic: {}
' + 'Detectors used for FAR calculation: {}.
' + 'Detectors used for localization: {}.
' + 'Detectors used only for localization: {}.' + ) + comment = comment.format( + ppdets(coinc_ifos), + args.ranking_statistic, + set(ifos) - self.skymap_only_ifos, + ppdets(followup_ifos), + ppdets(self.skymap_only_ifos) + ) ifar = coinc_results['foreground/ifar'] upload_checks = self.enable_gracedb_upload and self.ifar_upload_threshold < ifar @@ -604,7 +617,7 @@ class LiveEventManager(object): logging.info(f'Found {ifo} single with ifar {sifar}') followup_ifos = [i for i in active if i is not ifo] - followup_ifos = list(set(followup_ifos) | set(self.skymap_only_ifos)) + followup_ifos = list(set(followup_ifos) | self.skymap_only_ifos) # Don't recompute ifar considering other ifos sld = self.compute_followup_data( [ifo], @@ -636,10 +649,15 @@ class LiveEventManager(object): logging.info('Single-detector candidate! Saving as %s', fname) # Verbally explain some details not obvious from the other info - comment = ('Trigger produced as a {0} single. ' - 'FAR is based on {0} only.
' - 'Followup detectors: {1}') - comment = comment.format(ifo, ppdets(followup_ifos)) + comment = ( + 'Trigger produced as a {0} single.
' + 'Detectors used for FAR calculation: {0}.
' + 'Detectors used for localization: {1}.
' + 'Detectors used only for localization: {2}.' + ) + comment = comment.format( + ifo, ppdets(followup_ifos), ppdets(self.skymap_only_ifos) + ) # Has a coinc event at this time been uploaded recently? # If so, skip upload - Note that this means that we _always_ @@ -1017,8 +1035,6 @@ args = parser.parse_args() scheme.verify_processing_options(args, parser) fft.verify_fft_options(args, parser) Coincer.verify_args(args, parser) -ifos = set(args.channel_name.keys()) -analyze_singles = LiveSingle.verify_args(args, parser, ifos) if args.output_background is not None and len(args.output_background) != 2: parser.error('--output-background takes two parameters: period and path') @@ -1045,14 +1061,16 @@ if bank.min_f_lower < args.low_frequency_cutoff: 'minimum f_lower across all templates ' '({} Hz)'.format(args.low_frequency_cutoff, bank.min_f_lower)) -logging.info('Analyzing data from detectors %s', ppdets(ifos)) - evnt = LiveEventManager(args, bank) -logging.info('Detectors that only aid in the sky localization %s', ppdets(evnt.skymap_only_ifos)) + +logging.info('Analyzing data from detectors %s', ppdets(evnt.ifos)) +logging.info('Using %s for localization only', ppdets(evnt.skymap_only_ifos)) + +analyze_singles = LiveSingle.verify_args(args, parser, evnt.trigg_ifos) # include MPI rank and functional description into proctitle task_name = 'root' if evnt.rank == 0 else 'filtering' -setproctitle('PyCBC Live rank {:d} [{}]'.format(evnt.rank, task_name)) +setproctitle(f'PyCBC Live rank {evnt.rank:d} [{task_name}]') sg_chisq = SingleDetSGChisq.from_cli(args, bank, args.chisq_bins) @@ -1111,23 +1129,27 @@ with ctx: args.round_start_time logging.info('Starting from: %s', args.start_time) - # initialize the data readers for all detectors + # Initialize the data readers for all detectors. For rank 0, we need data + # from all detectors, including the localization-only ones. For higher + # ranks, we only need the detectors that can generate candidates. if args.max_length is not None: maxlen = args.max_length maxlen = int(maxlen) - data_reader = {ifo: StrainBuffer.from_cli(ifo, args, maxlen) - for ifo in ifos} + data_reader = { + ifo: StrainBuffer.from_cli(ifo, args, maxlen) + for ifo in (evnt.ifos if evnt.rank == 0 else evnt.trigg_ifos) + } evnt.data_readers = data_reader # create single-detector background "estimators" if analyze_singles and evnt.rank == 0: sngl_estimator = {ifo: LiveSingle.from_cli(args, ifo) - for ifo in ifos} + for ifo in evnt.trigg_ifos} - # Create double coincident background estimator for every combo + # Create double coincident background estimator + # for every pair of triggering interferometers if args.enable_background_estimation and evnt.rank == 0: - trigg_ifos = [ifo for ifo in ifos if ifo not in evnt.skymap_only_ifos] - ifo_combos = itertools.combinations(trigg_ifos, 2) + ifo_combos = itertools.combinations(evnt.trigg_ifos, 2) estimators = [] for combo in ifo_combos: logging.info('Will calculate %s background', ppdets(combo, "-")) @@ -1167,12 +1189,14 @@ with ctx: # main analysis loop data_end = lambda: data_reader[tuple(data_reader.keys())[0]].end_time last_bg_dump_time = int(data_end()) - psd_count = {ifo:0 for ifo in ifos} + psd_count = {ifo:0 for ifo in evnt.ifos} # Create dicts to track whether the psd has been recalculated and to hold # psd variation filters - psd_recalculated = {ifo: True for ifo in ifos} - psd_var_filts = {ifo: None for ifo in ifos} + psd_recalculated = { + ifo: True for ifo in (evnt.ifos if evnt.rank == 0 else evnt.trigg_ifos) + } + psd_var_filts = {ifo: None for ifo in evnt.trigg_ifos} while data_end() < args.end_time: t1 = pycbc.gps_now() @@ -1181,7 +1205,7 @@ with ctx: results = {} evnt.live_detectors = set() - for ifo in ifos: + for ifo in (evnt.ifos if evnt.rank == 0 else evnt.trigg_ifos): results[ifo] = False status = data_reader[ifo].advance( valid_pad, @@ -1259,23 +1283,27 @@ with ctx: # Calculate and add the psd variation for the results if args.psd_variation: - for ifo in results: - logging.info(f"Calculating PSD Variation Statistic for {ifo}") + logging.info("Calculating PSD Variation Statistic for %s", ifo) # A new filter is needed if the PSD has been recalculated if psd_recalculated[ifo] is True: - psd_var_filts[ifo] = variation.live_create_filter(data_reader[ifo].psd, - args.psd_segment_length, - int(args.sample_rate)) + psd_var_filts[ifo] = variation.live_create_filter( + data_reader[ifo].psd, + args.psd_segment_length, + int(args.sample_rate) + ) psd_recalculated[ifo] = False - psd_var_ts = variation.live_calc_psd_variation(data_reader[ifo].strain, - psd_var_filts[ifo], - args.increment) + psd_var_ts = variation.live_calc_psd_variation( + data_reader[ifo].strain, + psd_var_filts[ifo], + args.increment + ) - psd_var_vals = variation.live_find_var_value(results[ifo], - psd_var_ts) + psd_var_vals = variation.live_find_var_value( + results[ifo], psd_var_ts + ) results[ifo]['psd_var_val'] = psd_var_vals @@ -1295,7 +1323,7 @@ with ctx: gates = {ifo: data_reader[ifo].gate_params for ifo in data_reader} # map the results file to an hdf file - prefix = '{}-{}-{}-{}'.format(''.join(sorted(ifos)), + prefix = '{}-{}-{}-{}'.format(''.join(sorted(evnt.ifos)), args.file_prefix, data_end() - args.analysis_chunk, valid_pad) @@ -1310,8 +1338,9 @@ with ctx: data_end() - last_bg_dump_time > float(args.output_background[0]): last_bg_dump_time = int(data_end()) bg_dists = coinc_pool.broadcast(output_background, None) - bg_fn = '{}-LIVE_BACKGROUND-{}.hdf'.format(''.join(sorted(ifos)), - last_bg_dump_time) + bg_fn = '{}-LIVE_BACKGROUND-{}.hdf'.format( + ''.join(sorted(evnt.trigg_ifos)), last_bg_dump_time + ) bg_fn = os.path.join(args.output_background[1], bg_fn) with h5py.File(bg_fn, 'w') as bgf: for bg_ifos, bg_data, bg_time in bg_dists: From e8d671a3bbd2c109e65a5562384334b34ba50a8b Mon Sep 17 00:00:00 2001 From: Thomas Dent Date: Mon, 19 Feb 2024 17:19:59 +0100 Subject: [PATCH 14/29] Clean up options and help for fit over multiparam (#4640) * Clean up options and help for fit over multiparam Older code that this is copied from had the option of reading fit smoothing parameter values from an input fit file, rather than from the template bank (or calculating on the fly). This is no longer done so remove the (useless) option and references to it in help for other options * Update pycbc_fit_sngls_over_multiparam Few more fixes to options, in particular 'required' & defaults --- .../pycbc_fit_sngls_over_multiparam | 34 ++++++++----------- 1 file changed, 14 insertions(+), 20 deletions(-) diff --git a/bin/all_sky_search/pycbc_fit_sngls_over_multiparam b/bin/all_sky_search/pycbc_fit_sngls_over_multiparam index d4b2dac0f76..33e063203b3 100755 --- a/bin/all_sky_search/pycbc_fit_sngls_over_multiparam +++ b/bin/all_sky_search/pycbc_fit_sngls_over_multiparam @@ -174,18 +174,14 @@ parser = argparse.ArgumentParser(usage="", pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action=pycbc.version.Version) -parser.add_argument("--template-fit-file", +parser.add_argument("--template-fit-file", required=True, help="hdf5 file containing fit coefficients for each" " individual template. Required") -parser.add_argument("--bank-file", default=None, - help="hdf file containing template parameters. Required " - "unless reading param from template fit file") +parser.add_argument("--bank-file", required=True, + help="hdf file containing template parameters. Required") parser.add_argument("--output", required=True, help="Location for output file containing smoothed fit " - "coefficients. Required") -parser.add_argument("--use-template-fit-param", action="store_true", - help="Use parameter values stored in the template fit " - "file as template_param for smoothing.") + "coefficients. Required") parser.add_argument("--fit-param", nargs='+', help="Parameter(s) over which to regress the background " "fit coefficients. Required. Either read from " @@ -196,20 +192,19 @@ parser.add_argument("--fit-param", nargs='+', "multiple parameters, provide them as a list.") parser.add_argument("--approximant", default="SEOBNRv4", help="Approximant for template duration. Default SEOBNRv4") -parser.add_argument("--f-lower", type=float, default=0., - help="Starting frequency for calculating template " - "duration, if not reading from the template fit file") +parser.add_argument("--f-lower", type=float, + help="Start frequency for calculating template duration.") parser.add_argument("--min-duration", type=float, default=0., help="Fudge factor for templates with tiny or negative " "values of template_duration: add to duration values" " before fitting. Units seconds.") parser.add_argument("--log-param", nargs='+', - help="Take the log of the fit param before smoothing.") + help="Take the log of the fit param before smoothing. " + "Must be a list corresponding to fit params.") parser.add_argument("--smoothing-width", type=float, nargs='+', required=True, - help="Distance in the space of fit param values (or the " - "logs of them) to smooth over. Required. " - "This must be a list corresponding to the smoothing " - "parameters.") + help="Distance in the space of fit param values (or their" + " logs) to smooth over. Required. Must be a list " + "corresponding to fit params.") parser.add_argument("--smoothing-method", default="smooth_tophat", choices = _smooth_dist_func.keys(), help="Method used to smooth the fit parameters; " @@ -220,15 +215,14 @@ parser.add_argument("--smoothing-method", default="smooth_tophat", "the smoothing until 500 triggers are reached. " "'distance_weighted' weights the closest templates " "with a normal distribution of width smoothing-width " - "trucated at three smoothing-widths.") + "truncated at three smoothing-widths.") parser.add_argument("--smoothing-keywords", nargs='*', help="Keywords for the smoothing function, supplied " "as key:value pairs, e.g. total_trigs:500 to define " - "the number of templates in the n_closest smoothing " - "method") + "the number of templates for n_closest smoothing.") parser.add_argument("--output-fits-by-template", action='store_true', help="If given, will output the input file fits to " - "fit_by_template group") + "fit_by_template group.") args = parser.parse_args() if args.smoothing_keywords: From 42cd5256f5652c034b7ea7b5c2c23a9e9b39616e Mon Sep 17 00:00:00 2001 From: Gareth S Cabourn Davies Date: Tue, 20 Feb 2024 13:08:53 +0000 Subject: [PATCH 15/29] Bugfix: upload prep read psds function (#4644) --- bin/all_sky_search/pycbc_prepare_xml_for_gracedb | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/bin/all_sky_search/pycbc_prepare_xml_for_gracedb b/bin/all_sky_search/pycbc_prepare_xml_for_gracedb index 0baf6cb19dc..f915954884a 100755 --- a/bin/all_sky_search/pycbc_prepare_xml_for_gracedb +++ b/bin/all_sky_search/pycbc_prepare_xml_for_gracedb @@ -96,7 +96,7 @@ def read_psds(psd_files): psd = [group["psds"][str(i)] for i in range(len(group["psds"].keys()))] psds[ifo] = segmentlist(psd_segment(*segargs) for segargs in zip( psd, group["start_time"], group["end_time"])) - return psds + return psds psds = read_psds(args.psd_files) From e7e9666a1b91a2eed770b15f0dcad238f48b2c36 Mon Sep 17 00:00:00 2001 From: Tito Dal Canton Date: Tue, 20 Feb 2024 17:01:37 +0100 Subject: [PATCH 16/29] Fix comment bugs in #4635 (#4645) --- bin/pycbc_live | 22 ++++++++++++---------- 1 file changed, 12 insertions(+), 10 deletions(-) diff --git a/bin/pycbc_live b/bin/pycbc_live index 606f081cf9a..7cfb3af8290 100755 --- a/bin/pycbc_live +++ b/bin/pycbc_live @@ -542,6 +542,9 @@ class LiveEventManager(object): ) logging.info('Coincident candidate! Saving as %s', fname) + # Which IFOs were active? + live_ifos = [ifo for ifo in sld if 'snr_series' in sld[ifo]] + # Verbally explain some details not obvious from the other info comment = ( 'Trigger produced as a {} coincidence.
' @@ -553,9 +556,9 @@ class LiveEventManager(object): comment = comment.format( ppdets(coinc_ifos), args.ranking_statistic, - set(ifos) - self.skymap_only_ifos, - ppdets(followup_ifos), - ppdets(self.skymap_only_ifos) + ppdets(set(ifos) - self.skymap_only_ifos), + ppdets(live_ifos), + ppdets(set(live_ifos) & self.skymap_only_ifos) ) ifar = coinc_results['foreground/ifar'] @@ -577,9 +580,6 @@ class LiveEventManager(object): # even if not running it - do this before the thread so no # data buffers move on in a possible interim period - # Which IFOs were active? - live_ifos = [ifo for ifo in sld if 'snr_series' in sld[ifo]] - # Tell SNR optimized event about p_terr if hasattr(event, 'p_terr') and event.p_terr is not None: coinc_results['p_terr'] = event.p_terr @@ -648,6 +648,9 @@ class LiveEventManager(object): ) logging.info('Single-detector candidate! Saving as %s', fname) + # Which IFOs were active? + live_ifos = [ifo for ifo in sld if 'snr_series' in sld[ifo]] + # Verbally explain some details not obvious from the other info comment = ( 'Trigger produced as a {0} single.
' @@ -656,7 +659,9 @@ class LiveEventManager(object): 'Detectors used only for localization: {2}.' ) comment = comment.format( - ifo, ppdets(followup_ifos), ppdets(self.skymap_only_ifos) + ifo, + ppdets(live_ifos), + ppdets(set(live_ifos) & self.skymap_only_ifos) ) # Has a coinc event at this time been uploaded recently? @@ -694,9 +699,6 @@ class LiveEventManager(object): # where there is already a coinc continue - # Which IFOs were active? - live_ifos = [ifo for ifo in sld if 'snr_series' in sld[ifo]] - # Tell SNR optimized event about p_terr if hasattr(event, 'p_terr') and event.p_terr is not None: single['p_terr'] = event.p_terr From aa05af6a59cb06900f1da8e24387730245c8ba9a Mon Sep 17 00:00:00 2001 From: Ian Harry Date: Tue, 20 Feb 2024 16:21:05 +0000 Subject: [PATCH 17/29] Ensure approximant field is correct format (#4646) --- bin/pycbc_convertinjfiletohdf | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/bin/pycbc_convertinjfiletohdf b/bin/pycbc_convertinjfiletohdf index f0557fa03b2..c871d9a4189 100755 --- a/bin/pycbc_convertinjfiletohdf +++ b/bin/pycbc_convertinjfiletohdf @@ -151,7 +151,11 @@ class LVKNewStyleInjectionSet(object): elif lvk_name == 't_co_gps_add': continue else: - data[pycbc_name] = self.inj_file[f'{self.subdir}/{lvk_name}'][:] + lvk_file_dset = self.inj_file[f'{self.subdir}/{lvk_name}'] + if lvk_file_dset.dtype.char in ['U', 'O']: + data[pycbc_name] = lvk_file_dset[:].astype('S') + else: + data[pycbc_name] = lvk_file_dset[:] return data @@ -187,4 +191,4 @@ if injclass is not None: CBCHDFInjectionSet.write(args.output_file, samples) else: # Assume the injection is a HDF file (PyCBC format), only make a copy - shutil.copy(args.injection_file, args.output_file) \ No newline at end of file + shutil.copy(args.injection_file, args.output_file) From 232a4cef094d98e45809eda9d164cc0a6ab69a40 Mon Sep 17 00:00:00 2001 From: Tito Dal Canton Date: Tue, 20 Feb 2024 18:33:57 +0100 Subject: [PATCH 18/29] Various fixes to `EventManager`s (#4639) * Rename "snr_" datasets to just "snr" * Do not store template duration in trigger files * Cleanup and better document `write_events()` method * Fix how gating info is written * Fix imports * Detriplicate code * Codeclimate --- bin/pycbc_multi_inspiral | 15 +++- pycbc/events/eventmgr.py | 160 ++++++++++++++++----------------------- 2 files changed, 79 insertions(+), 96 deletions(-) diff --git a/bin/pycbc_multi_inspiral b/bin/pycbc_multi_inspiral index c0e8189c57e..ce58cf3cff4 100755 --- a/bin/pycbc_multi_inspiral +++ b/bin/pycbc_multi_inspiral @@ -352,10 +352,16 @@ with ctx: } network_names = sorted(network_out_vals.keys()) event_mgr = EventManagerCoherent( - args, args.instruments, ifo_names, - [ifo_out_types[n] for n in ifo_names], network_names, + args, + args.instruments, + ifo_names, + [ifo_out_types[n] for n in ifo_names], + network_names, [network_out_types[n] for n in network_names], - segments=segments, time_slides=time_slides) + segments=segments, + time_slides=time_slides, + gating_info={det: strain_dict[det].gating_info for det in strain_dict} + ) # Template bank: filtering and thinning logging.info("Read in template bank") @@ -697,6 +703,9 @@ with ctx: event_mgr.cluster_template_network_events( 'time_index', 'reweighted_snr', cluster_window, slide=slide) event_mgr.finalize_template_events() + +logging.info('Writing output') event_mgr.write_events(args.output) + logging.info("Finished") logging.info("Time to complete analysis: %d", int(time.time() - time_init)) diff --git a/pycbc/events/eventmgr.py b/pycbc/events/eventmgr.py index ec065900414..80e4d2c8a64 100644 --- a/pycbc/events/eventmgr.py +++ b/pycbc/events/eventmgr.py @@ -24,10 +24,13 @@ """This modules defines functions for clustering and thresholding timeseries to produces event triggers """ -import numpy, copy, os.path +import os.path +import copy +import itertools import logging -import h5py import pickle +import numpy +import h5py from pycbc.types import Array from pycbc.scheme import schemed @@ -177,6 +180,23 @@ def cluster_reduce(idx, snr, window_size): return idx.take(ind), snr.take(ind) +class H5FileSyntSugar(object): + """Convenience class that adds some syntactic sugar to h5py.File. + """ + def __init__(self, name, prefix=''): + self.f = h5py.File(name, 'w') + self.prefix = prefix + + def __setitem__(self, name, data): + self.f.create_dataset( + self.prefix + '/' + name, + data=data, + compression='gzip', + compression_opts=9, + shuffle=True + ) + + class EventManager(object): def __init__(self, opt, column, column_types, **kwds): self.opt = opt @@ -407,33 +427,21 @@ def save_performance(self, ncores, nfilters, ntemplates, run_time, self.write_performance = True def write_events(self, outname): - """ Write the found events to a sngl inspiral table + """Write the found events to a file. The only currently supported + format is HDF5, indicated by an .hdf or .h5 extension. """ self.make_output_dir(outname) - - if '.hdf' in outname: + if outname.endswith(('.hdf', '.h5')): self.write_to_hdf(outname) - else: - raise ValueError('Cannot write to this format') + return + raise ValueError('Unsupported event output file format') def write_to_hdf(self, outname): - class fw(object): - def __init__(self, name, prefix): - self.f = h5py.File(name, 'w') - self.prefix = prefix - - def __setitem__(self, name, data): - col = self.prefix + '/' + name - self.f.create_dataset(col, data=data, - compression='gzip', - compression_opts=9, - shuffle=True) - self.events.sort(order='template_id') th = numpy.array([p['tmplt'].template_hash for p in self.template_params]) tid = self.events['template_id'] - f = fw(outname, self.opt.channel_name[0:2]) + f = H5FileSyntSugar(outname, self.opt.channel_name[0:2]) if len(self.events): f['snr'] = abs(self.events['snr']) @@ -471,6 +479,10 @@ def __setitem__(self, name, data): # Not precessing f['sigmasq'] = self.events['sigmasq'] + # Template durations should ideally be stored in the bank file. + # At present, however, a few plotting/visualization codes + # downstream in the offline search workflow rely on durations being + # stored in the trigger files instead. template_durations = [p['tmplt'].template_duration for p in self.template_params] f['template_duration'] = numpy.array(template_durations, @@ -595,6 +607,31 @@ def add_template_events_to_ifo(self, ifo, columns, vectors): self.template_event_dict[ifo] = self.template_events self.template_events = None + def write_gating_info_to_hdf(self, hf): + """Write per-detector gating information to an h5py file object. + The information is laid out according to the following groups and + datasets: `//gating/{file, auto}/{time, width, pad}` where + "file" and "auto" indicate respectively externally-provided gates and + internally-generated gates (autogating), and "time", "width" and "pad" + indicate the gate center times, total durations and padding durations + in seconds respectively. + """ + if 'gating_info' not in self.global_params: + return + gates = self.global_params['gating_info'] + for ifo, gate_type in itertools.product(self.ifos, ['file', 'auto']): + if gate_type not in gates[ifo]: + continue + hf[f'{ifo}/gating/{gate_type}/time'] = numpy.array( + [float(g[0]) for g in gates[ifo][gate_type]] + ) + hf[f'{ifo}/gating/{gate_type}/width'] = numpy.array( + [g[1] for g in gates[ifo][gate_type]] + ) + hf[f'{ifo}/gating/{gate_type}/pad'] = numpy.array( + [g[2] for g in gates[ifo][gate_type]] + ) + class EventManagerCoherent(EventManagerMultiDetBase): def __init__(self, opt, ifos, column, column_types, network_column, @@ -656,7 +693,6 @@ def add_template_network_events(self, columns, vectors): """ Add a vector indexed """ # initialize with zeros - since vectors can be None, look for the # longest one that isn't - new_events = None new_events = numpy.zeros( max([len(v) for v in vectors if v is not None]), dtype=self.network_event_dtype @@ -681,19 +717,11 @@ def add_template_events_to_network(self, columns, vectors): self.template_events = None def write_to_hdf(self, outname): - class fw(object): - def __init__(self, name): - self.f = h5py.File(name, 'w') - - def __setitem__(self, name, data): - col = self.prefix + '/' + name - self.f.create_dataset( - col, data=data, compression='gzip', compression_opts=9, - shuffle=True) self.events.sort(order='template_id') th = numpy.array( [p['tmplt'].template_hash for p in self.template_params]) - f = fw(outname) + f = H5FileSyntSugar(outname) + self.write_gating_info_to_hdf(f) # Output network stuff f.prefix = 'network' network_events = numpy.array( @@ -721,8 +749,7 @@ def __setitem__(self, name, data): ifo_events = numpy.array([e for e in self.events if e['ifo'] == self.ifo_dict[ifo]], dtype=self.event_dtype) if len(ifo_events): - ifo_str = ifo.lower()[0] if ifo != 'H1' else ifo.lower() - f['snr_%s' % ifo_str] = abs(ifo_events['snr']) + f['snr'] = abs(ifo_events['snr']) f['event_id'] = ifo_events['event_id'] try: # Precessing @@ -736,8 +763,8 @@ def __setitem__(self, name, data): f['bank_chisq_dof'] = ifo_events['bank_chisq_dof'] f['cont_chisq'] = ifo_events['cont_chisq'] f['end_time'] = ifo_events['time_index'] / \ - float(self.opt.sample_rate[ifo_str]) + \ - self.opt.gps_start_time[ifo_str] + float(self.opt.sample_rate[ifo]) + \ + self.opt.gps_start_time[ifo] f['time_index'] = ifo_events['time_index'] f['slide_id'] = ifo_events['slide_id'] try: @@ -764,11 +791,6 @@ def __setitem__(self, name, data): dtype=numpy.float32) f['sigmasq'] = template_sigmasq[tid] - template_durations = [p['tmplt'].template_duration for p in - self.template_params] - f['template_duration'] = numpy.array(template_durations, - dtype=numpy.float32)[tid] - # FIXME: Can we get this value from the autochisq instance? # cont_dof = self.opt.autochi_number_points # if self.opt.autochi_onesided is None: @@ -820,17 +842,6 @@ def __setitem__(self, name, data): f['search/setup_time_fraction'] = \ numpy.array([float(self.setup_time) / float(self.run_time)]) - if 'gating_info' in self.global_params: - gating_info = self.global_params['gating_info'] - for gate_type in ['file', 'auto']: - if gate_type in gating_info: - f['gating/' + gate_type + '/time'] = numpy.array( - [float(g[0]) for g in gating_info[gate_type]]) - f['gating/' + gate_type + '/width'] = numpy.array( - [g[1] for g in gating_info[gate_type]]) - f['gating/' + gate_type + '/pad'] = numpy.array( - [g[2] for g in gating_info[gate_type]]) - def finalize_template_events(self): # Check that none of the template events have the same time index as an # existing event in events. I.e. don't list the same ifo event multiple @@ -965,41 +976,20 @@ def finalize_template_events(self, perform_coincidence=True, self.template_event_dict[ifo] = numpy.array([], dtype=self.event_dtype) - def write_events(self, outname): - """ Write the found events to a sngl inspiral table - """ - self.make_output_dir(outname) - - if '.hdf' in outname: - self.write_to_hdf(outname) - else: - raise ValueError('Cannot write to this format') - def write_to_hdf(self, outname): - class fw(object): - def __init__(self, name): - self.f = h5py.File(name, 'w') - - def __setitem__(self, name, data): - col = self.prefix + '/' + name - self.f.create_dataset(col, data=data, - compression='gzip', - compression_opts=9, - shuffle=True) - self.events.sort(order='template_id') th = numpy.array([p['tmplt'].template_hash for p in self.template_params]) tid = self.events['template_id'] - f = fw(outname) + f = H5FileSyntSugar(outname) + self.write_gating_info_to_hdf(f) for ifo in self.ifos: f.prefix = ifo ifo_events = numpy.array([e for e in self.events if e['ifo'] == self.ifo_dict[ifo]], dtype=self.event_dtype) if len(ifo_events): - ifo_str = ifo.lower()[0] if ifo != 'H1' else ifo.lower() - f['snr_%s' % ifo_str] = abs(ifo_events['snr']) + f['snr'] = abs(ifo_events['snr']) try: # Precessing f['u_vals'] = ifo_events['u_vals'] @@ -1012,8 +1002,8 @@ def __setitem__(self, name, data): f['bank_chisq_dof'] = ifo_events['bank_chisq_dof'] f['cont_chisq'] = ifo_events['cont_chisq'] f['end_time'] = ifo_events['time_index'] / \ - float(self.opt.sample_rate[ifo_str]) + \ - self.opt.gps_start_time[ifo_str] + float(self.opt.sample_rate[ifo]) + \ + self.opt.gps_start_time[ifo] try: # Precessing template_sigmasq_plus = numpy.array([t['sigmasq_plus'] for @@ -1034,11 +1024,6 @@ def __setitem__(self, name, data): dtype=numpy.float32) f['sigmasq'] = template_sigmasq[tid] - template_durations = [p['tmplt'].template_duration for p in - self.template_params] - f['template_duration'] = \ - numpy.array(template_durations, dtype=numpy.float32)[tid] - # FIXME: Can we get this value from the autochisq instance? cont_dof = self.opt.autochi_number_points if self.opt.autochi_onesided is None: @@ -1093,17 +1078,6 @@ def __setitem__(self, name, data): f['search/setup_time_fraction'] = \ numpy.array([float(self.setup_time) / float(self.run_time)]) - if 'gating_info' in self.global_params: - gating_info = self.global_params['gating_info'] - for gate_type in ['file', 'auto']: - if gate_type in gating_info: - f['gating/' + gate_type + '/time'] = numpy.array( - [float(g[0]) for g in gating_info[gate_type]]) - f['gating/' + gate_type + '/width'] = numpy.array( - [g[1] for g in gating_info[gate_type]]) - f['gating/' + gate_type + '/pad'] = numpy.array( - [g[2] for g in gating_info[gate_type]]) - __all__ = ['threshold_and_cluster', 'findchirp_cluster_over_window', 'threshold', 'cluster_reduce', 'ThresholdCluster', From 6b5fe834a81a071f73d253a00a0ef0932e49d1ef Mon Sep 17 00:00:00 2001 From: Thomas Dent Date: Fri, 23 Feb 2024 19:16:16 +0100 Subject: [PATCH 19/29] SNR optimizer option types (#4650) * SNR optimizer option types Fix issue in non-default settings for some options that Max complained about * Use standard verbosity argument --------- Co-authored-by: maxtrevor <65971534+maxtrevor@users.noreply.github.com> --- bin/pycbc_optimize_snr | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/bin/pycbc_optimize_snr b/bin/pycbc_optimize_snr index 5094a02dc62..513f1d2fac5 100755 --- a/bin/pycbc_optimize_snr +++ b/bin/pycbc_optimize_snr @@ -27,9 +27,9 @@ from pycbc.live import snr_optimizer parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=version.git_verbose_msg) -parser.add_argument('--verbose', action='store_true') parser.add_argument('--params-file', required=True, help='Location of the attributes file created by PyCBC ' 'Live') @@ -45,12 +45,12 @@ parser.add_argument('--psd-files', type=str, nargs='+', 'by PyCBC Live.') parser.add_argument('--approximant', required=True, help='Waveform approximant string.') -parser.add_argument('--snr-threshold', default=4.0, +parser.add_argument('--snr-threshold', type=float, default=4.0, help='If the SNR in ifo X is below this threshold do not ' 'consider it part of the coincidence. Not implemented') -parser.add_argument('--chirp-time-f-lower', default=20., +parser.add_argument('--chirp-time-f-lower', type=float, default=20., help='Starting frequency for chirp time window (Hz).') -parser.add_argument('--chirp-time-window', default=2., +parser.add_argument('--chirp-time-window', type=float, default=2., help='Chirp time window (s).') parser.add_argument('--gracedb-server', metavar='URL', help='URL of GraceDB server API for uploading events. ' From b4c1a9587676dcc033134f35269d4e990732985f Mon Sep 17 00:00:00 2001 From: Alex Nitz Date: Fri, 23 Feb 2024 15:08:19 -0500 Subject: [PATCH 20/29] Add pycbc support for the snowline sampler (#4587) * add support for snowline sampler * add to examples * add snowline to some docs, fixes * cc * fix * Update snowline.py * Update snowline.py * Update snowline.py * Update pool.py * Update snowline.py --- bin/inference/pycbc_inference_model_stats | 7 +- companion.txt | 1 + docs/inference/examples/sampler_platter.rst | 15 ++ examples/inference/samplers/run.sh | 3 +- examples/inference/samplers/snowline_stub.ini | 6 + pycbc/inference/io/__init__.py | 2 + pycbc/inference/io/snowline.py | 32 ++++ pycbc/inference/sampler/__init__.py | 2 + pycbc/inference/sampler/snowline.py | 168 ++++++++++++++++++ pycbc/pool.py | 6 + 10 files changed, 239 insertions(+), 3 deletions(-) create mode 100644 examples/inference/samplers/snowline_stub.ini create mode 100644 pycbc/inference/io/snowline.py create mode 100644 pycbc/inference/sampler/snowline.py diff --git a/bin/inference/pycbc_inference_model_stats b/bin/inference/pycbc_inference_model_stats index 263b6868cd7..9abf8c6dba8 100644 --- a/bin/inference/pycbc_inference_model_stats +++ b/bin/inference/pycbc_inference_model_stats @@ -27,6 +27,7 @@ import shutil import logging import time import numpy +import tqdm import pycbc from pycbc.pool import (use_mpi, choose_pool) @@ -116,7 +117,8 @@ logging.info('Getting samples') with loadfile(opts.input_file, 'r') as fp: # we'll need the shape; all the arrays in the samples group should have the # same shape - shape = fp[fp.samples_group]['loglikelihood'].shape + pick = list(fp[fp.samples_group].keys())[0] + shape = fp[fp.samples_group][pick].shape samples = {} for p in model.variable_params: samples[p] = fp[fp.samples_group][p][()].flatten() @@ -127,7 +129,8 @@ logging.info("Loaded %i samples", samples.size) # get the stats logging.info("Calculating stats") -data = list(pool.map(model_call, enumerate(samples))) +data = list(tqdm.tqdm(pool.imap(model_call, enumerate(samples)), + total=len(samples))) stats = [x[0] for x in data] rec = [x[1] for x in data] diff --git a/companion.txt b/companion.txt index 19ea6b0bb73..d49f4121e9e 100644 --- a/companion.txt +++ b/companion.txt @@ -18,6 +18,7 @@ https://github.com/willvousden/ptemcee/archive/master.tar.gz --extra-index-url https://download.pytorch.org/whl/cpu torch nessai>=0.11.0 +snowline # useful to look at PyCBC Live with htop setproctitle diff --git a/docs/inference/examples/sampler_platter.rst b/docs/inference/examples/sampler_platter.rst index 963232188f1..35f4158da99 100644 --- a/docs/inference/examples/sampler_platter.rst +++ b/docs/inference/examples/sampler_platter.rst @@ -79,6 +79,21 @@ or an external package (multinest). .. literalinclude:: ../../../examples/inference/samplers/multinest_stub.ini :language: ini +============================================================ +`Snowline `_ +============================================================ + +.. literalinclude:: ../../../examples/inference/samplers/snowline_stub.ini + :language: ini + + +============================================================ +`nessai `_ +============================================================ + +.. literalinclude:: ../../../examples/inference/samplers/nessai_stub.ini + :language: ini + If we run these samplers, we create the following plot: .. image:: ../../_include/sample.png diff --git a/examples/inference/samplers/run.sh b/examples/inference/samplers/run.sh index 91e41e7fbb7..bdeaefc8d6a 100755 --- a/examples/inference/samplers/run.sh +++ b/examples/inference/samplers/run.sh @@ -1,5 +1,5 @@ #!/bin/sh -for f in cpnest_stub.ini emcee_stub.ini emcee_pt_stub.ini dynesty_stub.ini ultranest_stub.ini epsie_stub.ini nessai_stub.ini; do +for f in cpnest_stub.ini emcee_stub.ini emcee_pt_stub.ini dynesty_stub.ini ultranest_stub.ini epsie_stub.ini nessai_stub.ini snowline_stub.ini; do echo $f pycbc_inference \ --config-files `dirname $0`/simp.ini `dirname $0`/$f \ @@ -17,6 +17,7 @@ ultranest_stub.ini.hdf:ultranest \ epsie_stub.ini.hdf:espie \ cpnest_stub.ini.hdf:cpnest \ nessai_stub.ini.hdf:nessai \ +snowline_stub.ini.hdf:snowline \ --output-file sample.png \ --plot-contours \ --plot-marginal \ diff --git a/examples/inference/samplers/snowline_stub.ini b/examples/inference/samplers/snowline_stub.ini new file mode 100644 index 00000000000..ac17879fe74 --- /dev/null +++ b/examples/inference/samplers/snowline_stub.ini @@ -0,0 +1,6 @@ +[sampler] +name = snowline + +num_gauss_samples = 500 +min_ess = 500 +max_improvement_loops = 3 diff --git a/pycbc/inference/io/__init__.py b/pycbc/inference/io/__init__.py index 4b3fd0ce909..3157f4408f6 100644 --- a/pycbc/inference/io/__init__.py +++ b/pycbc/inference/io/__init__.py @@ -37,6 +37,7 @@ from .multinest import MultinestFile from .dynesty import DynestyFile from .ultranest import UltranestFile +from .snowline import SnowlineFile from .nessai import NessaiFile from .posterior import PosteriorFile from .txt import InferenceTXTFile @@ -51,6 +52,7 @@ PosteriorFile.name: PosteriorFile, UltranestFile.name: UltranestFile, NessaiFile.name: NessaiFile, + SnowlineFile.name: SnowlineFile, } try: diff --git a/pycbc/inference/io/snowline.py b/pycbc/inference/io/snowline.py new file mode 100644 index 00000000000..4aa9edcbf13 --- /dev/null +++ b/pycbc/inference/io/snowline.py @@ -0,0 +1,32 @@ +# Copyright (C) 2024 Alex Nitz +# This program is free software; you can redistribute it and/or modify it +# under the terms of the GNU General Public License as published by the +# Free Software Foundation; either version 3 of the License, or (at your +# self.option) any later version. +# +# This program is distributed in the hope that it will be useful, but +# WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General +# Public License for more details. +# +# You should have received a copy of the GNU General Public License along +# with this program; if not, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + + +# +# ============================================================================= +# +# Preamble +# +# ============================================================================= +# +"""Provides IO for the snowline sampler. +""" +from .posterior import PosteriorFile + + +class SnowlineFile(PosteriorFile): + """Class to handle file IO for the ``snowline`` sampler.""" + + name = 'snowline_file' diff --git a/pycbc/inference/sampler/__init__.py b/pycbc/inference/sampler/__init__.py index 41da16f39c6..a25573b7ee2 100644 --- a/pycbc/inference/sampler/__init__.py +++ b/pycbc/inference/sampler/__init__.py @@ -25,6 +25,7 @@ from .ultranest import UltranestSampler from .dummy import DummySampler from .refine import RefineSampler +from .snowline import SnowlineSampler # list of available samplers samplers = {cls.name: cls for cls in ( @@ -32,6 +33,7 @@ UltranestSampler, DummySampler, RefineSampler, + SnowlineSampler, )} try: diff --git a/pycbc/inference/sampler/snowline.py b/pycbc/inference/sampler/snowline.py new file mode 100644 index 00000000000..6a54825879a --- /dev/null +++ b/pycbc/inference/sampler/snowline.py @@ -0,0 +1,168 @@ +# Copyright (C) 2023 Alex Nitz +# This program is free software; you can redistribute it and/or modify it +# under the terms of the GNU General Public License as published by the +# Free Software Foundation; either version 3 of the License, or (at your +# option) any later version. +# +# This program is distributed in the hope that it will be useful, but +# WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General +# Public License for more details. +# +# You should have received a copy of the GNU General Public License along +# with this program; if not, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + + +# +# ============================================================================= +# +# Preamble +# +# ============================================================================= +# +""" +This modules provides classes and functions for using the snowline sampler +packages for parameter estimation. +""" + +import sys +import logging + +from pycbc.inference.io.snowline import SnowlineFile +from pycbc.io.hdf import dump_state +from pycbc.pool import use_mpi +from .base import (BaseSampler, setup_output) +from .base_cube import setup_calls + + +# +# ============================================================================= +# +# Samplers +# +# ============================================================================= +# + +class SnowlineSampler(BaseSampler): + """This class is used to construct an Snowline sampler from the snowline + package. + + Parameters + ---------- + model : model + A model from ``pycbc.inference.models`` + """ + name = "snowline" + _io = SnowlineFile + + def __init__(self, model, **kwargs): + super().__init__(model) + + import snowline + log_likelihood_call, prior_call = setup_calls(model, copy_prior=True) + + self._sampler = snowline.ReactiveImportanceSampler( + list(self.model.variable_params), + log_likelihood_call, + transform=prior_call) + + do_mpi, _, rank = use_mpi() + self.main = (not do_mpi) or (rank == 0) + self.nlive = 0 + self.ndim = len(self.model.variable_params) + self.kwargs = kwargs + + def run(self): + self.result = self._sampler.run(**self.kwargs) + if not self.main: + sys.exit(0) + self._sampler.print_results() + + @property + def io(self): + return self._io + + @property + def niterations(self): + return self.result['niter'] + + @classmethod + def from_config(cls, cp, model, output_file=None, **kwds): + """ + Loads the sampler from the given config file. + """ + skeys = {} + opts = {'num_global_samples': int, + 'num_gauss_samples': int, + 'max_ncalls': int, + 'min_ess': int, + 'max_improvement_loops': int + } + for opt_name in opts: + if cp.has_option('sampler', opt_name): + value = cp.get('sampler', opt_name) + skeys[opt_name] = opts[opt_name](value) + inst = cls(model, **skeys) + + do_mpi, _, rank = use_mpi() + if not do_mpi or (rank == 0): + setup_output(inst, output_file) + return inst + + def checkpoint(self): + """ There is currently no checkpointing implemented""" + pass + + def resume_from_checkpoint(self): + """ There is currently no checkpointing implemented""" + pass + + def finalize(self): + logging.info("Writing samples to files") + for fn in [self.checkpoint_file, self.backup_file]: + self.write_results(fn) + + @property + def model_stats(self): + return {} + + @property + def samples(self): + samples = self.result['samples'] + params = list(self.model.variable_params) + samples_dict = {p: samples[:, i] for i, p in enumerate(params)} + return samples_dict + + def write_results(self, filename): + """Writes samples, model stats, acceptance fraction, and random state + to the given file. + + Parameters + ---------- + filename : str + The file to write to. The file is opened using the ``io`` class + in an an append state. + """ + with self.io(filename, 'a') as fp: + # write samples + fp.write_samples(self.samples, self.samples.keys()) + # write log evidence + fp.write_logevidence(self.logz, self.logz_err) + + # write full results + dump_state(self.result, fp, + path='sampler_info', + dsetname='presult') + + @property + def logz(self): + """Return bayesian evidence estimated by snowline sampler. + """ + return self.result['logz'] + + @property + def logz_err(self): + """Return error in bayesian evidence estimated by snowline sampler. + """ + return self.result['logzerr'] diff --git a/pycbc/pool.py b/pycbc/pool.py index 9750b4f4c91..8f2095181ee 100644 --- a/pycbc/pool.py +++ b/pycbc/pool.py @@ -121,6 +121,12 @@ def broadcast(self, fcn, args): def map(self, f, items): return [f(a) for a in items] + # This is single core, so imap and map + # would not behave differently. This is defined + # so that the general pool interfaces can use + # imap irrespective of the pool type. + imap = map + def use_mpi(require_mpi=False, log=True): """ Get whether MPI is enabled and if so the current size and rank """ From b4784f9eab741d24741c843b6987e243e987030e Mon Sep 17 00:00:00 2001 From: Ian Harry Date: Sat, 24 Feb 2024 20:54:25 +0000 Subject: [PATCH 21/29] Avoid pegasus 5.0.7 (#4651) * Avoid pegasus 5.0.7 * Needs to be here, not in API requirements * Revert "Avoid pegasus 5.0.7" This reverts commit c9c915f45c22b5de28118edbbca29665f2afd8cf. --- .github/workflows/inference-workflow.yml | 2 +- .github/workflows/search-workflow.yml | 2 +- .github/workflows/tmpltbank-workflow.yml | 2 +- .github/workflows/workflow-tests.yml | 2 +- 4 files changed, 4 insertions(+), 4 deletions(-) diff --git a/.github/workflows/inference-workflow.yml b/.github/workflows/inference-workflow.yml index 16e5b9cc65d..d95e1c5b714 100644 --- a/.github/workflows/inference-workflow.yml +++ b/.github/workflows/inference-workflow.yml @@ -25,7 +25,7 @@ jobs: wget -qO - https://download.pegasus.isi.edu/pegasus/gpg.txt | sudo apt-key add - echo "deb https://download.pegasus.isi.edu/pegasus/ubuntu bionic main" | sudo tee -a /etc/apt/sources.list sudo apt-get -o Acquire::Retries=3 update - sudo apt-get -o Acquire::Retries=3 install pegasus + sudo apt-get -o Acquire::Retries=3 install pegasus=5.0.6-1+ubuntu18 - run: sudo apt-get -o Acquire::Retries=3 install *fftw3* intel-mkl* - name: Install pycbc run: | diff --git a/.github/workflows/search-workflow.yml b/.github/workflows/search-workflow.yml index 18219411078..9383c5066cb 100644 --- a/.github/workflows/search-workflow.yml +++ b/.github/workflows/search-workflow.yml @@ -30,7 +30,7 @@ jobs: wget -qO - https://download.pegasus.isi.edu/pegasus/gpg.txt | sudo apt-key add - echo "deb https://download.pegasus.isi.edu/pegasus/ubuntu bionic main" | sudo tee -a /etc/apt/sources.list sudo apt-get -o Acquire::Retries=3 update - sudo apt-get -o Acquire::Retries=3 install pegasus + sudo apt-get -o Acquire::Retries=3 install pegasus=5.0.6-1+ubuntu18 - run: sudo apt-get -o Acquire::Retries=3 install *fftw3* intel-mkl* - name: Install pycbc run: | diff --git a/.github/workflows/tmpltbank-workflow.yml b/.github/workflows/tmpltbank-workflow.yml index bf29c9051d0..66a4f1ba3b2 100644 --- a/.github/workflows/tmpltbank-workflow.yml +++ b/.github/workflows/tmpltbank-workflow.yml @@ -29,7 +29,7 @@ jobs: wget -qO - https://download.pegasus.isi.edu/pegasus/gpg.txt | sudo apt-key add - echo "deb https://download.pegasus.isi.edu/pegasus/ubuntu bionic main" | sudo tee -a /etc/apt/sources.list sudo apt-get -o Acquire::Retries=3 update - sudo apt-get -o Acquire::Retries=3 install pegasus + sudo apt-get -o Acquire::Retries=3 install pegasus=5.0.6-1+ubuntu18 - run: sudo apt-get -o Acquire::Retries=3 install *fftw3* intel-mkl* - name: Install pycbc run: | diff --git a/.github/workflows/workflow-tests.yml b/.github/workflows/workflow-tests.yml index 754e25f8f11..86410066aff 100644 --- a/.github/workflows/workflow-tests.yml +++ b/.github/workflows/workflow-tests.yml @@ -34,7 +34,7 @@ jobs: wget -qO - https://download.pegasus.isi.edu/pegasus/gpg.txt | sudo apt-key add - echo "deb https://download.pegasus.isi.edu/pegasus/ubuntu bionic main" | sudo tee -a /etc/apt/sources.list sudo apt-get -o Acquire::Retries=3 update - sudo apt-get -o Acquire::Retries=3 install pegasus + sudo apt-get -o Acquire::Retries=3 install pegasus=5.0.6-1+ubuntu18 - run: sudo apt-get -o Acquire::Retries=3 install *fftw3* intel-mkl* - name: Install pycbc run: | From 6b1707887dc6507eefdab87206df10928d34c322 Mon Sep 17 00:00:00 2001 From: Rajiv Mayani Date: Tue, 27 Feb 2024 00:11:07 -0800 Subject: [PATCH 22/29] Fix Pegasus selector profile key (#4652) Change the selector profile to `execution.site` from `execution_site`. See: https://pegasus.isi.edu/documentation/reference-guide/configuration.html --- pycbc/workflow/pegasus_workflow.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pycbc/workflow/pegasus_workflow.py b/pycbc/workflow/pegasus_workflow.py index e442072e28a..9a80ff09fcd 100644 --- a/pycbc/workflow/pegasus_workflow.py +++ b/pycbc/workflow/pegasus_workflow.py @@ -75,7 +75,7 @@ def set_num_retries(self, number): self.add_profile("dagman", "retry", number) def set_execution_site(self, site): - self.add_profile("selector", "execution_site", site) + self.add_profile("selector", "execution.site", site) class Executable(ProfileShortcuts): From 98007f4bf7add107b06cf770609d15fa3784740d Mon Sep 17 00:00:00 2001 From: Gareth S Cabourn Davies Date: Tue, 27 Feb 2024 12:41:28 +0000 Subject: [PATCH 23/29] Standardize logging: bin, bin/bank and pycbc/tmpltbank (#4631) * Standardize bank logging setup * use loggers in tmpltbank * cc * change the verbose short option to lowercase, to avoid a clash with vary-fupper * add executables in bin * import properly in pycbc_multi_inspiral * mised import * CC * more missed imports * typo --- bin/bank/pycbc_aligned_bank_cat | 9 ++----- bin/bank/pycbc_aligned_stoch_bank | 11 +++----- bin/bank/pycbc_bank_verification | 17 ++++++------- bin/bank/pycbc_brute_bank | 13 +++++++--- bin/bank/pycbc_coinc_bank2hdf | 4 +-- bin/bank/pycbc_geom_aligned_2dstack | 10 +++----- bin/bank/pycbc_geom_aligned_bank | 12 +++------ bin/bank/pycbc_geom_nonspinbank | 10 +++----- bin/bank/pycbc_tmpltbank_to_chi_params | 11 +++----- bin/pycbc_banksim | 7 +++--- bin/pycbc_banksim_combine_banks | 6 +++-- bin/pycbc_banksim_match_combine | 7 ++++-- bin/pycbc_banksim_skymax | 9 ++++--- bin/pycbc_coinc_time | 10 ++++++-- bin/pycbc_compress_bank | 3 ++- bin/pycbc_convertinjfiletohdf | 3 +-- bin/pycbc_create_injections | 11 ++++---- bin/pycbc_data_store | 13 +++++++--- bin/pycbc_faithsim | 4 +-- bin/pycbc_faithsim_collect_results | 8 +++++- bin/pycbc_fit_sngl_trigs | 10 ++------ bin/pycbc_gwosc_segment_query | 12 ++++----- bin/pycbc_hdf5_splitbank | 8 +++--- bin/pycbc_hdf_splitinj | 6 ++++- bin/pycbc_inj_cut | 35 +++++++++++++++----------- bin/pycbc_inspiral | 18 +++++++------ bin/pycbc_inspiral_skymax | 9 ++++--- bin/pycbc_make_banksim | 4 ++- bin/pycbc_make_faithsim | 4 ++- bin/pycbc_make_html_page | 14 +++++++---- bin/pycbc_make_sky_grid | 32 +++++++++++++++++------ bin/pycbc_make_skymap | 8 ++++-- bin/pycbc_merge_inj_hdf | 9 +++---- bin/pycbc_multi_inspiral | 7 +++--- bin/pycbc_optimal_snr | 32 ++++++++++++----------- bin/pycbc_optimize_snr | 4 +-- bin/pycbc_process_sngls | 6 +++-- bin/pycbc_randomize_inj_dist_by_optsnr | 6 ++++- bin/pycbc_single_template | 14 ++++++++--- bin/pycbc_source_probability_offline | 3 ++- bin/pycbc_split_inspinj | 7 +++++- bin/pycbc_splitbank | 6 +++-- bin/pycbc_upload_xml_to_gracedb | 9 +++++-- pycbc/__init__.py | 6 ++--- pycbc/tmpltbank/bank_conversions.py | 4 +++ pycbc/tmpltbank/bank_output_utils.py | 5 ++++ pycbc/tmpltbank/brute_force_methods.py | 4 +++ pycbc/tmpltbank/calc_moments.py | 4 +++ pycbc/tmpltbank/coord_utils.py | 10 +++++--- pycbc/tmpltbank/lambda_mapping.py | 7 +++++- pycbc/tmpltbank/lattice_utils.py | 4 +++ pycbc/tmpltbank/option_utils.py | 20 +++++++++------ pycbc/tmpltbank/partitioned_bank.py | 5 +++- 53 files changed, 307 insertions(+), 203 deletions(-) diff --git a/bin/bank/pycbc_aligned_bank_cat b/bin/bank/pycbc_aligned_bank_cat index e09b6f274d9..f564c3b9bad 100644 --- a/bin/bank/pycbc_aligned_bank_cat +++ b/bin/bank/pycbc_aligned_bank_cat @@ -48,9 +48,8 @@ __program__ = "pycbc_aligned_bank_cat" parser = argparse.ArgumentParser(description=__doc__, formatter_class=tmpltbank.IndentedHelpFormatterWithNL) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=__version__) -parser.add_argument("--verbose", action="store_true", default=False, - help="verbose output") parser.add_argument("-i", "--input-glob", help="file glob the list of paramters") parser.add_argument("-I", "--input-files", nargs='+', @@ -76,11 +75,7 @@ tmpltbank.insert_base_bank_options(parser, match_req=False) options = parser.parse_args() -if options.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -logging.basicConfig(format='%(asctime)s %(message)s', level=log_level) +pycbc.init_logging(options.verbose) # Sanity check options if not options.output_file: diff --git a/bin/bank/pycbc_aligned_stoch_bank b/bin/bank/pycbc_aligned_stoch_bank index 9b9a473d3fa..e615fd75415 100644 --- a/bin/bank/pycbc_aligned_stoch_bank +++ b/bin/bank/pycbc_aligned_stoch_bank @@ -23,6 +23,7 @@ Stochastic aligned spin bank generator. import argparse import numpy import logging + import pycbc import pycbc.version from pycbc import tmpltbank @@ -45,8 +46,7 @@ parser = argparse.ArgumentParser(description=_desc, # Begin with code specific options parser.add_argument("--version", action="version", version=__version__) -parser.add_argument("--verbose", action="store_true", default=False, - help="verbose output") +pycbc.add_common_pycbc_options(parser) parser.add_argument("-V", "--vary-fupper", action="store_true", default=False, help="Use a variable upper frequency cutoff in laying " "out the bank. OPTIONAL.") @@ -96,12 +96,7 @@ tmpltbank.insert_ethinca_metric_options(parser) opts = parser.parse_args() -if opts.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -log_format='%(asctime)s %(message)s' -logging.basicConfig(format=log_format, level=log_level) +pycbc.init_logging(opts.verbose) # delete defaults for redundant options if not varying fupper if not opts.vary_fupper: diff --git a/bin/bank/pycbc_bank_verification b/bin/bank/pycbc_bank_verification index 9c76e28aaf5..6758daf9556 100644 --- a/bin/bank/pycbc_bank_verification +++ b/bin/bank/pycbc_bank_verification @@ -26,14 +26,16 @@ import argparse import logging import numpy import h5py +import matplotlib +matplotlib.use('Agg') +import pylab + from ligo.lw import lsctables, utils as ligolw_utils + import pycbc import pycbc.version from pycbc import tmpltbank, psd, strain from pycbc.io.ligolw import LIGOLWContentHandler -import matplotlib -matplotlib.use('Agg') -import pylab __author__ = "Ian Harry " @@ -48,8 +50,7 @@ parser = argparse.ArgumentParser(description=__doc__, # Begin with code specific options parser.add_argument("--version", action="version", version=__version__) -parser.add_argument("--verbose", action="store_true", default=False, - help="verbose output") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--histogram-output-file", action="store", default=None, help="Output a histogram of fitting factors to the " "supplied file. If not given no histogram is produced.") @@ -111,11 +112,7 @@ pycbc.strain.insert_strain_option_group(parser) opts = parser.parse_args() -if opts.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -logging.basicConfig(format='%(asctime)s %(message)s', level=log_level) +pycbc.init_logging(opts.verbose) # Sanity check options tmpltbank.verify_metric_calculation_options(opts, parser) diff --git a/bin/bank/pycbc_brute_bank b/bin/bank/pycbc_brute_bank index e887be710a1..583793b43c9 100644 --- a/bin/bank/pycbc_brute_bank +++ b/bin/bank/pycbc_brute_bank @@ -19,16 +19,21 @@ """Generate a bank of templates using a brute force stochastic method. """ -import numpy, h5py, logging, argparse, numpy.random +import numpy +import h5py +import logging +import argparse +import numpy.random +from scipy.stats import gaussian_kde + import pycbc.waveform, pycbc.filter, pycbc.types, pycbc.psd, pycbc.fft, pycbc.conversions from pycbc import transforms from pycbc.waveform.spa_tmplt import spa_length_in_time from pycbc.distributions import read_params_from_config from pycbc.distributions.utils import draw_samples_from_config, prior_from_config -from scipy.stats import gaussian_kde parser = argparse.ArgumentParser(description=__doc__) -parser.add_argument('--verbose', action='store_true') +pycbc.add_common_pycbc_options(parser) parser.add_argument('--output-file', required=True, help='Output file name for template bank.') parser.add_argument('--input-file', @@ -71,10 +76,10 @@ parser.add_argument('--tau0-end', type=float) parser.add_argument('--tau0-cutoff-frequency', type=float, default=15.0) pycbc.psd.insert_psd_option_group(parser) args = parser.parse_args() + pycbc.init_logging(args.verbose) numpy.random.seed(args.seed) - # Read the .ini file if it's in the input. if args.input_config is not None: config_parser = pycbc.types.config.InterpolatingConfigParser() diff --git a/bin/bank/pycbc_coinc_bank2hdf b/bin/bank/pycbc_coinc_bank2hdf index fb03caea604..db35f48d4be 100644 --- a/bin/bank/pycbc_coinc_bank2hdf +++ b/bin/bank/pycbc_coinc_bank2hdf @@ -23,6 +23,7 @@ with their template. import argparse import logging import numpy + import pycbc from pycbc.waveform import bank @@ -56,6 +57,7 @@ def parse_parameters(parameters): parser = argparse.ArgumentParser() parser.add_argument('--version', action='version', version=pycbc.version.git_verbose_msg) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--bank-file', required=True, help="The bank file to load. Must end in '.xml[.gz]' " "and must contain a SnglInspiral table or must end " @@ -82,8 +84,6 @@ bank.add_approximant_arg(parser, parser.add_argument("--force", action="store_true", default=False, help="Overwrite the given hdf file if it exists. " "Otherwise, an error is raised.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Be verbose.") args = parser.parse_args() pycbc.init_logging(args.verbose) diff --git a/bin/bank/pycbc_geom_aligned_2dstack b/bin/bank/pycbc_geom_aligned_2dstack index 179c74f8216..215440e76c1 100644 --- a/bin/bank/pycbc_geom_aligned_2dstack +++ b/bin/bank/pycbc_geom_aligned_2dstack @@ -31,6 +31,7 @@ import copy import numpy import logging import h5py + import pycbc.tmpltbank import pycbc.version from pycbc import pnutils @@ -47,10 +48,9 @@ usage = """usage: %prog [options]""" _desc = __doc__[1:] parser = argparse.ArgumentParser(usage, description=_desc, formatter_class=pycbc.tmpltbank.IndentedHelpFormatterWithNL) +pycbc.add_common_pycbc_options(parser) # Code specific options parser.add_argument('--version', action='version', version=__version__) -parser.add_argument("--verbose", action="store_true", default=False,\ - help="verbose output") parser.add_argument("--pn-order", action="store", type=str,\ default=None,\ help="Determines the PN order to use. Note that if you "+\ @@ -106,11 +106,7 @@ opts = parser.parse_args() opts.eval_vec4_depth = not opts.skip_vec4_depth opts.eval_vec5_depth = not opts.skip_vec5_depth -if opts.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -logging.basicConfig(format='%(asctime)s %(message)s', level=log_level) +pycbc.init_logging(opts.verbose) # Sanity check options if not opts.pn_order: diff --git a/bin/bank/pycbc_geom_aligned_bank b/bin/bank/pycbc_geom_aligned_bank index 97da6a5674c..d18c3490f91 100644 --- a/bin/bank/pycbc_geom_aligned_bank +++ b/bin/bank/pycbc_geom_aligned_bank @@ -31,8 +31,10 @@ import logging import h5py import configparser from scipy import spatial + from ligo.lw import ligolw from ligo.lw import utils as ligolw_utils + import pycbc import pycbc.psd import pycbc.strain @@ -175,9 +177,8 @@ parser = argparse.ArgumentParser(description=_desc, formatter_class=pycbc.tmpltbank.IndentedHelpFormatterWithNL) # Begin with code specific options +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=__version__) -parser.add_argument("-V", "--verbose", action="store_true", default=False, - help="verbose output") parser.add_argument("-s", "--stack-distance", action="store", type=float,\ default=0.2, help="Minimum metric spacing before we "+\ "stack.") @@ -249,12 +250,7 @@ outdoc.appendChild(ligolw.LIGO_LW()) create_process_table(outdoc, __program__, options=vars(opts)) ligolw_utils.write_filename(outdoc, opts.metadata_file) -if opts.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -log_format='%(asctime)s %(message)s' -logging.basicConfig(format=log_format, level=log_level) +pycbc.init_logging(opts.verbose) opts.max_mismatch = 1 - opts.min_match diff --git a/bin/bank/pycbc_geom_nonspinbank b/bin/bank/pycbc_geom_nonspinbank index 79298eff547..dcf5aa71890 100644 --- a/bin/bank/pycbc_geom_nonspinbank +++ b/bin/bank/pycbc_geom_nonspinbank @@ -25,6 +25,7 @@ import math import copy import numpy import logging + import pycbc import pycbc.version import pycbc.psd @@ -48,8 +49,7 @@ parser = argparse.ArgumentParser( # Begin with code specific options parser.add_argument('--version', action='version', version=__version__) -parser.add_argument("-V", "--verbose", action="store_true", default=False, - help="verbose output") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--random-seed", action="store", type=int, default=None, help="""Random seed to use when calling numpy.random @@ -78,11 +78,7 @@ tmpltbank.insert_ethinca_metric_options(parser) opts = parser.parse_args() -if opts.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -logging.basicConfig(format='%(asctime)s %(message)s', level=log_level) +pycbc.init_logging(opts.verbose) opts.max_mismatch = 1 - opts.min_match tmpltbank.verify_metric_calculation_options(opts, parser) diff --git a/bin/bank/pycbc_tmpltbank_to_chi_params b/bin/bank/pycbc_tmpltbank_to_chi_params index df0201ccc0f..ad7429b5f86 100644 --- a/bin/bank/pycbc_tmpltbank_to_chi_params +++ b/bin/bank/pycbc_tmpltbank_to_chi_params @@ -33,7 +33,9 @@ import argparse import logging import numpy import h5py + from ligo.lw import lsctables, utils as ligolw_utils + from pycbc import tmpltbank, psd, strain from pycbc.io.ligolw import LIGOLWContentHandler @@ -44,8 +46,7 @@ parser = argparse.ArgumentParser(description=__doc__, # Begin with code specific options parser.add_argument("--version", action="version", version=__version__) -parser.add_argument("--verbose", action="store_true", default=False, - help="verbose output") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--input-bank", action="store", required=True, help="The template bank to use an input.") parser.add_argument("--output-file", action="store", required=True, @@ -88,11 +89,7 @@ pycbc.strain.insert_strain_option_group(parser) opts = parser.parse_args() -if opts.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -logging.basicConfig(format='%(asctime)s %(message)s', level=log_level) +pycbc.init_logging(opts.verbose) # Sanity check options tmpltbank.verify_metric_calculation_options(opts, parser) diff --git a/bin/pycbc_banksim b/bin/pycbc_banksim index a7598f43226..66cf7e9f531 100644 --- a/bin/pycbc_banksim +++ b/bin/pycbc_banksim @@ -21,8 +21,11 @@ import logging from tqdm import tqdm from numpy import complex64, array from argparse import ArgumentParser +from math import ceil, log + from ligo.lw import utils as ligolw_utils from ligo.lw import lsctables + from pycbc.pnutils import mass1_mass2_to_mchirp_eta from pycbc.pnutils import mass1_mass2_to_tau0_tau3 from pycbc.pnutils import nearest_larger_binary_number @@ -32,7 +35,6 @@ from pycbc import DYN_RANGE_FAC from pycbc.types import FrequencySeries, TimeSeries, zeros, complex_same_precision_as from pycbc.filter import match, sigmasq from pycbc.io.ligolw import LIGOLWContentHandler -from math import ceil, log import pycbc.psd, pycbc.scheme, pycbc.fft, pycbc.strain, pycbc.version from pycbc.detector import overhead_antenna_pattern as generate_fplus_fcross from pycbc.waveform import TemplateBank @@ -144,8 +146,7 @@ aprs = sorted(list(set(td_approximants() + fd_approximants()))) parser = ArgumentParser(description=__doc__) parser.add_argument("--match-file", dest="out_file", metavar="FILE", required=True, help="File to output match results") -parser.add_argument("--verbose", action='store_true', default=False, - help="Print verbose statements") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) diff --git a/bin/pycbc_banksim_combine_banks b/bin/pycbc_banksim_combine_banks index 3c3570b27f4..56f382f19ee 100644 --- a/bin/pycbc_banksim_combine_banks +++ b/bin/pycbc_banksim_combine_banks @@ -26,6 +26,7 @@ from os.path import isfile import argparse import logging from numpy import * + import pycbc import pycbc.version @@ -37,8 +38,7 @@ _desc = __doc__[1:] parser = argparse.ArgumentParser(description=_desc) parser.add_argument('--version', action=pycbc.version.Version) -parser.add_argument("--verbose", action="store_true", default=False, - help="verbose output") +pycbc.add_common_pycbc_options(parser) parser.add_argument("-I", "--input-files", nargs='+', help="Explicit list of input files.") parser.add_argument("-o", "--output-file", required=True, @@ -46,6 +46,8 @@ parser.add_argument("-o", "--output-file", required=True, options = parser.parse_args() +pycbc.init_logging(options.verbose) + dtypef = dtype([('match', float64), ('bank', unicode_, 256), ('bank_i', int32), ('sim', unicode_, 256), ('sim_i', int32), ('sigmasq', float64)]) diff --git a/bin/pycbc_banksim_match_combine b/bin/pycbc_banksim_match_combine index f2008744887..1692504e4c1 100644 --- a/bin/pycbc_banksim_match_combine +++ b/bin/pycbc_banksim_match_combine @@ -26,7 +26,9 @@ import imp import argparse import numpy as np import h5py + from ligo.lw import utils, table + import pycbc from pycbc import pnutils from pycbc.waveform import TemplateBank @@ -43,8 +45,7 @@ __program__ = "pycbc_banksim_match_combine" parser = argparse.ArgumentParser(description=__doc__) parser.add_argument("--version", action="version", version=__version__) -parser.add_argument("--verbose", action="store_true", default=False, - help="verbose output") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--match-files", nargs='+', help="Explicit list of match files.") #parser.add_argument("--inj-files", nargs='+', @@ -62,6 +63,8 @@ parser.add_argument("--filter-func-file", default=None, "as numpy arrays.") options = parser.parse_args() +pycbc.init_logging(options.verbose) + dtypem = np.dtype([('match', np.float64), ('bank', np.unicode_, 256), ('bank_i', np.int32), ('sim', np.unicode_, 256), ('sim_i', np.int32), ('sigmasq', np.float64)]) diff --git a/bin/pycbc_banksim_skymax b/bin/pycbc_banksim_skymax index 40a208b9132..91e0973ca4f 100644 --- a/bin/pycbc_banksim_skymax +++ b/bin/pycbc_banksim_skymax @@ -21,8 +21,12 @@ import sys import logging from numpy import complex64, float32, sqrt, argmax, real, array from argparse import ArgumentParser +from math import ceil, log +from tqdm import tqdm + from ligo.lw import utils as ligolw_utils from ligo.lw import lsctables + from pycbc.pnutils import mass1_mass2_to_mchirp_eta, f_SchwarzISCO from pycbc.pnutils import nearest_larger_binary_number from pycbc.pnutils import mass1_mass2_to_tau0_tau3 @@ -35,11 +39,9 @@ from pycbc.filter import overlap_cplx, matched_filter from pycbc.filter import compute_max_snr_over_sky_loc_stat from pycbc.filter import compute_max_snr_over_sky_loc_stat_no_phase from pycbc.io.ligolw import LIGOLWContentHandler -from math import ceil, log import pycbc.psd, pycbc.scheme, pycbc.fft, pycbc.strain, pycbc.version from pycbc.detector import overhead_antenna_pattern as generate_fplus_fcross from pycbc.waveform import TemplateBank -from tqdm import tqdm ## Remove the need for these functions ######################################## @@ -150,8 +152,7 @@ aprs = sorted(list(set(td_approximants() + fd_approximants()))) parser = ArgumentParser(description=__doc__) parser.add_argument("--match-file", dest="out_file", metavar="FILE", required=True, help="File to output match results") -parser.add_argument("--verbose", action='store_true', default=False, - help="Print verbose statements") +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) diff --git a/bin/pycbc_coinc_time b/bin/pycbc_coinc_time index 7f0803fa5c4..34e4687ea63 100644 --- a/bin/pycbc_coinc_time +++ b/bin/pycbc_coinc_time @@ -1,7 +1,13 @@ #!/bin/env python -import argparse, logging, pycbc.version, ligo.segments, numpy +import argparse +import logging +import numpy from dqsegdb.apicalls import dqsegdbQueryTimes as query + +import ligo.segments + #from pycbc.workflow.segment import cat_to_veto_def_cat as convert_cat +import pycbc.version def sane(seg_list): """ Convert list of len two lists containing strs to segment list """ @@ -85,7 +91,7 @@ def get_vetoes(veto_def, ifo, server, veto_name, start_default, end_default): parser = argparse.ArgumentParser() parser.add_argument('--version', action='version', version=pycbc.version.git_verbose_msg) -parser.add_argument('--verbose', action='store_true') +pycbc.add_common_pycbc_options(parser) parser.add_argument('--gps-start-time', type=int, required=True, help="integer gps start time") diff --git a/bin/pycbc_compress_bank b/bin/pycbc_compress_bank index 851e6bd78b9..0d98bbffc75 100755 --- a/bin/pycbc_compress_bank +++ b/bin/pycbc_compress_bank @@ -30,6 +30,7 @@ import argparse import numpy import h5py import logging + import pycbc from pycbc import psd, DYN_RANGE_FAC from pycbc.waveform import compress @@ -40,6 +41,7 @@ from pycbc import filter parser = argparse.ArgumentParser(description=__description__) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--bank-file", type=str, required=True, help="Bank hdf file to load.") parser.add_argument("--output", type=str, required=True, @@ -96,7 +98,6 @@ parser.add_argument("--precision", type=str, choices=["double", "single"], parser.add_argument("--force", action="store_true", default=False, help="Overwrite the given hdf file if it exists. " "Otherwise, an error is raised.") -parser.add_argument("--verbose", action="store_true", default=False) # Insert the PSD options pycbc.psd.insert_psd_option_group(parser, include_data_options=False) diff --git a/bin/pycbc_convertinjfiletohdf b/bin/pycbc_convertinjfiletohdf index c871d9a4189..0998024f30a 100755 --- a/bin/pycbc_convertinjfiletohdf +++ b/bin/pycbc_convertinjfiletohdf @@ -160,6 +160,7 @@ class LVKNewStyleInjectionSet(object): parser = argparse.ArgumentParser() +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=pycbc.version.git_verbose_msg) parser.add_argument('--injection-file', required=True, @@ -167,8 +168,6 @@ parser.add_argument('--injection-file', required=True, "and must contain a SimInspiral table") parser.add_argument('--output-file', required=True, help="The ouput file name. Must end in '.hdf'.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Be verbose.") args = parser.parse_args() pycbc.init_logging(args.verbose) diff --git a/bin/pycbc_create_injections b/bin/pycbc_create_injections index a9616808a94..3fac44c3d66 100644 --- a/bin/pycbc_create_injections +++ b/bin/pycbc_create_injections @@ -108,10 +108,14 @@ along with any numpy ufunc. So, in the above example, mass ratio (q) is constrained to be <= 4 by using a function from the conversions module. """ -import os, sys +import os +import sys import argparse import logging import numpy +import h5py +from numpy.random import uniform + import pycbc import pycbc.version from pycbc.inject import InjectionSet @@ -122,12 +126,11 @@ from pycbc.io import record from pycbc.waveform import parameters from pycbc.workflow import configuration from pycbc.workflow import WorkflowConfigParser -from numpy.random import uniform -import h5py parser = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) configuration.add_workflow_command_line_group(parser) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=pycbc.version.git_verbose_msg, help='Prints version information.') @@ -163,8 +166,6 @@ parser.add_argument('--static-params-section', default="static_params", parser.add_argument("--force", action="store_true", default=False, help="If the output-file already exists, overwrite it. " "Otherwise, an OSError is raised.") -parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging messages.") opts = parser.parse_args() pycbc.init_logging(opts.verbose) diff --git a/bin/pycbc_data_store b/bin/pycbc_data_store index eb501e709a3..a6b17bac67f 100755 --- a/bin/pycbc_data_store +++ b/bin/pycbc_data_store @@ -1,16 +1,23 @@ #!/usr/bin/env python """ Create HDF strain cache file """ -import logging, argparse, numpy, h5py -import pycbc, pycbc.strain, pycbc.dq +import logging +import argparse +import numpy +import h5py + +import pycbc +import pycbc.strain +import pycbc.dq from pycbc.version import git_verbose_msg as version from pycbc.fft.fftw import set_measure_level from pycbc.events.veto import segments_to_start_end + set_measure_level(0) parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action='version', version=version) -parser.add_argument('--verbose', action="store_true") parser.add_argument("--science-name", help="Science flag definition") parser.add_argument("--segment-server") parser.add_argument("--veto-definer-file") diff --git a/bin/pycbc_faithsim b/bin/pycbc_faithsim index c1405316909..e1bf979a15a 100644 --- a/bin/pycbc_faithsim +++ b/bin/pycbc_faithsim @@ -27,6 +27,7 @@ import logging from numpy import complex64 import argparse import sys + from ligo.lw import utils as ligolw_utils from ligo.lw import lsctables @@ -84,8 +85,7 @@ parser = argparse.ArgumentParser(usage='', description="Calculate faithfulness for a set of waveforms.") parser.add_argument('--version', action='version', version=pycbc.version.git_verbose_msg) -parser.add_argument("-V", "--verbose", action="store_true", - help="print extra debugging information", default=False) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--param-file", dest="bank_file", metavar="FILE", help="Sngl or Sim Inspiral Table containing waveform " "parameters.") diff --git a/bin/pycbc_faithsim_collect_results b/bin/pycbc_faithsim_collect_results index 7d98ac02b22..a574e5e0415 100755 --- a/bin/pycbc_faithsim_collect_results +++ b/bin/pycbc_faithsim_collect_results @@ -5,12 +5,16 @@ Program for collecting the results of pycbc_faithsim comparing two approximants computing the match between them and creating a .dat file with the results. """ +import argparse import numpy as np + from ligo.lw import utils, table, lsctables + +from pycbc import add_common_pycbc_options, init_logging from pycbc.io.ligolw import LIGOLWContentHandler -import argparse parser = argparse.ArgumentParser(description=__doc__) +add_common_pycbc_options(parser) parser.add_argument( "--match-inputs", action="append", @@ -26,6 +30,8 @@ parser.add_argument( parser.add_argument("--output", required=True, help="name of the output .dat file") args = parser.parse_args() +init_logging(options.verbose) + mfields = ("match", "overlap", "time_offset", "sigma1", "sigma2") bfields = ( "match", diff --git a/bin/pycbc_fit_sngl_trigs b/bin/pycbc_fit_sngl_trigs index cef76488606..9f0c347656d 100644 --- a/bin/pycbc_fit_sngl_trigs +++ b/bin/pycbc_fit_sngl_trigs @@ -55,10 +55,8 @@ def get_bins(opt, pmin, pmax): parser = argparse.ArgumentParser(usage="", description="Perform maximum-likelihood fits of single inspiral trigger" "distributions to various functions") - +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action=pycbc.version.Version) -parser.add_argument("-V", "--verbose", action="store_true", - help="Print extra debugging information", default=False) parser.add_argument("--inputs", nargs="+", help="Input file or space-separated list of input files " "containing single triggers. Currently .xml(.gz) " @@ -134,11 +132,7 @@ binopt.add_argument("--irregular-bins", opt = parser.parse_args() -if opt.verbose: - log_level = logging.DEBUG -else: - log_level = logging.WARN -logging.basicConfig(format='%(asctime)s : %(message)s', level=log_level) +pycbc.init_logging(opt.verbose) statname = opt.sngl_stat.replace("_", " ") paramname = opt.fit_param.replace("_", " ") diff --git a/bin/pycbc_gwosc_segment_query b/bin/pycbc_gwosc_segment_query index 4403af29168..0ac34512a33 100644 --- a/bin/pycbc_gwosc_segment_query +++ b/bin/pycbc_gwosc_segment_query @@ -6,14 +6,11 @@ import json import argparse import shutil from urllib.request import urlopen -import ligo.segments -from pycbc.workflow import SegFile +import ligo.segments -# Logging formatting from pycbc optimal snr -log_fmt = '%(asctime)s %(message)s' -log_date_fmt = '%Y-%m-%d %H:%M:%S' -logging.basicConfig(level=logging.INFO, format=log_fmt, datefmt=log_date_fmt) +import pycbc +from pycbc.workflow import SegFile # Function to query json segment data from GWOSC @@ -64,6 +61,7 @@ def write_xml_file(ifo, summary_segment, segments, filename): parser = argparse.ArgumentParser() +pycbc.add_common_pycbc_options(parser) parser.add_argument('--gps-start-time', type=int, required=True) parser.add_argument('--gps-end-time', type=int, required=True) parser.add_argument('--include-segments', type=str, required=True) @@ -71,6 +69,8 @@ parser.add_argument('--output-file', type=str, required=True) parser.add_argument('--protract-hw-inj', type=int, default=0) args = parser.parse_args() +pycbc.init_logging(args.verbose) + gps_start_time = args.gps_start_time gps_end_time = args.gps_end_time duration = gps_end_time - gps_start_time diff --git a/bin/pycbc_hdf5_splitbank b/bin/pycbc_hdf5_splitbank index e86e649c28b..d080d5774bb 100755 --- a/bin/pycbc_hdf5_splitbank +++ b/bin/pycbc_hdf5_splitbank @@ -25,13 +25,15 @@ import argparse import numpy import h5py import logging +from numpy import random + import pycbc, pycbc.version from pycbc.waveform import bank -from numpy import random __author__ = "Soumi De " parser = argparse.ArgumentParser(description=__doc__[1:]) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) parser.add_argument("--bank-file", type=str, @@ -61,11 +63,11 @@ parser.add_argument("--random-seed", type=int, parser.add_argument("--force", action="store_true", default=False, help="Overwrite the given hdf file if it exists. " "Otherwise, an error is raised.") -parser.add_argument("--verbose", action="store_true", default=False) - args = parser.parse_args() +pycbc.init_logging(args.verbose) + # input checks if args.mchirp_sort and (args.random_seed is not None): raise RuntimeError("Can't supply a random seed if not sorting randomly!") diff --git a/bin/pycbc_hdf_splitinj b/bin/pycbc_hdf_splitinj index 0d255e3850a..6a6050f5abd 100644 --- a/bin/pycbc_hdf_splitinj +++ b/bin/pycbc_hdf_splitinj @@ -6,14 +6,16 @@ Split sets are organized to maximize time between injections. """ import argparse -from pycbc.inject import InjectionSet import h5py import numpy as np + +from pycbc.inject import InjectionSet import pycbc.version # Parse command line parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) parser.add_argument("-f", "--output-files", nargs='*', required=True, @@ -23,6 +25,8 @@ parser.add_argument("-i", "--input-file", required=True, args = parser.parse_args() +pycbc.init_logging(args.verbose) + # Read in input file as both an hdf file and an InjectionSet object inj_file = h5py.File(args.input_file, 'r') inj_set = InjectionSet(args.input_file) diff --git a/bin/pycbc_inj_cut b/bin/pycbc_inj_cut index 0622d21f0d2..1ca70bcb5ca 100644 --- a/bin/pycbc_inj_cut +++ b/bin/pycbc_inj_cut @@ -29,14 +29,17 @@ __email__ = "stephen.fairhurst@ligo.org" import argparse import logging import numpy + +from ligo.lw import utils as ligolw_utils +from ligo.lw import lsctables + import pycbc import pycbc.inject -from ligo.lw import utils as ligolw_utils from pycbc.types import MultiDetOptionAction -from ligo.lw import lsctables import pycbc.version parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action=pycbc.version.Version) parser.add_argument('--input', dest='inj_xml', required=True, help='Input LIGOLW injections file.') parser.add_argument('--output-missed', dest='output_missed', required=False, @@ -53,20 +56,13 @@ parser.add_argument('--snr-columns', nargs='+', action=MultiDetOptionAction, ' no useful data in it, good candidates are usually' \ ' alpha1, alpha2 etc.') parser.add_argument("-z", "--write-compress", action="store_true", help="Write compressed xml.gz files.") -parser.add_argument("-v", "--verbose", action="store_true", default=False, - help="Extended standard output.") opts = parser.parse_args() -log_fmt = '%(asctime)s %(message)s' -log_date_fmt = '%Y-%m-%d %H:%M:%S' -logging.basicConfig(level=logging.INFO, format=log_fmt, datefmt=log_date_fmt) +pycbc.init_logging(opts.verbose) -if opts.verbose: - logging.info("Loading injections") injections = pycbc.inject.InjectionSet(opts.inj_xml) - # injections that we should do: table and name of the file that will store it output_inject = opts.output_inject if opts.write_compress: @@ -89,11 +85,20 @@ for i, inj in enumerate(injections.table): else: out_sim_inspiral_inject.append(inj) -logging.info('Found %d/%d (potentially) found injections. Storing them to %s.', len(out_sim_inspiral_inject), len(injections.table), output_inject) -logging.info('Found %d/%d missed injections.', len(out_sim_inspiral_missed), len(injections.table)) - -if opts.verbose: - logging.info('Writing output') +logging.info( + 'Found %d/%d (potentially) found injections. Storing them to %s.', + len(out_sim_inspiral_inject), + len(injections.table), + output_inject +) + +logging.info( + 'Found %d/%d missed injections.', + len(out_sim_inspiral_missed), + len(injections.table) +) + +logging.info('Writing output') llw_doc = injections.indoc llw_root = llw_doc.childNodes[0] llw_root.removeChild(injections.table) diff --git a/bin/pycbc_inspiral b/bin/pycbc_inspiral index 4d5e6e9a585..3d4924c3e0b 100644 --- a/bin/pycbc_inspiral +++ b/bin/pycbc_inspiral @@ -16,19 +16,24 @@ # with this program; if not, write to the Free Software Foundation, Inc., # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. -import sys, os, copy -import logging, argparse, numpy, itertools +import sys +import os +import copy +import logging +import argparse +import numpy +import itertools +import time +from multiprocessing import Pool + import pycbc import pycbc.version from pycbc import vetoes, psd, waveform, strain, scheme, fft, DYN_RANGE_FAC, events from pycbc.vetoes.sgchisq import SingleDetSGChisq from pycbc.filter import MatchedFilterControl, make_frequency_series, qtransform from pycbc.types import TimeSeries, FrequencySeries, zeros, float32, complex64 -from multiprocessing import Pool -import pycbc.version import pycbc.opt import pycbc.inject -import time last_progress_update = -1.0 @@ -48,9 +53,8 @@ tstart = time.time() parser = argparse.ArgumentParser(usage='', description="Find single detector gravitational-wave triggers.") +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action=pycbc.version.Version) -parser.add_argument("-V", "--verbose", action="store_true", - help="print extra debugging information", default=False ) parser.add_argument("--update-progress", help="updates a file 'progress.txt' with a value 0 .. 1.0 when this amount of (filtering) progress was made", type=float, default=0) diff --git a/bin/pycbc_inspiral_skymax b/bin/pycbc_inspiral_skymax index f082f7d3841..2f2c477f72d 100755 --- a/bin/pycbc_inspiral_skymax +++ b/bin/pycbc_inspiral_skymax @@ -16,7 +16,11 @@ # with this program; if not, write to the Free Software Foundation, Inc., # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. -import logging, argparse, sys, numpy +import logging +import argparse +import sys +import numpy + from pycbc import vetoes, psd, waveform, events, strain, scheme, fft from pycbc.vetoes.sgchisq import SingleDetSGChisq from pycbc import DYN_RANGE_FAC @@ -28,8 +32,7 @@ import pycbc.opt parser = argparse.ArgumentParser(usage='', description="Find single detector precessing gravitational-wave triggers.") -parser.add_argument("-V", "--verbose", action="store_true", - help="print extra debugging information", default=False ) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--output", type=str, help="FIXME: ADD") parser.add_argument("--bank-file", type=str, help="FIXME: ADD") parser.add_argument("--snr-threshold", diff --git a/bin/pycbc_make_banksim b/bin/pycbc_make_banksim index cfbc8985bca..09b0efbad5c 100644 --- a/bin/pycbc_make_banksim +++ b/bin/pycbc_make_banksim @@ -9,6 +9,7 @@ import tempfile from optparse import OptionParser from glue.pipeline import CondorDAGJob, CondorDAGNode, CondorDAG, CondorJob +from pycbc import init_logging class BaseJob(CondorDAGJob, CondorJob): def __init__(self, log_dir, executable, cp, section, gpu=False, @@ -131,7 +132,8 @@ parser.add_option('--config', type=str) if options.config is None: raise ValueError("Config file is required") -logging.basicConfig(format='%(asctime)s %(message)s', level=logging.INFO) +# logging INFO +init_logging(1) confs = ConfigParser.ConfigParser() confs.read(options.config) diff --git a/bin/pycbc_make_faithsim b/bin/pycbc_make_faithsim index 8c1f33bbc3e..c27934745c3 100644 --- a/bin/pycbc_make_faithsim +++ b/bin/pycbc_make_faithsim @@ -9,6 +9,7 @@ import tempfile from optparse import OptionParser from glue.pipeline import CondorDAGJob, CondorDAGNode, CondorDAG, CondorJob +from pycbc import init_logging class BaseJob(CondorDAGJob, CondorJob): def __init__(self, log_dir, executable, cp, section, accounting_group=None): @@ -65,7 +66,8 @@ parser.add_option('--config', type=str) if options.config is None: raise ValueError("Config file is required") -logging.basicConfig(format='%(asctime)s %(message)s', level=logging.INFO) +# logging INFO +init_logging(1) confs = ConfigParser.ConfigParser() confs.read(options.config) diff --git a/bin/pycbc_make_html_page b/bin/pycbc_make_html_page index 1ed5612dd48..f4ea0f2e809 100644 --- a/bin/pycbc_make_html_page +++ b/bin/pycbc_make_html_page @@ -17,14 +17,17 @@ # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. import argparse -import os, stat -import pycbc.results +import os +import stat import random import shutil import zipfile import codecs -from ligo import segments from jinja2 import Environment, FileSystemLoader + +from ligo import segments + +import pycbc.results from pycbc.results.render import get_embedded_config, render_workflow_html_template, setup_template_render from pycbc.workflow import segment import pycbc.version @@ -166,6 +169,7 @@ parser = argparse.ArgumentParser(usage='pycbc_make_html_page \ description="Create static html pages of a filesystem's content.") parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) +pycbc.add_common_pycbc_options(parser) parser.add_argument('-f', '--template-file', type=str, help='Template file to use for skeleton html page.') parser.add_argument('-b', '--output-path', type=str, @@ -181,12 +185,12 @@ parser.add_argument('-s', '--analysis-subtitle', type=str, parser.add_argument('-l', '--logo', type=str, default=default_logo_location, help='File location of logo to include at top of page') -parser.add_argument('-v', '--verbose', action='store_true', - help='Print extra debugging information.', default=False) parser.add_argument('-P', '--protect-results', action='store_true', help='Make the output web pages read only.', default=False) opts = parser.parse_args() +pycbc.init_logging(opts.verbose) + # edit command line options analysis_title = opts.analysis_title.strip('"').rstrip('"') analysis_subtitle = opts.analysis_subtitle.strip('"').rstrip('"') diff --git a/bin/pycbc_make_sky_grid b/bin/pycbc_make_sky_grid index 8a1809441a4..e6e2bf85ddd 100644 --- a/bin/pycbc_make_sky_grid +++ b/bin/pycbc_make_sky_grid @@ -10,10 +10,13 @@ The grid is constructed following the method described in Section V of https://a import numpy as np import h5py import argparse -from pycbc.detector import Detector import itertools from scipy.spatial.transform import Rotation as R +import pycbc +from pycbc.detector import Detector + + def spher_to_cart(sky_points): """Convert spherical coordinates to cartesian coordinates. """ @@ -32,16 +35,29 @@ def cart_to_spher(sky_points): return spher parser = argparse.ArgumentParser(description=__doc__) -parser.add_argument('--ra', type=float, help="Right ascension (in rad) of the center of the external trigger error box") -parser.add_argument('--dec', type=float, help="Declination (in rad) of the center of the external trigger error box") -parser.add_argument('--instruments', nargs="+", type=str, required=True, help="List of instruments to analyze.") -parser.add_argument('--sky-error', type=float, required=True, help="3-sigma confidence radius (in rad) of the external trigger error box") -parser.add_argument('--trigger-time', type=int, required=True, help="Time (in s) of the external trigger") -parser.add_argument('--timing-uncertainty', type=float, default=0.0001, help="Timing uncertainty (in s) we are willing to accept") -parser.add_argument('--output', type=str, required=True, help="Name of the sky grid") +pycbc.add_common_pycbc_options(parser) +parser.add_argument('--ra', type=float, + help="Right ascension (in rad) of the center of the external trigger " + "error box") +parser.add_argument('--dec', type=float, + help="Declination (in rad) of the center of the external trigger " + "error box") +parser.add_argument('--instruments', nargs="+", type=str, required=True, + help="List of instruments to analyze.") +parser.add_argument('--sky-error', type=float, required=True, + help="3-sigma confidence radius (in rad) of the external trigger error " + "box") +parser.add_argument('--trigger-time', type=int, required=True, + help="Time (in s) of the external trigger") +parser.add_argument('--timing-uncertainty', type=float, default=0.0001, + help="Timing uncertainty (in s) we are willing to accept") +parser.add_argument('--output', type=str, required=True, + help="Name of the sky grid") args = parser.parse_args() +pycbc.init_logging(args.verbose) + if len(args.instruments) == 1: parser.error('Can not make a sky grid for only one detector.') diff --git a/bin/pycbc_make_skymap b/bin/pycbc_make_skymap index 4b011831011..89e5933db88 100755 --- a/bin/pycbc_make_skymap +++ b/bin/pycbc_make_skymap @@ -19,6 +19,9 @@ mpl_use_backend('agg') import numpy as np import h5py + +from ligo.gracedb.rest import GraceDb + import pycbc from pycbc.filter import sigmasq from pycbc.io import live, WaveformArray @@ -29,7 +32,6 @@ from pycbc.pnutils import nearest_larger_binary_number from pycbc.waveform.spa_tmplt import spa_length_in_time from pycbc import frame, waveform, dq from pycbc.psd import interpolate -from ligo.gracedb.rest import GraceDb def default_frame_type(time, ifo): @@ -476,9 +478,9 @@ def main(trig_time, mass1, mass2, spin1z, spin2z, f_low, f_upper, sample_rate, shutil.rmtree(tmpdir) if __name__ == '__main__': - pycbc.init_logging(True) parser = argparse.ArgumentParser(description=__doc__) + pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action=pycbc.version.Version) # note that I am not using a MultiDetOptionAction for --trig-time as I # explicitly want to handle cases like `--trig-time 1234` and @@ -563,6 +565,8 @@ if __name__ == '__main__': opt = parser.parse_args() + pycbc.init_logging(opt.verbose) + frame_type_dict = {f.split(':')[0]: f.split(':')[1] for f in opt.frame_type} \ if opt.frame_type is not None else None chan_name_dict = {f.split(':')[0]: f for f in opt.channel_name} \ diff --git a/bin/pycbc_merge_inj_hdf b/bin/pycbc_merge_inj_hdf index fe0d698e935..9b208738d47 100755 --- a/bin/pycbc_merge_inj_hdf +++ b/bin/pycbc_merge_inj_hdf @@ -24,6 +24,7 @@ import logging import argparse import numpy as np import h5py + import pycbc import pycbc.inject import pycbc.version @@ -43,20 +44,16 @@ def get_gc_end_time(injection): if __name__ == '__main__': parser = argparse.ArgumentParser(description=__doc__) + pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action=pycbc.version.Version) parser.add_argument('--injection-files', '-i', dest='injection_file', required=True, nargs='+', help='Input HDF5 files defining injections') parser.add_argument('--output-file', '-o', dest='out_file', required=True, help='Output HDF5 file') - parser.add_argument("--verbose", action="store_true", default=False, - help="Print logging messages") opts = parser.parse_args() - log_level = logging.INFO if opts.verbose else logging.WARN - logging.basicConfig(level=log_level, - format='%(asctime)s %(message)s', - datefmt='%Y-%m-%d %H:%M:%S') + pycbc.init_logging(opts.verbose) logging.info("Loading injections") diff --git a/bin/pycbc_multi_inspiral b/bin/pycbc_multi_inspiral index ce58cf3cff4..3d9967f42d0 100755 --- a/bin/pycbc_multi_inspiral +++ b/bin/pycbc_multi_inspiral @@ -31,9 +31,10 @@ from collections import defaultdict import argparse import numpy as np import h5py + from pycbc import ( detector, fft, init_logging, inject, opt, psd, scheme, strain, vetoes, - waveform, DYN_RANGE_FAC + waveform, DYN_RANGE_FAC, add_common_pycbc_options ) from pycbc.events import ranking, coherent as coh, EventManagerCoherent from pycbc.filter import MatchedFilterControl @@ -46,9 +47,7 @@ from pycbc.vetoes import sgchisq # pycbc_multi_inspiral executable. time_init = time.time() parser = argparse.ArgumentParser(description=__doc__) -parser.add_argument("-V", "--verbose", action="store_true", - help="prints extra debugging information during runtime", - default=False) +add_common_pycbc_options(parser) parser.add_argument("--output", type=str) parser.add_argument("--instruments", nargs="+", type=str, required=True, help="List of instruments to analyze.") diff --git a/bin/pycbc_optimal_snr b/bin/pycbc_optimal_snr index 7e99bed39f4..5511c06880f 100644 --- a/bin/pycbc_optimal_snr +++ b/bin/pycbc_optimal_snr @@ -27,12 +27,14 @@ import argparse import multiprocessing import numpy as np import h5py + +from ligo.lw import utils as ligolw_utils +from ligo.lw import lsctables + import pycbc import pycbc.inject import pycbc.psd import pycbc.version -from ligo.lw import utils as ligolw_utils -from ligo.lw import lsctables from pycbc.filter import sigma, make_frequency_series from pycbc.types import TimeSeries, FrequencySeries, zeros, float32, \ MultiDetOptionAction, load_frequencyseries @@ -117,6 +119,7 @@ def get_gc_end_time(injection): if __name__ == '__main__': parser = argparse.ArgumentParser(description=__doc__) + pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action=pycbc.version.Version) parser.add_argument('--input-file', '-i', dest='injection_file', required=True, @@ -159,11 +162,6 @@ if __name__ == '__main__': parser.add_argument('--ignore-waveform-errors', action='store_true', help='Ignore errors in waveform generation and keep ' 'the corresponding column unchanged') - parser.add_argument('--verbose', action='store_true', - help='Print which injection is going to be attempted ' - 'next and when it completes. Use this to debug ' - 'misbehaving approximants which crash or quit ' - 'the Python interpreter') parser.add_argument('--progress', action='store_true', help='Show a progress bar (requires tqdm)') psd_group = pycbc.psd.insert_psd_option_group_multi_ifo(parser) @@ -187,10 +185,7 @@ if __name__ == '__main__': if not opts.time_varying_psds: pycbc.psd.verify_psd_options_multi_ifo(opts, parser, detectors) - log_level = logging.INFO if opts.verbose else logging.WARN - logging.basicConfig(level=log_level, - format='%(asctime)s %(message)s', - datefmt='%Y-%m-%d %H:%M:%S') + pycbc.init_logging(opts.verbose) seg_len = opts.seg_length sample_rate = opts.sample_rate @@ -240,20 +235,27 @@ if __name__ == '__main__': psd = psds[det](injection_time) if psd is None: continue - logging.info('Trying injection %s at %s', inj.simulation_id, det) + logging.debug('Trying injection %s at %s', inj.simulation_id, det) try: wave = get_injection(injections, det, injection_time, simulation_id=inj.simulation_id) except Exception as e: if opts.ignore_waveform_errors: - logging.warn('%s: waveform generation failed, skipping (%s)', - inj.simulation_id, e) + logging.debug( + '%s: waveform generation failed, skipping (%s)', + inj.simulation_id, + e + ) continue else: logging.error('%s: waveform generation failed with the ' 'following exception', inj.simulation_id) raise - logging.info('Injection %s at %s completed', inj.simulation_id, det) + logging.debug( + 'Injection %s at %s completed', + inj.simulation_id, + det + ) sval = sigma(wave, psd=psd, low_frequency_cutoff=f_low) if ligolw: setattr(inj, column, sval) diff --git a/bin/pycbc_optimize_snr b/bin/pycbc_optimize_snr index 513f1d2fac5..fb36ba6c721 100755 --- a/bin/pycbc_optimize_snr +++ b/bin/pycbc_optimize_snr @@ -5,14 +5,14 @@ import os import argparse import logging +import numpy +import h5py # we will make plots on a likely headless machine, so make sure matplotlib's # backend is set appropriately from matplotlib import use as mpl_use_backend mpl_use_backend('agg') -import numpy -import h5py import pycbc from pycbc import ( fft, scheme, version diff --git a/bin/pycbc_process_sngls b/bin/pycbc_process_sngls index c740415ecdd..3e5909e449d 100644 --- a/bin/pycbc_process_sngls +++ b/bin/pycbc_process_sngls @@ -6,12 +6,14 @@ import argparse import logging import numpy import h5py + +import pycbc from pycbc.io import SingleDetTriggers from pycbc.events import stat, coinc, veto parser = argparse.ArgumentParser(description=__doc__) - +pycbc.add_common_pycbc_options(parser) parser.add_argument('--single-trig-file', required=True, help='Path to file containing single-detector triggers in ' 'HDF5 format. Required') @@ -41,7 +43,7 @@ stat.insert_statistic_option_group(parser) args = parser.parse_args() -logging.basicConfig(format='%(asctime)s %(message)s', level=logging.INFO) +pycbc.init_logging(args.verbose) # Munge together SNR cut and any other filter specified snr_filter = 'self.snr>%f' % (args.min_snr) if args.min_snr > 0. else None diff --git a/bin/pycbc_randomize_inj_dist_by_optsnr b/bin/pycbc_randomize_inj_dist_by_optsnr index 1bc4f507c27..c3b3a334209 100644 --- a/bin/pycbc_randomize_inj_dist_by_optsnr +++ b/bin/pycbc_randomize_inj_dist_by_optsnr @@ -14,6 +14,8 @@ from argparse import ArgumentParser from ligo.lw import lsctables from ligo.lw import utils as ligolw_utils from ligo.lw import table + +import pycbc from pycbc.io.ligolw import LIGOLWContentHandler @@ -50,6 +52,7 @@ def get_weights(distribution, r, r1, r2): return min_vol, weight parser = ArgumentParser(description = __doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--xml-file', required = True, help = 'Input LIGOLW file defining injections.') parser.add_argument('--output-file', required=True, @@ -90,10 +93,11 @@ parser.add_argument('--scale-by-chirp-dist', action='store_true', parser.add_argument('--seed', type = int, default = int((time.time()*100)% 1e6), help = 'Set the seed to use for the random number generator. ' \ 'If none specified, will use the current time.') -parser.add_argument('--verbose', action = 'store_true', help = 'Be verbose.') opts = parser.parse_args() +pycbc.init_logging(opts.verbose) + # parse the distributions # beta is the only special one if opts.distribution in ['volume', 'distance', 'logdist']: diff --git a/bin/pycbc_single_template b/bin/pycbc_single_template index 062f2c11815..ee78facdcc1 100755 --- a/bin/pycbc_single_template +++ b/bin/pycbc_single_template @@ -17,7 +17,13 @@ """ Calculate the SNR and CHISQ timeseries for either a chosen template, or a specific Nth loudest coincident event. """ -import sys, logging, argparse, numpy, pycbc, h5py +import sys +import logging +import argparse +import numpy +import h5py + +import pycbc from pycbc import vetoes, psd, waveform, strain, scheme, fft, filter from pycbc.io import WaveformArray from pycbc import events @@ -86,11 +92,10 @@ def select_segments(fname, anal_name, data_name, ifo, time, pad_data): parser = argparse.ArgumentParser(usage='', description="Single template gravitational-wave followup") +pycbc.add_common_pycbc_options(parser) parser.add_argument('--version', action=pycbc.version.Version) parser.add_argument('--output-file', required=True) parser.add_argument('--subtract-template', action='store_true') -parser.add_argument("-V", "--verbose", action="store_true", - help="print extra debugging information", default=False ) parser.add_argument("--low-frequency-cutoff", type=float, help="The low frequency cutoff to use for filtering (Hz)") parser.add_argument("--high-frequency-cutoff", type=float, @@ -182,6 +187,8 @@ scheme.insert_processing_option_group(parser) fft.insert_fft_option_group(parser) opt = parser.parse_args() +pycbc.init_logging(opt.verbose) + # Exclude invalid options if opt.window and opt.trigger_time is None: raise RuntimeError("Can't use --window option without a valid trigger time!") @@ -210,7 +217,6 @@ strain.verify_strain_options(opt, parser) strain.StrainSegments.verify_segment_options(opt, parser) scheme.verify_processing_options(opt, parser) fft.verify_fft_options(opt,parser) -pycbc.init_logging(opt.verbose) ctx = scheme.from_cli(opt) gwstrain = strain.from_cli(opt, pycbc.DYN_RANGE_FAC) diff --git a/bin/pycbc_source_probability_offline b/bin/pycbc_source_probability_offline index 2fc8eaceb84..f12447603ef 100755 --- a/bin/pycbc_source_probability_offline +++ b/bin/pycbc_source_probability_offline @@ -9,12 +9,14 @@ import tqdm import argparse import logging import numpy as np + import pycbc from pycbc.io import hdf from pycbc.pnutils import mass1_mass2_to_mchirp_eta from pycbc import mchirp_area parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument('--trigger-file', required=True) parser.add_argument('--bank-file', required=True) parser.add_argument('--single-detector-triggers', nargs='+', required=True) @@ -27,7 +29,6 @@ parser.add_argument('--ifar-threshold', type=float, default=None, 'above threshold.') parser.add_argument('--include-mass-gap', action='store_true', help='Option to include the Mass Gap region.') -parser.add_argument('--verbose', action='count') parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) mchirp_area.insert_args(parser) diff --git a/bin/pycbc_split_inspinj b/bin/pycbc_split_inspinj index ce38299eda1..28e5ee9c1ad 100644 --- a/bin/pycbc_split_inspinj +++ b/bin/pycbc_split_inspinj @@ -5,12 +5,14 @@ import argparse from ligo.lw import utils as ligolw_utils from ligo.lw import lsctables from itertools import cycle + import pycbc.version from pycbc.io.ligolw import LIGOLWContentHandler, get_table_columns # Parse command line parser = argparse.ArgumentParser() +pycbc.add_common_pycbc_options(parser) parser.add_argument("--version", action="version", version=pycbc.version.git_verbose_msg) group = parser.add_mutually_exclusive_group(required=True) @@ -30,7 +32,10 @@ if args.output_files and args.output_dir: # Read in input file xmldoc = ligolw_utils.load_filename( - args.input_file, verbose=True, contenthandler=LIGOLWContentHandler) + args.input_file, + verbose=args.verbose, + contenthandler=LIGOLWContentHandler +) tabletype = lsctables.SimInspiralTable allinjs = tabletype.get_table(xmldoc) diff --git a/bin/pycbc_splitbank b/bin/pycbc_splitbank index 42d2e6ba437..8919c31ac2f 100644 --- a/bin/pycbc_splitbank +++ b/bin/pycbc_splitbank @@ -29,9 +29,12 @@ import argparse from numpy import random, ceil + from ligo.lw import ligolw from ligo.lw import lsctables from ligo.lw import utils as ligolw_utils + +import pycbc from pycbc import version from pycbc.io.ligolw import LIGOLWContentHandler, create_process_table from pycbc.conversions import mchirp_from_mass1_mass2 @@ -46,6 +49,7 @@ __program__ = "pycbc_splitbank" parser = argparse.ArgumentParser(description=__doc__) parser.add_argument('--version', action='version', version=version.git_verbose_msg) +pycbc.add_common_pycbc_options(parser) group = parser.add_mutually_exclusive_group(required=True) group.add_argument('--templates-per-bank', metavar='SAMPLES', help='number of templates in the output banks', type=int) @@ -61,8 +65,6 @@ group.add_argument("-O", "--output-filenames", nargs='*', default=None, parser.add_argument("-o", "--output-prefix", default=None, help="Prefix to add to the template bank name (name becomes output#.xml[.gz])" ) -parser.add_argument("-V", "--verbose", action="store_true", - help="Print extra debugging information", default=False ) parser.add_argument("-t", "--bank-file", metavar='INPUT_FILE', help='Template bank to split', required=True) parser.add_argument("--sort-frequency-cutoff", diff --git a/bin/pycbc_upload_xml_to_gracedb b/bin/pycbc_upload_xml_to_gracedb index a922d0923cb..df9fc57848a 100755 --- a/bin/pycbc_upload_xml_to_gracedb +++ b/bin/pycbc_upload_xml_to_gracedb @@ -23,24 +23,26 @@ Take a coinc xml file containing multiple events and upload to gracedb. import os import argparse import logging -from ligo.gracedb.rest import GraceDb import numpy as np import h5py import matplotlib matplotlib.use('agg') -import pycbc +from ligo.gracedb.rest import GraceDb import lal import lal.series from ligo.lw import lsctables from ligo.lw import utils as ligolw_utils from ligo.segments import segment, segmentlist + +import pycbc from pycbc.io.live import gracedb_tag_with_version from pycbc.io.ligolw import LIGOLWContentHandler from pycbc.psd import interpolate from pycbc.types import FrequencySeries from pycbc.results import generate_asd_plot + def check_gracedb_for_event(gdb_handle, query, far): """ Check if there is an event in gracedb with the queried string @@ -68,6 +70,7 @@ logging.basicConfig(format='%(asctime)s:%(levelname)s : %(message)s', level=logging.INFO) parser = argparse.ArgumentParser(description=__doc__) +pycbc.add_common_pycbc_options(parser) parser.add_argument("--psd-files", nargs='+', required=True, help='HDF file(s) containing the PSDs to upload') parser.add_argument('--input-file', required=True, type=str, @@ -109,6 +112,8 @@ parser.add_argument('--generate-plots', action='store_true', args = parser.parse_args() +pycbc.init_logging(args.verbose) + if args.production_server: gracedb = GraceDb() else: diff --git a/pycbc/__init__.py b/pycbc/__init__.py index 97acf32e7a9..e716caf4a82 100644 --- a/pycbc/__init__.py +++ b/pycbc/__init__.py @@ -79,11 +79,11 @@ def add_common_pycbc_options(parser): title="PyCBC common options", description="Common options for PyCBC executables.", ) - group.add_argument('-V', '--verbose', action='count', default=0, + group.add_argument('-v', '--verbose', action='count', default=0, help='Add verbosity to logging. Adding the option ' 'multiple times makes logging progressively ' - 'more verbose, e.g. --verbose or -V provides ' - 'logging at the info level, but -VV or ' + 'more verbose, e.g. --verbose or -v provides ' + 'logging at the info level, but -vv or ' '--verbose --verbose provides debug logging.') def init_logging(verbose=False, diff --git a/pycbc/tmpltbank/bank_conversions.py b/pycbc/tmpltbank/bank_conversions.py index 9859d4e7a0a..f987d1df0eb 100644 --- a/pycbc/tmpltbank/bank_conversions.py +++ b/pycbc/tmpltbank/bank_conversions.py @@ -26,10 +26,14 @@ specific values from PyCBC template banks. """ +import logging import numpy as np + from pycbc import conversions as conv from pycbc import pnutils +logger = logging.getLogger('pycbc.tmpltbank.bank_conversions') + # Convert from parameter name to helper function # some multiple names are used for the same function conversion_options = ['mass1', 'mass2', 'spin1z', 'spin2z', 'duration', diff --git a/pycbc/tmpltbank/bank_output_utils.py b/pycbc/tmpltbank/bank_output_utils.py index 2fae646fb1e..4d5e281ad16 100644 --- a/pycbc/tmpltbank/bank_output_utils.py +++ b/pycbc/tmpltbank/bank_output_utils.py @@ -1,7 +1,10 @@ +import logging import numpy import h5py + from lal import PI, MTSUN_SI, TWOPI, GAMMA from ligo.lw import ligolw, lsctables, utils as ligolw_utils + from pycbc import pnutils from pycbc.tmpltbank.lambda_mapping import ethinca_order_from_string from pycbc.io.ligolw import ( @@ -10,6 +13,8 @@ from pycbc.waveform import get_waveform_filter_length_in_time as gwflit +logger = logging.getLogger('pycbc.tmpltbank.bank_output_utils') + def convert_to_sngl_inspiral_table(params, proc_id): ''' Convert a list of m1,m2,spin1z,spin2z values into a basic sngl_inspiral diff --git a/pycbc/tmpltbank/brute_force_methods.py b/pycbc/tmpltbank/brute_force_methods.py index c5ff02858c0..34ca473c5df 100644 --- a/pycbc/tmpltbank/brute_force_methods.py +++ b/pycbc/tmpltbank/brute_force_methods.py @@ -1,6 +1,10 @@ +import logging import numpy + from pycbc.tmpltbank.coord_utils import get_cov_params +logger = logging.getLogger('pycbc.tmpltbank.brute_force_methods') + def get_physical_covaried_masses(xis, bestMasses, bestXis, req_match, massRangeParams, metricParams, fUpper, diff --git a/pycbc/tmpltbank/calc_moments.py b/pycbc/tmpltbank/calc_moments.py index 305fb8834e6..c472dfacd91 100644 --- a/pycbc/tmpltbank/calc_moments.py +++ b/pycbc/tmpltbank/calc_moments.py @@ -14,9 +14,13 @@ # with this program; if not, write to the Free Software Foundation, Inc., # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. +import logging import numpy + from pycbc.tmpltbank.lambda_mapping import generate_mapping +logger = logging.getLogger('pycbc.tmpltbank.calc_moments') + def determine_eigen_directions(metricParams, preserveMoments=False, vary_fmax=False, vary_density=None): diff --git a/pycbc/tmpltbank/coord_utils.py b/pycbc/tmpltbank/coord_utils.py index c1bb38a7f62..0d0e4db509d 100644 --- a/pycbc/tmpltbank/coord_utils.py +++ b/pycbc/tmpltbank/coord_utils.py @@ -16,11 +16,15 @@ import logging import numpy + from pycbc.tmpltbank.lambda_mapping import get_chirp_params from pycbc import conversions from pycbc import pnutils from pycbc.neutron_stars import load_ns_sequence +logger = logging.getLogger('pycbc.tmpltbank.coord_utils') + + def estimate_mass_range(numPoints, massRangeParams, metricParams, fUpper,\ covary=True): """ @@ -235,7 +239,7 @@ def get_random_mass(numPoints, massRangeParams, eos='2H'): warn_msg += "(%s). " %(max_ns_g_mass) warn_msg += "The code will proceed using the latter value " warn_msg += "as the boundary mass." - logging.warn(warn_msg) + logger.warn(warn_msg) boundary_mass = max_ns_g_mass # Empty arrays to store points that pass all cuts @@ -711,7 +715,7 @@ def find_max_and_min_frequencies(name, mass_range_params, freqs): warn_msg += "for the metric: %s Hz. " %(freqs.min()) warn_msg += "Distances for these waveforms will be calculated at " warn_msg += "the lowest available metric frequency." - logging.warn(warn_msg) + logger.warn(warn_msg) if upper_f_cutoff > freqs.max(): warn_msg = "WARNING: " warn_msg += "Highest frequency cutoff is %s Hz " %(upper_f_cutoff,) @@ -719,7 +723,7 @@ def find_max_and_min_frequencies(name, mass_range_params, freqs): warn_msg += "for the metric: %s Hz. " %(freqs.max()) warn_msg += "Distances for these waveforms will be calculated at " warn_msg += "the largest available metric frequency." - logging.warn(warn_msg) + logger.warn(warn_msg) return find_closest_calculated_frequencies(cutoffs, freqs) diff --git a/pycbc/tmpltbank/lambda_mapping.py b/pycbc/tmpltbank/lambda_mapping.py index ed0121ed912..3e3516dbc96 100644 --- a/pycbc/tmpltbank/lambda_mapping.py +++ b/pycbc/tmpltbank/lambda_mapping.py @@ -14,12 +14,17 @@ # with this program; if not, write to the Free Software Foundation, Inc., # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. import re +import logging import numpy -import pycbc.libutils + from lal import MTSUN_SI, PI, CreateREAL8Vector +import pycbc.libutils + lalsimulation = pycbc.libutils.import_optional('lalsimulation') +logger = logging.getLogger('pycbc.tmpltbank.lambda_mapping') + # PLEASE ENSURE THESE ARE KEPT UP TO DATE WITH THE REST OF THIS FILE pycbcValidTmpltbankOrders = ['zeroPN','onePN','onePointFivePN','twoPN',\ 'twoPointFivePN','threePN','threePointFivePN'] diff --git a/pycbc/tmpltbank/lattice_utils.py b/pycbc/tmpltbank/lattice_utils.py index c938a8b3612..d6db4f98c86 100644 --- a/pycbc/tmpltbank/lattice_utils.py +++ b/pycbc/tmpltbank/lattice_utils.py @@ -14,10 +14,14 @@ # with this program; if not, write to the Free Software Foundation, Inc., # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. +import logging import copy import numpy + import lal +logger = logging.getLogger('pycbc.tmpltbank.lattice_utils') + def generate_hexagonal_lattice(maxv1, minv1, maxv2, minv2, mindist): """ This function generates a 2-dimensional lattice of points using a hexagonal diff --git a/pycbc/tmpltbank/option_utils.py b/pycbc/tmpltbank/option_utils.py index 722c250b205..e9451a7e716 100644 --- a/pycbc/tmpltbank/option_utils.py +++ b/pycbc/tmpltbank/option_utils.py @@ -19,11 +19,15 @@ import textwrap import numpy import os + from pycbc.tmpltbank.lambda_mapping import get_ethinca_orders, pycbcValidOrdersHelpDescriptions from pycbc import pnutils from pycbc.neutron_stars import load_ns_sequence from pycbc.types import positive_float, nonnegative_float +logger = logging.getLogger('pycbc.tmpltbank.option_utils') + + class IndentedHelpFormatterWithNL(argparse.ArgumentDefaultsHelpFormatter): """ This class taken from @@ -589,21 +593,21 @@ def verify_mass_range_options(opts, parser, nonSpin=False): # inform him/her about the caveats involved in this. if hasattr(opts, 'remnant_mass_threshold') \ and opts.remnant_mass_threshold is not None: - logging.info("""You have asked to exclude EM dim NS-BH systems from the - target parameter space. The script will assume that m1 is - the BH and m2 is the NS: make sure that your settings - respect this convention. The script will also treat the - NS as non-spinning: use NS spins in the template bank - at your own risk!""") + logger.info("""You have asked to exclude EM dim NS-BH systems from the + target parameter space. The script will assume that m1 + is the BH and m2 is the NS: make sure that your settings + respect this convention. The script will also treat the + NS as non-spinning: use NS spins in the template bank + at your own risk!""") if opts.use_eos_max_ns_mass: - logging.info("""You have asked to take into account the maximum NS + logger.info("""You have asked to take into account the maximum NS mass value for the EOS in use.""") # Find out if the EM constraint surface data already exists or not # and inform user whether this will be read from file or generated. # This is the minumum eta as a function of BH spin and NS mass # required to produce an EM counterpart if os.path.isfile('constraint_em_bright.npz'): - logging.info("""The constraint surface for EM bright binaries + logger.info("""The constraint surface for EM bright binaries will be read in from constraint_em_bright.npz.""") # Assign min/max total mass from mass1, mass2 if not specified diff --git a/pycbc/tmpltbank/partitioned_bank.py b/pycbc/tmpltbank/partitioned_bank.py index 29505b2a03e..9b22658aaa8 100644 --- a/pycbc/tmpltbank/partitioned_bank.py +++ b/pycbc/tmpltbank/partitioned_bank.py @@ -17,8 +17,11 @@ import copy import numpy import logging + from pycbc.tmpltbank import coord_utils +logger = logging.getLogger('pycbc.tmpltbank.partitioned_bank') + class PartitionedTmpltbank(object): """ This class is used to hold a template bank partitioned into numerous bins @@ -524,7 +527,7 @@ def add_point_by_masses(self, mass1, mass2, spin1z, spin2z, warn_msg += "convention is that mass1 > mass2. Swapping mass1 " warn_msg += "and mass2 and adding point to bank. This message " warn_msg += "will not be repeated." - logging.warn(warn_msg) + logger.warn(warn_msg) self.spin_warning_given = True # These that masses obey the restrictions of mass_range_params From 7f9cb87e0606ca8f8eba17094765f105078d4a9a Mon Sep 17 00:00:00 2001 From: Gareth S Cabourn Davies Date: Tue, 27 Feb 2024 15:43:17 +0000 Subject: [PATCH 24/29] Add a template for pull requests (#4585) * Create standard_template.md * Update .github/PULL_REQUEST_TEMPLATE/standard_template.md * Update .github/PULL_REQUEST_TEMPLATE/standard_template.md * Update standard_template.md * Update standard_template.md Use some 'delete as appropriate' lists to make using the template a little less daunting --- .../standard_template.md | 52 +++++++++++++++++++ 1 file changed, 52 insertions(+) create mode 100644 .github/PULL_REQUEST_TEMPLATE/standard_template.md diff --git a/.github/PULL_REQUEST_TEMPLATE/standard_template.md b/.github/PULL_REQUEST_TEMPLATE/standard_template.md new file mode 100644 index 00000000000..79a06d213f5 --- /dev/null +++ b/.github/PULL_REQUEST_TEMPLATE/standard_template.md @@ -0,0 +1,52 @@ + + + + +- [ ] The author of this pull request confirms they will adhere to the [code of conduct](https://github.com/gwastro/pycbc/blob/master/CODE_OF_CONDUCT.md) + +## Standard information about the request + + +This is a: bug fix, new feature, efficiency update, other (please describe) + + +This change affects: the offline search, the live search, inference, PyGRB + + +This change changes: documentation, result presentation / plotting, scientific output + + +This change: has appropriate unit tests, follows style guidelines (See e.g. [PEP8](https://peps.python.org/pep-0008/)), has been proposed using the [contribution guidelines](https://github.com/gwastro/pycbc/blob/master/CONTRIBUTING.md) + + +This change will: break current functionality, require additional dependencies, require a new release, other (please describe) + +## Motivation + + +## Contents + + +## Links to any issues or associated PRs + + +## Testing performed + + +## Additional notes + From a3f5e92a42c3d7d7ece0ef4cf10b84b1e7b5b72f Mon Sep 17 00:00:00 2001 From: kkacanja <123669569+kkacanja@users.noreply.github.com> Date: Tue, 27 Feb 2024 10:44:16 -0500 Subject: [PATCH 25/29] Allow option for workflow to end at inspiral jobs (#4612) * Updates to stopping after data inspiral jobs and statmap * Updates to stopping after data inspiral jobs and statmap * Changed typing of stop-after * Allow pycbc_make_offline_search_workflow to stop after inspiral jobs. * Fixed comment --- .../pycbc_make_offline_search_workflow | 164 ++++++++++-------- 1 file changed, 95 insertions(+), 69 deletions(-) diff --git a/bin/workflows/pycbc_make_offline_search_workflow b/bin/workflows/pycbc_make_offline_search_workflow index ac96e5a9815..25fc1e627ba 100755 --- a/bin/workflows/pycbc_make_offline_search_workflow +++ b/bin/workflows/pycbc_make_offline_search_workflow @@ -57,6 +57,92 @@ def ifo_combos(ifos): for ifocomb in combinations: yield ifocomb +def finalize(container, workflow, finalize_workflow): + # Create the final log file + log_file_html = wf.File(workflow.ifos, 'WORKFLOW-LOG', workflow.analysis_time, + extension='.html', directory=rdir['workflow']) + + gen_file_html = wf.File(workflow.ifos, 'WORKFLOW-GEN', workflow.analysis_time, + extension='.html', directory=rdir['workflow']) + + # Create a page to contain a dashboard link + dashboard_file = wf.File(workflow.ifos, 'DASHBOARD', workflow.analysis_time, + extension='.html', directory=rdir['workflow']) + dashboard_str = """

Pegasus Dashboard Page

""" + kwds = {'title': "Pegasus Dashboard", + 'caption': "Link to Pegasus Dashboard", + 'cmd': "PYCBC_SUBMIT_DAX_ARGV", } + save_fig_with_metadata(dashboard_str, dashboard_file.storage_path, **kwds) + + # Create pages for the submission script to write data + wf.makedir(rdir['workflow/dax']) + wf.makedir(rdir['workflow/input_map']) + wf.makedir(rdir['workflow/output_map']) + wf.makedir(rdir['workflow/planning']) + + wf.make_results_web_page(finalize_workflow, os.path.join(os.getcwd(), + rdir.base)) + + container += workflow + container += finalize_workflow + + container.add_subworkflow_dependancy(workflow, finalize_workflow) + + container.save() + + open_box_result_path = rdir['open_box_result'] + if os.path.exists(open_box_result_path): + os.chmod(open_box_result_path, 0o0700) + else: + pass + + logging.info("Written dax.") + + # Close the log and flush to the html file + logging.shutdown() + with open(wf_log_file.storage_path, "r") as logfile: + logdata = logfile.read() + log_str = """ +

Workflow generation script created workflow in output directory: %s

+

Workflow name is: %s

+

Workflow generation script run on host: %s

+ 
%s
+ """ % (os.getcwd(), args.workflow_name, socket.gethostname(), logdata) + kwds = {'title': 'Workflow Generation Log', + 'caption': "Log of the workflow script %s" % sys.argv[0], + 'cmd': ' '.join(sys.argv), } + save_fig_with_metadata(log_str, log_file_html.storage_path, **kwds) + + # Add the command line used to a specific file + args_to_output = [sys.argv[0]] + for arg in sys.argv[1:]: + if arg.startswith('--'): + # This is an option, add tab + args_to_output.append(' ' + arg) + else: + # This is a parameter, add two tabs + args_to_output.append(' ' + arg) + + gen_str = '
' + ' \\\n'.join(args_to_output) + '
' + kwds = {'title': 'Workflow Generation Command', + 'caption': "Command used to generate the workflow.", + 'cmd': ' '.join(sys.argv), } + save_fig_with_metadata(gen_str, gen_file_html.storage_path, **kwds) + layout.single_layout(rdir['workflow'], ([dashboard_file, gen_file_html, log_file_html])) + sys.exit(0) + +def check_stop(job_name, container, workflow, finalize_workflow): + #This function will finalize the workflow and stop it at the job specified. + #Under [workflow] supply the option stop-after = and the job name you want the workflow to stop at. Current options are [inspiral, hdf_trigger_merge, statmap] + if workflow.cp.has_option('workflow', 'stop-after'): + stop_after = workflow.cp.get('workflow', 'stop-after') + if stop_after == job_name: + logging.info("Search worflow will be stopped after " + str(job_name)) + finalize(container, workflow, finalize_workflow) + else: + pass + + parser = argparse.ArgumentParser(description=__doc__[1:]) parser.add_argument('--version', action='version', version=__version__) parser.add_argument('--verbose', action='count', @@ -183,10 +269,15 @@ ind_insps = insps = wf.setup_matchedfltr_workflow(workflow, analyzable_segs, datafind_files, splitbank_files_fd, output_dir, tags=['full_data']) +#check to see if workflow should stop at inspiral jobs. +check_stop('inspiral', container, workflow, finalize_workflow) + insps = wf.merge_single_detector_hdf_files(workflow, hdfbank, insps, output_dir, tags=['full_data']) +check_stop('hdf_trigger_merge', container, workflow, finalize_workflow) + # setup sngl trigger distribution fitting jobs # 'statfiles' is list of files used in calculating statistic # 'dqfiles' is the subset of files containing data quality information @@ -264,6 +355,9 @@ for ifo in ifo_ids.keys(): workflow, hdfbank, inspcomb, statfiles, final_veto_file, final_veto_name, output_dir, tags=ctagsngl) +#check to see if workflow should stop at statmap jobs. +check_stop('statmap', container, workflow, finalize_workflow) + ifo_sets = list(ifo_combos(ifo_ids.keys())) if analyze_singles: ifo_sets += [(ifo,) for ifo in ifo_ids.keys()] @@ -836,72 +930,4 @@ wf.make_versioning_page( ) ############################ Finalization #################################### - -# Create the final log file -log_file_html = wf.File(workflow.ifos, 'WORKFLOW-LOG', workflow.analysis_time, - extension='.html', directory=rdir['workflow']) - -gen_file_html = wf.File(workflow.ifos, 'WORKFLOW-GEN', workflow.analysis_time, - extension='.html', directory=rdir['workflow']) - -# Create a page to contain a dashboard link -dashboard_file = wf.File(workflow.ifos, 'DASHBOARD', workflow.analysis_time, - extension='.html', directory=rdir['workflow']) -dashboard_str = """

Pegasus Dashboard Page

""" -kwds = { 'title' : "Pegasus Dashboard", - 'caption' : "Link to Pegasus Dashboard", - 'cmd' : "PYCBC_SUBMIT_DAX_ARGV", } -save_fig_with_metadata(dashboard_str, dashboard_file.storage_path, **kwds) - -# Create pages for the submission script to write data -wf.makedir(rdir['workflow/dax']) -wf.makedir(rdir['workflow/input_map']) -wf.makedir(rdir['workflow/output_map']) -wf.makedir(rdir['workflow/planning']) - -wf.make_results_web_page(finalize_workflow, os.path.join(os.getcwd(), - rdir.base)) - -container += workflow -container += finalize_workflow - -container.add_subworkflow_dependancy(workflow, finalize_workflow) - -container.save() - -# Protect the open box results folder -os.chmod(rdir['open_box_result'], 0o0700) - -logging.info("Written dax.") - -# Close the log and flush to the html file -logging.shutdown() -with open (wf_log_file.storage_path, "r") as logfile: - logdata=logfile.read() -log_str = """ -

Workflow generation script created workflow in output directory: %s

-

Workflow name is: %s

-

Workflow generation script run on host: %s

-
%s
-""" % (os.getcwd(), args.workflow_name, socket.gethostname(), logdata) -kwds = { 'title' : 'Workflow Generation Log', - 'caption' : "Log of the workflow script %s" % sys.argv[0], - 'cmd' :' '.join(sys.argv), } -save_fig_with_metadata(log_str, log_file_html.storage_path, **kwds) - -# Add the command line used to a specific file -args_to_output = [sys.argv[0]] -for arg in sys.argv[1:]: - if arg.startswith('--'): - # This is an option, add tab - args_to_output.append(' ' + arg) - else: - # This is a parameter, add two tabs - args_to_output.append(' ' + arg) - -gen_str = '
' + ' \\\n'.join(args_to_output) + '
' -kwds = { 'title' : 'Workflow Generation Command', - 'caption' : "Command used to generate the workflow.", - 'cmd' :' '.join(sys.argv), } -save_fig_with_metadata(gen_str, gen_file_html.storage_path, **kwds) -layout.single_layout(rdir['workflow'], ([dashboard_file, gen_file_html, log_file_html])) +finalize(container, workflow, finalize_workflow) From fd2681702e452739a3e08c9755a0f7f5c3a724b7 Mon Sep 17 00:00:00 2001 From: Ian Harry Date: Tue, 27 Feb 2024 22:30:31 +0000 Subject: [PATCH 26/29] Need to use atol for isclose here (#4655) --- pycbc/types/timeseries.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pycbc/types/timeseries.py b/pycbc/types/timeseries.py index 3bf94065473..70d8802cb23 100644 --- a/pycbc/types/timeseries.py +++ b/pycbc/types/timeseries.py @@ -188,10 +188,10 @@ def time_slice(self, start, end, mode='floor'): start_idx = float(start - self.start_time) * self.sample_rate end_idx = float(end - self.start_time) * self.sample_rate - if _numpy.isclose(start_idx, round(start_idx)): + if _numpy.isclose(start_idx, round(start_idx), rtol=0, atol=1E-3): start_idx = round(start_idx) - if _numpy.isclose(end_idx, round(end_idx)): + if _numpy.isclose(end_idx, round(end_idx), rtol=0, atol=1E-3): end_idx = round(end_idx) if mode == 'floor': From 56b544fdb462899a569a96ac4018761652987a1d Mon Sep 17 00:00:00 2001 From: Shichao Wu Date: Wed, 28 Feb 2024 22:44:24 +0100 Subject: [PATCH 27/29] add support for multi-returning-value functions in transform (#4301) * add support for multi-value functions in transform * fix cc issue * Update transforms.py * fix cc issue * fix * Update transforms.py * Update transforms.py * fix cc issue * Update transforms.py * add 7 classes * fix cc issues * move LISA stuff to another PR * Update transforms.py --- pycbc/transforms.py | 93 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 93 insertions(+) diff --git a/pycbc/transforms.py b/pycbc/transforms.py index c8fa44a9a9b..487de5ab103 100644 --- a/pycbc/transforms.py +++ b/pycbc/transforms.py @@ -318,6 +318,98 @@ def from_config(cls, cp, section, outputs): return cls(inputs, outputs, transform_functions, jacobian=jacobian) +class CustomTransformMultiOutputs(CustomTransform): + """Allows for any transform to be defined. Based on CustomTransform, + but also supports multi-returning value functions. + + Parameters + ---------- + input_args : (list of) str + The names of the input parameters. + output_args : (list of) str + The names of the output parameters. + transform_functions : dict + Dictionary mapping input args to a string giving a function call; + e.g., ``{'q': 'q_from_mass1_mass2(mass1, mass2)'}``. + jacobian : str, optional + String giving a jacobian function. The function must be in terms of + the input arguments. + """ + + name = "custom_multi" + + def __init__(self, input_args, output_args, transform_functions, + jacobian=None): + super(CustomTransformMultiOutputs, self).__init__( + input_args, output_args, transform_functions, jacobian) + + def transform(self, maps): + """Applies the transform functions to the given maps object. + Parameters + ---------- + maps : dict, or FieldArray + Returns + ------- + dict or FieldArray + A map object containing the transformed variables, along with the + original variables. The type of the output will be the same as the + input. + """ + if self.transform_functions is None: + raise NotImplementedError("no transform function(s) provided") + # copy values to scratch + self._copytoscratch(maps) + # ensure that we return the same data type in each dict + getslice = self._getslice(maps) + # evaluate the functions + # func[0] is the function itself, func[1] is the index, + # this supports multiple returning values function + out = { + p: self._scratch[func[0]][func[1]][getslice] if + len(self._scratch[func[0]]) > 1 else + self._scratch[func[0]][getslice] + for p, func in self.transform_functions.items() + } + return self.format_output(maps, out) + + @classmethod + def from_config(cls, cp, section, outputs): + """Loads a CustomTransformMultiOutputs from the given config file. + + Example section: + + .. code-block:: ini + + [{section}-outvar1+outvar2] + name = custom_multi + inputs = inputvar1, inputvar2 + outvar1, outvar2 = func1(inputs) + jacobian = func2(inputs) + """ + tag = outputs + outputs = list(outputs.split(VARARGS_DELIM)) + all_vars = ", ".join(outputs) + inputs = map(str.strip, + cp.get_opt_tag(section, "inputs", tag).split(",")) + # get the functions for each output + transform_functions = {} + output_index = slice(None, None, None) + for var in outputs: + # check if option can be cast as a float + try: + func = cp.get_opt_tag(section, var, tag) + except Exception: + func = cp.get_opt_tag(section, all_vars, tag) + output_index = slice(outputs.index(var), outputs.index(var)+1) + transform_functions[var] = [func, output_index] + s = "-".join([section, tag]) + if cp.has_option(s, "jacobian"): + jacobian = cp.get_opt_tag(section, "jacobian", tag) + else: + jacobian = None + return cls(inputs, outputs, transform_functions, jacobian=jacobian) + + # # ============================================================================= # @@ -2725,6 +2817,7 @@ def from_config(cls, cp, section, outputs, # dictionary of all transforms transforms = { CustomTransform.name: CustomTransform, + CustomTransformMultiOutputs.name: CustomTransformMultiOutputs, MchirpQToMass1Mass2.name: MchirpQToMass1Mass2, Mass1Mass2ToMchirpQ.name: Mass1Mass2ToMchirpQ, MchirpEtaToMass1Mass2.name: MchirpEtaToMass1Mass2, From f3c887eba66a012f0e47314ec7f38c5453c02824 Mon Sep 17 00:00:00 2001 From: Han Wang <68005354+HumphreyWang@users.noreply.github.com> Date: Sun, 3 Mar 2024 23:06:33 +0000 Subject: [PATCH 28/29] avoid using for loop in space psds; fix typo (#4657) * avoid using for loop in space psds; fix typo * change back to astropy constants --- pycbc/psd/analytical_space.py | 166 ++++++++++++++-------------------- 1 file changed, 68 insertions(+), 98 deletions(-) diff --git a/pycbc/psd/analytical_space.py b/pycbc/psd/analytical_space.py index 1f388ae59bd..4237a724124 100644 --- a/pycbc/psd/analytical_space.py +++ b/pycbc/psd/analytical_space.py @@ -29,7 +29,7 @@ """ import numpy as np -from astropy import constants +from astropy.constants import c from pycbc.psd.read import from_numpy_arrays @@ -49,11 +49,11 @@ def psd_lisa_acc_noise(f, acc_noise_level=3e-15): The PSD value or array for acceleration noise. Notes ----- - Pease see Eq.(11-13) in for more details. + Please see Eq.(11-13) in for more details. """ s_acc = acc_noise_level**2 * (1+(4e-4/f)**2)*(1+(f/8e-3)**4) s_acc_d = s_acc * (2*np.pi*f)**(-4) - s_acc_nu = (2*np.pi*f/constants.c.value)**2 * s_acc_d + s_acc_nu = (2*np.pi*f/c.value)**2 * s_acc_d return s_acc_nu @@ -74,10 +74,10 @@ def psd_lisa_oms_noise(f, oms_noise_level=15e-12): The PSD value or array for OMS noise. Notes ----- - Pease see Eq.(9-10) in for more details. + Please see Eq.(9-10) in for more details. """ s_oms_d = oms_noise_level**2 * (1+(2e-3/f)**4) - s_oms_nu = s_oms_d * (2*np.pi*f/constants.c.value)**2 + s_oms_nu = s_oms_d * (2*np.pi*f/c.value)**2 return s_oms_nu @@ -96,13 +96,13 @@ def lisa_psd_components(f, acc_noise_level=3e-15, oms_noise_level=15e-12): Returns ------- - [low_freq_component, high_freq_component] : list + low_freq_component, high_freq_component : The PSD value or array for acceleration and OMS noise. """ low_freq_component = psd_lisa_acc_noise(f, acc_noise_level) high_freq_component = psd_lisa_oms_noise(f, oms_noise_level) - return [low_freq_component, high_freq_component] + return low_freq_component, high_freq_component def omega_length(f, len_arm=2.5e9): @@ -120,7 +120,7 @@ def omega_length(f, len_arm=2.5e9): omega_len : float or numpy.array The value of 2*pi*f*LISA_arm_length. """ - omega_len = 2*np.pi*f * len_arm/constants.c.value + omega_len = 2*np.pi*f * len_arm/c.value return omega_len @@ -151,23 +151,18 @@ def analytical_psd_lisa_tdi_1p5_XYZ(length, delta_f, low_freq_cutoff, The TDI-1.5 PSD (X,Y,Z channel) for LISA. Notes ----- - Pease see Eq.(19) in for more details. + Please see Eq.(19) in for more details. """ len_arm = np.float64(len_arm) acc_noise_level = np.float64(acc_noise_level) oms_noise_level = np.float64(oms_noise_level) - psd = [] fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - for f in fr: - [s_acc_nu, s_oms_nu] = lisa_psd_components( - f, acc_noise_level, oms_noise_level) - omega_len = omega_length(f, len_arm) - psd.append(16*(np.sin(omega_len))**2 * - (s_oms_nu+s_acc_nu*(3+np.cos(omega_len)))) - fseries = from_numpy_arrays(fr, np.array(psd), - length, delta_f, low_freq_cutoff) + s_acc_nu, s_oms_nu = lisa_psd_components( + fr, acc_noise_level, oms_noise_level) + omega_len = omega_length(fr, len_arm) + psd = 16*(np.sin(omega_len))**2 * (s_oms_nu+s_acc_nu*(3+np.cos(omega_len))) - return fseries + return from_numpy_arrays(fr, psd, length, delta_f, low_freq_cutoff) def analytical_psd_lisa_tdi_2p0_XYZ(length, delta_f, low_freq_cutoff, @@ -196,23 +191,19 @@ def analytical_psd_lisa_tdi_2p0_XYZ(length, delta_f, low_freq_cutoff, The TDI-2.0 PSD (X,Y,Z channel) for LISA. Notes ----- - Pease see Eq.(20) in for more details. + Please see Eq.(20) in for more details. """ len_arm = np.float64(len_arm) acc_noise_level = np.float64(acc_noise_level) oms_noise_level = np.float64(oms_noise_level) - psd = [] fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - for f in fr: - [s_acc_nu, s_oms_nu] = lisa_psd_components( - f, acc_noise_level, oms_noise_level) - omega_len = omega_length(f, len_arm) - psd.append(64*(np.sin(omega_len))**2 * (np.sin(2*omega_len))**2 * - (s_oms_nu+s_acc_nu*(3+np.cos(2*omega_len)))) - fseries = from_numpy_arrays(fr, np.array(psd), - length, delta_f, low_freq_cutoff) + s_acc_nu, s_oms_nu = lisa_psd_components( + fr, acc_noise_level, oms_noise_level) + omega_len = omega_length(fr, len_arm) + psd = (64*(np.sin(omega_len))**2 * (np.sin(2*omega_len))**2 * + (s_oms_nu+s_acc_nu*(3+np.cos(2*omega_len)))) - return fseries + return from_numpy_arrays(fr, psd, length, delta_f, low_freq_cutoff) def analytical_csd_lisa_tdi_1p5_XY(length, delta_f, low_freq_cutoff, @@ -241,22 +232,19 @@ def analytical_csd_lisa_tdi_1p5_XY(length, delta_f, low_freq_cutoff, The CSD between LISA's TDI-1.5 channel X and Y. Notes ----- - Pease see Eq.(56) in for more details. + Please see Eq.(56) in for more details. """ len_arm = np.float64(len_arm) acc_noise_level = np.float64(acc_noise_level) oms_noise_level = np.float64(oms_noise_level) - csd = [] fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - for f in fr: - omega_len = omega_length(f, len_arm) - [s_acc_nu, s_oms_nu] = lisa_psd_components( - f, acc_noise_level, oms_noise_level) - csd.append(-8*np.sin(omega_len)**2 * np.cos(omega_len) * - (s_oms_nu+4*s_acc_nu)) - fseries = from_numpy_arrays(fr, np.array(csd), - length, delta_f, low_freq_cutoff) - return fseries + s_acc_nu, s_oms_nu = lisa_psd_components( + fr, acc_noise_level, oms_noise_level) + omega_len = omega_length(fr, len_arm) + csd = (-8*np.sin(omega_len)**2 * np.cos(omega_len) * + (s_oms_nu+4*s_acc_nu)) + + return from_numpy_arrays(fr, csd, length, delta_f, low_freq_cutoff) def analytical_psd_lisa_tdi_1p5_AE(length, delta_f, low_freq_cutoff, @@ -285,24 +273,20 @@ def analytical_psd_lisa_tdi_1p5_AE(length, delta_f, low_freq_cutoff, The PSD of LISA's TDI-1.5 channel A and E. Notes ----- - Pease see Eq.(58) in for more details. + Please see Eq.(58) in for more details. """ len_arm = np.float64(len_arm) acc_noise_level = np.float64(acc_noise_level) oms_noise_level = np.float64(oms_noise_level) - psd = [] fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - for f in fr: - [s_acc_nu, s_oms_nu] = lisa_psd_components( - f, acc_noise_level, oms_noise_level) - omega_len = omega_length(f, len_arm) - psd.append(8*(np.sin(omega_len))**2 * - (4*(1+np.cos(omega_len)+np.cos(omega_len)**2)*s_acc_nu + - (2+np.cos(omega_len))*s_oms_nu)) - fseries = from_numpy_arrays(fr, np.array(psd), - length, delta_f, low_freq_cutoff) + s_acc_nu, s_oms_nu = lisa_psd_components( + fr, acc_noise_level, oms_noise_level) + omega_len = omega_length(fr, len_arm) + psd = (8*(np.sin(omega_len))**2 * + (4*(1+np.cos(omega_len)+np.cos(omega_len)**2)*s_acc_nu + + (2+np.cos(omega_len))*s_oms_nu)) - return fseries + return from_numpy_arrays(fr, psd, length, delta_f, low_freq_cutoff) def analytical_psd_lisa_tdi_1p5_T(length, delta_f, low_freq_cutoff, @@ -331,23 +315,19 @@ def analytical_psd_lisa_tdi_1p5_T(length, delta_f, low_freq_cutoff, The PSD of LISA's TDI-1.5 channel T. Notes ----- - Pease see Eq.(59) in for more details. + Please see Eq.(59) in for more details. """ len_arm = np.float64(len_arm) acc_noise_level = np.float64(acc_noise_level) oms_noise_level = np.float64(oms_noise_level) - psd = [] fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - for f in fr: - [s_acc_nu, s_oms_nu] = lisa_psd_components( - f, acc_noise_level, oms_noise_level) - omega_len = omega_length(f, len_arm) - psd.append(32*np.sin(omega_len)**2 * np.sin(omega_len/2)**2 * - (4*s_acc_nu*np.sin(omega_len/2)**2 + s_oms_nu)) - fseries = from_numpy_arrays(fr, np.array(psd), - length, delta_f, low_freq_cutoff) + s_acc_nu, s_oms_nu = lisa_psd_components( + fr, acc_noise_level, oms_noise_level) + omega_len = omega_length(fr, len_arm) + psd = (32*np.sin(omega_len)**2 * np.sin(omega_len/2)**2 * + (4*s_acc_nu*np.sin(omega_len/2)**2 + s_oms_nu)) - return fseries + return from_numpy_arrays(fr, psd, length, delta_f, low_freq_cutoff) def averaged_lisa_fplus_sq_approx(f, len_arm=2.5e9): @@ -367,7 +347,7 @@ def averaged_lisa_fplus_sq_approx(f, len_arm=2.5e9): The sky and polarization angle averaged squared antenna response. Notes ----- - Pease see Eq.(36) in for more details. + Please see Eq.(36) in for more details. """ from scipy.interpolate import interp1d from astropy.utils.data import download_file @@ -405,7 +385,7 @@ def averaged_response_lisa_tdi_1p5(f, len_arm=2.5e9): The sky and polarization angle averaged TDI-1.5 response to GW. Notes ----- - Pease see Eq.(39) in for more details. + Please see Eq.(39) in for more details. """ omega_len = omega_length(f, len_arm) ave_fp2 = averaged_lisa_fplus_sq_approx(f, len_arm) @@ -431,7 +411,7 @@ def averaged_response_lisa_tdi_2p0(f, len_arm=2.5e9): The sky and polarization angle averaged TDI-2.0 response to GW. Notes ----- - Pease see Eq.(40) in for more details. + Please see Eq.(40) in for more details. """ omega_len = omega_length(f, len_arm) response_tdi_1p5 = averaged_response_lisa_tdi_1p5(f, len_arm) @@ -469,21 +449,19 @@ def sensitivity_curve_lisa_semi_analytical(length, delta_f, low_freq_cutoff, LISA's sensitivity curve (6-links). Notes ----- - Pease see Eq.(42-43) in for more details. + Please see Eq.(42-43) in for more details. """ - sense_curve = [] len_arm = np.float64(len_arm) acc_noise_level = np.float64(acc_noise_level) oms_noise_level = np.float64(oms_noise_level) fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) fp_sq = averaged_lisa_fplus_sq_approx(fr, len_arm) - for i in range(len(fr)): - [s_acc_nu, s_oms_nu] = lisa_psd_components( - fr[i], acc_noise_level, oms_noise_level) - omega_len = 2*np.pi*fr[i] * len_arm/constants.c.value - sense_curve.append((s_oms_nu + s_acc_nu*(3+np.cos(2*omega_len))) / - (omega_len**2*fp_sq[i])) - fseries = from_numpy_arrays(fr, np.array(sense_curve)/2, + s_acc_nu, s_oms_nu = lisa_psd_components( + fr, acc_noise_level, oms_noise_level) + omega_len = 2*np.pi*fr * len_arm/c.value + sense_curve = ((s_oms_nu + s_acc_nu*(3+np.cos(2*omega_len))) / + (omega_len**2*fp_sq)) + fseries = from_numpy_arrays(fr, sense_curve/2, length, delta_f, low_freq_cutoff) return fseries @@ -509,19 +487,15 @@ def sensitivity_curve_lisa_SciRD(length, delta_f, low_freq_cutoff): LISA's sensitivity curve in SciRD. Notes ----- - Pease see Eq.(114) in for more details. + Please see Eq.(114) in for more details. """ - sense_curve = [] fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - for f in fr: - s_I = 5.76e-48 * (1+(4e-4/f)**2) - s_II = 3.6e-41 - R = 1 + (f/2.5e-2)**2 - sense_curve.append(10/3 * (s_I/(2*np.pi*f)**4+s_II) * R) - fseries = from_numpy_arrays(fr, sense_curve, - length, delta_f, low_freq_cutoff) + s_I = 5.76e-48 * (1+(4e-4/fr)**2) + s_II = 3.6e-41 + R = 1 + (fr/2.5e-2)**2 + sense_curve = 10/3 * (s_I/(2*np.pi*fr)**4+s_II) * R - return fseries + return from_numpy_arrays(fr, sense_curve, length, delta_f, low_freq_cutoff) def sensitivity_curve_lisa_confusion(length, delta_f, low_freq_cutoff, @@ -557,7 +531,7 @@ def sensitivity_curve_lisa_confusion(length, delta_f, low_freq_cutoff, LISA's sensitivity curve with Galactic confusion noise. Notes ----- - Pease see Eq.(85-86) in for more details. + Please see Eq.(85-86) in for more details. """ if base_model == "semi": base_curve = sensitivity_curve_lisa_semi_analytical( @@ -571,13 +545,11 @@ def sensitivity_curve_lisa_confusion(length, delta_f, low_freq_cutoff, if duration < 0 or duration > 10: raise Exception("Must between 0 and 10.") fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) - sh_confusion = [] f1 = 10**(-0.25*np.log10(duration)-2.7) fk = 10**(-0.27*np.log10(duration)-2.47) - for f in fr: - sh_confusion.append(0.5*1.14e-44*f**(-7/3)*np.exp(-(f/f1)**1.8) * - (1.0+np.tanh((fk-f)/(0.31e-3)))) - fseries_confusion = from_numpy_arrays(fr, np.array(sh_confusion), + sh_confusion = (0.5*1.14e-44*fr**(-7/3)*np.exp(-(fr/f1)**1.8) * + (1.0+np.tanh((fk-fr)/0.31e-3))) + fseries_confusion = from_numpy_arrays(fr, sh_confusion, length, delta_f, low_freq_cutoff) fseries = from_numpy_arrays(base_curve.sample_frequencies, base_curve+fseries_confusion, @@ -622,7 +594,7 @@ def sh_transformed_psd_lisa_tdi_XYZ(length, delta_f, low_freq_cutoff, noise, transformed from LISA sensitivity curve. Notes ----- - Pease see Eq.(7,41-43) in for more details. + Please see Eq.(7,41-43) in for more details. """ fr = np.linspace(low_freq_cutoff, (length-1)*2*delta_f, length) if tdi == "1.5": @@ -680,13 +652,11 @@ def semi_analytical_psd_lisa_confusion_noise(length, delta_f, low_freq_cutoff, raise Exception("The version of TDI, currently only for 1.5 or 2.0.") fseries_response = from_numpy_arrays(fr, np.array(response), length, delta_f, low_freq_cutoff) - sh_confusion = [] f1 = 10**(-0.25*np.log10(duration)-2.7) fk = 10**(-0.27*np.log10(duration)-2.47) - for f in fr: - sh_confusion.append(0.5*1.14e-44*f**(-7/3)*np.exp(-(f/f1)**1.8) * - (1.0+np.tanh((fk-f)/(0.31e-3)))) - fseries_confusion = from_numpy_arrays(fr, np.array(sh_confusion), + sh_confusion = (0.5*1.14e-44*fr**(-7/3)*np.exp(-(fr/f1)**1.8) * + (1.0+np.tanh((fk-fr)/0.31e-3))) + fseries_confusion = from_numpy_arrays(fr, sh_confusion, length, delta_f, low_freq_cutoff) psd_confusion = 2*fseries_confusion.data * fseries_response.data fseries = from_numpy_arrays(fseries_confusion.sample_frequencies, From f35d04bb218c004eb59fe96034c4f61aa30a3be2 Mon Sep 17 00:00:00 2001 From: Tito Dal Canton Date: Fri, 8 Mar 2024 17:58:39 +0100 Subject: [PATCH 29/29] Fix import broken by change in lscsoft-glue (#4662) * Drop lscsoft-glue dependency * Glue still needed, do not drop it! --- pycbc/results/render.py | 2 -- pycbc/results/versioning.py | 20 +++----------------- pycbc/workflow/core.py | 3 ++- pycbc/workflow/datafind.py | 7 ++++++- 4 files changed, 11 insertions(+), 21 deletions(-) diff --git a/pycbc/results/render.py b/pycbc/results/render.py index d10c4989e6d..3f58c4bf6fb 100644 --- a/pycbc/results/render.py +++ b/pycbc/results/render.py @@ -1,5 +1,3 @@ -#!/usr/bin/python - # Copyright (C) 2015 Christopher M. Biwer # # This program is free software; you can redistribute it and/or modify it diff --git a/pycbc/results/versioning.py b/pycbc/results/versioning.py index 64536377d22..529b39eeb91 100644 --- a/pycbc/results/versioning.py +++ b/pycbc/results/versioning.py @@ -1,5 +1,3 @@ -#!/usr/bin/python - # Copyright (C) 2015 Ian Harry # # This program is free software; you can redistribute it and/or modify it @@ -20,8 +18,9 @@ import subprocess import urllib.parse -import lal, lalframe -import pycbc.version, glue.git_version +import lal +import lalframe +import pycbc.version def get_library_version_info(): """This will return a list of dictionaries containing versioning @@ -88,19 +87,6 @@ def add_info_new_version(info_dct, curr_module, extra_str): pass library_list.append(lalsimulationinfo) - glueinfo = {} - glueinfo['Name'] = 'LSCSoft-Glue' - glueinfo['ID'] = glue.git_version.id - glueinfo['Status'] = glue.git_version.status - glueinfo['Version'] = glue.git_version.version - glueinfo['Tag'] = glue.git_version.tag - glueinfo['Author'] = glue.git_version.author - glueinfo['Builder'] = glue.git_version.builder - glueinfo['Branch'] = glue.git_version.branch - glueinfo['Committer'] = glue.git_version.committer - glueinfo['Date'] = glue.git_version.date - library_list.append(glueinfo) - pycbcinfo = {} pycbcinfo['Name'] = 'PyCBC' pycbcinfo['ID'] = pycbc.version.version diff --git a/pycbc/workflow/core.py b/pycbc/workflow/core.py index a455b37c81e..2853e50d528 100644 --- a/pycbc/workflow/core.py +++ b/pycbc/workflow/core.py @@ -36,7 +36,6 @@ import lal import lal.utils import Pegasus.api # Try and move this into pegasus_workflow -from glue import lal as gluelal from ligo import segments from ligo.lw import lsctables, ligolw from ligo.lw import utils as ligolw_utils @@ -1572,6 +1571,8 @@ def convert_to_lal_cache(self): """ Return all files in this object as a glue.lal.Cache object """ + from glue import lal as gluelal + lal_cache = gluelal.Cache([]) for entry in self: try: diff --git a/pycbc/workflow/datafind.py b/pycbc/workflow/datafind.py index 9dc898effcd..42055d25949 100644 --- a/pycbc/workflow/datafind.py +++ b/pycbc/workflow/datafind.py @@ -34,7 +34,6 @@ import urllib.parse from ligo import segments from ligo.lw import utils, table -from glue import lal from gwdatafind import find_urls as find_frame_urls from pycbc.workflow.core import SegFile, File, FileList, make_analysis_dir from pycbc.io.ligolw import LIGOLWContentHandler @@ -694,6 +693,8 @@ def setup_datafind_from_pregenerated_lcf_files(cp, ifos, outputDir, tags=None): datafindOuts : pycbc.workflow.core.FileList List of all the datafind output files for use later in the pipeline. """ + from glue import lal + if tags is None: tags = [] @@ -821,6 +822,8 @@ def get_missing_segs_from_frame_file_cache(datafindcaches): missingFrames: Dict. of ifo keyed lal.Cache instances The list of missing frames """ + from glue import lal + missingFrameSegs = {} missingFrames = {} for cache in datafindcaches: @@ -947,6 +950,8 @@ def run_datafind_instance(cp, outputDir, observatory, frameType, Cache file listing all of the datafind output files for use later in the pipeline. """ + from glue import lal + if tags is None: tags = []