Kuwahara Filter

Kuwahara filter는 에지를 보존하면서 노이즈를 제거하는 대표적인 비선형 필터이다. 먼저 각 픽셀을 중심으로 한 윈도를 설정한 뒤, 중앙 픽셀이 각각의 한 꼭짓점이 되도록 정의된 네 개의 겹치는 subregion으로 나눈다. 각 부분영역의 평균과 분산을 계산한 후, 분산이 가장 작은 영역의 평균값으로 중앙 픽셀의 값을 대체한다.

균일한 영역에서는 주변 픽셀들의 분산이 작으므로 spike noise는 주변 영역의 평균값으로 대체되어 효과적으로 제거된다. 반면 에지 부근에서는 에지를 가로지르는 부분영역의 분산은 크고, 에지를 포함하지 않는 부분영역의 분산은 상대적으로 작다. 따라서 분산이 가장 작은 부분영역이 선택되어 에지를 넘는 평균화가 일어나지 않으므로, 에지의 위치와 형태를 비교적 잘 보존할 수 있다.

예를 들어 윈도의 크기가 \(5×5\)인 경우에는 네 개의 \(3×3\) 부분영역에서 각각 평균과 분산을 계산하여 중앙 픽셀의 값을 결정한다. 이때 integral image를 이용하면 각 부분영역의 픽셀 합을 매우 빠르게 계산할 수 있으며, 픽셀 제곱값에 대한 integral image를 함께 사용하면 분산 역시 반복 계산 없이 효율적으로 구할 수 있다. 따라서 Kuwahara filter의 계산량을 크게 줄일 수 있다.

center_pixel = mean[argmin[variance[i]]]

원본 이미지, 필터링된 이미지(경계부분은 미처리됨)

Kuwahara filter보다 더 향상된 filter로는 symmetric nearest neighbour (SNN) filter가 있다. 이 filter는 기준 픽셀에서 대칭 지점에 있는 두 지점의 픽셀 값 중에서 중앙 픽셀 값에 가까운 것을 취해서 평균값을 취한다.

 

아래 코드는 gray 이미지에 대해 Kuwahara filter을 구현한 것으로 컬러 이미지의 경우는 RGB 채널별로 따로 적용하면 된다. 그리고 경계 영역에서는 처리를 하지 않도록 만들었다.

#define GET_SUM(x,y)  sum[(y)*width+(x)]
#define GET_SUM2(x,y) sum2[(y)*width+(x)]
#define BLK_SUM(x,y) ((GET_SUM((x),(y)) + GET_SUM((x)-hs1,(y)-hs1) \
                        - GET_SUM((x)-hs1,(y)) - GET_SUM((x),(y)-hs1)))
#define BLK_SUM2(x,y) ((GET_SUM2((x),(y)) + GET_SUM2((x)-hs1,(y)-hs1) \
                        - GET_SUM2((x)-hs1,(y)) - GET_SUM2((x),(y)-hs1)))
void integral_image(BYTE *image, int width, int height, int *sum);
void integral_image_sq(BYTE *image, int width, int height, int *sum2);
// output -> input;
void kuwahara_filter (std::vector<BYTE> &image, int width, int height, int size) 
{
    // integral image; for window mean;
    std::vector<int> sum(width * height, 0);
    // integral image of square; for window variance;
    stt::vector<int> sum2(width * height, 0);
    
    for (int y = 0, yw = 0; y < height; y++, yw += width ) {
        int s = 0, ss = 0;
        for (int x = 0; x < width; x++ ) {
            int a = image[yw + x];
            s += a; ss += a * a;
            sum[yw + x] = s;
            sum2[yw + x] = ss;
        }
    }
    for (int x = 0; x < width; x++ )
        for (int y = 1, yw = width; y < height; y++, yw += width ) {
            sum[yw + x]  += sum[yw - width + x];
            sum2[yw + x] += sum2[yw - width + x];
        }
    
    // do filtering;
    size = ((size >> 1) << 1) + 1; // make window size be an odd number;
    const int offset = ( size - 1 ) >> 1;
    const int hs1 = ( size + 1 ) >> 1;
    const int nn = hs1 * hs1;    // # of pixels within the  window;
    const int dx[] = {0, offset, 0, offset};
    const int dy[] = {0, 0, offset, offset};
    int var[4], mean[4];
    for (int y = height - hs1; y-->hs1;)
        for (int x = width - hs1; x-->hs1;) {
            int vmin = INT_MAX, winner = 0;
            for (int k = 0; k < 4; ++k) {
                mean[k] = BLK_SUM(x+dx[k], y+dy[k])/nn;
                var[k] = BLK_SUM2(x+dx[k], y+dy[k])/nn - mean[k]*mean[k];
                if (vmin > var[k])
                    vmin = var[k], winner = k;
            }
            int a = mean[winner];
            image[y*width + x] = a > 255 ? 255 : a < 0 ? 0 : a;
        }
};