def test_new_img_like_side_effect():
    img1 = Nifti1Image(np.ones((2, 2, 2, 2)), affine=np.eye(4))
    hash1 = joblib.hash(img1)
    new_img_like(img1, np.ones((2, 2, 2, 2)), img1.affine.copy(),
                 copy_header=True)
    hash2 = joblib.hash(img1)
    assert_equal(hash1, hash2)

def test_new_img_like_mgz():
    """Check that new images can be generated with bool MGZ type
    This is usually when computing masks using MGZ inputs, e.g.
    when using plot_stap_map
    """

    ref_img = nibabel.load(os.path.join(datadir, 'test.mgz'))
    data = np.ones(ref_img.get_data().shape, dtype=np.bool)
    affine = ref_img.affine
    new_img_like(ref_img, data, affine, copy_header=False)

def test_new_img_like():
    # Give a list to new_img_like
    data = np.zeros((5, 6, 7))
    data[2:4, 1:5, 3:6] = 1
    affine = np.diag((4, 3, 2, 1))
    img = nibabel.Nifti1Image(data, affine=affine)
    img2 = new_img_like([img, ], data)
    np.testing.assert_array_equal(img.get_data(), img2.get_data())

    # test_new_img_like_with_nifti2image_copy_header
    img_nifti2 = nibabel.Nifti2Image(data, affine=affine)
    img2_nifti2 = new_img_like([img_nifti2, ], data, copy_header=True)
    np.testing.assert_array_equal(img_nifti2.get_data(), img2_nifti2.get_data())

def test_new_img_like_mgz():
    """Check that new images can be generated with bool MGZ type
    This is usually when computing masks using MGZ inputs, e.g.
    when using plot_stap_map
    """

    if not LooseVersion(nibabel.__version__) >= LooseVersion('1.2.0'):
        # Old nibabel do not support MGZ files
        raise SkipTest

    ref_img = nibabel.load(os.path.join(datadir, 'test.mgz'))
    data = np.ones(ref_img.get_data().shape, dtype=np.bool)
    affine = ref_img.get_affine()
    new_img_like(ref_img, data, affine, copy_header=False)

def test_encode_nii():
    mni = datasets.load_mni152_template()
    encoded = html_stat_map._encode_nii(mni)
    decoded = html_stat_map._decode_nii(encoded)
    assert np.allclose(mni.get_data(), decoded.get_data())

    mni = image.new_img_like(mni, np.asarray(mni.get_data(), dtype='>f8'))
    encoded = html_stat_map._encode_nii(mni)
    decoded = html_stat_map._decode_nii(encoded)
    assert np.allclose(mni.get_data(), decoded.get_data())

    mni = image.new_img_like(mni, np.asarray(mni.get_data(), dtype='<i4'))
    encoded = html_stat_map._encode_nii(mni)
    decoded = html_stat_map._decode_nii(encoded)
    assert np.allclose(mni.get_data(), decoded.get_data())

def flip_img_lr(img):
    """ Convenience function to flip image on X axis"""
    # This won't work for all image formats! But
    # does work for those that we're working with...
    assert isinstance(img, nib.nifti1.Nifti1Image)
    img = new_img_like(img, data=img.get_data()[::-1], copy_header=True)
    return img

def plot_tbss(img, mean_FA_skeleton, start, end, row_l=6, step=1, title='',
    axis='z', pngfile=None):
    ''' Inspired from plot_two_maps. Plots a TBSS contrast map over the
    skeleton of a mean FA map'''

    # Dilate tbss map
    import numpy as np
    from skimage.morphology import cube, dilation
    from nilearn import image
    d = np.array(image.load_img(img).dataobj)
    dil_tbss = dilation(d, cube(2))
    dil_tbss_img = image.new_img_like(img, dil_tbss)

    slice_nb = int(abs(((end - start) / float(step))))
    images = []

    for line in range(int(slice_nb/float(row_l) + 1)):
        opt = {'title':{True:title,
                        False:None}[line==0],
               'colorbar':False,
               'black_bg':True,
               'display_mode':axis,
               'threshold':0.2,
               'cmap': cm.Greens,
               'cut_coords':range(start + line * row_l * step,
                                       start + (line+1) * row_l * step,
                                       step)}
        method = 'plot_stat_map'
        opt.update({'stat_map_img': mean_FA_skeleton})

        t = getattr(plotting, method).__call__(**opt)


        try:
            # Add overlay
            t.add_overlay(dil_tbss_img, cmap=cm.hot, threshold=0.95, colorbar=True)
        except TypeError:
            print img, 'probably empty tbss map'
            pass

        # Converting to PIL and appending it to the list
        buf = io.BytesIO()
        t.savefig(buf)
        buf.seek(0)
        im = Image.open(buf)
        images.append(im)

    # Joining the images
    imsize = images[0].size
    out = Image.new('RGBA', size=(imsize[0], len(images)*imsize[1]))
    for i, im in enumerate(images):
        box = (0, i * imsize[1], imsize[0], (i+1) * imsize[1])
        out.paste(im, box)

    if pngfile is None:
        import tempfile
        pngfile = tempfile.mkstemp(suffix='.png')[1]
    print 'Saving to...', pngfile, '(%s)'%title

    out.save(pngfile)

def clean_img(img):
    """ Remove nan/inf entries."""
    img = check_niimg(img)
    img_data = img.get_data()
    img_data[np.isnan(img_data)] = 0
    img_data[np.isinf(img_data)] = 0
    return new_img_like(img, img_data, copy_header=True)

def test_view_img():
    mni = datasets.load_mni152_template()
    with warnings.catch_warnings(record=True) as w:
        # Create a fake functional image by resample the template
        img = image.resample_img(mni, target_affine=3 * np.eye(3))
        html_view = html_stat_map.view_img(img)
        _check_html(html_view)
        html_view = html_stat_map.view_img(img, threshold='95%')
        _check_html(html_view)
        html_view = html_stat_map.view_img(img, bg_img=mni)
        _check_html(html_view)
        html_view = html_stat_map.view_img(img, bg_img=None)
        _check_html(html_view)
        html_view = html_stat_map.view_img(img, threshold=2., vmax=4.)
        _check_html(html_view)
        html_view = html_stat_map.view_img(img, symmetric_cmap=False)
        img_4d = image.new_img_like(img, img.get_data()[:, :, :, np.newaxis])
        assert len(img_4d.shape) == 4
        html_view = html_stat_map.view_img(img_4d, threshold=2., vmax=4.)
        _check_html(html_view)

    # Check that all warnings were expected
    warnings_set = set(warning_.category for warning_ in w)
    expected_set = set([FutureWarning, UserWarning,
                       DeprecationWarning])
    assert warnings_set.issubset(expected_set), (
        "the following warnings were not expected: {}").format(
        warnings_set.difference(expected_set))

def test_resample_stat_map():

    # Start with simple simulated data
    bg_img, data = _simulate_img()

    # Now double the voxel size and mess with the affine
    affine = 2 * np.eye(4)
    affine[3, 3] = 1
    affine[0, 1] = 0.1
    stat_map_img = Nifti1Image(data, affine)

    # Make a mask for the stat image
    mask_img = new_img_like(stat_map_img, data > 0, stat_map_img.affine)

    # Now run the resampling
    stat_map_img, mask_img = html_stat_map._resample_stat_map(
        stat_map_img, bg_img, mask_img, resampling_interpolation='nearest')

    # Check positive isotropic, near-diagonal affine
    _check_affine(stat_map_img.affine)
    _check_affine(mask_img.affine)

    # Check voxel size matches bg_img
    assert stat_map_img.affine[0, 0] == bg_img.affine[0, 0], (
        "stat_map_img was not resampled at the resolution of background")
    assert mask_img.affine[0, 0] == bg_img.affine[0, 0], (
        "mask_img was not resampled at the resolution of background")

def new_nii_like(ref_img, data, affine=None, copy_header=True):
    """
    Coerces `data` into NiftiImage format like `ref_img`

    Parameters
    ----------
    ref_img : str or img_like
        Reference image
    data : (S [x T]) array_like
        Data to be saved
    affine : (4 x 4) array_like, optional
        Transformation matrix to be used. Default: `ref_img.affine`
    copy_header : bool, optional
        Whether to copy header from `ref_img` to new image. Default: True

    Returns
    -------
    nii : :obj:`nibabel.nifti1.Nifti1Image`
        NiftiImage
    """

    ref_img = check_niimg(ref_img)
    nii = new_img_like(ref_img,
                       data.reshape(ref_img.shape[:3] + data.shape[1:]),
                       affine=affine,
                       copy_header=copy_header)
    nii.set_data_dtype(data.dtype)

    return nii

def test_new_img_like():
    # Give a list to new_img_like
    data = np.zeros((5, 6, 7))
    data[2:4, 1:5, 3:6] = 1
    affine = np.diag((4, 3, 2, 1))
    img = nibabel.Nifti1Image(data, affine=affine)
    img2 = new_img_like([img, ], data)
    np.testing.assert_array_equal(img.get_data(), img2.get_data())

    def join_networks(self, network_indices):
        """Return a NiftiImage containing a binarised version of the sum of
        the RSN images of each of the `network_indices`."""
        oimg = self._get_img(network_indices[0]).get_data()
        for idx in network_indices[1:]:
            oimg += self._get_img(idx).get_data()

        return niimg.new_img_like(self._get_img(network_indices[0]),
                                  oimg.astype(bool))

    def _run_interface(self, runtime):
        t1_img = nli.load_img(self.inputs.in_file)
        t1_data = t1_img.get_data()
        epi_data = nli.load_img(self.inputs.ref_file).get_data()

        # We assume the image is already masked
        mask = t1_data > 0

        t1_min, t1_max = np.unique(t1_data)[[1, -1]]
        epi_min, epi_max = np.unique(epi_data)[[1, -1]]
        scale_factor = (epi_max - epi_min) / (t1_max - t1_min)

        inv_data = mask * ((t1_max - t1_data) * scale_factor + epi_min)

        out_file = fname_presuffix(self.inputs.in_file, suffix='_inv', newpath=runtime.cwd)
        nli.new_img_like(t1_img, inv_data, copy_header=True).to_filename(out_file)
        self._results['out_file'] = out_file
        return runtime

def inject_skullstripped(subjects_dir, subject_id, skullstripped):
    mridir = op.join(subjects_dir, subject_id, 'mri')
    t1 = op.join(mridir, 'T1.mgz')
    bm_auto = op.join(mridir, 'brainmask.auto.mgz')
    bm = op.join(mridir, 'brainmask.mgz')

    if not op.exists(bm_auto):
        img = nb.load(t1)
        mask = nb.load(skullstripped)
        bmask = new_img_like(mask, mask.get_data() > 0)
        resampled_mask = resample_to_img(bmask, img, 'nearest')
        masked_image = new_img_like(img, img.get_data() * resampled_mask.get_data())
        masked_image.to_filename(bm_auto)

    if not op.exists(bm):
        copyfile(bm_auto, bm, copy=True, use_hardlink=True)

    return subjects_dir, subject_id

def rescale(source, target):
    import nibabel as nib
    import numpy as np
    from nilearn import image

    n = nib.load(source)
    d = np.array(n.dataobj)
    s = n.dataobj.slope
    i = image.new_img_like(n, d/s)
    i.to_filename(target)
    log.info('Rescaling done: %s rescaled to %s'%(source, target))

def compute_confounds(imgs, mask_img, n_confounds=5, get_randomized_svd=False,
                      compute_not_mask=False):
    """
    """
    confounds = []
    if not isinstance(imgs, collections.Iterable) or \
            isinstance(imgs, _basestring):
        imgs = [imgs, ]

    img = _utils.check_niimg_4d(imgs[0])
    shape = img.shape[:3]
    affine = get_affine(img)

    if isinstance(mask_img, _basestring):
        mask_img = _utils.check_niimg_3d(mask_img)

    if not _check_same_fov(img, mask_img):
        mask_img = resample_img(
            mask_img, target_shape=shape, target_affine=affine,
            interpolation='nearest')

    if compute_not_mask:
        print("Non mask based confounds extraction")
        not_mask_data = np.logical_not(mask_img.get_data().astype(np.int))
        whole_brain_mask = masking.compute_multi_epi_mask(imgs)
        not_mask = np.logical_and(not_mask_data, whole_brain_mask.get_data())
        mask_img = new_img_like(img, not_mask.astype(np.int), affine)

    for img in imgs:
        print("[Confounds Extraction] {0}".format(img))
        img = _utils.check_niimg_4d(img)
        print("[Confounds Extraction] high ariance confounds computation]")
        high_variance = high_variance_confounds(img, mask_img=mask_img,
                                                n_confounds=n_confounds)
        if compute_not_mask and get_randomized_svd:
            signals = masking.apply_mask(img, mask_img)
            non_constant = np.any(np.diff(signals, axis=0) != 0, axis=0)
            signals = signals[:, non_constant]
            signals = signal.clean(signals, detrend=True)
            print("[Confounds Extraction] Randomized SVD computation")
            U, s, V = randomized_svd(signals, n_components=n_confounds,
                                     random_state=0)
            if high_variance is not None:
                confound_ = np.hstack((U, high_variance))
            else:
                confound_ = U
        else:
            confound_ = high_variance
        confounds.append(confound_)

    return confounds

def get_largest_blobs(ic_maps):
    """ Generator for the largest blobs in each IC spatial map.
    These should be masked and thresholded.

    Parameters
    ----------
    ic_maps: sequence of niimg-like

    Returns
    -------
    blobs: generator of niimg-like
    """
    # store the average value of the blob in a list
    for i, icimg in enumerate(iter_img(ic_maps)):
        yield niimg.new_img_like(icimg, largest_connected_component(icimg.get_data()))

    def fit(self, X=None, y=None):  # noqa
        super(HemisphereMasker, self).fit(X, y)

        # x, y, z
        hemi_mask_data = reorder_img(self.mask_img_).get_data().astype(np.bool)

        xvals = hemi_mask_data.shape[0]
        midpt = np.ceil(xvals / 2.0)
        if self.hemi == "r":
            other_hemi_slice = slice(midpt, xvals)
        else:
            other_hemi_slice = slice(0, midpt)

        hemi_mask_data[other_hemi_slice] = False
        mask_data = self.mask_img_.get_data() * hemi_mask_data
        self.mask_img_ = new_img_like(self.mask_img_, data=mask_data)

        return self

def split_bilateral_rois(maps_img):
    """Convenience function for splitting bilateral ROIs
    into two unilateral ROIs"""

    new_rois = []

    for map_img in iter_img(maps_img):
        for hemi in ["L", "R"]:
            hemi_mask = HemisphereMasker(hemisphere=hemi)
            hemi_mask.fit(map_img)
            if hemi_mask.mask_img_.get_data().sum() > 0:
                hemi_vectors = hemi_mask.transform(map_img)
                hemi_img = hemi_mask.inverse_transform(hemi_vectors)
                new_rois.append(hemi_img.get_data())

    new_maps_data = np.concatenate(new_rois, axis=3)
    new_maps_img = new_img_like(maps_img, data=new_maps_data, copy_header=True)
    print("Changed from %d ROIs to %d ROIs" % (maps_img.shape[-1], new_maps_img.shape[-1]))
    return new_maps_img