동일 직선상에 놓이지 않은 세점 A, B, C에 의해서 만들어지는 삼각형의 외접원 (circumcircle)에 한점 D가 포함되는가를 안되는가를 판단하는 것은 평면상의 점들을 삼각화화는 알고리즘에서 매우 중요하다. 물론 직접 외접원의 중심과 반지름을 구하고 이것을 가지고 중심과 점 D와의 거리를 계산하여서 판단을 할 수 있지만 이 경우에 수치계산의 오류가 들어갈 여지가 매우 많다(직접 중심과 반지름을 구하는 식을 만들어 보면 쉽게 알 수 있다). 여기서는 이것을 회피하는 방법을 고안하여서 보다 robust한 판단을 주는 식을 유도하려고 한다.

이것은 공간상의 네점 A(ax, ay, az), B(bx, by, bz), C(cx, cy, cz), D(dx, dy, dz)로 만들어지는 사변형의 체적을 구하는 공식이 아래의 행렬식으로 주어진다는 사실에서 출발한다.

사용자 삽입 이미지

이 사실을 이용하면 평면상의 3점에 의해서 형성이 된 원에 주어진 한점이 포함이 되느냐가를 판별하는 식 역시 3차원 공간에서의 체적으로 표현이 가능하다.

z= x^2 + y^2으로 표현이 되는 포물체(paraboloid)의 한 점 (a, b, a^2 + b^2)에서의 접평면의 방정식을 구하면,

         z=2ax + 2bx - (a^2 + b^2) ;

으로 주어진다. 이 평면을 z방향으로 r^2만큼 이동시키면

         z=2ax + 2by -(a^2 +b^2) + r^2;

으로 주어지고, 이 평면과 포물체가 만나는 경계선이 바닥에 정사영하면 아래의 방정식을 얻는다.

        (x-a)^2 + (y-b)^2 = r^2 ;

즉, 단면의 프로젝션은 중심이 (a,b)이고 반지름이 r^2인 원이 된다.

따라서, 원주상의 점들은 포물체상의 곡선상에 존재하고, 원의 내부점은 평면아래의 포물체면에 존재하게되고, 원주 밖의 점은 평면위의 포물체면에 존재하게 된다. 세점이 그리는 원인 경우에 원 내부와 외부점에 의한 사면체의 꼭지점이 평면을 기준으로 반대로 향하게 된다. 따라서 이 사면체의 체적도 부호를 달리하게 된다. 즉, 체적의 부호를 판별함으로써, 외접원내에 있는지 밖에 있는지를 판별할 수 있는 것이다.



다시 정리하면, 세점 a(ax, ay), b(bx, by) c(cx, cy) d(dx, dy)가 주어진 경우에 a, b, c(반시계방향)에 의해서 형성이 되는 원에 점 d가 포함이 되는지의 여부는 공간상의 네점

        A = (ax, ay, ax^2+ay^2),
        B = (bx, by, bx^2+by^2),
        C = (cx, cy, cx^2+cy^2),
        D = (dx, dy, dx^2+dy^2)

가 형성하는 체적의 부호에 의해서 결정이 된다.

사용자 삽입 이미지


위의 행렬식은 좀 더 단순화 시킬 수 있다. 1, 2, 3행에서 4행을 행단위로 빼도 값이 변화지 않으므로 오른쪽 식을 얻고(4열에 대해서 전개 후), 1열에 2*dx를 곱한것과 2열에 2*dy를 곱한것을 각각 3열에 빼면 아래줄의 식을 얻는다.

사용자 삽입 이미지

이 3x3 행렬식은 쉽게 계산이 된다. ccw(a,b,c)를 계산하여서 전체의 부호를 고정시킨다.

//A,B,C가 만드는 외접원에 D가 들어가는가?;
// det=0인 경우는 세점이 일직선상이나, D가 원주상에 놓인 경우다;
// 마지막 판단을 수정하여서 원주상에 놓인 경우도 해결할 수 있다.
int incircle(POINT A, POINT B, POINT C, POINT D){
    double ccw=CCW(A,B,C) ;
    double adx=A.x-D.x, ady=A.y-D.y,
               bdx=B.x-D.x, bdy=B.y-D.y,
               cdx=C.x-D.x, cdy=C.y-D.y,
               bdxcdy = bdx * cdy, cdxbdy = cdx * bdy,
               cdxady = cdx * ady, adxcdy = adx * cdy,
               adxbdy = adx * bdy, bdxady = bdx * ady,
               alift = adx * adx + ady * ady,
               blift = bdx * bdx + bdy * bdy,
               clift = cdx * cdx + cdy * cdy;
    double det = alift * (bdxcdy - cdxbdy)
                    + blift * (cdxady - adxcdy)
                    + clift * (adxbdy - bdxady);

    if(ccw>0) return det>0; //CCW인 경우에..
    else      return det<0;    //CW인  경우에..    
}

'Computational Geometry' 카테고리의 다른 글

Circle Drawing Algorithm  (1) 2008.06.03
Wu's Line Algorithm  (1) 2008.06.02
삼각형 외접원의 Inclusion Test  (0) 2008.05.29
Brute Force Triangulation  (0) 2008.05.28
Polygon Triangulation (II)  (0) 2008.05.26
Polygon Triangulation  (4) 2008.05.25
Posted by helloktk