平面点云的轮廓线计算-alpha shapes算法原理和实现

alpha shape算法又称为滚球法，是一种提取边界点的算法。跟凸壳提取相比，alpha shape算法能够了凹包情形，且对多个点云时 能勾勒出多个边界线，这是他的优势。

研究alpha shape算法的博文虽然不多，但也有相当数量了。但给出的算法大多有错误，或者只是部分实现。现将算法原理再梳理一遍。

如上图所示，alpha shape算法就好像在散点上设置一个半径为alpha的球在上面滚动，最后出来的线就是轮廓线。因此这个算法的滚球半径不能太小，否则可能不能将散点全部包起来。

alpha shape算法的具体流程：

1. 输入点云坐标向量p，确定滚球半径r（即前文的alpha）
2. 在p中任选一点p0（可以按p向量的排序来选，遍历所有点，也可以按照一定规则优选边界附近的）
3. 以p0为中心2r为半径画圆，圆内的剩余点记为pr_cir
4. 在pr_cir内任选一点（仍然可以按照序号来选，方便写成循环，尽管这样效率不高，但也是一种简单的实现方法）作为p1，然后计算出过p0和p1两点的圆心（当然应该有2个，圆心分别记为p2和p3）具体计算公式下面给出。
5. 将去除p0和p1的剩余点记为pr_cirr，然后计算pr_cirr所有点和p2 p3两个圆心的距离，可以用一个循环实现，将所有计算出的距离放入一个数组d中，由于有两个圆心，d可以分为两个d1和d2.
6. 判定是否为边界点的条件：只要pr_cirr内所有点到某一个圆的距离都大于r，就算是边界点。都大于r的意思是如果d1中最小值也大于r，且两个圆只要有一个就行，是“或”的关系。写成代码就是min(d1)>= r或min(d2)>=r。如果满足该要求，p1就是边界点，就可以定义新的p1开始下一轮循环了。
7. 如果不能满足上述条件，进入步骤4继续计算，即遍历pr_cir，如果pr_cir所有点都不能满足要求，则不为边界点。

上面就是alpha shape算法的具体流程，有些博文对判定条件的描述有些模糊，可能会造成误解。注意，必须是两个圆有一个满足要求，就算边界点。有些点对两个圆是不能同时满足大于r的要求的，这样就找不到边界点。

计算圆心坐标的公式：

 也可以用伪代码描述上述算法：

% 读入数据和半径
read data p
read r


% 遍历p
outline = []
for i = 1: length(p)
    % 定义p0
    p0 = p(0)
    % 定义pr
    pr = p(remove p0)
    % 定义pr_cir
    pr_cir = pr(inner 2*r)
    
    % 遍历pr_cir
    for j = 1: length(pr_cir)
        % 定义p1
        p1 = pr_cir(0)
        % 定义pr_cirr
        pr_cirr = pr_cir(remove p1)
        % 计算圆心
        p2 = cir1(p0,p1)
        p3 = cir2(p0,p1)
        % 计算pr_cirr到圆心的距离
        for k = 1: length(pr_cirr)
            pk = pr_cirr(k)
            d1(k) = distance(p2,pk)
            d2(k) = distance(p3,pk)
        end
        mind1 = min(d1)
        mind2 = min(d2)
        if mind1>=r || mind2 >= r
            outline = [outline p1]
            break
        end
    end
end
obtain outline

按照这个算法，就能编出没有bug的alpha shape边界点算法了。

还参考了斯坦福这个教程，不过有点复杂，先放在这里吧。

https://graphics.stanford.edu/courses/cs268-11-spring/handouts/AlphaShapes/as_fisher.pdf

来源：Intimes