code: def __compute_edges_nd(label_image):
      
      Er = set()
     
      def append(v1, v2):
          if v1 != v2:
              Er.update([(min(v1, v2), max(v1, v2))])
          
      vappend = scipy.vectorize(append)
     
      for dim in range(label_image.ndim):
          slices_x = []
          slices_y = []
          for di in range(label_image.ndim):
              slices_x.append(slice(None, -1 if di == dim else None))
              slices_y.append(slice(1 if di == dim else None, None))
          vappend(label_image[slices_x], label_image[slices_y])
          
      return Er
docs: Computes the region neighbourhood defined by a star shaped n-dimensional structuring
      element (as returned by scipy.ndimage.generate_binary_structure(ndim, 1)) for the
      supplied region/label image.
      Note The returned set contains neither duplicates, nor self-references
      (i.e. (id_1, id_1)), nor reversed references (e.g. (id_1, id_2) and (id_2, id_1).
      
      @param label_image An image with labeled regions (nD).
      @param return A set with tuples denoting the edge neighbourhood.