728x90

// 수퍼 삼각형은 입력 점집합을 모두 포함하도록 충분히 큰 삼각형을 선택한다;
// 삼각화한 후 최외각이 convex가 아니면 수퍼 삼각형을 더 크게 다시 잡아야 함.
// 마지막 수정: 2021.4.7, 코드 간결화;
#define REMOVE_FLAG    (-1)
//
// max num of triangles = 2 * (N - 2) - 1;
// max num of edegs     = 3 * (N - 2);
int inc_delaunay_triangulation(POINT P[], int N, TRIANGLE *V) {
    // 메모리는 입력 점집합의 크기가 (N+3)인 경우를 기준으로 할당되어야 한다;
    int emax = 3 * (N + 1) ;
    EDGE *E = (EDGE *)malloc(sizeof(EDGE) * emax) ;
    set_supertriangle(P, N) ;
    // First element of triangle index is supertriangle vertex;
    push_tri(&V[0], N, N + 1, N + 2);
    int nt = 1;
    for (int i = 0; i < N; i++) {
        int ne = 0;
        /* 추가점(P[i])를 외접원에 포함시키는 삼각형을 모두 찾아 제거한다.
        ** 이들 삼각형의 공유되는 내부 edge를 제외한 나머지 외부 edge와
        ** 추가점으로 새로운 삼각형을 만들어 리스트에 추가한다.
        */ 
        for (int j = 0; j < nt; j++) {
            TRIANGLE *t = &V[j];
            if (incircle(P[t->a], P[t->b], P[t->c], P[i])) {
                // 제거되는 삼각형의 edge는 백업한다; 내부 edge는 중복된다.
                push_edge(&E[ne++], t->a, t->b);
                push_edge(&E[ne++], t->b, t->c);
                push_edge(&E[ne++], t->c, t->a);
                // P[i]를 포함하는 외접원 삼각형은 제거;
                // swap(j, nt-1) to remove triangle at j
                V[j--] = V[--nt];  
            }
        }
        // 중복되어 있는 내부 edge를 제거;
        remove_double_edges(E, ne);
        
        // 외부 edge와 P[i]가 만드는 새 삼각형을 리스트에 추가;
        for (int j = 0; j < ne; j++) {
            if (is_removed_edge(&E[j])) continue;
            push_tri(&V[nt++], E[j].a, E[j].b, i);
        }
    }
    //done with super_triangle;
    free(E);
    return remove_supertriangle_relics(N, V, nt) ;
}
void remove_double_edges(EDGE *E, int ne) {
    for (int j = 0; j < ne - 1; j++)
    	for (int k = j + 1; k < ne; k++)
            if (is_double_edge(&E[j], &E[k])) {
            	pop_edge(&E[j]); pop_edge(&E[k]);
            }
}
void push_tri(TRIANGLE *V, int a, int b, int c) {
    V->a = a; V->b = b; V->c = c;
}
int is_removed_edge(const EDGE *E) {
    return (E->a == REMOVE_FLAG || E->b == REMOVE_FLAG);
};
int is_double_edge(const EDGE* E1, const EDGE* E2) {
    return ((E1->a == E2->a && E1->b == E2->b) ||
            (E1->a == E2->b && E1->b == E2->a)) ;
}
void push_edge(EDGE *E, int a, int b) {
    E->a = a; E->b = b;
}
void pop_edge(EDGE* E) {
    E->a = E->b = REMOVE_FLAG;
}

// 삼각형 중에서 수퍼 삼각형의 꼭지점(N,N+1,N+2)을 가지고 있는 것을 제거;
int remove_supertriangle_relics(int N, TRIANGLE *V, int nt) {
    for (int i = 0; i < nt; i++)
        if (V[i].a >= N || V[i].b >= N || V[i].c >= N)
            // swap(i, nt - 1) to remove triangle at j
            V[i--] = V[--nt];
    return nt;
}

// A,B,C가 만드는 외접원에 D가 들어가는가?;
int incircle(POINT A, POINT B, POINT C, POINT D) {
	POINT AD = A - D, BD = B - D, CD = C - D;
    double ccw = CCW(A, B, C);
    double det = norm2(AD) * cross(BD, CD)
               + norm2(BD) * cross(CD, AD)
               + norm2(CD) * cross(AD, BD)
    if (ccw > 0) return det > 0; //CCW인 경우에..
    else         return det < 0; //CW인  경우에..    
}

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

삼각형 외접원의 Inclusion Test  (0) 2020.12.30
Point in Polygon  (2) 2020.12.14
Incremental Delaunay Triangulation  (1) 2020.12.01
Chain Hull  (2) 2012.09.16
Quick Hull  (2) 2012.09.16
Monotone Cubic Interpolation  (0) 2012.09.09
Posted by helloktk

댓글을 달아 주세요

  1. hgmhc 2020.12.02 20:40 신고  댓글주소  수정/삭제  댓글쓰기

    그 개어렵다는 델로네 삼각분할이군요...
    나중에 배울때 참고하겠습니다^^

728x90

실행화일.zip
다운로드

마우스 클릭으로 생긴 이미지 상의 입력 점을 가지고 triangulation을 한 후에 각각의 삼각형 내의 영역을 그 영역의 평균 컬러 값으로 채우는 효과이다. voronoi 다이어그램을 이용해서 tiling 하는 효과와 유사하다.

그림: Lena사진(일부)을 이용한 예:

사용 algorithm
  - incremental delaunay triangulation;
  - polygon fill;



Voronoi Tiling:

 

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

한 점에서 선분까지 거리  (1) 2010.01.16
삼각형 채우기  (0) 2009.12.19
Triangulation을 이용한 Imge Effect  (0) 2009.12.18
Ellipse Drawing Algorithm  (0) 2008.06.04
Circle Drawing Algorithm  (1) 2008.06.03
Wu's Line Algorithm  (1) 2008.06.02
Posted by helloktk

댓글을 달아 주세요

728x90

주어진 점집합에서 세 점을 뽑아서 만든 삼각형의 외접원에 다른 점이 하나도 포함하지 않으면 triangulation의 기본 삼각형 cell이 된다. 주어진 점으로 만들 수 있는 삼각형의 총개수가 ${}_nC_3$ 이므로, 기본 삼각형을 찾기 위해서는 이들 각각의 삼각형에 대서 나머지 점을 가지고 incircle 테스트를 수행해야 한다. 따라서 이 알고리즘은 ${\cal O}(n^4)$의 스텝이 필요하게 된다.

/*brute force attack*/
foreach point p1
      foreach point p2 
             foreach point p3
                   foreach point p4 
                          if(incircle(p1,p2,p3,p4))
                                 iscell=false;
                                 break ;
                   endfor;
                   if(iscell) 
                           add(triangle(p1,p2,p3));
              endfor;
       endfor;
endfor;

세 점에 의해서 형성이 되는 외접원은 대수적으로 쉽게 구할 수 있다. 여기서는 좀 더 기하학적인 접근을 쓰면, 평면의 점은 

$$(x, y)\rightarrow (x, y, z=x^2 + y^2)$$

의 mapping에 의해서 3차원 paraboloid 곡면의 점으로 옮길 수 있다. paraboloid 위의 세 점이 형성하는 3차원에서 평면이 paraboloid를 절단하는 곡선을 $x-y$ 평면으로 정사영하면 원이 된다는 것을 쉽게 알 수 있다.(incircle 포스팅 참조). 따라서 주어진 점이 세 점의 외접원에 포함되는가를 테스트하는 작업은 이 점을 paraboloid로 올렸을 때의 점과 (paraboloid로 올려진) 외접원의 3점이 형성하는 3차에서의 평면과 관계를 조사하는 것으로 바꿀 수 있다.

주어진 점이 외접원에 포함되면 paraboloid로 변환된 점은 평면의 아래에 놓이고, 외접원 밖의 점이면 평면 위에 놓이게 된다. 물론 외접원 위의 점은 평면에 놓인다. 따라서 평면의 법선 벡터 구하고, 삼각형의 한 꼭짓점을 기준한 주어진 점의 변위 벡터와 내적을 구하면 내적의 값은 평면 위인지, 아래인지 또는 평면에 놓인 점인가에 따라서 부호가 달라진다. 평면의 수직 벡터를 고정하면(예제는 아래 방향: $n_z < 0$), 평면 위에 놓인 점과의 내적은 음수, 평면 아래에 놓인 점과의 내적은 양수가 되고, 평면의 점과의 내적은 0이다. 

주어진 세 점이 만드는 외접원 내부(and 경계)에 들어가는 점이 없으면 이 삼각형을 선택한다.

** 참고 : Computational Geometry in C(2nd Edition) by Joseph O'Rourke

// input  x[0], x[1],....,x[n-1],
// input  y[0], y[1],....,y[n-1];
// calculate z[0]=x[0]^2+y[0]^2, z[1]=x[1]^2+y[1]^2,...;
int dt4(int n, double x[], double y[]) {
    double *z = new double [n] ;
    for(int i = 0; i < n; i++) 
        z[i] = x[i] * x[i] + y[i] * y[i] ;

    int ntri = 0;
    /* For each triple (i,j,k) */
    for (int i = 0; i < n - 2; i++ )
        for (int j = i + 1; j < n; j++ )
            for (int k = i + 1; k < n; k++ )
                if ( j != k ) {
                    /* Compute normal to triangle (i,j,k)::  outter_product(j-i, k-i)*/
                    double nx = (y[j] - y[i]) * (z[k] - z[i]) - (y[k] - y[i]) * (z[j] - z[i]); 
                    double ny = (x[k] - x[i]) * (z[j] - z[i]) - (x[j] - x[i]) * (z[k] - z[i]);
                    double nz = (x[j] - x[i]) * (y[k] - y[i]) - (x[k] - x[i]) * (y[j] - y[i]);
                    
                    /* Only examine faces on bottom of paraboloid: nz < 0. */
                    int flag = (nz < 0);
                    if (flag) {
                        /* For each other point m */
                        for (int m = 0; m < n; m++) {
                            /* Check if m above (i,j,k)::i점을 기준으로 m 과 
                            ** normal 간의 내적으로 체크(내적이 양수이면 
                            ** m이 원의 내부에 완전히 들어 있는 경우가 된다. 
                            ** 0인 경우는 원주상에 놓인 경우이므로 배제한다
                            */
                            flag &= ((x[m]-x[i])*nx + (y[m]-y[i])*ny + (z[m]-z[i])*nz <= 0);
                        }
                    }
                    if (flag) {
                        ntri++;
                        // (i, j, k)의 외접원이 다른 점을 포함하지 않으므로 이 삼각형은 
                        // 삼각분할의 한 면을 형성하게 된다.
                        // addTriangle(tri(i, j, k));
                    }
                }
    delete[] z ;
    return ntri;
}

사용자 삽입 이미지

 

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

Circle Drawing Algorithm  (1) 2008.06.03
Wu's Line Algorithm  (1) 2008.06.02
Brute Force Triangulation  (1) 2008.05.28
Polygon Triangulation (II)  (0) 2008.05.26
Polygon Triangulation  (4) 2008.05.25
Polygon Fill  (0) 2008.05.22
Posted by helloktk

댓글을 달아 주세요

  1. hgmhc 2020.12.30 15:18 신고  댓글주소  수정/삭제  댓글쓰기

    계속 보게 되네요 ㅋㅋ

728x90

monotone 폴리곤은 한 축 방향으로 꼭짓점의 좌표값이 단조 증가하는 형태를 갖는 단순 폴리곤이다. 단순 폴리곤을 monotone 폴리곤으로 분할하는 방법은 잘 알려져 있고 O(n log n)의 복잡도를 갖는다(나중에 이 부분을 자세히 리뷰할 예정이다). monotone 폴리곤은 한 방향에 대해서 이미 정렬이 되어 있기 때문이 그 방향으로 sweep-line 기법을 쓰면 선형 시간 내에 triangulation이 가능하다. 이러한 알고리즘은 이미지 처리에 적용하면 매우 효율적으로 이용이 가능하다. raster 상의 데이터는 scanline 방향으로 항상 정렬이 되어 있기 때문이다.


참고 : http://www.personal.kent.edu/~rmuhamma/Compgeometry/MyCG/PolyPart/polyPartition.htm

Algorithm Monotone_Triangulation
Input: Monotone polygon P (ccw-ness needed!)
Output: Set of triangles
(1) 주어진 폴리곤의 꼭지점을 x-좌표를 기준으로 정렬하여 V에 저장 (O(n) due to monotonicity).
(2) sweep-line 리스트(L)을 초기화 :: V의 처음 두 포인트(이벤트)를 넣음.
(3) while(V is not empty) do
       get next vertex (p) from V ;
       if (p is opposite to points in L) then
            while(size(L)>1) do
                out triangle { first(L), second(L), p }
                remove(first(L))
            end while
            add p to L;
       else 
            while ( size(L)>1 && ccw(last(L), previous(last(L)), p)) do
                 out triangle { last(L), previous(last(L)), p }
                 remove(last(L))
            end while
            add p to L
       end if
     end while/*입력이 ccw이면 ccw-test도 ccw로, cw면 test도 cw로 해야한다.*/ 

/**
**  http://blog.naver.com/helloktk/80051120809 에서 옮김.
*/

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

Polygon Fill  (0) 2008.05.22
Fortune's Sweep Algorithm  (0) 2008.05.22
Triangulating Monotone Polygon  (0) 2008.05.22
Trapezoidalization  (0) 2008.05.22
Optimizing Polygon Triangulation  (0) 2008.05.22
Sweep Line Algorithm  (0) 2008.05.22
Posted by helloktk

댓글을 달아 주세요