Basic Usage


All the healpix maps must have a proper header defining their :

  • frame using the COORDSYS keyword,
  • order using the ORDERING keyword.

However you can correct the headers in the construction of the list of maps

maps = [ {'HFI 100': {'filename': 'data/HFI_SkyMap_100_2048_R2.00_full.fits', 'COORDSYS': 'C'}}]

There is two main way to use hpproj, the first way is to use the standalone program on the command line, this will efficiently produce cuts for similar maps, or use it programmatically from within a python script or program which will offer an additional speed-up on high memory system.

From the command line - cutsky

The command line program is called cutsky and takes 3 arguments at minimum, the longitude and latitude of the desired projection (by default in galactic coordinate, but see below) and a list of healpix map to cut from :

$ cutsky 0.0 0.0 --mapfilenames data/HFI_SkyMap_100_2048_R2.00_full.fits  data/HFI_SkyMap_857_2048_R2.00_full.fits

by default this will produce two png files centered on galactic longitude and latitude (0,0). Fits images of central photometries can be obtain using the --fits or --phot options. Help on cutsky can be obtain by

$ cutsky -h

usage: cutsky [-h] [--npix NPIX | --radius RADIUS] [--pixsize PIXSIZE]
              [--coordframe {galactic,fk5}]
              [--mapfilenames MAPFILENAMES [MAPFILENAMES ...]] [--fits]
              [--png] [--votable aperture [aperture ...]] [--outdir OUTDIR] [-v | -q] [--conf CONF]
              lon lat

Reproject the spherical sky onto a plane.

positional arguments:
  lon                   longitude of the projection [deg]
  lat                   latitude of the projection [deg]

optional arguments:
  -h, --help            show this help message and exit
  --npix NPIX           number of pixels (default 256)
  --radius RADIUS       radius of the requested region [deg]
  --pixsize PIXSIZE     pixel size [arcmin] (default 1)
  --coordframe {galactic,fk5}
                        coordinate frame of the lon. and lat. of the
                        projection and the projected map (default: galactic)
                        any projection code supported by wcslib (default:TAN)

input maps:
  one of the two options must be present
  --mapfilenames MAPFILENAMES [MAPFILENAMES ...]
                        absolute path to the healpix maps
  --conf CONF           absolute path to a config file

  --fits                output fits file
  --png                 output png file (Default: True if nothing else)
  --votable aperture [aperture ...]
                        list of aperture [arcmin] to make circular aperture photometry
  --outdir OUTDIR       output directory (default:.)

  -v, --verbose         verbose mode
  -q, --quiet           quiet mode

It takes two float arguments, the latitude and longitude center of the requested projection, either in galactic or equatorial coordinate frame (controled by the --coordframe option) and a list of healpix maps, either on the command line with the --mapfilenames argument or describe in a config file (with the --conf option). Several other optional arguments can also be set like --npix the number of pixels, their size (--pixsize) or the projection type --ctype.

The cutted maps can be saved as fits (--fits) or png (--png) and central circular aperture photometry can be performed and saved as a votable (--votable aperture). The output products directory can be tune using the --outdir option. All theses options can also be provided by the config file.

The config file follows a simple ini syntax with a global section [cutsky] to gather all previous options. The rest of the sections is used to describe the healpix maps used by cutsky. The section name [test] will be used as a legend and index by cutsky.

npix = 256
pixsize = 2
coordframe = galactic
png = True

filename = hpproj/data/CMB_I_SMICA_128_R2.00.fits
docut = False

[HFI 100]
filename = hpproj/data/HFI_SkyMap_100_128_R2.00_RING.fits

[HFI 857]
filename = hpproj/data/HFI_SkyMap_857_128_R2.00_NESTED.fits
docut = True
docontour = True

As a function call - cutsky()

It is also possible to call cutsky from a python program or script, as a function. You first need to define a list of maps on which to perform the cuts, as list of tuple with at minimum (filename.fits, {'legend': "legend")) given the full path to the healpix map, and a dictionnary with a least the key legend

from hpproj import cutsky

maps = [('data/HFI_SkyMap_100_2048_R2.00_full.fits', {'legend': 'HFI 100', 'aperture': [1, 2, 3]}),
        ('data/HFI_SkyMap_857_2048_R2.00_full.fits', {'legend': 'HFI 857', 'docontour': True})]

result = cutsky([0.0, 0.0], maps=maps)

The first argument is the latitude and longitude of the requested maps, by default in galactic frame (see the coordframe keyword), and the maps list define the healpix maps.

This will produce a list of dictionnaries containing 4 keys:

  • legend,
  • fits an,
  • png, a b61encoded png image of the fits
  • phot, the corresponding photometry

Additionnal parameters can by passed to the function :

  • patch=[256,1] : the size of the patch in pixel, and the size of the pixels in arcmin
  • ctype='TAN' : the desired type of projection

As an object - CutSky

It is however more efficient to use cutsky as an object :

from hpproj import CutSky, to_coord

maps = [('data/HFI_SkyMap_100_2048_R2.00_full.fits', {'legend': 'HFI 100', 'aperture': [1, 2, 3]}),
        ('data/HFI_SkyMap_857_2048_R2.00_full.fits', {'legend': 'HFI 857', 'docontour': True})]

cutsky = CutSky(maps, low_mem=False)

coord = to_coord([0.0,0.0])
result = cutsky.cut_fits(coord) # Will only produce the 'fits' key
result = cutsky.cut_png(coord)  # Will only produce the 'png' key (and 'fits' if absent)
result = cutsky.cut_phot(coord) # Will only produce the 'phot' key (and fits' if absent)

The result product should be similar to the cutsky() function. However with the low_mem keyword the healpix maps will be read only once in memory, for all cut_* calls. Similar to cutsky() several keyword parameters can be passed to CutSky() :

  • npix=256 : the size of the patch in pixels
  • pixsize=1 : the size of the pixels in arcmin
  • ctype='TAN': the desired type of projection

As internal calls - hp_helper

Alternatively if you simply want to get a projected array, you can use the hp_project() function

from import fits
from astropy.coordinates import SkyCoord
import healpy as hp
from hpproj import hp_project

hp_data, hp_header = hp.read_map('data/HFI_SkyMap_100_2048_R2.00_full.fits', h=True)
hp_header = fits.Header(hp_header)

hdu = hp_project(hp_data, hp_header, SkyCoord(0, 0, unit='deg'))

Or, if you prefer to get full control, you can also use the internal functions like build_wcs() and hp_to_wcs()

from import fits
import healpy as hp
import hpproj as hpp

hp_data, hp_header = hp.read_map('data/HFI_SkyMap_100_2048_R2.00_full.fits', h=True)
hp_header = fits.Header(hp_header)
hp_hdu = fits.ImageHDU(hp_data, hp_header)

w = hpp.build_wcs(0, 0)

proj_map = hpp.hp_to_wcs(hp_data, hp_header, w)

Note that both hp_project and hp_to_wcs accept either an, or both hp_data, hp_header