Pastro

Pastro is a package to estimate the probabilities of astrophysical origins of the gravitational wave candidates detected by the GstLAL pipeline, for the purpose of assisting in real-time multi-messenger follow-up observations. Events like the Binary Neutron Star Coalescence (BNS) also emit EM follow-up signals, such as Gamma-Ray bursts, while producing gravitational wave signals from its collision. We develop a software package that can estimate the probability of the sources of these GW events.

To generate pastro models with pastro, one needs to prepare a mass model with mass cuts implemented for the NS boundaries. In O4, GstLAL uses manifold to generate the mass models and mass models with mass cuts.

GitLab Repo: https://git.ligo.org/gstlal/pastro.git

Method: FGMC

Ray, Anarya. et. al. (2023), dcc-link: https://dcc.ligo.org/P2300141

The FGMC method implements Bayesain rate estimation from gravitational wave trigger data by constructing a Poisson mixture model for the various astrophysical and terrestrial source categories. The inferred rates are then used to compute the probability of astrophysical origin for every trigger found above threshold [1, 2]. After installing a writable GstLAL container with pastro package complied, you can generate the pastro model following the example code below (on CIT). You can find the instruction to install the container and pastro in the README.md of the pastro repo above. There is also an example Makefile in the repo: /path/to/pastro/examples/fgmc/Makefile. (citation: [1] Kapadia, S. Caudill, S. Creighton, J. Farr, W. et al. arxiv: 1903.06881; [2] Farr, W. Gair, J. Mandel, I. Cutler, C. arxiv 1302.5341)

If you want to run the code on the ICDS, simply change the environment variable from $TMPDIR to /ligo.

If you already have the template-weights files and estimated rates, please skip the part (a) & (b), directly generate the pastro model from part (c).

(a) Rate Estimation (O4)

Rate estimation is important to a low-latency pastro calculation, which approximate the expected rates of detectable GW sources in the absence of the full list of the triggers.

To calculate the expected rate of detectable triggers, we use the astrophysical rates estimated using FGMC through O1+O2+O3 (BBH, BNS, NSBH). These numbers go into the --rates-astro option in the example given below. We also need the sensitive hypervolumes for the past runs (BBH,BNS,NSBH). These comprise the --VT-old. We scale the hypervolumes to O4 sensitivity using the projected BNS range of O4 and the BNS ranges of O1+O2+O3, the run-times of the past observing runs, all taken from the observing scenarios doc. We also need the FAR threshold used in the past observing run (--threshold-old) FGMC estimates as well as their combined terrestrial count whichis close to 11541 ie the total number of gstlal triggers found above threshold in O1O2O3. For details of the scaling calculation and its derivation, see T2200165.

The calculation is implemented using the script calculate_rates which lives in the gwsci-pastro module. This script needs cosmology_utils.py which needs to be copied to the run directory. To estimate such rate for O4 (using O1+O2+O3 triggers):

singularity exec -B $TMPDIR gstlal-dev/ ./calculate_rates --runtime-old '{"O2"
:0.75,"O1":0.333333333333,"O3":1.0}' --dbns-old '{"O1":0.08,"O2":0.1,"O3":0.12}' --threshold-old 8760.12 --dbns-new
0.160  --rates-astro '{"BBH": 40.0, "NSBH":68.0,"BNS":163.0}' --VT-old '{"BBH":1.63,"NSBH":0.048,"BNS":0.014}'  --terrestrial-old 11541 --T-inj-start 1262304018 --T-inj-stop 1265760018  --N-inj '{"BBH":118036.67,"BNS":118036.67, "NSBH":118036.67}' --z-max '{"BBH":1.9, "BNS":0.15, "NSBH":0.25}' --output rates_V_O4.json 
{"rates_astrophysical": {"BBH": 40.0, "NSBH": 68.0, "BNS": 163.0}, "rates_injections": {"BBH": 4045.100324417181, "BNS": 1086074.9231302636, "NSBH": 269801.54621791927}, "V": {"BBH": 3.9860764449083943, "NSBH": 0.1173813922427012, "BNS": 0.03423623940412118}}

Note, the same script (calculate_rates) also calculates the rate of detectable injections using the total number of injections generated, the simulated run time and the redshift distributions of the injections. The other options options (T_inj,N_inj,and z_max) are for calculating the injection rates, taken from this repo . These numbers do not affect the generation of astrophysical rates. Even though this functionality of calculatign injection rates is there and was used in some MDCs, it will not be used in O4, and hence both the injection options and the injection rates output can be ignored. We did not split up the script as we did not want to make changes to the code this late unless absolutely necessary.

This script generates a file called rates_V.json which contains the astrophysical and injection rates as well as the scaled hyper volumes. The detectable rates required for pastro are obtained by multiplying these two, which the pastro module does by itself.

(b) Compute Template Weights (O4)

The template weights are the Bayes factors, or marginalized likelihood ratio, in the calculation of GW trigger probability.

First, create and enter template weight working directory:

$ mkdir -p /path/to/working/directory/template-weights
$ cd /path/to/working/directory/template-weights

using BBH as an example:

$ mkdir BBH
$ cd BBH
$ export GWSCI=/path/to/container/gstlal-dev
$ scp $GWSCI/src/manifold/examples/O4-bank-dag.yaml .

Edit the mass-model-cuts section of the yaml file to the following:

mass-model-cuts:
    m1: [2.5,500]
    m2: [2.5,500]
    x1: [-0.95,0.95]
    x2: [-0.95,0.95]

So as to demarcate the region of parameter space used for BBH. Copy the bank:

$ scp -r /home/anarya.ray/manifold-template-weights-official/new-bank/bank /path/to/working/directory/template-weights

Then make and run the dag, and plot template weights for various SNR’s(5,20,50):

$ singularity exec ${GWSCI} manifold_cbc_bank_mass_model_dag --output-h5 H1L1V1-O4_MANIFOLD_BANK_SALPETER_BBH_-0-2000000000.h5 --bank-file /path/to/working/directory/template-weights/bank/H1L1V1-O4_MANIFOLD_BANK-0-2000000000.h5 --singularity-image ${GWSCI} /path/to/working/directory/template-weights/BBH/O4-bank-dag.yaml
$ condor_submit_dag manifold_cbc_bank_mass_model.dag
$ singularity exec ${GWSCI} manifold_cbc_bank_plot_mass_model H1L1V1-O4_MANIFOLD_BANK_SALPETER_BBH_-0-2000000000.h5 5
$ singularity exec ${GWSCI} manifold_cbc_bank_plot_mass_model H1L1V1-O4_MANIFOLD_BANK_SALPETER_BBH_-0-2000000000.h5 20
$ singularity exec ${GWSCI} manifold_cbc_bank_plot_mass_model H1L1V1-O4_MANIFOLD_BANK_SALPETER_BBH_-0-2000000000.h5 50

For other categories, changing the name from BBH to BNS or NSBH. You will need the mass model of all 3 categories. Put them into the same directory for the $PATH to your template-weight files, and name them as ceff_<category>.h5. The $PATH will be the input of --weights-dir when generating the p-astro model below.

(c) Generate P(astro) Model

To generate a pastro model file for the final posterior probability (on CIT), you need to have the previously generated template-weights files, a ranking-statistic file, and a picked FAR threshold:

singularity exec -B $TMPDIR /path/to/gstlal-dev/ initialize_pastro_fgmc_model --rankstat-filename <your-rank-stat-pdf-file> --weights-dir <path-to-your-template-weights> --V-new <your-sensitive-volume V> --far-threshold <far-threshold> --output <output-model-h5file-name> --rates <rates> --rates-inj <injection-rate>

Once you have the pastro model generated, you can estimate the events by:

singularity exec -B $TMPDIR /path/to/gstlal-dev/ calculate_pastro_fgmc --fgmc-model pastro_model_fgmc.h5 --templateid {template_id} --snr {snr-thresh}  --lnlr {likelihood-thresh} {"BBH": <estimated-BBH-rate>, "BNS": <estimated-BNS-rate>, "NSBH": <estimated-NSBH-rate>, "Terr": <noise-rate>}

For example, you can estimate the p-astro of injection events:

singularity exec -B  /path/to/gstlal-dev/ calculate_pastro_fgmc_test --output-file pastro_inj.txt --fgmc-model pastro_model_fgmc.h5 --injection-file ${injdb}

Optional: Create Mass Model file for p-astro from manifold

To generate mass models with manifold, one needs a hdf5 format template bank. This can be generated with manifold as described here.

First, you need to have the manifold installed within the singularity. After you have the manifold installed, run the following command:

singularity exec -B $TMPDIR /path/to/gstlal-dev/ manifold_cbc_bank_mass_model_dag --singularity-image /path/to/gstlal-dev/ --output-h5 /path/to/output-h5-filename/ --bank-file /path/to/template-bank-h5-filename/ /path/to/bank-dag-yaml-file/

for O4, we use H1L1V1-O4_MANIFOLD_BANK-0-2000000000.h5 for the template-bank-h5, and a config-O4-bank for the yaml file.

Then, submit the dag with

condor_submit_dag manifold_cbc_bank_mass_model_dag.dag

After generating the mass_model files for all three categories (BBH, NSBH, BNS) from manifold, put 3 .h5 files into the same directory, which will be the path to --weights-dir to generate the pastro model as above.

Debug

If you run into the following error when running p-astro jobs:

Traceback (most recent call last):
  File "/usr/local/bin/your-pastro-job-name", line 4, in <module>
    __import__('pkg_resources').run_script('pastro==0.0.3', 'your-pastro-job-name')
  File "/usr/lib/python3.6/site-packages/pkg_resources/__init__.py", line 654, in run_script
    self.require(requires)[0].run_script(script_name, ns)
  File "/usr/lib/python3.6/site-packages/pkg_resources/__init__.py", line 1425, in run_script
    .format(**locals()),
pkg_resources.ResolutionError: Script 'scripts/your-pastro-job-name' not found in metadata at '/usr/local/lib/python3.6/site-packages/gwsci_pastro-0.0.2-py3.6.egg-info'

especially when you make changes and re-install the p-astro. Within the singularity, you need to delete the previously compiled program from --old-and-unmanageable option by:

rm -r /path/to/gstlal-dev/usr/local/lib/python3.6/site-packages/gwsci_pastro-0.0.2-py3.6.egg-info

and re-install the package without the --old-and-unmanageable option:

Singularity> python3 setup.py install

Plotting and Analysis

(add docs here once the plotting code is merged into repo) We majorly use Sankey diagram and gg-style plot to demonstrate and compare our result.