Minimum Volume Box

주어진 점집합을 감싸는 최소 면적의 bounding box(외접 직사각형)는 먼저 점집합의 Convex Hull을 구한 후 Rotating Calipers 알고리즘을 이용하여 계산할 수 있다. 다음 알고리즘은 반시계 방향으로 정렬된 Convex Hull을 입력으로 받아 최소 면적의 bounding box를 구하는 방법이다.

먼저 좌표축에 평행한 초기 bounding box를 구성한 후, Convex Hull의 네 방향에 접하는 네 개의 caliper를 생각한다. 각 caliper는 현재 접하고 있는 Convex Hull의 변과 함께 회전하며, 다음 변과 이루는 회전각을 계산한다. 네 개의 caliper에 대해 계산된 회전각 가운데 가장 작은 각만큼 네 개의 caliper를 동시에 반시계 방향으로 회전시킨다. 회전이 끝나면 새로운 방향에서의 bounding box를 계산하고 그 면적을 구한다. 이 과정을 반복하면서 면적이 가장 작은 bounding box를 선택하면 된다.

이 알고리즘은 Convex Hull을 바닥 위에서 굴리면서 각 방향에서 물체를 감싸는 bounding box를 구하는 과정으로 생각하면 직관적으로 이해하기 쉽다. Convex Hull을 한 바퀴 굴리는 동안 모든 변은 정확히 한 번씩 bounding box의 한 변과 평행하게 되므로, 최소 면적의 bounding box가 될 수 있는 모든 후보를 빠짐없이 검사하게 된다.

Convex Hull의 꼭짓점 수를 (N)이라 하면, 각 caliper는 Convex Hull의 각 변을 최대 한 번씩만 따라 이동하므로 Rotating Calipers 알고리즘의 시간 복잡도는 \(O(N)\)이다. 따라서 Convex Hull이 이미 계산되어 있다면 최소 면적의 bounding box는 선형 시간에 구할 수 있다. Convex Hull 계산까지 포함한 전체 알고리즘의 시간 복잡도는 \(O(n\log n)+O(N)\)이며, 일반적으로 \(N \le n\)이므로 전체 시간 복잡도는 \(O(n\log n)\)이 된다.

std::vector<CfPt> MinVolumeBox(const std::vector<CfPt>& chull) {
    if (chull.size() < 3) return std::vector<CfPt> ();
    std::vector<CfPt> edges(chull.size());
    std::vector<int>  visited(chull.size(), 0);  // initialized = 0;
    enum {BB_NONE, BB_LEFT, BB_RIGHT, BB_BOTTOM, BB_TOP};
    for (int i = 0; i < chull.size(); i++) {
        const CfPt e = chull[(i+1) % chull.size()] - chull[i];
        edges[i] = e / e.Norm();// make a unit vector;
    }
    // 먼저 좌표의 최소-최대 값을 구한다;
    // index을 순환할 때, 항상 최대-최소가 (index+chull.size)%chull.size가 되게
    double xmin = chull[0].x, xmax = xmin;
    double ymin = chull[0].y, ymax = curr;
    int left = 0, bot = 0, right = 0, top = 0 ;
    for (int i = 1; i < chull.size(); i++) {
        if (chull[i].x <= xmin)      { xmin = chull[i].x; left = i;}
        else if (chull[i].x >= xmax) { xmax = chull[i].x; right = i;}
        if (chull[i].y <= ymin)      { ymin = chull[i].y; bot = i;}
        else if (chull[i].y >= ymax) { ymax = chull[i].y; top = i;}
    }
    if (chull[0].x <= xmin)       { xmin = chull[0].x; left = 0;}
    else if ( chull[0].x >= xmax) { xmax = chull[0].x; right = 0;}
    if (chull[0].y <= ymin)       { ymin = chull[0].y; bot = 0;}
    else if (chull[0].y >= ymax)  { ymax = chull[0].y; top = 0;}

    /*____Initial mvbox Setting____*/
    CfPt center = CfPt(xmin + xmax, ymin + ymax)/2;
    CfPt ax0 = CfPt(1,0), ax1 = CfPt(0,1);	// mvbox의 unit vectors;
    CfPt eh = ax0, ev = ax1;				// working unit vectors;
    double hscale  = (xmax - xmin) / 2; 
    double vscale  = (ymax - ymin) / 2;
    double minarea = hscale * vscale;		// area / 4;
    int done = 0;
    while (!done) {
        // edges[bot]는 chull[bot+1]-chull[bot];
        // edges[bot],edges[right],edges[top],edges[left] 중에서
        // 현재의 mvbox의 변(eh, ev,-eh,-ev)과 가장 작은 각(max_cosine)을 이루는 
        // edge를 찾음;
        int iflag = BB_NONE;
        double	max_cosine = 0;
        double	cosine = eh * edges[bot];
        if (cosine > max_cosine) { max_cosine = cosine; iflag = BB_BOTTOM;}

        cosine = ev * edges[right]; ;
        if (cosine > max_cosine) { max_cosine = cosine; iflag = BB_RIGHT;}
        
        cosine = -eh * edges[top]; 
        if (cosine > max_cosine) { max_cosine = cosine; iflag = BB_TOP;}
        
        cosine = -ev * edges[left]; 
        if (cosine > max_cosine) { max_cosine = cosine; iflag = BB_LEFT;}
        
        switch (iflag) {
        case BB_BOTTOM:
            if (!visited[bot]) {// edges[bot]가 mvbox의 한변에 포함됨;
                eh = edges[bot]; ev = eh.Rot90(); 
                // 다음 변으로 회전;
                visited[bot] = 1;
                if (++bot == chull.size()) bot = 0;
            } else done = 1;
            break;
        case BB_RIGHT:
            if (!visited[right]) {// edges[right]가 mvbox의 한변에 포함됨;
                ev = edges[right]; eh = -ev.Rot90(); 
                // 다음 변으로 회전;
                visited[right] = 1;
                if (++right == chull.size()) right = 0;
            } else done = 1;
            break;
        case BB_TOP:
            if (!visited[top]) { // edges[top]이 mvbox의 한변에 포함됨;
                eh = -edges[top]; ev = eh.Rot90(); 
                // 다음 변으로 회전;
                visited[top] = 1;
                if (++top == chull.size()) top = 0;
            } else done = 1;
            break;
        case BB_LEFT:
            if (!visited[left]) {// edges[left]가 mvbox의 한변에 포함됨;
                ev = -edges[left]; eh = -ev.Rot90(); 
                // 다음 변으로 회전;
                visited[left] = 1;
                if (++left == chull.size()) left = 0;
            } else done = 1;
            break;
        case BB_NONE: // polygon이 직사각형임;
            done = 1;
            break;
        }
        // check area of mvbox;
        double len0 = eh * (chull[right] - chull[left]) / 2;
        double len1 = ev * (chull[top] - chull[bot]) / 2;
        double area = len0 * len1; //면적의 1/4 임;
        if (area < minarea) {
            minarea	= area;
            ax0	= eh; ax1 = ev;
            hscale = len0; vscale = len1;
            // box center: tt = left 기준 top-bot을 연결하는 vector/2;
            CfPt tt = (chull[top] + chull[bot]) / 2 - chull[left];
            double zz = ev * tt; 
            center = eh * len0 + ev * zz + chull[left];
        }
    }
    eh = ax0 * hscale; // recover original scale;
    ev = ax1 * vscale; // recover original scale;
    std::vector<CfPt> mvbox(4);
    mvbox[0] = center + eh + ev;
    mvbox[1] = center - eh + ev;
    mvbox[2] = center - eh - ev;
    mvbox[3] = center + eh - ev;
    return mvbox;
}