pythonhealpy

How should I plot healpix binning


I was plotting celestial sphere with equal area binning. My first intention was making a binning based on the basic calculus knowledge which I learned from 1st year University calculus class. Following is what I made.

import numpy as np
import matplotlib.pyplot as plt

def Drawpixel(lm,lM,bm,bM,co,w):
    l0=[lm,lM]
    b0=[bM,bM]
    l1=[lm,lM]
    b1=[bm,bm]
    l2=[lm,lm]
    b2=[bm,bM]
    l3=[lM,lM]
    b3=[bm,bM]
    plt.plot(l0,b0,c=co,linewidth=w)
    plt.plot(l1,b1,c=co,linewidth=w)
    plt.plot(l2,b2,c=co,linewidth=w)
    plt.plot(l3,b3,c=co,linewidth=w)

def Paint(lm,lM,bm,bM,co):
    l=[lm,lM]
    b1=[bM,bM]
    b2=[bm,bm]
    plt.fill_between(l,b1,b2,color=co)
    

plt.subplot(111,projection="aitoff")
plt.title("Equal Binning", fontsize = 30)
co='k'
w=1
for i in range(16):
    bM=np.pi/2-np.arccos((8-i)/8)
    bm=np.pi/2-np.arccos((7-i)/8)
    for j in range(32):
        lm=(j-16)/16*np.pi*-1
        lM=(j-15)/16*np.pi*-1
        Drawpixel(lm,lM,bm,bM,co,w)

ax=plt.gca()
plt.pause(0.1)
ax.set_xticklabels(ax.get_xticklabels()[::-1])
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlabel("L",fontsize=25)
plt.ylabel("B",fontsize=25)
plt.show()

With this I would get something like thisenter image description here

After plotting the plot I made, I learned that there is a binning method named healpix and it would provide equal area binning with different shape. Also I learned that there was a python software named healpy, however, I could not find how could I generate similar plot like above from the healpy instruction. Is there a simple method I can make a plot using healpix? (If there is a simple function supported from healpy, that would be also appreciated.)


Solution

  • Healpy is commonly used in astrophysics for all sky data. You can use the functions in healpy to make such a map. Healpy typically uses fits files for input and output.

    import healpy as hp 
    import numpy as np 
    import pyplot as plt 
    nside = 128 # parameter that defines the pixel size 
    npix = hp.nside2npix(nside) 
    image = np.zeros(npix) # data needs to be an array of size npix or use hp.read_map
    hp.mollview(image) # you can use another projection from hp.visufunc
    hp.graticules()
    plt.show()