code: def indices(image, voxelspacing = None, mask = slice(None)):
      r
      if type(image) == tuple or type(image) == list:
          image = image[0]
          
      if not type(mask) is slice:
          mask = numpy.array(mask, copy=False, dtype=numpy.bool)
          
      if voxelspacing is None:
          voxelspacing = [1.] * image.ndim
  
      return join(*[a[mask].ravel() * vs for a, vs in zip(numpy.indices(image.shape), voxelspacing)])
docs: r"""
      Takes an image and returns the voxels ndim-indices as voxel-wise feature. The voxel
      spacing is taken into account, i.e. the indices are not array indices, but millimeter
      indices.
      
      This is a multi-element feature where each element corresponds to one of the images
      axes, e.g. x, y, z, ...
      
      Parameters
      ----------
      image : array_like or list/tuple of array_like 
          A single image or a list/tuple of images (for multi-spectral case).
      voxelspacing : sequence of floats
          The side-length of each voxel.
      mask : array_like
          A binary mask for the image. 
      
      Returns
      -------
      indices : ndarray
          Each voxels ndimensional index.
  
      Notes
      -----
      This feature is independent of the actual image content, but depends
      solely on its shape. Therefore always a one-dimensional feature is returned, even if
      a multi-spectral image has been supplied.