I'm trying to implement 2D cross-correlation to acquire the displacement of these two small images (interrogation windows). The issue is that the result from Cross-Correlation (CC) is different compared to the FFT approach; therefore, the displacement/fitting would be wrong. The ffwtools package in R uses the C library: Fastest Fourier Transform in the West [(FFTW)](https://www.fftw.org/) to do the conversions. The R code and the issue: ```{r} # Assume Window1 at t1 is: Window1 <- structure(c( 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ), dim = c(8L, 8L)) # Assume Window2 at t2 is: Window2 <- structure(c( 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1 ), dim = c(8L, 8L)) # Getting the cross-correlation using gsignal package CC_gsignal <- gsignal::xcorr2(Window1, Window2) ``` Now using FFT approach (i.e. fftwtools in R), I have generally followed the example/function [here.](https://de.mathworks.com/matlabcentral/fileexchange/53570-xcorr2_fft-a-b) ```{r} # dimensions of the CC_fft matrix cc_row <- nrow(Window1) + nrow(Window2) - 1 cc_col <- ncol(Window1) + ncol(Window2) - 1 # padding arrays padded_Window1 <- matrix(0, nrow = cc_row, ncol = cc_col) padded_Window1[1:nrow(Window1), 1:ncol(Window1)] <- Window1 # For padding Window2, I have tried multiple approaches; # option 3 produces the correct values but not the locations padded_Window2 <- matrix(0, nrow = cc_row, ncol = cc_col) # # Option 1 # padded_Window2[1:nrow(Window2), # 1:ncol(Window2)] <- Window2 # # Option 2 # padded_Window2[1:nrow(Window2), # 1:ncol(Window2)] <- Window2 [nrow(Window2):1, # ncol(Window2):1] # option 3 padded_Window2[nrow(Window2):cc_row, ncol(Window2):cc_col] <- Window2 # # Option 4 # padded_Window2[nrow(Window2):cc_row, # ncol(Window2):cc_col] <- Window2 [nrow(Window2):1, # ncol(Window2):1] # fft fftWindow1 <- matrix(fftwtools::fftw2d(padded_Window1), nrow(padded_Window1)) fftWindow2 <- matrix(fftwtools::fftw2d(padded_Window2), nrow(padded_Window2)) # CC_fftwtools CC_fftwtools <- round( matrix( Re(fftwtools::fftw2d(Conj(fftWindow1) * fftWindow2, inverse = 1)), nrow(padded_Window1) ) / length(padded_Window1), digits = 3 ) # plotting library(plot.matrix) par(mfrow = c(2, 2)) plot(Window1, col = topo.colors, main = "Window1") plot(Window2, col = topo.colors, main = "Window2") plot(CC_gsignal, col = topo.colors, main = "CC_gsignal") plot(CC_fftwtools, col = topo.colors, main = "CC_fftwtools") ``` [![enter image description here][1]][1] From the two bottom images, you can see the difference between CC from gsginal and fftwtools. I have tried to solve the issue by changing how Window2 is padded (see code above). I have also tried to shift the `CC_fftwtools` array using `gsignal::fftshift` without success. One dirty solution would be to flip the CC_fftwtools result vertically & horizontaly. However, I believe I'm doing something wrong, and there is a better explanation for my error. Any hints or tips are much appreciated! New Info: When exchanging Window1 with Window2, the results of CC exchange as well. Meaning the left image CC_gsignal will look like the current cc_fftwtools and vice versa. Edit: As per @OverLordGoldDragon's request, I add below a Python version using Scipy as well as FFT approach: ```python import numpy as np import matplotlib.pyplot as plt from numpy.fft import fft, ifft, fft2, ifft2, ifftshift from PIL import Image from scipy import signal Window1 =array([[0., 0., 0., 0., 0., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 0., 0.], [0., 0., 1., 1., 1., 1., 0., 0.], [0., 0., 1., 1., 1., 0., 0., 0.], [0., 0., 1., 1., 1., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 0., 0.]]) Window2 = array([[0., 0., 0., 0., 0., 0., 0., 0.], [0., 1., 1., 1., 1., 0., 0., 0.], [0., 1., 1., 1., 0., 0., 0., 0.], [0., 1., 1., 1., 0., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 0., 0.], [0., 0., 0., 0., 0., 0., 1., 1.], [0., 0., 0., 0., 0., 0., 1., 1.]]) corr_scipy = signal.correlate2d(Window1, Window2, boundary='fill', mode='full') def cross_correlate_2d(x, h): h = ifftshift(ifftshift(h, axes=0), axes=1) return ifft2(fft2(x) * np.conj(fft2(h))).real # no zero padding corr_fft = cross_correlate_2d(Window1, Window2) fig = plt.figure(figsize=(8, 8)) fig.add_subplot(1,2, 1); plt.imshow(corr_scipy, cmap='turbo') fig.add_subplot(1,2, 2); plt.imshow(corr_fft, cmap='turbo') plt.show() ``` [![enter image description here][2]][2] [1]: https://i.sstatic.net/HXG7Z.png [2]: https://i.sstatic.net/M0yme.png