Monday, October 18, 2010

Learning About HEALPix

I've been playing with Shirley's LRG catalog and mask. The mask uses something called HEALPix which I've never used before. It basically breaks up the spherical sky into equal area pixels:

The map is in the form of weights for each pixel. Shirley says that if the weight is greater than 0.2 she considers that pixel in the mask.

So to create a random catalog, I just go through the weights, and if it is greater than 0.2 I place objects in that pixel.

The key is to now convert from a Healpix mask to ra/dec.
There are two routines that do this that I have found:


The problem is that these routines return the ra/dec of the center of the pixel. And to avoid putting extra power at the separations of the pixels, we want to select random points within the pixel. However, as you can see below, the ra range of these pixels changes dramatically as a function of cos (theta) (i.e. large at the poles and small at the equator) and so it isn't a simple procedure to fill in the area of a pixel evenly with random points:

I've searched to see if there is an existing package that does this, and I can't seem to find anything. In the HEALPix documentation there is a comment on random number generation, it doesn't seem to be what I am interested in here. I also have a email out to Shirley about this.

There are many functions that deal with HEALPix on their web page and well as in IDLUTILS. I'm currently combing through them to see if there is a function that gives me the ra/dec ranges for a given pixel, or even better, automatically picks a random point given a pixel.

Some questions:
  1. How do I assign a random place in a HEALPix pixel and convert this to a ra/dec position?
  2. The mask contains a "weight" for each HEALPix pixel, how to do factor in the weight into the random catalog. My intuition is that pixel's with higher weight should have more objects in them, but will this add clustering signal to the random catalog that we don't want?
  3. How do we distinguish what is signal because of observation effects (we only looked at that patch of the sky once, or we didn't target as many objects there) and what is signal due to clustering of the data (there is a huge galaxy cluster in that part of the sky, and that is why there is more data there)?


  1. Hi Jessica,
    I'm new to healpix.I would really be glad if u could tell me about the Healpix coordinates.How is every pixel identified or named, like a point on earth being identified by latitude and longitude.How do u convert ra and dec to healpix coordinates?

  2. Dear Anonymous,
    Healpix breaks up a sphere into equal area pixel regions. Each pixel is assigned a number. In IDL, using IDLUtils you can translate between ra and dec and pixel number using ang2pix_nest or ang2ring. You can also get the center position of each pixel in ra/dec or x,y,z with pix2ang_nest and pix2vec_nest

    See for more details. There is also information on the JPL web page about utilities for different programming languages if you don't want to use IDL. Hope this is helpful.


  3. Dear Jessica,
    Can you please help me in plotting a HEALPix grid. How do you form a HEALPix grid with Nside=16 and how do u convert a normal fits image with RA/DEC into a fits image with Healpix map.
    Thank You.

  4. Hi.
    As far as i'm aware i'm not able to use RA and DEC input arrays in ang2pix_ring. I am trying to convert the co-ordinates (given in RA and DEC) of known LAT GRBs into positions on my HEALPix map. My pixel numbers are usually given by inputting theta and phi angles, along with Nside. Is there a simple way of me plotting these GRBs on my HEALPix maps? Thanks.

  5. Anonymous,
    I've never used ang2pix_ring. What I do is the following.

    1) I convert ra/dec to theta/phi (healpix expected coordinates)
    phi = pra*!pi/180.
    theta = (90.0-pdec)*!pi/180.

    2) Find out which pixel these coordinates correspond to:
    ang2pix_nest, nside_in, theta, phi, ipnest

    3) Read in my mask:
    read_fits_map, mask, maskweight, nside=nside_in, ordering=order

    4) See if they are inside my mask or not:
    dataweight = maskweight[ipnest]
    inMask = where(dataweight GT 0)

    Hope this helps.

  6. Rather than generating random points on a pixel by pixel basis, you could try generating a set of random points (i.e. ra/dec) that uniformly cover the whole sphere (this should be fairly straightforward) and then using ang2pix to determine whether a given point is in the mask or not. Of course, this requires more computational overhead than doing it your way but it is also very easy to implement.

  7. Hey Jessica,

    I'm having the same problem as you. Did you ever find a script/code that solved your problem or did you go with what the latest Anonymous comment suggested.



  8. Just saw the date on this post... can't image you're actually still working on this.

    1. Christopher,
      I suggest you look at some of my other posts regarding healpix. It might answer your questions:


  9. This comment has been removed by the author.

  10. What happens if you increase the number of pixels to (let's say) a trillion or so, and then use their centers?