Phase Correlation in OpenCV

Nov 4, 2008 | Tags: Correlation | del.icio.us del.icio.us | digg Digg

Phase Correlation is a method to check the similarity of two images with equal size. It can be used for template matching, object tracking, motion estimation, etc. In this article, I will explain the basics and the code to perform Phase Correlation in OpenCV. The source code used in this tutorial is available for you to download. The package contains all necessary files and some installation notes.

Table of Contents:

  1. The Basics
  2. Main Code
  3. Phase Correlation Function
  4. Results
  5. Summary

1. The Basics

To obtain the Phase Correlation of two images, perform these steps:

  1. Load two images, f and g.
  2. Perform FFT on each image, resulting in F and G
  3. Obtain the cross power spectrum using this formula:
    Cross power spectrum formula
    where G bar is the complex conjugate of G.
  4. Obtain the phase correlation by performing IFFT on R

The result is a 2D array with each element has a value between 0 to 1. The location of the highest value corresponds with the object translation movement from image 1 to image 2.

2. Main Code

The main code performs a very simple task. It loads the reference and template images, call phase_correlation() function and display the result.

Note: In the code listings below I removed all error checking lines to make things simple. Remember that you should always have some error checking in your programs.

First, load the template and the reference image. Since Phase Correlation computes 2 equal images, the width and the height of both images need to be checked.

Listing 1: Load reference and template image

  1. /* load reference and template image */
  2. IplImage *ref = cvLoadImage( argv[1], CV_LOAD_IMAGE_GRAYSCALE );
  3. IplImage *tpl = cvLoadImage( argv[2], CV_LOAD_IMAGE_GRAYSCALE );
  4.  
  5. /* both images' size should be equal */
  6. if( ( tpl->width != ref->width ) || ( tpl->height != ref->height ) ) {
  7.     fprintf( stderr, "Both images must have equal width and height!\n" );
  8.     return 1;
  9. }
  10.  

Next create a new image, poc, for the Phase Correlation result. Then call phase_correlation() function with tpl, ref and poc passed as parameters.

Listing 2: Perform Phase Correlation

  1. /* create a new image, to store phase correlation result */
  2. IplImage *poc = cvCreateImage( cvSize( tpl->width, tpl->height ),
  3.                                IPL_DEPTH_64F, 1 );
  4.  
  5. /* get phase correlation of input images */
  6. phase_correlation( ref, tpl, poc );
  7.  

Now that poc contains the result from Phase Correlation, find its maximum value and the location using cvMinMaxLoc. The result then displayed in the command line.

Listing 3: Get maximum value

  1. /* find the maximum value and its location */
  2. CvPoint  minloc, maxloc;
  3. double   minval, maxval;
  4. cvMinMaxLoc( poc, &minval, &maxval, &minloc, &maxloc, 0 );
  5.  
  6. /* print it */
  7. fprintf( stdout,
  8.          "Maxval at (%d, %d) = %2.4f\n", maxloc.x, maxloc.y, maxval );
  9.  

The last step, display the images and free memory.

Listing 4: Display images and free memory

  1. /* display images and free memory */
  2. cvNamedWindow( "tpl", CV_WINDOW_AUTOSIZE );
  3. cvNamedWindow( "ref", CV_WINDOW_AUTOSIZE );    
  4. cvShowImage( "tpl", tpl );
  5. cvShowImage( "ref", ref );
  6.  
  7. cvWaitKey( 0 );
  8.  
  9. cvDestroyWindow( "tpl" );
  10. cvDestroyWindow( "ref" );    
  11. cvReleaseImage( &tpl );
  12. cvReleaseImage( &ref );
  13. cvReleaseImage( &poc );

In the next section I will explain the phase_correlation() function, line by line.

3. Phase Correlation Function

This function performs Phase Correlation operation. It takes 3 parameters:

  1. IplImage *ref: The reference image
  2. IplImage *tpl: The template image
  3. IplImage *poc: Empty image to store the result

To obtain the DFT of the images, I use the FFTW library rather than OpenCV's DFT functions. This will increase the speed since FFTW computes FFT very fast.

First get the width and the height of the image. These properties are needed for the rest of the function. Also setup the pointers to copy the images to FFTW input and vice versa.

Listing 5: Phase Correlation Function

  1. void phase_correlation( IplImage *ref, IplImage *tpl, IplImage *poc )
  2. {
  3.     int     i, j, k;
  4.     double  tmp;
  5.    
  6.     /* get image properties */
  7.     int width    = ref->width;
  8.     int height   = ref->height;
  9.     int step     = ref->widthStep;
  10.     int fft_size = width * height;
  11.  
  12.     /* setup pointers to images */
  13.     uchar   *ref_data = ( uchar* ) ref->imageData;
  14.     uchar   *tpl_data = ( uchar* ) tpl->imageData;
  15.     double  *poc_data = ( double* )poc->imageData;
  16.  

Note that the function above assumes that the three images have equal width and height. Next initialize some FFTW input/output arrays and plans, load the reference image to img1 and the template image to img2.

Listing 6: Initialize FFTW input and output arrays

  1.     /* allocate FFTW input and output arrays */
  2.     fftw_complex *img1 = ( fftw_complex* )fftw_malloc( sizeof( fftw_complex ) * width * height );
  3.     fftw_complex *img2 = ( fftw_complex* )fftw_malloc( sizeof( fftw_complex ) * width * height );
  4.     fftw_complex *res  = ( fftw_complex* )fftw_malloc( sizeof( fftw_complex ) * width * height );    
  5.    
  6.     /* setup FFTW plans */
  7.     fftw_plan fft_img1 = fftw_plan_dft_1d( width * height, img1, img1, FFTW_FORWARD,  FFTW_ESTIMATE );
  8.     fftw_plan fft_img2 = fftw_plan_dft_1d( width * height, img2, img2, FFTW_FORWARD,  FFTW_ESTIMATE );
  9.     fftw_plan ifft_res = fftw_plan_dft_1d( width * height, res,  res,  FFTW_BACKWARD, FFTW_ESTIMATE );
  10.    
  11.     /* load images' data to FFTW input */
  12.     for( i = 0, k = 0 ; i < height ; i++ ) {
  13.         for( j = 0 ; j < width ; j++, k++ ) {
  14.             img1[k][0] = ( double )ref_data[i * step + j];
  15.             img1[k][1] = 0.0;
  16.            
  17.             img2[k][0] = ( double )tpl_data[i * step + j];
  18.             img2[k][1] = 0.0;
  19.         }
  20.     }
  21.  

In the listing above res will hold the end result of the phase correlation operation. Next is the heart of this function. It performs Phase Correlation and store the result to poc.

Listing 7: Compute phase correlation

  1.     /* obtain the FFT of img1 */
  2.     fftw_execute( fft_img1 );
  3.    
  4.     /* obtain the FFT of img2 */
  5.     fftw_execute( fft_img2 );
  6.    
  7.     /* obtain the cross power spectrum */
  8.     for( i = 0; i < fft_size ; i++ ) {
  9.         res[i][0] = ( img2[i][0] * img1[i][0] ) - ( img2[i][1] * ( -img1[i][1] ) );
  10.         res[i][1] = ( img2[i][0] * ( -img1[i][1] ) ) + ( img2[i][1] * img1[i][0] );
  11.  
  12.         tmp = sqrt( pow( res[i][0], 2.0 ) + pow( res[i][1], 2.0 ) );
  13.  
  14.         res[i][0] /= tmp;
  15.         res[i][1] /= tmp;
  16.     }
  17.  
  18.     /* obtain the phase correlation array */
  19.     fftw_execute(ifft_res);
  20.  
  21.     /* normalize and copy to result image */
  22.     for( i = 0 ; i < fft_size ; i++ ) {
  23.         poc_data[i] = res[i][0] / ( double )fft_size;
  24.     }
  25.  

Note that we have to normalize the IFFT result since FFTW computes an unnormalized DFT. Now all we have to do is just free memory and return to main function.

Listing 8: Free Memory

  1.     /* deallocate FFTW arrays and plans */
  2.     fftw_destroy_plan( fft_img1 );
  3.     fftw_destroy_plan( fft_img2 );
  4.     fftw_destroy_plan( ifft_res );
  5.     fftw_free( img1 );
  6.     fftw_free( img2 );
  7.     fftw_free( res );
  8. }

The Phase Correlation result is stored in poc. The main code should find the maximum value and its location from this result.

4. Results

Here are some results I've got when computing phase correlation using this code.

Table 1. Phase Correlation Results

# Image 1 Image 2 Result
1. image 1 image 2
(Image 2 is identical
with image 1)		 Max. value = 1.0
in (0,0)
2. image 1 image 2
(image 1 with some noise)		 Max. value = 0.2729
in (0,0)
3. image 1 image 2
(object moved 6 pixels to
the right and 15 pixels to
the bottom)		 Max. value = 0.6344
in (6,15)

The first result shows the maximum value is 1.0 at (0, 0). This means both images are identical.

The second result shows the maximum value is 0.2729 at (0, 0). This means there is no translation from image 1 to image 2, but image 2 has some noise in it.

The last result shows the maximum value is 0.6344 at (6, 15). This means image 1 has a translative movement 6 pixels to the right and 15 pixels to the bottom.

5. Summary

This article covered Phase Correlation function in OpenCV. The code makes use FFTW library to obtain the DFTs rather than using OpenCV's DFT functions. This will increase the speed since FFTW computes FFT very fast.

The source code in this tutorial is available for you to download. The package contains all necessary files and some installation notes.

Send your comments, questions and bug reports to me [at] nashruddin.com.

Related Articles

Recommended Book

book cover Learning OpenCV: Computer Vision with the OpenCV Library
By: Gary Bradski, Adrian Kaehler
This book is the "de facto" OpenCV User's Manual. It provides a practical, pragmatic, accessible book on computer vision, with algorithmic explanation and concrete example code snippets. Written by the creators of OpenCV, no doubt you should obtain a copy.

The Downloads

22 Comments

diyana on Sep 10, 2008:

very cool website! I like it...actually, I'm going to study openCV..so this blog is helpful.

Marcos Rodrigues on Jan 14, 2009:

Congratulations on your web site, very useful indeed. I need to search for a template in a larger image but it needs to be invariant to scale, rotation and translation (your example shows invariance to translation). I wonder if using fft one would be able to achieve scale and rotation invariance? Have you tried this?

Thanks in anticipation.

Nash on Jan 14, 2009:

Check this one out:
http://www.lps.usp.br/~hae/software/cirateg/index.html

ABHAY SHANKAR on Feb 26, 2009:

I have some queries ..........

1.Is FFTW similar to OpenCV or is it different from it?

2.What is the procedure to use rather install FFTW and will it effect the OpenCV already installed ?

3.After installing do we need to change ant settings or provide a path to use FFTW functions or will they be taken automatically?

4.Are the functions included in FFTW also included in OpenCV ?

Nash on Feb 26, 2009:

Hi Abhay,

1. FFTW is a library of fast C routines for computing Discrete Fourier Transforms (DFT). It's not same with OpenCV. Note that OpenCV also has functions to compute DFTs.

2. Just extract the package into a directory and add the directory to your system path. It doesn't affect OpenCV installation.

3. Yes, add the directory to your system path. See the Makefile included in the download above for an example on how compiling fftw programs.

4. OpenCV also has functions to compute DFTs, but FFTW performs alot faster.

Asad on Feb 27, 2009:

Hi,
Asad here i need ur help.
i sent u a mail.
plz chk it.

jane_doe on Feb 27, 2009:

how about more objects... how to?

ABHAY SHANKAR on Mar 2, 2009:

Hi Nash,
I tried to compile the phase correlation code in my Linux machine.
The code runs fine without any error but it is not showing the correct output. When Iam comparing two similar images the maxvalue output comes 0.00 instead of 1.00. Same is the case with other images also. Everytime the maxvalue output is coming 0.00.
What could be the problem????
Can u plzz explain where am I making mistake??

Nash on Mar 2, 2009:

It's very hard for me to say what's wrong, since everything works fine in my system. Note that the Makefile in the package above for use in Windows. I'd like to see your Makefile, with the compilation messages if any.

I also suggest you to try compiling and running some simple programs first:
OpenCV Examples 1
FFTW Example
FFTW with OpenCV

If you have successfully compiling and running those samples, I guarantee 99% you'll have no problem running my phase correlation code

Chandra Praksh on Mar 2, 2009:

Hi Nash,

I am compiling phase correlation function on my linux machine. I dont know how to generate a makefile so i combine both the phase_correlation.c and main.c togather. Then i compile it by..

gcc `pkg-config --cflags opencv` -o fftwmerge fftwmerge.c `pkg-config --libs opencv` -lfftw3

and execute it by the below..

./fftwmerge1 0201_1.jpg 0201_1.jpg

and it returns..
Maxval at (0, 0) = 0.0000

it is giving the same output for all the images which u have tested.

how i can solve it.

ABHAY SHANKAR on Mar 2, 2009:

Nash,
Thankx for ur reply but i have'nt used any makefile. I have merged both the main code and the phase correlation code together. After doing so I have compiled the program.
I think that makefile is only needed when both the codes are separately used . But I have used them together so there is no need of makefile.
If I am doing it wrong or if I have to follow some other steps after combining the two codes together and then compiling then plz tell me.

Nash on Mar 3, 2009:

Ok guys, I'll try to run it in Linux and figure it out. Meanwhile, just compile and run the code in Windows. It works.

Deb on Apr 8, 2009:

Hi, when I run mingw32-make I get this error:

C:\Phase_Correlation>mingw32-make
gcc main.o phase_correlation.o -L"C:\OpenCV\lib" -L"C:\OpenCV\fftw" -lcxcore -lcv -lcvaux -lhighgui -lml -lcvcam -lfftw3-3 -o phase_correlation
C:\MinGW\bin\..\lib\gcc\mingw32\3.4.5\..\..\..\..\mingw32\bin\ld.exe: cannot find -lcvcam
collect2: ld returned 1 exit status
mingw32-make: *** [all] Error 1

lcvcam isn't where it's supposed to be... where can I get it?If I erase lcvcam from Makefile doesn't work

Nash on Apr 9, 2009:

@Chandra
@Abhay

I've merge the files to make it simpler. See the full listing of phase_correlation.c. I've compile and test the code in my FreeBSD box and it runs very well.

$ gcc phase_correlation.c -o phase_correlation `pkg-config --cflags opencv` `pkg-config --cflags fftw` `pkg-config --libs opencv` `pkg-config --libs fftw` -lcxcore -lcv -lhighgui -lfftw3

$ ./phase_correlation 0201_1.jpg 0201_1.jpg
Maxval at (0, 0) = 1.0000

$ ./phase_correlation 0201_4.jpg 0201_5.jpg
Maxval at (0, 0) = 0.2729

$ ./phase_correlation 0201_2.jpg 0201_3.jpg
Maxval at (6, 15) = 0.6344

Nash on Apr 9, 2009:

@Deb

Yes, lcvcam should be removed since the code doesn't use it. I've rewrite the code to make it simpler. Please download the package above and read how to compile in the README file.

mclish on Apr 14, 2009:

what the location of peak is always be positive .
what about the image with diplcement to left and up ???

Guenter Sandner on Apr 23, 2009:

just add the following to lines of code:

if ( maxloc.x >= width/2 )
    maxloc.x -= width;
if ( maxloc.y >= height/2 )
    maxloc.y -= height;

Ed Elston on Aug 6, 2009:

Nash, very nice explanation and sample code. I have a couple of questions before starting with this:

What would be the expected result if you have more than 1 (similar or dissimilar) object in the image? Could this be dealt with using sub-images to avoid the problem?

Thanks, Ed

Elmond on Aug 13, 2009:

Hello there,

When will we actually see a real 2D FFT? or by theory is this how it is being done? How can FFTW be used to show the spectrum of an image using opencv?

Kind regards, and thank you very much!!!
Elmond

Bhavik on Apr 1, 2010:

very nice explanation and sample code.................
but,i need to find the next step i.e. Co-efficient of correlation
based on this value , i will be using threshold and then to decide whether to reject or accept two images

stn on May 4, 2010:

hi
1. why do you use 1D Fourier Transform instead of 2D Fourier Transform of an image?

thanks

Nash on May 4, 2010:

Simply because it is much simpler.

Leave a comment

Name (required)
Email (will not be published) (required)
Website

Characters left = 1000

Tags

Correlation  HTTP  Image  Networking  OpenCV  PHP  Projects  Regex  SEO  Security  WebService 

Recent Posts

  1. OpenCV Utility: Reading Image Pixels Value
  2. OpenCV Circular ROI
  3. OpenCV 2.0 Installation on Windows XP and Visual Studio 2008
  4. Runtime ROI Selection using Mouse
  5. Real Time Eye Tracking and Blink Detection
View Archives

About the Author

avatar Cool PHP programmer writing cool PHP scripts. Feel free to contact
Tel. +62 31 8662872
+62 856 338 6017
ICQ 489571630
Skype dede_bl4ckheart
Yahoo dede_bl4ckheart
Google nashruddin.amin

Recommended Sites:

Hacker's HTTP Client
HTML and CSS Tutorials
Stop Dreaming Start Action
Online Quran and Translation