极坐标转换为直角坐标的快速算法

Question

我有一个在极坐标网格上的图像。这个图像应该被转换成笛卡尔网格，但我知道的唯一算法对此非常慢。现在我使用笛卡尔网格，对于每个点，我找到r和theta值，然后查找两个向量中定义的最小误差：

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

这会在外部嵌套循环中产生一个嵌套的for循环，因此我的复杂度为O(N^4)。一个512x512的图像需要整整一分钟才能完成。当然，这样的复杂度不能使用，所以我想知道是否有任何更快的算法来做到这一点？

我有这个图像和这两个向量。图像的X轴是角度，而图像的Y轴是距离中心的长度。角度始终为0-2pi，范围从0到r_max。

谢谢你们的帮助。

编辑：范围从0到r_max，而不是从-r_max到r_max，如之前所述。我看到有一些误解。我使用正常的反向转换，即：

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

问题是我必须先将x和y值转换为x'和y'值，因为网格在结果图像中从-r_max到r_max，但在数据中是以像素为单位的。所以我有一个512x512的图像，但r_max可能是3.512之类的东西。因此，我必须将每个像素值转换为网格值，然后找到r和theta值。当我找到r和theta值时，我必须遍历两个向量，range和th_vec，以查找与原始图像匹配的像素：

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

这使得我的复杂度为O(n^4)，因为th_vec和range向量的大小与图像相同。因此，如果我有一个512x512元素的正方形矩阵，我必须遍历68 719 476 736个元素，这非常慢。所以我想知道是否有更快的算法？我无法更改输入数据，所以据我所知，如果您不从三角测量等开始，这是唯一的方法，但在记忆时间上这太昂贵了。

- martiert

---

这是用来做什么的？另外，为什么你没有从0到pi的角度或从0到r_max的范围？2*pi可以表示整个圆，那么为什么需要负距离呢？ - Jonathan Graehl

你的极坐标网格是否相对于极坐标均匀分割？ - Richard Dunlap

如果你从 x,y 中找到了 r_0 和 th_0 的某些浮点值，那么你只需要查看极坐标图像中的四个对应点 (r,th)，即 (r_0,th_0) 的四个最近邻居，它们是 floor(r_0),ceil(r_0) 和 floor(th_0),ceil(th_0) 的四种组合，其中 floor() 和 ceil() 产生的结果会被舍入到你的极坐标网格。 - jilles de wit

Answer 1

如何呢？

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

这是从极坐标转换为直角坐标的标准方法，除非你要使用某种表查找，否则没有更快的选项。

编辑：wrang wrang有一个很好的观点。如果您想将在极坐标下的图像I(angle, r)转换为直角坐标下的图像I_new(x, y)，那么最好使用反向变换，具体如下：

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

通常情况下，angle和r不会是整数，所以您需要在图像 I 中进行某种插值。最简单的方法是将angle和r四舍五入；这将给您提供最近邻插值。如果需要更好的质量，请尝试更复杂类型的插值，例如双线性或双三次插值。

- Martin B

---

你必须描述如何填充整个图像，因此除非你使用一些聪明的插值方法，否则反向变换更有用。 - Jonathan Graehl

Answer 2

你可以循环遍历极坐标图像中的每个像素，然后在笛卡尔图像平面上呈现出所得到的弧形部分：

极坐标转换为笛卡尔坐标系 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);
    }
}

许多绘图库都内置了绘制弧形的函数，但您也可以用简单的多边形来近似绘制它：

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);
}

- e.James

Answer 3

如果您不关心平滑处理，为什么不计算每个目标笛卡尔像素坐标的极坐标并读取颜色值呢？如果需要帮助，请参见http://en.wikipedia.org/wiki/Polar_coordinate_system#Converting_between_polar_and_Cartesian_coordinates。

- Jonathan Graehl

---

如果您希望最终图像被填充，最好进行反向变换并在源图像中进行插值。 - Martin Beckett

Answer 4

如果您的网格相对于极坐标均匀分割，则如果利用最接近点到（r，theta）将是其所包含的网格元素的四个角之一的事实，您的算法可以简化为O（N ^ 2）。

在更一般的情况下，其中网格是r和theta维度的任意分区的乘积，如果您必须在每个分区中搜索点的位置，则可能会增长到O（（N log N）^ 2）。但是，如果分区是系统地构建的，则应该能够回到O（N ^ 2）。

- Richard Dunlap

Answer 5

如果你所有的图片都是512x512，那么我会使用查找表将极坐标图像中一组加权像素映射到直角坐标图像。虽然这一开始需要很多工作，但最终计算的复杂度只会是O(n^2)。如果查找表不可行，我会使用：

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

对于极坐标图像中的每个像素，将其映射到笛卡尔图像中的“a”像素，其中输出像素是落在其上的所有输入像素的平均值。然后重复执行膨胀操作，直到没有未初始化的像素为止。对于膨胀操作，使用3x3结构元素，并且仅当输出像素先前没有值时，才用中心像素的值替换输出像素的值。最后，对整个图像应用高斯滤波器以平滑硬边缘。这是我能想到的最快的方法，在合理的时间内生成一个令人愉悦的图像。

- jilles de wit

Answer 6

记忆有时会出错，但可能存在一种涉及FFT的快速版本算法。曾经我上过一门医学成像课程，当反变换/光栅化CT扫描时，似乎会出现这种情况。一些搜索关键词包括Radon变换、滤波反投影算法和CT扫描。我在维基百科上简要查阅了这些内容，但没有找到什么有用的信息，也许更彻底的研究会发现一些宝藏。

- RAC

Answer 7

O(N²log(N)) 算法：

• 数组S将用于最接近源（极坐标）坐标的笛卡尔坐标。
• S开始填充一个“尚未初始化”的值。（Python：None，Haskell：Nothing等）
• O(N²) - 迭代所有的极坐标。
  • 转换为笛卡尔坐标
  • 在目标图像中找到最接近的笛卡尔坐标。（通过四舍五入和应用边界）
  • 使用此坐标填充S中的相关单元格
• O(N²log(N)) - 执行如下所述的修改后的Dijkstra算法：
  • 我们搜索算法的“图”如下：
    • S的所有单元格都是节点
    • 单元格的邻居是国王可以从其移动到的那些单元格
  • 如果单元格未初始化，则其“分数”为无限大，并且为指向它的极坐标的未触及的笛卡尔坐标的距离
  • 更新单元格N的邻居时，我们将N单元格中的值放入其中（但仅在它的分数比当前分数更好时才这样做，就像在Dijkstra中一样）
  • 起点是上面描述的数组S初始化

- yairchu