Maximum of cross-correlation not moving

Question

I already asked this question here yesterday, but it was very poorly worded I think. I made a much more detailled post explaining my problem of stackoverflow, as it might also be a code problem.

Here is a copy-paste of this post :

I am working on a C# project involving calculating cross-correlations between arrays. The goal is to precisely track the shift of the pattern present on the images (with sub-pixel accuracy). I am using the MathNet library.

Here is the formula used for the cross correlation :

corr(a, b) = ifft(fft(a_and_zeros) * conj(fft(b_and_zeros)))

The double[1024]arrays I am working on are derived from 1024*768 grayscale images by calculating the mean of the signal alongside each 768 pixels row. Here is an exemple of an image and the corresponding array :

(the green curve is irrelevant here, don't mind it ; note that these are just an exemple, for my actual experiments I use un-saturated images)

Here is the "algorithm" I am using :

1. Isolating one 1024 * 768 image from the camera
2. Shifting the camera by a known distance (e.g. 3 pixels)
3. Isolating another image from the camera
4. Obtaining data from each image by averaging the pixels color of the 768 pixels on each of the 1024 row (I obtain two double[1024] array)
5. Zero-padding both signals : placing each double[1024] array in the center of a double[2048] array
6. Computing the cross correlation of the two signals
7. Looking for the Max of the cross-corrleation
8. Applying a parabolic fit function for sub-sample accuracy

My problem is the following : when there is a strong signal on the edge of the initial array (as you can see on the exemple), and if the shift of the image is small (less than 10 pixels I'd say, but I can't be sure this is the actual limit), the cross correlation will indicate a shift less than 0.1 pixels.

I understand that this is probably coming from the use of FFT/IFFT, but shouldn't the use of a zero-padding allow to avoid this behaviour ? What could I do to avoid it ?

As my goal is to measure very small shifts, it is completely preventing me from progressing.

Here is the code I wrote for the cross-correlation (ExperimentInit and ExperimentFinal are my 2 initials arrays) :

ExperimentInitWindowed = new double[ExperimentInit.Length];
ExperimentFinalWindowed = new double[ExperimentFinal.Length];
ExpInitLarge = new Complex[2 * ExperimentInit.Length];
ExpFinalLarge = new Complex[2 * ExperimentFinal.Length];
CrossCorrExpArrayReShifted = new double[RefInitLarge.Length];

//Zero-padding
for (int i = 0; i < ExpFinalLarge.Length / 2; i++)
{
    ExpInitLarge[i + ExpInitLarge.Length / 4] = ExperimentInit[i];
    ExpFinalLarge[i + ExpFinalLarge.Length / 4] = ExperimentFinal[i];
}

//FFT
Accord.Math.Transforms.FourierTransform2.FFT(ExpInitLarge, FourierTransform.Direction.Forward);
Accord.Math.Transforms.FourierTransform2.FFT(ExpFinalLarge, FourierTransform.Direction.Forward);

//Conjugating ExpFinal
Complex[] CompConjExpFinal = new Complex[RefInitLarge.Length];
for (int l = 0; l < ExpInitLarge.Length; l++)
{
    CompConjExpFinal[l] = Complex.Conjugate(ExpFinalLarge[l]);
    ExpFinalLarge[l] = CompConjExpFinal[l];
}

//Element-wise multiplication of the complex arrays
Complex[] ExpMultipliedArray = new Complex[ExpFinalLarge.Length];
for (int l = 0; l < RefFinalLarge.Length; l++)
{
    ExpMultipliedArray[l] = Complex.Multiply(ExpInitLarge[l], ExpFinalLarge[l]);
}

//InverseFFT
Accord.Math.Transforms.FourierTransform2.FFT(ExpMultipliedArray, FourierTransform.Direction.Backward);
double[] CrossCorrExpArrayRe = ExpMultipliedArray.Re();

//Reorganizing the data inside the array to center the zero frequency component
double[] LeftHalfExp = new double[ExpFinalLarge.Length / 2];
double[] RightHalfExp = new double[ExpFinalLarge.Length / 2];
double[] LeftHalfRef = new double[ExpFinalLarge.Length / 2];
double[] RightHalfRef = new double[ExpFinalLarge.Length / 2];

for (int l = 0; l < ExpFinalLarge.Length / 2; l++)
{
    LeftHalfExp[l] = CrossCorrExpArrayRe[l];
    RightHalfExp[l] = CrossCorrExpArrayRe[l + (ExpFinalLarge.Length / 2)];
    LeftHalfRef[l] = CrossCorrRefArrayRe[l];
    RightHalfRef[l] = CrossCorrRefArrayRe[l + (ExpFinalLarge.Length / 2)];
    CrossCorrExpArrayReShifted[l] = RightHalfExp[l];
    CrossCorrExpArrayReShifted[l + (ExpFinalLarge.Length / 2)] = LeftHalfExp[l];
}

And here is the code for finding the max of the crosscorrelation and the shift :

//Finding max of the CrossCorr, and index of it
CrossCorrExpArrayReShiftedMax = CrossCorrExpArrayReShifted.Max();
CrossCorrExpArrayReShiftedMaxIndex = Array.IndexOf(CrossCorrExpArrayReShifted, CrossCorrExpArrayReShiftedMax);

//Performing LeastSquare Parabolic fitting
double[] xExp = { CrossCorrExpArrayReShiftedMaxIndex - 1, CrossCorrExpArrayReShiftedMaxIndex, CrossCorrExpArrayReShiftedMaxIndex + 1 };
double[] yExp = { CrossCorrExpArrayReShifted[CrossCorrExpArrayReShiftedMaxIndex - 1], CrossCorrExpArrayReShiftedMax, CrossCorrExpArrayReShifted[CrossCorrExpArrayReShiftedMaxIndex + 1] };
double[] CoefficientsExp = Fit.Polynomial(xExp, yExp, 2);
AtermExp = CoefficientsExp[2];
BtermExp = CoefficientsExp[1];
CtermExp = CoefficientsExp[0];

//Calculating shift
shiftLeastSquareExp = (-BtermExp / (2 * AtermExp)) - (CrossCorrExpArrayReShifted.Length / 2);

I will gladly provide any additional informations you could need. Thanks !

EDIT : I just tried using the cross-correlation method described here, and I get the exact same results than with my homemade code. This indicates the problem does not come from the code, but rather from my understanding of the way cross-correlation works, I think. Anybody can shed some light on this ?

Answer 1

I suspect that you need to pre-whiten the images.

Here's an example with and without prewhitening.

I've found a simple way to do pre-whitening (though in no way optimal) is to just do a column-wise diff on the images first.

Some more text here.