Minimum Enclosing Circle

평면 위에서 주어진 점집합을 포함하는 가장 작은 원을 찾는 문제는 가장 단순하게는 임의의 두 점이 지름이 되는 원과 임의의 세 점이 만드는 외접원을 모두 조사하여 찾으면 된다. 따라서 $O(N^4)$의 복잡도를 가질 것이다. 그러나 이 문제는 점집합의 개수에 비례하는 시간 복잡도($O(N)$)를 가지는 알고리즘이 알려져 있고, 여기서는 재귀적 방법을 이용한 Welzl's algorithm을 구현한다.

 

Welzl's algorithm은 평면에 주어진 여러 점을 모두 포함하는 최소 크기의 원(Smallest Enclosing Circle)을 구하는 대표적인 무작위 재귀 알고리즘이다. 2차원에서는 최소 원이 경계에 있는 최대 3개의 점에 의해 결정된다는 성질을 이용한다. 알고리즘의 핵심 과정을 정리하면 

1. 입력 점들의 순서를 무작위로 섞고 한 점 $p$를 선택한다.
2. $p$를 제외한 점들에 대해 최소 포함 원을 재귀적으로 구한다.
3. $p$가 그 원의 내부에 있으면 현재 원을 그대로 사용한다.
4. $p$가 원 밖에 있다면, 최종 최소 원은 반드시 $p$를 경계에 포함해야 한다. 따라서 $p$를 경계점 집합 $R$에 추가하고 다시 계산한다.
5. $R$에 점이 3개가 되면 세 점을 지나는 원이 결정되므로 재귀를 종료한다. $R$이 2개라면 두 점을 지름의 양 끝으로 하는 원, 1개라면 그 점을 중심으로 하는 반지름 0의 원에서 시작할 수 있다.

즉, 모든 점으로 가능한 원을 조사하는 대신 원 밖으로 벗어나는 점만 최소 원의 경계를 결정하는 후보로 추가하는 것이 핵심이다. 무작위화 덕분에 평균적으로 경계점이 변경되는 경우가 드물어 기대 $O(N)$의 효율적인 계산이 가능하다. 이러한 이유로 Welzl 알고리즘은 최소 포함 원이나 3차원의 최소 포함 구(bounding sphere)를 계산하는 데 널리 알려진 방법이다.

int MinEncCir ( CfPt *P, int n, CfPt* boundary, int b, CfPt& center, double& rad ) {
    // exiting cases
    if ( b == 3 ) CalcCircle3 ( boundary, center, rad ); // a circle passing 3 points;
    else if ( ( n == 1 ) && ( b == 0 ) ) {
        rad = 0; b = 1;
        center = boundary[0] = P[0];
    } 
    else if ( ( n == 2 ) && ( b == 0 ) ) {
        boundary[0] = P[0]; boundary[1] = P[1]; b = 2;
        CalcCircle2 (boundary, center, rad );  // a circle with diagonal consisting of 2 points;
    }
    else if ( ( n == 0 ) && ( b == 2 ) ) {
        CalcCircle2 ( boundary, center, rad );
    }
    else if ( ( n == 1 ) && ( b == 1 ) ) {
        boundary[1] = P[0]; b = 2;
        CalcCircle2 ( boundary, center, rad );
    }
    else {// general case; ( b < 3 ) && ( n + b > 2 )
        // choose a random pivot;
        int k = rand() % n;
        if ( k != 0 ) SWAP( P[0], P[k] );
        int b1 = MinEncCir ( &P[1], n - 1, boundary, b, center, rad );
        if (!InCircle(P[0], center, rad) ) {
            // Now, P[0] belongs to the boundary.
            boundary[b++] = P[0];
            return MinEncCir ( &P[1], n - 1, boundary, b, center, rad );
        } else return b1;
    }
    return b;
}

minimum_enclosing_ball_lecture13.pdf

0.16MB