Search code examples
pythonphysicsastronomyhealpy

Coordinate system used for healpy.query_disc()


The healpy.query_disc() function takes an argument vec which is a three-component unit vector defining the center of the disc. What coordinate system is being used here - why is there a third dimension for a 2-d projection? What point is the "tail" of this vector llocated at?


Solution

  • Very good you found the solution yourself, for later reference here is a full working code example:

    import healpy as hp
    import numpy as np
    
    
    # `lonlat=True` switches `ang2vec` from requiring colatitude $\theta$ and longitude $\phi$ in radians to longitude and latitude in degrees (notice that also the order changes)
    
    # in degrees
    lon = 60
    lat = 30
    vec = hp.ang2vec(lon, lat, lonlat=True)
    
    
    nside = 256
    large_disc = hp.query_disc(nside, vec, radius=np.radians(20))
    small_disc = hp.query_disc(nside, vec, radius=np.radians(8))
    tiny_disc = hp.query_disc(nside, vec, radius=np.radians(2))
    
    
    # `query_disc` returns a list of pixels, by default in RING ordering, let's check their length:
    
    list(map(len, [large_disc, small_disc, tiny_disc]))
    
    
    # ## Create a map and plot it in Mollweide projection
    
    m = np.zeros(hp.nside2npix(nside))
    
    
    m[large_disc] = 1
    m[small_disc] = 2
    m[tiny_disc] = 3
    
    
    hp.mollview(m)
    hp.graticule()
    

    mollweide projection

    See the notebook with plots here: https://zonca.dev/2020/10/example-healpy-query_disc.html