docs: def labelHealpix(pixels, values, nside, threshold=0, xsize=1000):
          """
          Label contiguous regions of a (sparse) HEALPix map. Works by mapping 
          HEALPix array to a Mollweide projection and applying scipy.ndimage.label
       
          Assumes non-nested HEALPix map.
          
          Parameters:
          pixels    : Pixel values associated to (sparse) HEALPix array
          values    : (Sparse) HEALPix array of data values
          nside     : HEALPix dimensionality
          threshold : Threshold value for object detection
          xsize     : Size of Mollweide projection
          
          Returns:
          labels, nlabels
          """
          proj = healpy.projector.MollweideProj(xsize=xsize)
          vec = healpy.pix2vec(nside,pixels)
          xy = proj.vec2xy(vec)
          ij = proj.xy2ij(xy)
          xx,yy = proj.ij2xy()
       
          # Convert to Mollweide
          searchims = []
          if values.ndim < 2: iterate = [values]
          else:               iterate = values.T
          for i,value in enumerate(iterate):
              logger.debug("Labeling slice %i...")
              searchim = numpy.zeros(xx.shape,dtype=bool)
              select = (value > threshold)
              yidx = ij[0][select]; xidx = ij[1][select]
              searchim[yidx,xidx] |= True
              searchims.append( searchim )
          searchims = numpy.array(searchims)
  
          # Full binary structure
          s = ndimage.generate_binary_structure(searchims.ndim,searchims.ndim)
       
          ### # Dilate in the z-direction
          logger.info("  Dilating image...")
          searchims = ndimage.binary_dilation(searchims,s,1)
          
          # Do the labeling
          logger.info("  Labeling image...")
          labels,nlabels = ndimage.label(searchims,structure=s)
  
          # Convert back to healpix
          pix_labels = labels[:,ij[0],ij[1]].T
          pix_labels = pix_labels.reshape(values.shape)
          pix_labels *= (values > threshold) # re-trim
  
          return pix_labels, nlabels