CN102034241A - 一种球面距离快速变换方法 - Google Patents

Abstract

本发明公开一种球面距离快速变换的方法，属于图像处理技术领域。本发明方法包括：在经度-纬度参数域上将球面参数化成一个二维网格，行表示纬线，列表示经线；在二维网格上标注边界像素和非边界像素；沿表示纬线的行扫描，对于二维网格的每个像素，得到在相同纬度上与其具有最小距离的边界像素的索引并记录于该像素中；沿表示经线的列扫描，将每个网格像素保存的索引指向的边界像素依次加入到一个队列中，并在每次加入一个新的边界像素时做判断，以删除不满足一定条件的边界像素；沿表示经度的列再次扫描，根据上述队列依次得到列中的每个像素的球面距离变换。

Description

一种球面距离快速变换方法

技术领域

本发明涉及一种球面距离快速变换方法，属于图像处理技术领域。

背景技术

在二维空间中，假设图像是由边界像素和非边界像素构成，距离变换（距离场）是指将每一个像素映射成以距离函数定义的从该像素到最近边界像素的距离值。球面距离变换是指定义在球面上的球面测地距离变换。球面距离变换在球面网格剖分、气象建模等领域有着重要的应用。

直接的球面距离变换方法是通过比较球面上任一非边界点到所有边界点的距离进行的，这种方法的原理和实现简单，但是时间复杂度却非常高，因此几乎无法实际应用。采用基于顺序向量扩散的方法一定程度上缓解了球面距离变换的时间消耗（相关技术可参考Michikawa T, Suzuki H. Spherical Distance Transforms[C]，Proceedings of the 2008 International Conference on Computational Sciences and Its Applications，p 405--412），但是由于该方法采用了扩散的方法，因此只能单线程执行球面距离变换，效率仍然受限。在二维图像空间中，有一种独立扫描算法（Rosenfeld A, Pfaltz J L. Sequential Operations in Digital Picture Processing[J]. J. ACM. 1966, 13(4): 471-494）充分利用了二维距离变换过程中行、列之间的局部特性，使得距离变换能够并行执行。然而球面距离是指球面上两点之间的测地距离，与一般的欧氏距离存在本质上的区别，因此二维距离变换的算法无法被直接应用到球面距离变换。

发明内容

本发明的目的是提供一种球面距离快速变换方法。

为实现上述目的，本发明所采取的技术方案是：

（1）在经度-纬度参数域上，将球面参数化成二维网格；其中所述二维网格的行表示纬线，列表示经线；

（2）根据球面上边界点的分布要求，在所述二维网格中将与球面上边界点对应的网格像素标注为边界像素，其余网格像素标注为非边界像素；

（3）先并行地对所述二维网格的每一行从行首向行尾做第一次扫描，当扫描至行尾后，再从行首向行尾继续扫描，在该继续扫描过程中，若在不超过行的长度的一半时遇到边界像素，则结束第一次扫描，否则在达到行的长度的一半时结束第一次扫描；并且，在整个第一次扫描过程中跟踪每个边界像素并将每个边界像素的有效索引保存到在相应边界像素之后被扫描到的网格像素内；如果某一行中不存在边界像素，则在第一次扫描过程之后该行的所有网格像素保存无效索引；

当第一次扫描结束后，再并行地对所述二维网格的每一行从行尾到行首做第二次扫描，当扫描至行首后，再从行尾向行首继续扫描，在该继续扫描过程中，若在不超过行的长度的一半时遇到边界像素，则结束第二次扫描，否则在达到行的长度的一半时结束第二次扫描；并且，在整个第二次扫描过程中跟踪每个边界像素并将每个边界像素到当前网格像素的距离与第一次扫描过程中当前网格像素所保存的有效索引所指向的边界像素到当前网格像素的距离进行比较：如果在第二次扫描过程中跟踪的边界像素到当前网格像素的距离小于第一次扫描过程中当前网格像素所保存的有效索引所指向的边界像素到当前网格像素的距离，则将第二次扫描过程中所跟踪的边界像素的有效索引保存到当前网格像素内；

（4）对所述二维网格的每一列分别构造一个用于保存边界像素的先进后出的候选队列，且只能在所述候选队列的尾端插入或取出元素；接着并行地对所述二维网格的每一列从头到尾进行扫描，在该扫描过程中，依次将各列中的网格像素所保存的有效索引指向的边界像素从所在列的候选队列的尾端插入到该候选队列中，并且当出现候选队列中的边界像素的数量≥3时，按以下步骤进行：

（a）顺序取出该候选队列中的三个边界像素并分别依次记作q₁、q₂、q₃；在对应的正在扫描的列的中心线所在的经度上找两个点P1和P2，使得P1到q1、q2的距离相等且P2到q₂、q₃的距离相等，并进行以下判断：如果P1的纬度值小于P2的纬度值，或者P1和P2同时处在以两个相邻的网格像素的中心点为端点的球面线段上，则按照q₃、q₁的次序将这两个边界像素放回候选队列并执行步骤（b），否则，按照q₃、q₂、q₁的次序将这三个边界像素放回对应的候选队列并继续扫描该列的余下网格像素；

（b）判断该候选队列中的边界像素的数量是否≥3，若是，则返回步骤（a）；若不是，则继续扫描该列的余下网格像素；

（5）再次并行地沿着所述二维网格的每一列从头到尾进行扫描，顺序从各列对应的候选队列的首端取出两个边界像素，并分别依次记作q₄、q₅；从当前列的第一个网格像素开始，对于扫描过程中的每一个网格像素，分别计算它的中心点到q₄、q₅的球面距离，并记作d₁和d₂；如果d₁<d₂，则将d₁作为该网格像素的球面距离变换结果保存；否则，将d₂作为该网格像素的球面距离变换结果保存，然后将q₅赋值给q₄，并从当前列对应的候选队列的首端取出下一个边界像素赋值给q₅；当扫描完所有列的最后一个网格像素后，球面距离的快速变换即完成。

与现有技术相比，本发明的有益效果是：本发明将球面在经度-纬度参数域上参数化成一个二维网格后，在给定的球面边界点分布的前提下，提供了一种通过独立并行的扫描式的方法，快速执行球面上每个点的球面距离变换。如果令经度和纬度上的参数化粒度均为N，则使用本发明的方法进行球面距离快速变换的时间复杂度为O(N²)，即是关于参数化后二维网格像素数量的线性复杂度；同时由于本发明的方法中针对每行和每列的扫描是相互独立的，因此本发明在并行系统中的时间复杂度将于并行系统计算核心的数量成反比，执行效率将得到进一步的加速；例如在使用最新的GPU加速后，本发明提供的方法相比目前最新提出的基于顺序向量扩散的方法（相关技术可参考Michikawa T, Suzuki H. Spherical Distance Transforms[C]，Proceedings of the 2008 International Conference on Computational Sciences and Its Applications，p 405--412）有最高达300多倍的效率提升。

附图说明

图1-图5是本发明第一次扫描过程的示意图；

图6-图8是本发明第二次扫描过程中用于判断是否删除当前列的边界像素候选队列中的每三个边界像素的中间那个边界像素时三种情形的示意图；其中，图6和图7为不放回位于中间的边界像素的情形，图8为三个边界像素全部放回候选队列的情形；

图9-图10是本发明根据当前列的边界像素候选队列进行球面距离变换时两种情形的示意图，其中，图9中d1＜d2,图10中d1≥d2。

图11是一个应用本发明的球面距离快速变换的实例的效果图。

具体实施方式

下面结合附图，对本发明的球面距离快速变换的方法作进一步说明，具体步骤如下：

步骤1）：在球面经度-纬度参数域上，将单位球（半径为1个单位）的球面参数化成一个按照经度和纬度分别等间距划分的二维网格，其中二维网格的行表示纬线，列表示经线；令经度和纬度的参数化粒度均为N，即沿着球面的经线和纬线分别作N等分（在经度和纬度上可以使用不同的参数化粒度）。

步骤2）：根据给定的在球面上边界点的分布要求（分布要求可由用户根据实际的特定要求指定，或者仅由程序随机设定），在所述二维网格中将对应的网格像素标注为边界像素，其余的则是非边界像素；

步骤3）：并行地沿着所述二维网格的每一行扫描，得到在当前行的纬度下，距离所述二维网格每一行的每个像素最近的边界像素，并将该边界像素的索引记录在所述二维网格的每个像素中；如果某一网格像素本身是一个边界像素，则记录其自身的索引；如果某一行中不存在边界像素，则该行的所有网格像素保存一个无效的索引；此步骤将按照如下两次扫描完成：

3a）如图1所示为所述二维网格的某一行，A、B、C、D为在步骤2）中标注的4个边界像素，二维网格的每一行表示一条纬线，因此行首和行尾是相连的；第一次扫描将沿该行的行首（左）至行尾（右）的方向扫描最多3/2个圆周，即扫描至行尾后，再从行首继续扫描直至遇到边界像素或超过一半的行的长度；扫描过程中使用一个临时的索引T1，用以指向一个边界像素；扫描开始时，该索引T1不指向任何边界像素；当扫描进度从行首开始后，如果当前网格像素不是边界像素时，将所述的临时索引T1记录在当前的网格像素中；如果当前网格像素是一个边界像素时（如标注为A、B、C、D的网格像素），则将所述的临时索引T1改写为指向当前的边界网格像素，并将该索引T1记录在当前的网格像素中；从行首至行尾的扫描完成后，本实施例中，其结果将如图2所示，相应的小写字母（a、b、c、d）指向其对应大写字母标注的边界网格像素；继续回到行首扫描，如果在不超过行的长度的一半时遇到边界像素或者扫描的路程已经达到行的长度的一半，则结束扫描；扫描完成后，本实施例中，其结果如图3所示。

3b）第二次扫描将沿该行的行尾（右）至行首（左）的方向扫描最多3/2个圆周，即扫描至行首后，再从行尾继续扫描最多一半的行的长度；扫描过程中使用一个临时的索引T2，用以指向一个边界像素；扫描开始时，该索引T2不指向任何边界像素，即为一个无效索引；当扫描进度从行尾开始后，如果当前网格像素不是边界像素且临时索引T2指向一个有效的边界像素时，计算临时索引T2指向的边界像素到当前网格像素的经度距离，记作d₁；同时计算当前网格像素所记录的索引指向的边界像素到当前网格像素的经度距离，记作d₂，其中经度距离定义为两个网格像素在球面上的经度之差的绝对值，其范围为                                                

；如果d₁<d₂，则将临时索引T2记录在当前的网格像素中；如果当前网格像素本身就是一个边界像素时（如标注为A、B、C、D的网格像素），则将所述的临时索引T2改写为指向当前的网格像素的索引即可。在本实施例中，从行尾至行首的扫描完成后，其结果如图4所示；继续回到行尾扫描，如果在不超过行的长度的一半时遇到边界像素或者扫描的路程已经达到行的长度的一半，则结束扫描；在本实施例中，扫描完成后，结果如图5所示；其中图4和图5中所示的浅灰色背景的网格像素，表示其保存的边界像素的索引在第二次扫描过程中被更新。

步骤4）：并行地沿着表示经度的每一列扫描，得到所述二维网格的与每一列对应的一个存有用于每列的网格像素进行球面距离变换的边界像素的队列；经过步骤3）后，对于所述二维网格的每一列中的每一个网格像素，如果该网格像素所在的行（纬度）存在边界网格像素，则其保存了指向其中距离它最近的边界网格像素的索引；如果该网格像素所在的行（纬度）不存在边界网格像素，则其保存的是一个不指向任何地方的索引（即无效索引）；该步骤按如下方法扫描：

并行地沿着所述二维网格的每一列从低纬度（即纬度值为0）到高纬度（即纬度值为

）的方向（也可以从高纬度到低纬度扫描）扫描，令Q为一个保存与当前列相关的边界像素的先进后出的队列，即队列中的元素保持其插入到该队列中时的相对顺序，且按照其插入到队列中的次序的逆序取出元素，称Q为候选队列；在扫描过程中，每当当前被扫描的网格像素保存的索引为有效索引，则将其指向的边界像素插入到队列Q中；每次有新的边界元素插入到队列Q中后，记此时队列Q中的元素个数为n；如果n≥3，则进行如下操作：

（a）从队列Q的尾端依次取出三个边界像素并分别记作q1、q2、q3，在对应的正在扫描的列的中心线所在的经度

上找两个点P1和P2（如图6、图7、图8），使得P1到q1、q2的距离相等，P2到q₂、q₃的距离相等，并进行以下判断：如果点P1的纬度值小于P2的纬度值（如图6），或者P1和P2同时处在以两个相邻的网格像素的中心点O1和O2为端点的球面线段上（如图7），则按照q₃、q₁的次序将这三个边界像素放回候选队列Q并执行步骤（b）；否则（如图8），按照q₃、q₂、q₁的次序将这两个边界像素放回对应的候选队列Q中并继续扫描该列的余下网格像素；

（b）判断此时候选队列Q中的边界像素的数量n是否仍满足n≥3，若是，则返回步骤（a）；若不是，则继续扫描该列的余下网格像素；

步骤5）：再次并行地沿着所述二维网格的每一列从低纬度（即纬度值为0）到高纬度（即纬度值为

）的方向（也可以从高纬度到低纬度扫描）扫描。此时将步骤（4）得到的各列对应的候选队列的存取性质改为先进先出，即从队列首端取出元素，而从队列尾端插入元素；顺序取出各列对应的候选队列中的两个边界像素，并分别依次记作q₄、q₅；从当前列的第一个网格像素开始，对于扫描过程中的每一个网格像素，分别计算它的像素中心点P（如图9和图10）到q₄、q₅的球面距离，并记作d₁和d₂；如果d₁<d₂（如图9），则将d₁作为该网格像素的球面距离变换结果保存；否则（如图10），将d₂作为该网格像素的球面距离变换结果保存，然后将q₅赋值给q₄，并从当前列对应的候选队列中取出下一个边界像素赋值给q₅。当所有列都扫描完它们的最后一个网格像素后，球面距离变换即完成，即最终的二维网格的每个网格像素保存的距离为该网格像素到最近的边界像素的距离。

使用本发明的实施例，对于参数化为2048×2048的二维网格，并随机选取100个网格像素作为边界像素后进行球面距离的快速变换，使用一般方法需要44秒得到正确的球面距离变换，而本实施例仅使用0.013秒即可得到正确的球面距离变换。在不同参数化的二维网格下并随机选取1000个网格像素作为边界像素的条件下，与目前最新的方法（相关技术可参考Michikawa T, Suzuki H. Spherical Distance Transforms[C]，Proceedings of the 2008 International Conference on Computational Sciences and Its Applications，p 405--412）相比较的过程中，本实施例达到了最高近300倍的加速比，并且在8192×8192的参数化二维网格测试中，本实施例能达到15fps的球面距离变换性能。图11为一个应用球面距离快速变换的实例效果图，用于建模海洋中每个点到其最近的海岸线的距离变换，图11中，位于海洋中的条带为等距线，其余部分为大陆的示意图。

Claims (1)

1.一种球面距离快速变换方法，其特征在于包括如下步骤：

（1）在经度-纬度参数域上，将球面参数化成二维网格；其中所述二维网格的行表示纬线，列表示经线；

（2）根据球面上边界点的分布要求，在所述二维网格中将与球面上边界点对应的网格像素标注为边界像素，其余网格像素标注为非边界像素；

（3）先并行地对所述二维网格的每一行从行首向行尾做第一次扫描，当扫描至行尾后，再从行首向行尾继续扫描，在该继续扫描过程中，若在不超过行的长度的一半时遇到边界像素，则结束第一次扫描，否则在达到行的长度的一半时结束第一次扫描；并且，在整个第一次扫描过程中跟踪每个边界像素并将每个边界像素的有效索引保存到在相应边界像素之后被扫描到的网格像素内；如果某一行中不存在边界像素，则在第一次扫描过程之后该行的所有网格像素保存无效索引；

当第一次扫描结束后，再并行地对所述二维网格的每一行从行尾到行首做第二次扫描，当扫描至行首后，再从行尾向行首继续扫描，在该继续扫描过程中，若在不超过行的长度的一半时遇到边界像素，则结束第二次扫描，否则在达到行的长度的一半时结束第二次扫描；并且，在整个第二次扫描过程中跟踪每个边界像素并将每个边界像素到当前网格像素的距离与第一次扫描过程中当前网格像素所保存的有效索引所指向的边界像素到当前网格像素的距离进行比较：如果在第二次扫描过程中跟踪的边界像素到当前网格像素的距离小于第一次扫描过程中当前网格像素所保存的有效索引所指向的边界像素到当前网格像素的距离，则将第二次扫描过程中所跟踪的边界像素的有效索引保存到当前网格像素内；

（4）对所述二维网格的每一列分别构造一个用于保存边界像素的先进后出的候选队列，且只能在所述候选队列的尾端插入或取出元素；接着并行地对所述二维网格的每一列从头到尾进行扫描，在该扫描过程中，依次将各列中的网格像素所保存的有效索引指向的边界像素从所在列的候选队列的尾端插入到该候选队列中，并且当出现候选队列中的边界像素的数量≥3时，按以下步骤进行：

（a）顺序取出该候选队列中的三个边界像素并分别依次记作q₁、q₂、q₃；在对应的正在扫描的列的中心线所在的经度上找两个点P1和P2，使得P1到q1、q2的距离相等且P2到q₂、q₃的距离相等，并进行以下判断：如果P1的纬度值小于P2的纬度值，或者P1和P2同时处在以两个相邻的网格像素的中心点为端点的球面线段上，则按照q₃、q₁的次序将这两个边界像素放回候选队列并执行步骤（b），否则，按照q₃、q₂、q₁的次序将这三个边界像素放回对应的候选队列并继续扫描该列的余下网格像素；

（b）判断该候选队列中的边界像素的数量是否≥3，若是，则返回步骤（a）；若不是，则继续扫描该列的余下网格像素；

（5）再次并行地沿着所述二维网格的每一列从头到尾进行扫描，顺序从各列对应的候选队列的首端取出两个边界像素，并分别依次记作q₄、q₅；从当前列的第一个网格像素开始，对于扫描过程中的每一个网格像素，分别计算它的中心点到q₄、q₅的球面距离，并记作d₁和d₂；如果d₁<d₂，则将d₁作为该网格像素的球面距离变换结果保存；否则，将d₂作为该网格像素的球面距离变换结果保存，然后将q₅赋值给q₄，并从当前列对应的候选队列的首端取出下一个边界像素赋值给q₅；当扫描完所有列的最后一个网格像素后，球面距离的快速变换即完成。

Images