Phase correlation

Question

How can rotation angle be determined by phase correlation(using fft) of 2 images? The algorithm given in http://en.wikipedia.org/wiki/Phase_correlation returns linear shift, not angular. It also mentions images have to be converted to log-polar coordinates to compute rotation. How is this conversion achieved in python? And post conversion do the same steps of the algorithm hold?

Answer 1

Log polar transformation is actually rotation and scale invariant.. Rotation corresponds to shift in y axis and scaling corresponds to shift in x axis in log polar transformation

So simple steps are as follows for finding an image x in image y:

1. Find image x in image y (use phase correlation in cartesian coordinates)

2. Compute log polar transforms of both x and y (this is a whole other problem, see references below), make sure to center on the same feature in both images.

3. Find FFT of x and y, say F(X) and F(y)

4. Find phase correlation of F(x) and F(y), call it R

5. Find the IFFT (inverse FFT) of R. The peak value of R corresponds to the rotation deviation in the Y Axis and to the Scaling deviation in the X Axis from the original Image.

References:

1. http://etd.lsu.edu/docs/available/etd-07072005-113808/unrestricted/Thunuguntla_thesis.pdf

Answer 2

I have been working on the same problem for a while. I took the week-end to write this. It's not the cleanest code there is, but I'm only a physicist, not a programmer...

The phase correlation itself is simple: use your favorite convolution algorithm to convolve two images. The peak position gives you the rotation/scaling difference. It's well explained on Wikipedia (in the link mentioned in the question).

My problem was I couldn't find a good log-polar converter, so I wrote one. It's not fool-proof, but it gets the job done. Anyone willing to rewrite it to make it clearer, please do so!

import scipy as sp
from scipy import ndimage
from math import *

def logpolar(input,silent=False):
    # This takes a numpy array and returns it in Log-Polar coordinates.

    if not silent: print("Creating log-polar coordinates...")
    # Create a cartesian array which will be used to compute log-polar coordinates.
    coordinates = sp.mgrid[0:max(input.shape)*2,0:360]
    # Compute a normalized logarithmic gradient
    log_r = 10**(coordinates[0,:]/(input.shape[0]*2.)*log10(input.shape[1]))
    # Create a linear gradient going from 0 to 2*Pi
    angle = 2.*pi*(coordinates[1,:]/360.)

    # Using scipy's map_coordinates(), we map the input array on the log-polar 
    # coordinate. Do not forget to center the coordinates!
    if not silent: print("Interpolation...")
    lpinput = ndimage.interpolation.map_coordinates(input,
                                            (log_r*sp.cos(angle)+input.shape[0]/2.,
                                             log_r*sp.sin(angle)+input.shape[1]/2.),
                                            order=3,mode='constant')

    # Returning log-normal...
    return lpinput

Warning: This code is designed for grayscale images. It can easily be adapter to work on color images by looping the line with map_coordinates() on each separate color frame.

EDIT: Now, the code to do the correlation is simple. After your script has imported both images as image and target, do the following:

# Conversion to log-polar coordinates
lpimage = logpolar(image)
lptarget = logpolar(target)

# Correlation through FFTs    
Fcorr = np.fft.fft2(lpimage)*np.conj(np.fft.fft2(lptarget))
correlation = np.fft.ifft2(Fcorr)

The array correlation should contain a peak which coordinates are the size difference and the angle difference. Also, instead of using FFTs, you could simply use numpy's np.correlate() function:

# Conversion to log-polar coordinates
lpimage = logpolar(image)
lptarget = logpolar(target)

# Correlation
correlation = np.correlate(lpimage,lptarget)