Posterior Samples to FITS Files (ligo-skymap-from-samples)


This example shows how to extract the conditional distance posterior for a given sky location from the kernel density estimator.

First, run ligo-skymap-from-samples to create the KDE and the 3D sky map from the posterior sample chain:

$ ligo-skymap-from-samples --maxpts 1000 posterior_samples.dat

Then you can load the pickle in a Python interpreter and evaluate the KDE at any point:

>>> from astropy.coordinates import SkyCoord
>>> from matplotlib import pyplot as plt
>>> import numpy as np
>>> import pickle
>>> with open('skypost.obj', 'rb') as f:  
...     skypost = pickle.load(f)  
>>> coord = SkyCoord.from_name('NGC 4993')
>>> distance = np.arange(1, 100)
>>> coords = np.column_stack((np.tile(coord.ra.rad, len(distance)),
...                           np.tile(coord.dec.rad, len(distance)),
...                           distance))  
>>> post = skypost.posterior_spherical(coords)  
>>> plt.plot(distance, post)  

Generate a FITS sky map file from posterior samples using clustering and kernel density estimation.

The input file should be an HDF5 file with the following columns:

  • ra, rightascension, or right_ascension: J2000 right ascension in


  • dec or declination: J200 declination in radians

  • dist, distance, or luminosity_distance: luminosity distance in Mpc (optional)

The output consist of two files:

  • skypost.obj, a pickle representation of the kernel density estimator

  • skymap.fits.gz, a 3D localization in HEALPix/FITS format

usage: ligo-skymap-from-samples [-h] [--seed SEED] [--version]
                                [-l CRITICAL|ERROR|WARNING|INFO|DEBUG|NOTSET]
                                [--outdir OUTDIR]
                                [--fitsoutname SKYMAP.fits[.gz]]
                                [--loadpost SKYPOST.obj] [--maxpts MAXPTS]
                                [--trials TRIALS] [--enable-distance-map]
                                [--top-nside TOP_NSIDE] [-j [JOBS]]
                                [--instruments H1|L1|V1|... [H1|L1|V1|... ...]]
                                [--objid OBJID] [--path PATH]
                                [--tablename TABLENAME]

Positional Arguments


posterior samples file

Named Arguments


show program’s version number and exit

-l, --loglevel

Default: INFO

--outdir, -o

output directory

Default: .


filename for the FITS file

Default: “skymap.fits”


filename for pickled posterior state


maximum number of posterior points to use; if omitted or greater than or equal to the number of posterior samples, then use all samples


number of trials at each clustering number

Default: 5

--enable-distance-map, --disable-distance-map

generate HEALPix map of distance estimates

Default: True

--enable-multiresolution, --disable-multiresolution

generate a multiresolution HEALPix map

Default: True


choose a start nside before HEALPix refinement steps (must be a valid nside)

Default: 16

-j, --jobs

Number of threads

Default: 1


instruments to store in FITS header


event ID to store in FITS header


The path of the dataset within the HDF5 file


The name of the table to search for recursively within the HDF5 file. By default, search for posterior_samples

Default: “posterior_samples”

random number generator options

Options that affect the Numpy pseudo-random number genrator


Pseudo-random number generator seed [default: initialized from /dev/urandom or clock]