Find the shift between two images using the DFT method

    Parameters
    ----------
    im0: 2d array
        First image
    im1: 2d array
        Second image
    isccs: Boolean, default false
        Set to True if the images are alredy DFT and in CCS representation
    subpix: boolean, default True
        Set to True (default) if you want subpixel precision

    Returns
    -------
    [y, x]: 2 numbers
        The offset

    Notes
    -----
    This algorithm detect a shift using the global phase difference of the DFTs

    If the images are already DFT and in the CCS format, set isccs to true.
    In that case the images should have the same size.

    If subpix is True, a gaussian fit is used for subpix precision

def find_shift_dft(im0, im1, isccs=False, subpix=True):
    """Find the shift between two images using the DFT method

    Parameters
    ----------
    im0: 2d array
        First image
    im1: 2d array
        Second image
    isccs: Boolean, default false
        Set to True if the images are alredy DFT and in CCS representation
    subpix: boolean, default True
        Set to True (default) if you want subpixel precision

    Returns
    -------
    [y, x]: 2 numbers
        The offset

    Notes
    -----
    This algorithm detect a shift using the global phase difference of the DFTs

    If the images are already DFT and in the CCS format, set isccs to true.
    In that case the images should have the same size.

    If subpix is True, a gaussian fit is used for subpix precision
    """
    # sanitize input
    im0 = np.asarray(im0, dtype=np.float32)
    im1 = np.asarray(im1, dtype=np.float32)

    # check input
    if not isccs:
        im0, im1 = dft_optsize_same(im0, im1)
    else:
        # Work only if the shapes are the same
        assert(im0.shape == im1.shape)

    # f0*conj(f1)
    mulSpec = cv2.mulSpectrums(im0, im1, flags=0, conjB=True)
    # norm(f0)*norm(f1)
    normccs = cv2.sqrt(cv2.mulSpectrums(im0, im0, flags=0, conjB=True) *
                       cv2.mulSpectrums(im1, im1, flags=0, conjB=True))
    # compute the inverse DFT
    xc = cv2.dft(ccs_normalize(mulSpec, normccs),
                 flags=cv2.DFT_REAL_OUTPUT | cv2.DFT_INVERSE)
    # Blur xc to remove some noise and improve the subpixel detection
    # workaround as GaussianBlur doesn't work with BORDER_WRAP
    blurRadii = 2
    xc = cv2.copyMakeBorder(xc, blurRadii, blurRadii, blurRadii, blurRadii,
                            borderType=cv2.BORDER_WRAP)
    xc = cv2.GaussianBlur(xc, (2 * blurRadii + 1, 2 * blurRadii + 1), 1.5)
    xc = xc[blurRadii:-blurRadii, blurRadii:-blurRadii]
    # save shape
    shape = np.asarray(xc.shape)
    # find max
    idx = np.asarray(np.unravel_index(np.argmax(xc), shape))

    """
    from matplotlib import pyplot as plt
    from numpy.fft import fftshift
    plt.figure()
    plt.imshow(np.log(np.abs(fftshift(im0))))
    plt.figure()
    plt.imshow(np.log(np.abs(fftshift(im1))))
    plt.figure()
    plt.imshow(fftshift(ccs_normalize(mulSpec,normccs)))
    plt.figure()
    extent= (-np.shape(xc)[1]/2, np.shape(xc)[1]/2, -np.shape(xc)[0]/2, np.shape(xc)[0]/2  )
    plt.imshow(np.log(np.abs(fftshift(xc))),extent = extent)

    #"""
    # plt.imshow(fftshift(xc))
    # print(idx)
    # plt.figure()
#    if toremove:
#        plt.figure(1)
#        l=len(xc[:,0])
#        plt.plot(np.arange(l)/l,xc[:,0])
#        print(l,xc[-1,0])
#        plt.figure(2)

    #"""

    if subpix:
        # update idx
        idx = np.asarray([get_peak_pos(xc[:, idx[1]], wrap=True),
                          get_peak_pos(xc[idx[0], :], wrap=True)])
    else:
        # restrics to reasonable values
        idx[idx > shape // 2] -= shape[idx > shape // 2]
    return idx