Fast algorithm for polar -> cartesian conversion

Question

I have an image on a polar grid. This image should be transformed into a cartesian grid, but the only algorithm I know of is really slow for this. Now I use the cartesian grid, for each point I find the r and theta values, and then I look in two vectors to find the smallest error defined by:

min{(th_vec - theta)^2 + (range - r)^2}

This gives a nested for-loop inside of the outer nested for-loop, so I have a complexity of O(N^4). A 512x512 image uses a whole minute to complete. Of course, a complexity like that can not be used, so I'm wondering if anyone know of any faster algorithms to do this?

I have the image, and the two vectors. The X-axis of the image is the angle, while the Y-axis of the image is the length from the center. The angle is always from 0-2pi, and the range goes from 0 to r_max.

Thank you in advance.

EDIT: The range goes from 0 to r_max, not -r_max to r_max as it stood before. I see that there have been some missunderstandings. I have used the normal, inverse, conversion with;

r=sqrt(x^2 + y^2);
theta=atan2(y,x);

The problem is that I have to first convert the x and y values to x' and y' values, since the grid is from -r_max to r_max in the resulting image, but in pixels in the data. So I have a 512x512 image, but r_max can be something like 3.512. So I have to convert each pixel value into the grid value, then find the r and theta values. When I have found the r and theta values I have to run trough two vectors, range and th_vec, to find the pixel in the original image that matches:

min{(range - r)^2 + (th_vec - theta)^2}

This gives me a complexity of O(n^4), since the th_vec and range vectors are the same size as the image. So if I have a square matrix of 512x512 elements, I have to run trough 68 719 476 736 elements, which is way slow. So I'm wondering if there is a faster algorithm? I can't change the input data, so as far as I know, this is the only way to do it if you don't start with triangulation and stuff, but this is to expensive in times of memory.

Answer 1

How about

x=r*cos(angle)
y=r*sin(angle)

This is the standard way of converting from polar to Cartesian, and unless you're going to use some kind of table lookup, there's not really a faster option.

Edit: wrang wrang has a good point. If you're trying to transform an image in polar coordinates I(angle, r) to an image in Cartesian coordinates I_new(x, y), you are definitely better off using the inverse transformation, as follows:

for x=1,...,width
    for y=1,...,height
        angle=atan2(y, x)
        r=sqrt(x^2+y^2)
        I_new(x, y)=I(angle, r)
    end
end

As a rule, angle and r will not be integer, so you have to do some kind of interpolation in the image I. The easiest way to do this is simply to round angle and r; this will give you nearest-neighbour interpolation. If you need better quality, try more sophisticated types of interpolation such as bilinear or bicubic interpolation.

Answer 2

You could loop over each pixel in the polar image map and then render the resulting arc section in the cartesian image plane:

polar to cartesian conversion http://img24.imageshack.us/img24/4635/polartocartesian.png

const float dR = 2*r_max / polar_image_height;
const float dA = 2*pi / polar_image_width;

float angle;
float radius;
for (int polar_x = 0; polar_x < polar_image_width; polar_x++)
{
    for (int polar_y = 0; polar_y < polar_image_height; polar_y++)
    {
        angle = polar_x * dA;
        radius = polar_y * dR - r_max;
        DrawArcSection(radius, radius+dR, angle, angle+dA);
    }
}

Many drawing libraries have built-in functions for drawing that arc section, but you could always just approximate it with a simple polygon:

void DrawArcSection(float minRadius, float maxRadius,
                    float minAngle, float maxAngle)
{
    point P1 = MakePoint(minRadius * cos(minAngle) + image_width/2,
                         minRadius * sin(minAngle) + image_height/2);
    point P2 = MakePoint(minRadius * cos(maxAngle) + image_width/2,
                         minRadius * sin(maxAngle) + image_height/2);
    point P3 = MakePoint(maxRadius * cos(minAngle) + image_width/2,
                         maxRadius * sin(minAngle) + image_height/2);
    point P3 = MakePoint(maxRadius * cos(maxAngle) + image_width/2,
                         maxRadius * sin(maxAngle) + image_height/2);

    DrawPolygon(P1, P2, P3, P4);
}