docs: def pix2ang(nside, pix):
      r"""Convert pixel to angle :math:`\theta` :math:`\phi`.
  
          nside and pix are broadcast with numpy rules.
  
          Returns: theta, phi
  
          This is translated from chealpix.c; but refer to Section 4.1 of
          http://adsabs.harvard.edu/abs/2005ApJ...622..759G
      """
      nside, pix = numpy.lib.stride_tricks.broadcast_arrays(nside, pix)
  
      ncap = nside * (nside - 1) * 2
      npix = 12 * nside * nside
  
      def northpole(pix, npix):
          iring = (1 + ((1 + 2 * pix) ** 0.5)).astype('i8') // 2
          iphi = (pix + 1) - 2 * iring * (iring - 1)
          z = 1.0 - (iring*iring) * 4. / npix
          phi = (iphi - 0.5) * 0.5 * numpy.pi / iring
          return z, phi
  
      def equatorial(pix, nside, npix, ncap):
          ip = pix - ncap
          iring = ip // (4 * nside) + nside
          iphi = ip % (4 * nside) + 1
          fodd = (((iring + nside) &1) + 1.) * 0.5
          z = (2 * nside - iring) * nside * 8.0 / npix
          phi = (iphi - fodd) * (0.5 * numpy.pi) / nside
          return z, phi
  
      def southpole(pix, npix):
          ip = npix - pix
          iring = (1 + ((2 * ip - 1)**0.5).astype('i8')) // 2
          iphi = 4 * iring + 1 - (ip - 2 * iring * (iring - 1))
          z = -1 + (iring * iring) * 4. / npix
          phi = (iphi - 0.5 ) * 0.5 * numpy.pi / iring
          return z, phi
  
      mask1 = pix < ncap
  
      mask2 = (~mask1) & (pix < npix - ncap)
      mask3 = pix >= npix - ncap
  
      z = numpy.zeros(pix.shape, dtype='f8')
      phi = numpy.zeros(pix.shape, dtype='f8')
      z[mask1], phi[mask1] = northpole(pix[mask1], npix[mask1])
      z[mask2], phi[mask2] = equatorial(pix[mask2], nside[mask2], npix[mask2], ncap[mask2])
      z[mask3], phi[mask3] = southpole(pix[mask3], npix[mask3])
      return numpy.arccos(z), phi