如何计算二维向量形状的中轴线?

13

我有一个将2D形状作为路径元素存储在SVG中。这些形状由贝塞尔曲线和线段组成。

我还使用弧长参数化在形状上生成一组等间距点。

如何使用SVG或这些点来确定形状的中轴线?

我使用Python,但任何伪代码或算法建议都将不胜感激。


以下是我正在处理的形状类型示例,红点是沿曲线采样的点。

example


你能为我们提供一个示例形状让我们玩吗? - will
当然,这是上面示例中使用的链接:http://hastebin.com/esiyojehik.svg 这是沿着这个形状的一组点:http://hastebin.com/usosaruyup.py - flutillie
你所说的中轴线是指什么?是指经典的主惯性轴(http://www.mukimuki.fr/flashblog/2009/05/20/magic-moments/)还是指一条曲线,始终保持在形状的中心,并试图沿着其边界最终呈现出曲线的特性? - Diego Mazzaro
1
我在Github上几乎完成了一个类似的项目(https://github.com/FlorisSteenkamp/MAT)。你已经可以查看演示(http://mat-demo.appspot.com)。一旦项目完成,我会尽快发布完整的答案。 - Floris
2个回答

24

有点晚,但是让我们开始吧:

Medial and Scale Axis Transform Scale Axis Transform Medial Axis Transform 2-Prong 3-Prong

上面的图片显示: 1. 叠加的中轴变换(MATSAT)。 2. 仅缩放轴变换。 3. 仅中轴变换。 4. 许多2叉(见下文)之一。 5. SAT中的单个3叉(MAT中有很多)。
为了找到中轴(Media Axis)中轴变换(Media Axis Transform)(上图紫色曲线),可以使用以下算法(基于Choi、Choi、Moon和Wee的paper,也可以看到一个演示实现here,也可以处理相交的形状和带孔的形状)。还存在其他算法。
这个算法比找到二进制图像(例如位图)骨架(也称为草火或离散变换)难得多,但具有几个优点(如是解析的)。为简化事情,下面的讨论只处理没有孔的简单(非相交)形状。
一些定义
  • Curve - 参数为t ∈ [0,1]的参数贝塞尔曲线。换句话说,是矢量图形中使用的典型曲线。
  • Shape (Ω) - (定义直接取自论文) - “由有限数量互不相交的简单封闭曲线界定的 ℝ² 中连通有界开子集的闭包。” 为了回答这个问题,还可以进一步假设shape没有孔。简单来说,它是由许多曲线给出的边界所围成的区域 - 上图中填充的灰色部分。
  • Boundary - shape的边界。换句话说,是一系列索引的curve,使得任何曲线在t = 0处的起始点对应于前一个曲线在t = 1处的端点。假定boundary是正向定向的(即沿着边界逆时针走,形状始终在您的左侧)。
  • Boundary point - 边界上的任何点,完全由标识curvet curve参数的索引确定。
  • Maximal disk - 形状内没有被形状内其他圆盘遮挡的圆盘。每个maximal disk都可以由1个或多个(通常是2个)boundary point和其中心点进行识别。
  • n-prong - 在n个点处切向接触形状边界的maximal disk
  • Branch point - 具有n>= 3的n-prong maximal disk center。换句话说,是MA上形成平面树中的一个顶点的点。
  • Medial Axis (MA) - 所有n-prongs的中心的联合。
  • Medial Axis Transform (MAT) - 所有n-prongs的联合。换句话说,maximal disks的半径包括在定义中。
  • Contact point (cp) - 链接到已找到的唯一maximal diskboundary point
  • Sharp corner - 通常在曲线参数t = 0或t = 1处的boundary point,不可微分(不平滑),并且形成内角 < 180°。
  • Dull corner - 与尖角相同,只是角度> 180°。
  • Cp graph - 具有contact point作为顶点的图形(不一定是平面图),因此还包含有关所有找到的maximal disk以及整个MAT的信息。图形中的每个顶点(cp)都具有以下边缘:
    • next - 通过沿着boundary逆时针走来自当前点的第一个contact point
    • prev - 通过沿着boundary顺时针走来自当前点的第一个contact point一些注意事项
      • Medial Axis Transform(MAT)形成一个无根平面树,每个1-prong中心(包括sharp corner)都有一片树叶;如果形状有孔,则MAT是一个平面图。
      • 您可以通过对MAT应用Scale Axis Transform(SAT)(上图中的红色曲线)来消除形状边界上的噪声。demo实现了这样的转换,基于this研究。然而,该演示改进了算法,确保SATMAT的子集并保证了拓扑保持。

      算法概述

      1. 找到所有的1-prongs,它们恰好是Medial Axis树的叶子节点。这些是尖角和1-prongs,与最大局部内曲率的边界点相接触。将1-prongs的接触点添加到cp-graph中。请注意,由于1-prong只有一个接触点,因此next-around和prev-around边将再次返回找到的接触点。
      2. 对于代表性的边界点集(例如由OP给出的那些红色点),计算它们的最大圆盘。它们通常是2-prongs。(见下面的“查找2-prong”)。通过连接适当的边将它们的接触点添加到cp-graph中。请注意,在这种情况下,由于我们发现了2-prongs,如果连续两次执行next-around(或prev-around),则会回到相同的接触点。
      3. 最后找到n-prongs,其中n>= 3。从每个接触点开始:

        1. 通过从接触点(例如cp1)开始并重复应用cp1.next.next-around来遍历cp-graph,直到回到cp1。
        2. 如果上述中的1或2次迭代中的任何一次发生,请移动到边界上的下一个接触点(使用next)。但是,如果需要3次或更多次迭代,则可以证明存在3-prong分支点,并且应在cp1和cp1.next之间的每个边界段上插入其接触点。

          在这种情况下,请插入3-prong(见下面如何找到这些)并重新从步骤1开始,从cp1再次开始。

        3. 重复步骤1和2,直到找到所有的3-prongs。
      4. MAT的结构现在已经隐含地找到了。 要提取MAT,需要遍历cp-graph。选择任何一个接触点(称其为cp1)作为根(为简单起见,假设它链接到2-prong),然后应用next(称其为cp2)。现在可以从cp1的链接最大圆盘中心到cp2的最大圆盘中心绘制一条线。这条线是MA的子集的近似值。通过插值计算点之间的曲线,可以得到更好的近似值,因为我们知道接触点处的边界切线。让两个线端点成为二次贝塞尔曲线的端点。由于MA在两个最大圆盘中心的方向完全相同(通过2个边界切线的平均值),因此我们也有端点处二次贝塞尔曲线的导数,并且我们可以在点之间构造最佳近似的二次贝塞尔曲线。如果cp2链接到3-prong,则意味着由MA形成的平面树分叉(即分支或叉),我们需要跟随两个MA边(例如使用堆栈或递归)。另一方面,如果cp2链接到1-prong,则已到达MA的叶子节点,并且该遍历对于该边缘停止。实际上,由于MA形成了平面树,因此遍历算法本质上与任何其他树数据结构的遍历算法相同。

      寻找2针插头

      我将在这里进行总结,但是Choi等人的论文非常清晰易懂,很容易理解。选择边界上的一个点作为我们2-prong的第一个接触点,并称其为bp1。从边界点(即要找到的第一个“prong”)向内绘制一条法线射线。现在迭代执行以下步骤: 1. 选择一个点p1(沿着该射线距离bp1为d1),使得以射线为圆心并切于bp1的圆至少会在边界上再相交一个点。 2. 找到最靠近p1的边界点,称之为bp2。最后,找到距离bp1和bp2等距离的射线上的点,并称其为p2。 3. 使用p2作为新的p1返回步骤2。 4. 重复步骤2和3,直到从p1到两个边界点的距离在某种容差范围内相等。现在已经找到了2-prong,其中bp1和bp2是其两个接触点,p1是其中心。

      寻找三叉戟

      这里我们与论文稍有不同,使用构建外接圆方法来代替使用势函数。

      1. 选择一个边界点(cp1)的链接至多的最大圆心作为要查找的三叉戟中心的初始猜测。称其为p
      2. 分别对每个3+边界片段找到距离p最近的3个点。
      3. 构造这些点的外接圆并将其圆心用作新的点p
      4. 重复步骤2和3,直到从p到3个边界片段的距离在一定公差内相等。
      5. 将找到的这3个边界点标记为三叉戟的接触点,点p是三叉戟的中心。

      如果有任何不清楚的地方,请随时提问。


我已经不再参与这个项目了,但是非常感谢你提供如此详尽清晰的答案。我很喜欢阅读它! - flutillie
谢谢您的解释。您知道这个算法的计算复杂度吗?如果所有边界线段都是线性的,是否可以进行其他简化? - headdab
该算法的平均运行时间复杂度为(对于指定的精度)O(n log n),最坏情况下为O(n^2),其中n是贝塞尔曲线的数量。但这只是一个很大的简化,因为它还取决于输入贝塞尔曲线的位数等因素。是的,如果所有边界线段都是线性的,那么该算法会简单得多,如果我没记错的话,已知有一种O(n log n)的平均情况算法。 - Floris
算法是否容易返回一组离散的最大球?据我所知,目前返回了MAT的路径,因此可以从中计算出最大球的位置。但是如何获取每个球的直径呢?顺便说一句,这是一个很棒的库,我还在努力理解论文! - fweth
@fweth 谢谢 :) 是的,它返回作为 MAT 一部分的最大球。MA(Medial Axis)是形状“中心”的路径。MAT(Medial Axis Transform)是路径加上路径上每个点的球半径(包含所有重建原始形状所需的信息)。例如,在库中迭代路径时(请参见例如 https://stackblitz.com/edit/typescript-xuruhe?file=src%2Fget-thinned-path.ts 的第 69-83 行),路径上的点作为 CpNode 对象返回。此对象包含一个具有 ...(继续) - Floris
反过来,它有一个circle属性,这正是您要寻找的。它实际上存储了额外的(在技术上数学上是冗余的信息),例如切线相交最大球体边界上的点等。 - Floris

2
您可以从skimage(scikit-image)中查看代码。您将找到骨架化和介轴的代码(skimage.morphology.medial_axis)。
源代码在此地址可用:https://github.com/scikit-image/scikit-image/blob/v0.12.2/skimage/morphology/_skeletonize.py#L103 该算法计算图像的介轴转换,即其距离变换的脊线。
该算法的不同步骤如下所示。
A lookup table is used, that assigns 0 or 1 to each configuration of
   the 3x3 binary square, whether the central pixel should be removed
   or kept. 

We want a point to be removed if it has more than one neighbor
   and if removing it does not change the number of connected components.


The distance transform to the background is computed, as well as
   the cornerness of the pixel.

The foreground (value of 1) points are ordered by
   the distance transform, then the cornerness.

A cython function is called to reduce the image to its skeleton. It
   processes pixels in the order determined at the previous step, and
   removes or maintains a pixel according to the lookup table. 

Because
   of the ordering, it is possible to process all pixels in only one
   pass.

我希望这能对你有所帮助


网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接