Resamples an image from a conventional grid to a log-polar grid, and back.
Images sampled on a log-polar grid have interesting and useful properties. Rotation and scaling are trivial operations, and this can lead to efficient algorithms for straight line detection and optic flow estimation. Such sampling schemes also provide a rough model of biological foveal vision, of interest in the development of active computer vision systems.
The functions in this file carry out resampling from a conventional image to a log-polar image, and back. There is usually considerable information loss in each direction, but resampling to log-polar can still be useful for computational experiments.
The toolbox function imtransform does the main work; these functions are a front end to it.
David Young (view profile)
Hi Bagus, Yes, you can do template matching on log-polar images. It's a trade off - you get matching across scale and rotation, but you lose translation invariance, so you need to use it as part of a system which can control the fixation point. See for example, http://homepages.inf.ed.ac.uk/rbf/CVonline/LOCAL_COPIES/YOUNG2/
Bagus Adhi Kusuma (view profile)
Hi David, whether polar able to measure the similarity of the two images as template matching?
David Young (view profile)
Hi Angel, thank you for your comment. I'd need more information to work out what is wrong - at least the full error message including the traceback, and also the type of the array you are passing.
Angel Lindo (view profile)
Hi David, the code is great, but it throws me an error:
'XData' and 'YData' must be real, double, 1-by-2 vectors.
I have tried to fix it, but I can't.
---Angel Lindo---
David Young (view profile)
dong yun: Thank you for your comment. The "circular pixels" condition is the condition for the radial separation between the centres of neighbouring pixels to be the same as the tangential separation. That is, a neighbouring pixel along a wedge is the same distance away as a neighbouring pixel round a ring.
If the condition is not met, the local aspect ratio of a pixel is not unity - that is, the pixel may represent a region which is elongated in either the ring or wedge direction.
Whether this matters or not depends on your application.
dong yun (view profile)
Hi David,mant thanks your code.But I have a probelm: what's the 'Circular pixels condition'? can you give me a detailed description?if we should not explicitly specify those parameters rmin,rmax,nr,nw? if these paramers doesn't statisfy the 'circular pixels condition',what effects it result in?
----Yunyun Dong
AMAL DUROU (view profile)
Phillip (view profile)
Works as advertised, good job. Am using this to recover rotation/scale for image registration.
David Young (view profile)
Hi Daniel, The differences from the original are almost inevitable, because the log-polar transform is usually lossy. Unless you use a huge output array, the pixels in the outer rings will be further apart than the original image pixels. When you transform back, there's no way to recover the lost information, and so the resolution suffers. Antialiasing won't be able to help fundamentally.
Both functions use *imtransform* with default interpolation. (I ought to update to *imwarp* sometime.) It might be possible to get some improvements in quality by changing the interpolation rule - you'd need to learn about *makeresampler* to generate a resampler to give to imtransform. It still will not restore the original image though.
Daniel (view profile)
This is a great start for me. I'm trying to model foveal vision, but when I use the logsampback function the output image will have significant variations from the original that look like aliasing. Do you know of a way to add anti-aliasing into the logsample function?
David Young (view profile)
ANCY: One of your comments is repeated 4 times and does not seem to relate to this submission. Please would you delete it; you could try putting it on Matlab Answers.
On your most recent comment, I am not sure what you mean by "full sampling". To do log-polar sampling and its inverse, you just call these functions like any other functions. They are fully described in the help information. An example of their use is at http://www.sussex.ac.uk/Users/davidy/log_lines/index.html
ANCY (view profile)
what is the function of logsample function.please share the code for full sampling
ANCY (view profile)
If any one know the code for segmentation based sampling,grid based sampling,keypoint based sampling .then please share it
ANCY (view profile)
If any one know the code for segmentation based sampling,grid based sampling,keypoint based sampling .then please share it
ANCY (view profile)
If any one know the code for segmentation based sampling,grid based sampling,keypoint based sampling .then please share it
ANCY (view profile)
If any one know the code for segmentation based sampling,grid based sampling,keypoint based sampling .then please share it
Gomathy (view profile)
Hi Dave, a quick update. The log-sample routine works fine for phase correlation too. Thanks,G
David Young (view profile)
Hi Gomathy, I'm sorry, but I don't know very much about phase correlation, and not enough to help with your code I'm afraid.
Gomathy (view profile)
Hi David, the correct LP parameters for the earlier post is
rmin = 5;
rmax = 125;
nw = 500;
I tried various ways to measure theta using Phase Correlation. It doesn't work. Can you please look at the test code ?
Gomathy (view profile)
Hi David, This is Gomathy again. I am trying to use Phase correlation to measure displacements after converting the image into LP. It fails to detect displacement, whereas the same LP converted image works fine for normxcorr2. I am not sure, where am i making mistakes. Here is the test code. I would greatly appreciate, if you can you can please help.
Iref = imread('cameraman.tif');
deg = 73;
I = imrotate(Iref, deg, 'bicubic', 'crop');
rmin = 5;
rmax = 300;
nw = 1500;
xc = size(Iref, 2)/2;
yc = size(Iref, 1)/2;
lp_Iref = logsample(Iref, rmin, rmax, xc, yc, [], nw);
lp_I = logsample(I, rmin, rmax, xc, yc, [], nw);
THETA_F1 =fftshift(fft2(lp_Iref));
THETA_F2 =fftshift(fft2(lp_I));
a1 = angle(THETA_F1); a2 = angle(THETA_F2);
THETA_CROSS = exp(-1* (a1 - a2)); % Phase Difference
THETA_PHASE = abs(fftshift(ifft2(THETA_CROSS)));
[peak, idx] = max(THETA_PHASE (:));
[ypeak,xpeak] = ind2sub(size(THETA_PHASE ),idx);
theta=xpeak-(size(lp_I,2)+1)/2;
Michael Meyers (view profile)
Hello David
How would your code be modified to include calculating the scale from the xpeak? I can't find any examples of that calculation on the web. I understand the angle as the y axis covers 2pi but what does the xaxis cover in terms of scale?
Bharath Ramesh (view profile)
Hello David,
Do you know how to find the image point corresponding a point in the log-polar transformed space? Is there an easy way to do it?
Janaka (view profile)
Hi David
Thanks for your submission. I have been following your answers and submissions and they all are very useful. Greatly appreciate your work.
I am trying to use imrotate for checking the rotation between two image frames (1024x1024 int16) coming from HH and HV polarized antennas of a marine radar. They are 16-bit data (comes as Matlab matrices) and the HV is rotated w.r.t HH. I tried to use the above piece of code but first I tried with the rotation of HH w.r.t. itself but drot always comes as 0. Is it due to the data? Should this valid for int8 only? Here is the modified code I used. Greatly appreciate your reply.
HHdble=double(HH);
HVdble=double(HV);
Iref = HHdble;
deg = 10;
I = imrotate(Iref, deg, 'bilinear', 'crop');
rmin = 5;
rmax = 512; %%125;
nw = 512;
xc = size(Iref, 2)/2;
yc = size(Iref, 1)/2;
lp_Iref = logsample(Iref, rmin, rmax, xc, yc, [], nw);
lp_I = logsample(I, rmin, rmax, xc, yc, [], nw);
cc = normxcorr2(lp_I, lp_Iref);
[max_cc, imax] = max(cc(:));
[ypeak, xpeak] = ind2sub(size(cc), imax)
drot = 360 * (ypeak - nw) / nw %% comes as '0' .
Gomathy (view profile)
Thanks David.
David Young (view profile)
Hi Gomathy. Measuring small rotations is difficult, but you could try my affine optic flow estimator at http://www.mathworks.co.uk/matlabcentral/fileexchange/27093-affine-optic-flow If you have difficulty with it, please let me know.
Gomathy (view profile)
Hi David, Thanks for your reply and explaining about estimation of deg from ypeak. The medical image sequences, what i am dealing with, have rotations less than 1degree. Even, 0.05,.. Though it is a small value, it turn up as a swing in registered video seq. I tried couple of other methods. But, Nothing works for less than a degree. Now, Log-polar also fails. Do you have any other recommendation for < 1deg? Thanks agian.
David Young (view profile)
Hi Gomathy. Your code is almost OK. One problem is that you have chosen very small rotations to test, so it is hard to see if the results are right. Instead of fractions of a degree, try using many degrees, so that you can tell if there's a big error or just a small inaccuracy.
The second problem is that your conversion from the y-coordinate of the correlation peak to an angle is not correct. Remember that the number of wedges is just a full circle, so in your case a change of 500 in the y-coordinate corresponds to 360 degrees or 2*pi radians. If you use this, then the code works.
Here is my test code, based on yours, but simplified to make it easier to analyse.
Iref = imread('cameraman.tif');
deg = 73;
I = imrotate(Iref, deg, 'bilinear', 'crop');
rmin = 5;
rmax = 125;
nw = 500;
xc = size(Iref, 2)/2;
yc = size(Iref, 1)/2;
lp_Iref = logsample(Iref, rmin, rmax, xc, yc, [], nw);
lp_I = logsample(I, rmin, rmax, xc, yc, [], nw);
cc = normxcorr2(lp_I, lp_Iref);
[max_cc, imax] = max(cc(:));
[ypeak, xpeak] = ind2sub(size(cc), imax)
drot = 360 * (ypeak - nw) / nw % prints 72.72 for initial deg=73
The choice of the parameters is up to you - you seem to have chosen good ones for this problem.
Gomathy (view profile)
Hi, i am using Log-Polar for measuring and correcting the image rotation. For test purpose, cameraman.tif is rotated to known value. And using log-polar , we like to measure the rotation and de-rotate to end. It is too difficult to find the suitable parameters for rmin/rmax/nr/nw. Below code doesn't measure the rotation. Can you please help to fix this issue.
***************
I(:,:,1)=imread('cameraman.tif');
deg=[0 0.1,-0.2,0.3,-0.4,0.5,-0.6,0.7,-0.8,0.9];
rad=deg*(pi/180);
for i=2:10 I(:,:,i)=imrotate(I(:,:,1),-deg(i),'bilinear','crop');
end
[h w]=size(I(:,:,1));
lp_Iref = logsample(I(:,:,1), 5, 125, w/2, h/2, [], 500); %mri
for i=2:10
clear lp_I ;
lp_I = logsample(I(:,:,i),5, 125, w/2, h/2, [], 500);
clear cc;
cc = normxcorr2(lp_I,lp_Iref);
[max_cc, imax] = max((cc(:)));
[ypeak(i),xpeak(i)] = ind2sub(size(cc), imax(1));
xpeak(i)=xpeak(i)+xx1;
ypeak(i)=ypeak(i)+yy1;
drot_r(i)=size(lp_Iref,1)-ypeak(i);
dscale(i)=size(lp_Iref,2)-xpeak(i);
end
drot_d=drot_r*(180/pi);
*************
The measured rotation in drot_d is not same as deg. It would be great if you can help to fix this.
Regards
David Young (view profile)
Yeoh cs: There are no general rules for setting those parameters. It depends on the original image and the purpose of doing the transform.
rmax determines the "field of view" of the transform. That is, it gives the radius of the region of the original image that is sampled.
nw determines the resolution at the outside ring. If you want to capture the full detail of the original image, with one transform pixel for each original pixel, you need to make nw equal to 2*pi*rmax. However, this results in massive oversampling nearer the centre, so you may want a smaller value.
I normally set only one of rmin and nr, and allow the other to default, so that the pixels in the sampled image have an aspect ratio of 1. There is always a "hole" at the centre of a log-polar image, and rmin determines how big this hole is. Near the centre of the log-polar image the original image is grossly oversampled, and the smaller rmin is, the worse this is.
In the end though, you have to look at what you are using the log-polar image for, and use either theoretical arguments or empirical tests to decide the parameters.
Yeoh cs (view profile)
Young: thx for ur advice, 1 more question .... how to defined the rmin /rmax /nr /nw value of different pattern ?
David Young (view profile)
Yeoh cs: I'm sorry, but I have not written a tutorial on the use of this transform, and there isn't room to write one here. There are many papers on applications of log-polar transforms, and it would be best to start with those and then ask specific questions on MATLAB Answers.
Yeoh cs (view profile)
can this code use to do the recognition? how do i compare the 2 image using function of logsample? i m new to log polar transform. i wish to use log polar trans to do the object recognition.
can u giv me an example that how to compare 2 logimage after going through the logsample function.. Thanks
David Young (view profile)
Masthan Zayid: Sorry, I don't know what the log polar gabor transform is.
masthan zayid (view profile)
hi...i need code for 2D log polar gabor transform for my project....
David Young (view profile)
Aravind: The words and equations in the documentation give the meanings of the arguments as clearly and precisely as I can - but it's true there isn't a tutorial there on the log-polar transform. Perhaps you need to look at the literature first to make sure you're familiar with the main ideas - perhaps start with Weiman & Chaikin's 1979 paper, which Google will find for you.
Here's an example of applying the function - you can change parameters to see what they do:
im = imread('cameraman.tif');
imshow(im);
logim = logsample(im, 5, 200, 150, 150, [], 300);
figure; imshow(logim);
There is a much more complex example in the code for my FEX contribution 27093, but unless you're already familiar with the basic ideas this is unlikely to shed a lot of light.
Aravind Natarajan (view profile)
Since I'm new to log polar transform I'm not able to understand the arguments for the function logsample.Can you please post an example of how to use the function.Thank you.
David Young (view profile)
Thanks Vlatko. Yes, the tilde as a marker for an usused result only came in with release 2009b. However, if you just omit it, the function will no longer work properly for RGB images. A better change is to use
[U, V, unused] = size(arr);
which is compatible with older versions and preserves the functionality.
Vlatko Becanovic (view profile)
At line 65 in the function logsample.m it says:
[U, V, ~] = size(arr);
If I change this to be:
[U, V] = size(arr);
it all works great.
Thanks for a nice program!