Sky Map I/O (ligo.skymap.io.fits)¶
Reading and writing HEALPix FITS files.
An example FITS header looks like this:
$ fitsheader test.fits.gz
# HDU 0 in test.fits.gz
SIMPLE = T / conforms to FITS standard
BITPIX = 8 / array data type
NAXIS = 0 / number of array dimensions
EXTEND = T
# HDU 1 in test.fits.gz
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / array data type
NAXIS = 2 / number of array dimensions
NAXIS1 = 4096 / length of dimension 1
NAXIS2 = 192 / length of dimension 2
PCOUNT = 0 / number of group parameters
GCOUNT = 1 / number of groups
TFIELDS = 1 / number of table fields
TTYPE1 = 'PROB '
TFORM1 = '1024E '
TUNIT1 = 'pix-1 '
PIXTYPE = 'HEALPIX ' / HEALPIX pixelisation
ORDERING= 'RING ' / Pixel ordering scheme, either RING or NESTED
COORDSYS= 'C ' / Ecliptic, Galactic or Celestial (equatorial)
EXTNAME = 'xtension' / name of this binary table extension
NSIDE = 128 / Resolution parameter of HEALPIX
FIRSTPIX= 0 / First pixel # (0 based)
LASTPIX = 196607 / Last pixel # (0 based)
INDXSCHM= 'IMPLICIT' / Indexing: IMPLICIT or EXPLICIT
OBJECT = 'FOOBAR 12345' / Unique identifier for this event
REFERENC= 'http://www.youtube.com/watch?v=0ccKPSVQcFk' / URL of this event
DATE-OBS= '2013-04-08T21:37:32.25' / UTC date of the observation
MJD-OBS = 56391.151064815 / modified Julian date of the observation
DATE = '2013-04-08T21:50:32' / UTC date of file creation
CREATOR = 'fits.py ' / Program that created this file
RUNTIME = 21.5 / Runtime in seconds of the CREATOR program
- ligo.skymap.io.fits.read_sky_map(filename, nest=False, distances=False, moc=False, **kwargs)[source] [edit on github]¶
Read a LIGO/Virgo/KAGRA-type sky map and return a tuple of the HEALPix array and a dictionary of metadata from the header.
- Parameters:
- filename: string
Path to the optionally gzip-compressed FITS file.
- nest: bool, optional
If omitted or False, then detect the pixel ordering in the FITS file and rearrange if necessary to RING indexing before returning.
If True, then detect the pixel ordering and rearrange if necessary to NESTED indexing before returning.
If None, then preserve the ordering from the FITS file.
Regardless of the value of this option, the ordering used in the FITS file is indicated as the value of the ‘nest’ key in the metadata dictionary.
- distances: bool, optional
If true, then read also read the additional HEALPix layers representing the conditional mean and standard deviation of distance as a function of sky location.
- moc: bool, optional
If true, then preserve multi-order structure if present.
Examples
Test that we can read a legacy IDL-compatible file (https://bugs.ligo.org/redmine/issues/5168):
>>> import tempfile >>> with tempfile.NamedTemporaryFile(suffix='.fits') as f: ... nside = 512 ... npix = ah.nside_to_npix(nside) ... ipix_nest = np.arange(npix) ... hp.write_map(f.name, ipix_nest, nest=True, column_names=['PROB']) ... m, meta = read_sky_map(f.name) ... np.testing.assert_array_equal(m, hp.ring2nest(nside, ipix_nest))
- ligo.skymap.io.fits.write_sky_map(filename, m, **kwargs)[source] [edit on github]¶
Write a gravitational-wave sky map to a file, populating the header with optional metadata.
- Parameters:
- filename: str
Path to the optionally gzip-compressed FITS file.
- m
astropy.table.Table,numpy.array If a Numpy record array or astorpy.table.Table instance, and has a column named ‘UNIQ’, then interpret the input as NUNIQ-style multi-order map [1]. Otherwise, interpret as a NESTED or RING ordered map.
- **kwargs
Additional metadata to add to FITS header. If m is an
astropy.table.Tableinstance, then the header is initialized from bothm.metaandkwargs.
References
[1]Górski, K.M., Wandelt, B.D., Hivon, E., Hansen, F.K., & Banday, A.J. 2017. The HEALPix Primer. The Unique Identifier scheme. http://healpix.sourceforge.net/html/intronode4.htm#SECTION00042000000000000000
Examples
Test header contents:
>>> order = 9 >>> nside = 2 ** order >>> npix = ah.nside_to_npix(nside) >>> prob = np.ones(npix, dtype=float) / npix
>>> import tempfile >>> from ligo.skymap import version >>> with tempfile.NamedTemporaryFile(suffix='.fits') as f: ... write_sky_map(f.name, prob, nest=True, ... vcs_version='foo 1.0', vcs_revision='bar', ... build_date='2018-01-01T00:00:00') ... for card in fits.getheader(f.name, 1).cards: ... print(str(card).rstrip()) XTENSION= 'BINTABLE' / binary table extension BITPIX = 8 / array data type NAXIS = 2 / number of array dimensions NAXIS1 = 8 / length of dimension 1 NAXIS2 = 3145728 / length of dimension 2 PCOUNT = 0 / number of group parameters GCOUNT = 1 / number of groups TFIELDS = 1 / number of table fields TTYPE1 = 'PROB ' TFORM1 = 'D ' TUNIT1 = 'pix-1 ' PIXTYPE = 'HEALPIX ' / HEALPIX pixelisation ORDERING= 'NESTED ' / Pixel ordering scheme: RING, NESTED, or NUNIQ COORDSYS= 'C ' / Ecliptic, Galactic or Celestial (equatorial) NSIDE = 512 / Resolution parameter of HEALPIX INDXSCHM= 'IMPLICIT' / Indexing: IMPLICIT or EXPLICIT VCSVERS = 'foo 1.0 ' / Software version VCSREV = 'bar ' / Software revision (Git) DATE-BLD= '2018-01-01T00:00:00' / Software build date
>>> uniq = moc.nest2uniq(np.uint8(order), np.arange(npix)) >>> probdensity = prob / hp.nside2pixarea(nside) >>> moc_data = np.rec.fromarrays( ... [uniq, probdensity], names=['UNIQ', 'PROBDENSITY']) >>> with tempfile.NamedTemporaryFile(suffix='.fits') as f: ... write_sky_map(f.name, moc_data, ... vcs_version='foo 1.0', vcs_revision='bar', ... build_date='2018-01-01T00:00:00') ... for card in fits.getheader(f.name, 1).cards: ... print(str(card).rstrip()) XTENSION= 'BINTABLE' / binary table extension BITPIX = 8 / array data type NAXIS = 2 / number of array dimensions NAXIS1 = 16 / length of dimension 1 NAXIS2 = 3145728 / length of dimension 2 PCOUNT = 0 / number of group parameters GCOUNT = 1 / number of groups TFIELDS = 2 / number of table fields TTYPE1 = 'UNIQ ' TFORM1 = 'K ' TTYPE2 = 'PROBDENSITY' TFORM2 = 'D ' TUNIT2 = 'sr-1 ' PIXTYPE = 'HEALPIX ' / HEALPIX pixelisation ORDERING= 'NUNIQ ' / Pixel ordering scheme: RING, NESTED, or NUNIQ COORDSYS= 'C ' / Ecliptic, Galactic or Celestial (equatorial) MOCORDER= 9 / MOC resolution (best order) INDXSCHM= 'EXPLICIT' / Indexing: IMPLICIT or EXPLICIT VCSVERS = 'foo 1.0 ' / Software version VCSREV = 'bar ' / Software revision (Git) DATE-BLD= '2018-01-01T00:00:00' / Software build date