Thread Subject: Register one image with another that has been both translated and

Hello,

I am having some difficulties with image registration, and was hoping
to get some help. My application involves a camera that takes images
every so often, the objects within them move uniformly (you can think
of this as the camera is moving over stationary objects). The
difference between the images will be BOTH a translation and rotation
of the entire image, with potentially large portions of one image
being outside that of the other (a significant overlap, but still
should be plenty to properly register the images). I need a method
that can extract the translation (in pixels) and rotation angle
required to register the images quickly.

I have found methods that use the FFT, and want to use them. For
example, the below code can successfully find the pixel translation
required to register two grayscale images that only contain objects
that are translated (no rotation from one image to the next):

% Take FFT of each image
F1 = fft2(im1);
F2 = fft2(im2);
% Create phase difference matrix using phase angles from the two FFTs
(i.e., phase correlation)
pdm = exp(i*(angle(F1)-angle(F2)));
% Solve for phase correlation function. The peak of this function
% corresponds to the registration point for translation.
pcf = real(ifft2(pdm));
figure; mesh(pcf);
% Find the linear index of where the maximum occurs, then convert that
into
% the 2D image subscripts
linear_index = find((pcf == (max(max(pcf)))));
[X, Y] = ind2sub(size(pcf),linear_index);
% Adjust subcripts to be correct offsets from (0,0) in pixels (should
return 0,0 if images are the same)
X = min(size(pcf,1)-X, X) - 1;
Y = -1*(min(size(pcf,2)-Y, Y) - 1);

I know that to find both the translation and the rotation, that I
first need to find the rotation. After finding the rotation, I rotate
the first image by that rotation and then use the above code to find
the translation between the images. I would also like to use an FFT-
based approach here too. I have found papers stating that the polar
transform (or log-polar transform) of each image's magnitude spectrum
can lead to finding out the rotation. I have a desire to use the FFT-
based approaches because they are faster than other regression or
least-squares approaches I have found.

Could someone please help me with the Matlab code for finding the
rotation angle between a pair of images, where the second image has
been translated and rotated (likely with large portions lying outside
of the image boundary after registration)? Any help is appreciated.
Thanks in advance.

> I have found papers stating that the polar
> transform (or log-polar transform) of each image's magnitude spectrum
> can lead to finding out the rotation

I wouldn't mind taking a look at those papers. It seems pretty counter-intuitive that analyzing the Fourier domain would help you determine rotation.

Rotation in the image causes a corresponding rotation of the spectrum, so you end up with the same analysis problem for the spectrum as you did for the original image.

Anyway, if the papers prove these methods work, what help do you really need? Why not just code them up?

On Sep 23, 9:51=A0am, "Matt " <mjacobson.removet...@xorantech.com>
wrote:
> > I have found papers stating that the polar
> > transform (or log-polar transform) of each image's magnitude spectrum
> > can lead to finding out the rotation
>
> I wouldn't mind taking a look at those papers. It seems pretty counter-in=
tuitive that analyzing the Fourier domain would help you determine rotation=
.
>
> Rotation in the image causes a corresponding rotation of the spectrum, so=
 you end up with the same analysis problem for the spectrum as you did for =
the original image.
>
> Anyway, if the papers prove these methods work, what help do you really n=
eed? Why not just code them up?

I have attempted to implement what the papers suggest, but they are
mostly specialized for their application, and they do not cover some
of the details that would be required to replicate the behavior. I
wanted to see if anyone has had any experience doing what I am trying
to do to help suggest how to best implement it or provide a code
segment. Or, if something else would work better or faster, to suggest
that. I was hoping to find Matlab code similar to what I posted for
determining translation because I can follow it and intuitively know
why it works. Thus, I am looking to see if anyone is willing to
explain or provide code segments on how to determine the amount of
rotation from one image to another that also has translation. Any help
is appreciated. Thanks!

On Sep 23, 9:51=A0am, "Matt " <mjacobson.removet...@xorantech.com>
wrote:
> > I have found papers stating that the polar
> > transform (or log-polar transform) of each image's magnitude spectrum
> > can lead to finding out the rotation
>
> I wouldn't mind taking a look at those papers. It seems pretty counter-in=
tuitive that analyzing the Fourier domain would help you determine rotation=
.
>
> Rotation in the image causes a corresponding rotation of the spectrum, so=
 you end up with the same analysis problem for the spectrum as you did for =
the original image.
>
> Anyway, if the papers prove these methods work, what help do you really n=
eed? Why not just code them up?

I have attempted to implement what the papers suggest, but they are
mostly specialized solutions and the details to fully reproduce the
behavior are not provided. I was thus hoping that someone may have
experience doing what I am trying to do, and would be willing to
provide code samples or examples. The code I listed in my original
message is great because I can follow it and understand why it works.
Something like this is what I was wanting to determine the rotation
angle between two images where one is both translated and rotated.
Once I find out the rotation, I can find out the translation using the
code I listed in the original post.

In particular, the following two example images are those which I am
trying to register to make sure things work before moving to more
complex images. im2.png has a 36 degree rotation and a translation
applied. What I want is something that can automatically extract the
36 degrees for me.
* https://www.cresis.ku.edu/~cgifford/access/im1.png
* https://www.cresis.ku.edu/~cgifford/access/im2.png

Any help is appreciated. Thanks!

If you have the image processing toolbox, you might try something like

theta=linspace(0,180,N);
Y1=fft(radon(im1),theta);
Y2=fft(radon(im2),theta);

the columns of Y1 and Y2 will be polar samples of the Fourier spectrum of im1 and im2 at angles theta. (Is this what you mean by "polar transform"?)

Anyway, information about the rotation is definitely encapsulated in Y1 and Y2, so comparing them in some way is probably a good idea.

On Sep 23, 3:40=A0pm, "Matt " <mjacobson.removet...@xorantech.com>
wrote:
> If you have the image processing toolbox, you might try something like
>
> theta=3Dlinspace(0,180,N);
> Y1=3Dfft(radon(im1),theta);
> Y2=3Dfft(radon(im2),theta);
>
> the columns of Y1 and Y2 will be polar samples of the Fourier spectrum of=
 im1 and im2 at angles theta. (Is this what you mean by "polar transform"?)
>
> Anyway, information about the rotation is definitely encapsulated in Y1 a=
nd Y2, so comparing them in some way is probably a good idea.

Matt,

Thank you for that tip. I do have the image processing toolbox. I will
play around with Y1 and Y2 a bit and see what I can come up with in
terms of how to compare them to extract the rotation angle. Any other
ideas are also welcome. I will report back. Thanks!

> I know that to find both the translation and the rotation, that I
> first need to find the rotation. After finding the rotation, I rotate
> the first image by that rotation and then use the above code to find
> the translation between the images.

That may be the order that appeals to you for some reason, but it's not obviously the necessary one. Why not find the translation first, for example, by comparing the positions of the centroids of the two images? Then, you could unshift first and determine rotation second

>I have a desire to use the FFT-
> based approaches because they are faster than other regression or
> least-squares approaches I have found.

That's understandable, but my suspicion is that the FFT approaches will be more sensitive to the field of view truncation than those other, slower methods.

It might be worth doing a fast and coarse FFT-based registration first and then clean it up with something like Sum of Squared Differences (SSD).

On Sep 23, 4:39=A0pm, "Matt " <mjacobson.removet...@xorantech.com>
wrote:
> >I have a desire to use the FFT-
> > based approaches because they are faster than other regression or
> > least-squares approaches I have found.
>
> That's understandable, but my suspicion is that the FFT approaches will b=
e more sensitive to the field of view truncation than those other, slower m=
ethods.
>
> It might be worth doing a fast and coarse FFT-basedregistrationfirst and =
then clean it up with something like Sum of Squared Differences (SSD).

That is probably true, but the application requires speed
unfortunately.

I played around with the code for a while and cannot seem to get it to
spit out the correct rotation angle for the two images I linked to in
a previous post. The code I used is below:

% First image
im1 =3D imread('im1.png');
im1c =3D im1;
im1 =3D double(im1);
im1 =3D im1/255;
im1 =3D im1(:,:,1);
% Second image
im2 =3D imread('im2.png');
im2 =3D double(im2);
im2 =3D im2/255;
im2 =3D im2(:,:,1);
% Columns contain polar samples from Fourier spectrum for thetas in
range [0,180]
theta=3Dlinspace(0,180,181);
Y1=3Dfft(radon(im1,theta));
Y2=3Dfft(radon(im2,theta));
figure; imshow(log(1+real(Y1)));
figure; imshow(log(1+real(Y2)));
% Correlate the two images to find the maximum
c =3D sum(fft2(Y1).*fft2(Y2));
p =3D real(ifft(c));
figure; plot(p);
[max_val, rot_angle] =3D max(p);

I am expecting that the rotations that will be experienced between
frames will be in the range of -90 to 90 degrees. The above code
returns 125 for rot_angle. I can visually see a shift between Y1 and
Y2, so how can I use them to properly correlate to get the rotation
angle back out in the -90 to 90 range? What is the most efficient
Matlab way of doing this? Thanks again for your help. It is much
appreciated.

On Sep 23, 8:48=A0pm, ChrisGifford <chris.m.giff...@gmail.com> wrote:
> On Sep 23, 4:39=A0pm, "Matt " <mjacobson.removet...@xorantech.com>
> wrote:
>
> > >I have a desire to use the FFT-
> > > based approaches because they are faster than other regression or
> > > least-squares approaches I have found.
>
> > That's understandable, but my suspicion is that the FFT approaches will=
 be more sensitive to the field of view truncation than those other, slower=
 methods.
>
> > It might be worth doing a fast and coarse FFT-basedregistrationfirst an=
d then clean it up with something like Sum of Squared Differences (SSD).
>
> That is probably true, but the application requires speed
> unfortunately.
>
> I played around with the code for a while and cannot seem to get it to
> spit out the correct rotation angle for the two images I linked to in
> a previous post. The code I used is below:
>
> % First image
> im1 =3D imread('im1.png');
> im1c =3D im1;
> im1 =3D double(im1);
> im1 =3D im1/255;
> im1 =3D im1(:,:,1);
> % Second image
> im2 =3D imread('im2.png');
> im2 =3D double(im2);
> im2 =3D im2/255;
> im2 =3D im2(:,:,1);
> % Columns contain polar samples from Fourier spectrum for thetas in
> range [0,180]
> theta=3Dlinspace(0,180,181);
> Y1=3Dfft(radon(im1,theta));
> Y2=3Dfft(radon(im2,theta));
> figure; imshow(log(1+real(Y1)));
> figure; imshow(log(1+real(Y2)));
> % Correlate the two images to find the maximum
> c =3D sum(fft2(Y1).*fft2(Y2));
> p =3D real(ifft(c));
> figure; plot(p);
> [max_val, rot_angle] =3D max(p);
>
> I am expecting that the rotations that will be experienced between
> frames will be in the range of -90 to 90 degrees. The above code
> returns 125 for rot_angle. I can visually see a shift between Y1 and
> Y2, so how can I use them to properly correlate to get the rotation
> angle back out in the -90 to 90 range? What is the most efficient
> Matlab way of doing this? Thanks again for your help. It is much
> appreciated.

Turns out that I forgot that the magnitude of the of the polar samples
of the Fourier spectrum is important for determining rotation. Having
the phase present causes the correlation to be higher in places where
a match really isn't present in the spatial domain. If I do that, it
appears to be working for the most part. Thanks for the discussion,
and definitely let me know if anyone dreams up an even more efficient
way to do this.

> > % Correlate the two images to find the maximum
> > c =3D sum(fft2(Y1).*fft2(Y2));
> > p =3D real(ifft(c));
> > figure; plot(p);
> > [max_val, rot_angle] =3D max(p);


Just one more thing. The shift you see between abs(Y1) and abs(Y2) should be a 1 dimensional shift along the row axis. So, you really should only need to correlate row-wise as opposed to a fully 2D correlation operation.

So, instead of fft2(abs(Y)), this can probably be done using fft(abs(Y),[],2) to take the FFT's row-wise, which should save you some CPU time.

ChrisGifford <chris.m.gifford@gmail.com> wrote in message <322bf070-5261-4e98-a321-7922681d8e8e@a1g2000hsb.googlegroups.com>...
> On Sep 23, 3:40=A0pm, "Matt " <mjacobson.removet...@xorantech.com>
> wrote:
> > If you have the image processing toolbox, you might try something like
> >
> > theta=3Dlinspace(0,180,N);
> > Y1=3Dfft(radon(im1),theta);
> > Y2=3Dfft(radon(im2),theta);
> >
> > the columns of Y1 and Y2 will be polar samples of the Fourier spectrum of=
> im1 and im2 at angles theta. (Is this what you mean by "polar transform"?)
> >
> > Anyway, information about the rotation is definitely encapsulated in Y1 a=
> nd Y2, so comparing them in some way is probably a good idea.
>
> Matt,
>
> Thank you for that tip. I do have the image processing toolbox. I will
> play around with Y1 and Y2 a bit and see what I can come up with in
> terms of how to compare them to extract the rotation angle. Any other
> ideas are also welcome. I will report back. Thanks!

I had a similar problem many years ago, when I was inspecting 'credit' cards coming through a vision system on a conveyor belt, and had to do a template match to pick up printing flaws. The method outlined here worked well for me, however without seeing examples of your images I couldn't say whether the approach would be suitable for your application.

On both your images, you need to convert it from an intensity image into an angle image (i.e. where the intensity of each pixel is a measure of the gradient direction at that point). I believe that there is an option in Matlab for calculating such an image directly, but as I was using C++ at the time, I did it using the ATAN2 function applied to the x and y gradients at each pixel position. You will find that the initial image will look meaningless, as the gradient in uniform regions is solely dictated by the noise. My solution to this was to mask the resulting angle image by a dilated edge image, found using something like the Sobel operator. The resulting image is then effectively the gradient angle of the image edges.

I then made a histogramme of the masked angle images, and the angle of rotation could be clearly seen in the shift in the histogramme shape, that can easily be quantified by undertaking a circular correlation of the two histogrammes. Once you know the relative rotation angle, you can apply a compensating twist to one of the images using the rotate() function.

Now you have two images with no rotational shift, just translational shift, which can be determined using a 2D correlation process.

Sorry I don't have any Matlab code, as I don't have the Image Processing toolbox on this machine.

Hope this helps

Regards

Dave Robinson

One method that I've had success with for rotation and translation
image registration is to use the Fourier-Mellin transform. If you
search around for it there is some Matlab code out there to
demonstrate the effectiveness of the algorithm. The demo code had a
few bugs but they should be pretty easy to figure out.

Hello again. I thought I would keep the below questions in this thread
since they involve this application. Below is the code that properly
finds the rotation angle between two frames.

% Only pay attention to 0 to 180 degrees of FFT
thetas = linspace(0, 180, 180);

% Polar transform the images and get the magnitude information, which
should
% only be different by a translation
F1 = abs(fft(radon(frame1, thetas)));
F2 = abs(fft(radon(frame2, thetas)));

% Correlate the magnitude components by doing the FFT2 of each image
and
% point-wise multiplying them and summing the columns. The IFFT can
then be
% used to find the maximum index, which represents the
% translation-invariant angle of rotation between the frames.
correlation = sum(fft2(F1) .* fft2(F2));
peaks = real(ifft(correlation));
angle = max(peaks);

Could someone please provide a more in-depth explanation of what the
code is doing and why it works? For example, how does the 1D FFT and
radon transform combine to do a polar transform? The correlation
section appears to be something like the Fourier Slice Theorem, but
how and why does that code provide the peak information for finding
the proper translation. Any help on the explanation is greatly
appreciated. If you notice any instances where the process could be
done more efficiently, please let me know. Thanks ahead of time!