Making skymaps from GWOSC .fits files with healpy

Data Analysis
1

I am trying to make skymaps from GWOSC data using healpy. I am using the following few lines of code

import healpy as hp
import numpy as np
from glob import glob

skymap_loc = ${SKYMAP_DIR/*}

for skymap in glob(skymap_loc):
    hpx = hp.read_map(skymap)

where I filled the skymap directory by downloading data from the “Skymap for ${event_name}” link on GWOSC.

This code results in the error:

ValueError: Wrong pixel number (it is not 12*nside**2)

Looking at the file header for one of these files with fitsheader (IGWN-GWTC3p0-v1-GW191103_012549_PEDataRelease_cosmo_reweight_C01_IMRPhenomXPHM.fits), I see that this may be due to the size of the data. In particular, the parameter

NAXIS2  =                16896 / length of dimension 2

seems incorrect. Comparing these lines to an example here, Working with Sky Maps — LIGO/Virgo Public Alerts User Guide 16 documentation, the value of NAXIS2 is equal to 12*2048^2. Later on the page, we see that NSIDE for the example fits file is 2048. sqrt(16896/12), the NAXIS2 for a GWOSC fits file, isn’t an integer. This is my suspicion of what’s happening, but I’m not very familiar with healpy. How can I correctly extract the data from the .fits files available on GWOSC?

1 Like
2

Hi Adrian,

Many of the skymaps uploaded on GWOSC are created with the package ligo.skymap and can be read using its I/O routines. These skymaps also contain layers representing the distance posterior conditioned on sky direction.

Below’s a small example of how to read a GWOSC skymap

from ligo.skymap.io.fits import read_sky_map
read_sky_map('GW190814_C01:IMRPhenomPv3HM.fits)

the read_sky_map skymap will return an array (corresponding to healpy indices) with the skymap and additional information. If you add the keywork distance=True it will also return information on the distance.

1 Like