Median-Cut 컬러 양자화

컬러 영상을 표현하는 가장 일반적인 방법은 RGB 컬러 모델이다. RGB 컬러에서는 빨강(\(R\)), 초록(\(G\)), 파랑(\(B\)) 각각에 8비트를 할당하여 256단계의 밝기를 표현한다. 따라서 전체적으로 \(256\times256\times256=16,777,216\)가지의 컬러를 나타낼 수 있다. 그러나 실제 영상에서는 이처럼 많은 컬러가 모두 사용되는 경우는 드물며, 대부분의 픽셀은 비교적 적은 수의 색상에 집중되어 분포한다. 따라서 많은 응용에서는 24비트 컬러를 그대로 사용하는 것보다 대표적인 컬러만을 이용하여 영상을 표현하면 메모리 사용량과 연산량을 크게 줄일 수 있다.

RGB 컬러 영상은 각 픽셀을 \((R,G,B)\)로 이루어진 하나의 벡터로 생각하여 3차원 \(RGB\) 공간의 한 점으로 표현할 수 있다. 실제 영상에서는 픽셀들이 \(RGB\) 공간 전체에 균일하게 분포하는 경우는 드물며, 대부분 비슷한 색을 갖는 픽셀들이 여러 개의 군집(cluster)을 이루어 분포한다. 따라서 같은 군집에 속하는 픽셀들을 하나의 대표 컬러로 대체하면, 원래는 수많은 점으로 이루어져 있던 컬러 공간이 소수의 대표 컬러만으로 표현된다. 이 과정을 컬러 양자화(color quantization)라고 하며, 대표 컬러들을 Lookup Table(또는 colormap)에 저장하고 각 픽셀은 해당 컬러의 인덱스만 저장하면 적은 메모리로도 원본과 유사한 영상을 표현할 수 있다.

이제 문제는 \(RGB\) 공간의 컬러를 어떻게 효율적으로 군집화할 것인가로 바뀐다. 가장 단순한 방법은 컬러 공간을 일정한 크기의 격자로 나누는 것이지만, 실제 컬러 분포를 반영하지 못하므로 비효율적이다. 보다 널리 사용되는 방법이 Median-Cut 알고리즘이다.

Median-Cut 알고리즘에서는 먼저 영상에 존재하는 모든 컬러를 포함하는 하나의 box를 만든다. 그런 다음 현재 박스에서 가장 길이가 긴 축(\(R\), \(G\), 또는 \(B\)축)을 선택하고, 그 축 방향으로 투영한 컬러 분포의 median을 기준으로 박스를 둘로 분할한다. 이렇게 하면 두 박스가 대략 비슷한 수의 픽셀을 포함하게 된다(한 축에 대한 메디안이므로 정확히 절반으로 나뉘지는 않는다). 이후에는 일반적으로 가장 많은 픽셀을 포함하거나 우선순위가 가장 높은 박스를 다시 같은 방법으로 분할한다. 원하는 개수의 박스가 만들어질 때까지 이 과정을 반복하면 최종적으로 원하는 수의 대표 컬러를 얻을 수 있다.

\(RGB\) 공간 전체를 그대로 사용하면 히스토그램의 크기가 너무 커져 메모리 사용량과 계산량이 증가한다. 따라서 실제 구현에서는 \(RGB\)의 하위 비트를 일부 버려 컬러 공간을 먼저 축소한다. 예를 들어 각 채널의 하위 3비트를 제거하면 \(R\), \(G\), \(B\)가 각각 5비트가 되어 가능한 컬러 수는 \(32\times32\times32=32768\) 개로 줄어든다. 이렇게 하면 히스토그램을 훨씬 작은 메모리로 관리할 수 있어 알고리즘의 수행 속도가 크게 향상된다.

최종적으로 각 박스에 포함된 픽셀들의 평균 \(RGB\) 값(centroid)을 대표 컬러로 선택한 뒤, 원본 영상의 각 픽셀은 \(RGB\) 공간에서 가장 가까운 대표 컬러(유클리드 거리가 가장 작은 컬러)로 대체된다. 이 과정을 통해 연속적으로 분포하던 컬러 공간은 소수의 대표 컬러로 양자화된다.

그러나 대표 컬러의 개수가 너무 적으면 인접한 픽셀 사이의 색상 변화가 계단 형태로 나타나는 색 띠(banding) 현상이 발생할 수 있다. 이러한 시각적인 왜곡을 줄이기 위해서는 Floyd-Steinberg dithering과 같은 dithering 기법을 적용하여 양자화 오차를 주변 픽셀로 분산시키는 후처리를 함께 사용하는 경우가 많다.

 

**네이버 블로그에서 이전;

source code: dev.w3.org/Amaya/libjpeg/jquant2.c

 

int MedCutQuantizer(CRaster& raster, CRaster& out8) {
    const int MAX_CUBES = 256;
    const int HSIZE =32 * 32 * 32; // # of 15-bit colors;
    int hist[HSIZE] = {0};
    CSize sz = raster.GetSize();
    for (int y = 0; y < sz.cy; y++) {
        BYTE *p = (BYTE *)raster.GetLinePtr(y);
        for (int x = 0; x < sz.cx; x++) {
            int b = *p++; int g = *p++;	int r = *p++;
            hist[bgr15(b, g, r)]++;
        }
    }
    std::vector<int> hist_ptr;
    for (int color15 = 0; color15 < HSIZE; color15++)
        if (hist[color15] != 0) { 
            // 0이 아닌 히스토그램의 bin에 해당하는 15비트 컬러를 줌;
            hist_ptr.push_back(color15);
        }
    
    std::vector<Cube> cubeList;
    Cube cube(0, hist_ptr.size()-1, sz.cx * sz.cy, 0);
    shrinkCube(cube, &hist_ptr[0]);
    cubeList.push_back(cube);
    while (cubeList.size() < MAX_CUBES) {
        int level = 255, splitpos = -1; 
        for (int k = cubeList.size(); k-->0;) {
            if (cubeList[k].lower == cubeList[k].upper) ;	// single color; cannot be split
            else if (cubeList[k].level < level) {
                level = cubeList[k].level;
                splitpos = k;
            }
        }
        if (splitpos == -1)	// no more cubes to split
            break;
        // Find longest dimension of this cube
        Cube cube = cubeList[splitpos];
        int dr = cube.rmax - cube.rmin;
        int dg = cube.gmax - cube.gmin;
        int db = cube.bmax - cube.bmin;
        int longdim = 0;
        if (dr >= dg && dr >= db) longdim = 0;
        if (dg >= dr && dg >= db) longdim = 1;
        if (db >= dr && db >= dg) longdim = 2;
        // 가장 넓게 분포한 color를 colorTab의 15비트 컬러의 상위비트 이동;
        reorderColors(&hist_ptr[0], cube.lower, cube.upper, longdim);
        // hist_ptr값을 증가 순서로 바꿈->가장 폭이 넓은 컬러 방향으로 정렬됨;
        quickSort(&hist_ptr[0], cube.lower, cube.upper);
        // hist_ptr을 다시 원컬러로 복구;
        restoreColorOrder(&hist_ptr[0], cube.lower, cube.upper, longdim);
        // Find median
        int count = 0;
        int med = 0;
        for (med = cube.lower; med < cube.upper; med++) {
            if (count >= cube.count/2) break;
            count += hist[hist_ptr[med]];
        }
        // cube에 들어있는 color의 median을 기준으로 두개 cube로 쪼김;
        // 낮은 쪽은 원래 cube를 대체하고 높은 쪽은 새로 list에 추가;
        Cube cubeLow(cube.lower, med - 1, count, cube.level + 1);
        shrinkCube(cubeLow, &hist_ptr[0]);
        cubeList[splitpos] = cubeLow;				// add in old slot
        // 
        Cube cubeHigh(med, cube.upper, cube.count - count, cube.level + 1);
        shrinkCube(cubeHigh, &hist_ptr[0]);
        cubeList.push_back(cubeHigh);	
    }
    // 각 cube에 포함된 컬러를 하나의 대표값으로 설정; 
    // 대표값은 cube내 컬러의 평균이다;
    BYTE rLUT[256], gLUT[256], bLUT[256];
    for (int k = cubeList.size(); k-->0;) {
        int rsum = 0, gsum = 0, bsum = 0;
        Cube &cube = cubeList[k];
        for (int i = cube.lower; i <= cube.upper; i++) {
            int color15 = hist_ptr[i];
            rsum += getRed(color15) * hist[color15];
            gsum += getGreen(color15) * hist[color15];
            bsum += getBlue(color15) * hist[color15];
        }
        rLUT[k] = int(rsum / (double)cube.count + .5);
        gLUT[k] = int(gsum / (double)cube.count + .5);
        bLUT[k] = int(bsum / (double)cube.count + .5);
    }
        
    // 각 cube에 포함된 color를 대표하는 컬러값을 쉽게 찾도록 lookup table
    // 설정함;
    BYTE inverseMap[HSIZE]; 
    for (int k = cubeList.size(); k-->0;) {
        // 각 cube 내의 컬러는 하나의 대표컬러로...;
        Cube& cube = cubeList[k];
        for (int i = cube.lower; i <= cube.upper; i++)
            inverseMap[hist_ptr[i]] = (BYTE)k;
    }	

    // 8-bit indexed color 이미지 생성;
    out8.SetDimensions(sz, 8);
    for (int y = 0; y < sz.cy; y++) {
        BYTE *p = (BYTE *)raster.GetLinePtr(y);
        BYTE *q = (BYTE *)out8.GetLinePtr(y);
        for (int x = 0; x < sz.cx; x++) {
            int b = *p++; int g = *p++;	int r = *p++;
            *q++ = inverseMap[bgr15(b, g, r)];
        }
    }
    // 8비트 이미지 팔레트 설정;
    for (int i = 0; i < 256; i++) {
        out8.m_palette[i].rgbBlue = bLUT[i];
        out8.m_palette[i].rgbGreen = gLUT[i];
        out8.m_palette[i].rgbRed = rLUT[i];
    }
    // SaveRaster(out8, "res_medcut.bmp");
    return 1;
}