def _extract_hemispheric_difference(image, mask = slice(None), sigma_active = 7, sigma_reference = 7, cut_plane = 0, voxelspacing = None):
    """
    Internal, single-image version of `hemispheric_difference`.
    """
    # constants
    INTERPOLATION_RANGE = int(10) # how many neighbouring values to take into account when interpolating the medial longitudinal fissure slice
    
    # check arguments
    if cut_plane >= image.ndim:
        raise ArgumentError('The suppliedc cut-plane ({}) is invalid, the image has only {} dimensions.'.format(cut_plane, image.ndim))
    
    # set voxel spacing
    if voxelspacing is None:
        voxelspacing = [1.] * image.ndim
    
    # compute the (presumed) location of the medial longitudinal fissure, treating also the special of an odd number of slices, in which case a cut into two equal halves is not possible
    medial_longitudinal_fissure = int(image.shape[cut_plane] / 2)
    medial_longitudinal_fissure_excluded = image.shape[cut_plane] % 2
    
    # split the head into a dexter and sinister half along the saggital plane
    # this is assumed to be consistent with a cut of the brain along the medial longitudinal fissure, thus separating it into its hemispheres
    slicer = [slice(None)] * image.ndim
    slicer[cut_plane] = slice(None, medial_longitudinal_fissure)
    left_hemisphere = image[slicer]

    slicer[cut_plane] = slice(medial_longitudinal_fissure + medial_longitudinal_fissure_excluded, None)
    right_hemisphere = image[slicer]
    
    # flip right hemisphere image along cut plane
    slicer[cut_plane] = slice(None, None, -1)
    right_hemisphere = right_hemisphere[slicer]

    # substract once left from right and once right from left hemisphere, including smoothing steps
    right_hemisphere_difference = _substract_hemispheres(right_hemisphere, left_hemisphere, sigma_active, sigma_reference, voxelspacing)
    left_hemisphere_difference = _substract_hemispheres(left_hemisphere, right_hemisphere, sigma_active, sigma_reference, voxelspacing)
    
    # re-flip right hemisphere image to original orientation
    right_hemisphere_difference = right_hemisphere_difference[slicer]
    
    # estimate the medial longitudinal fissure if required
    if 1 == medial_longitudinal_fissure_excluded:
        left_slicer = [slice(None)] * image.ndim
        right_slicer = [slice(None)] * image.ndim
        left_slicer[cut_plane] = slice(-1 * INTERPOLATION_RANGE, None)
        right_slicer[cut_plane] = slice(None, INTERPOLATION_RANGE)
        interp_data_left = left_hemisphere_difference[left_slicer]
        interp_data_right = right_hemisphere_difference[right_slicer]
        interp_indices_left = list(range(-1 * interp_data_left.shape[cut_plane], 0))
        interp_indices_right = list(range(1, interp_data_right.shape[cut_plane] + 1))
        interp_data = numpy.concatenate((left_hemisphere_difference[left_slicer], right_hemisphere_difference[right_slicer]), cut_plane)
        interp_indices = numpy.concatenate((interp_indices_left, interp_indices_right), 0)
        medial_longitudinal_fissure_estimated = interp1d(interp_indices, interp_data, kind='cubic', axis=cut_plane)(0)
        # add singleton dimension
        slicer[cut_plane] = numpy.newaxis
        medial_longitudinal_fissure_estimated = medial_longitudinal_fissure_estimated[slicer]

    # stich images back together
    if 1 == medial_longitudinal_fissure_excluded:
        hemisphere_difference = numpy.concatenate((left_hemisphere_difference, medial_longitudinal_fissure_estimated, right_hemisphere_difference), cut_plane)
    else:
        hemisphere_difference = numpy.concatenate((left_hemisphere_difference, right_hemisphere_difference), cut_plane)

    # extract intensities and return
    return _extract_intensities(hemisphere_difference, mask)
def _extract_hemispheric_difference ( image , mask = slice ( None ) , sigma_active = 7 , sigma_reference = 7 , cut_plane = 0 , voxelspacing = None ) : INTERPOLATION_RANGE = int ( 10 ) if cut_plane >= image . ndim : raise ArgumentError ( 'The suppliedc cut-plane ({}) is invalid, the image has only {} dimensions.' . format ( cut_plane , image . ndim ) ) if voxelspacing is None : voxelspacing = [ 1. ] * image . ndim medial_longitudinal_fissure = int ( image . shape [ cut_plane ] / 2 ) medial_longitudinal_fissure_excluded = image . shape [ cut_plane ] % 2 slicer = [ slice ( None ) ] * image . ndim slicer [ cut_plane ] = slice ( None , medial_longitudinal_fissure ) left_hemisphere = image [ slicer ] slicer [ cut_plane ] = slice ( medial_longitudinal_fissure + medial_longitudinal_fissure_excluded , None ) right_hemisphere = image [ slicer ] slicer [ cut_plane ] = slice ( None , None , - 1 ) right_hemisphere = right_hemisphere [ slicer ] right_hemisphere_difference = _substract_hemispheres ( right_hemisphere , left_hemisphere , sigma_active , sigma_reference , voxelspacing ) left_hemisphere_difference = _substract_hemispheres ( left_hemisphere , right_hemisphere , sigma_active , sigma_reference , voxelspacing ) right_hemisphere_difference = right_hemisphere_difference [ slicer ] if 1 == medial_longitudinal_fissure_excluded : left_slicer = [ slice ( None ) ] * image . ndim right_slicer = [ slice ( None ) ] * image . ndim left_slicer [ cut_plane ] = slice ( - 1 * INTERPOLATION_RANGE , None ) right_slicer [ cut_plane ] = slice ( None , INTERPOLATION_RANGE ) interp_data_left = left_hemisphere_difference [ left_slicer ] interp_data_right = right_hemisphere_difference [ right_slicer ] interp_indices_left = list ( range ( - 1 * interp_data_left . shape [ cut_plane ] , 0 ) ) interp_indices_right = list ( range ( 1 , interp_data_right . shape [ cut_plane ] + 1 ) ) interp_data = numpy . concatenate ( ( left_hemisphere_difference [ left_slicer ] , right_hemisphere_difference [ right_slicer ] ) , cut_plane ) interp_indices = numpy . concatenate ( ( interp_indices_left , interp_indices_right ) , 0 ) medial_longitudinal_fissure_estimated = interp1d ( interp_indices , interp_data , kind = 'cubic' , axis = cut_plane ) ( 0 ) slicer [ cut_plane ] = numpy . newaxis medial_longitudinal_fissure_estimated = medial_longitudinal_fissure_estimated [ slicer ] if 1 == medial_longitudinal_fissure_excluded : hemisphere_difference = numpy . concatenate ( ( left_hemisphere_difference , medial_longitudinal_fissure_estimated , right_hemisphere_difference ) , cut_plane ) else : hemisphere_difference = numpy . concatenate ( ( left_hemisphere_difference , right_hemisphere_difference ) , cut_plane ) return _extract_intensities ( hemisphere_difference , mask )
Internal, single-image version of `hemispheric_difference`.