Cubic Spline Kernel

일정한 간격 $h$마다 샘플링된 데이터 $\{ (x_k, f_k) \}$를 이용해서 이들 데이터를 표현하는 spline를 구해보자. 일반적으로 spline은 샘플링 데이터와 kernel이라고 불리는 함수의 convolution sum 형태로 표현할 수 있다.

$$ g(x) = \sum_k  f_k K \left( \frac{x-x_k}{h}\right)$$ kernel의 종류에 따라 재구성된 spline은 원래 샘플을 모두 통과하는 interpolation이거나, 샘플을 반드시 통과하지 않는 approximation이 될 수도 있다. 이미지의 resampling 과정에서는 다양한 spline kernel이 사용되며, 각각 서로 다른 특성을 갖는다.

 

3차 spline kernel은 중심을 기준으로 support가 \([-2,2]\)인 구간에서만 0이 아닌 piecewise cubic polynomial로 정의되며, 우함수의 성질을 갖는다. 따라서 일반적인 형태는

$$ K(s) = \left\{ \begin{matrix} A_1|s|^3 + B_1 |s|^2 +C_1 |s| + D_1    &  |s| <1 \\ A_2 |s|^3 + B_2 |s|^2 + C_2 |s| + D_2 & 1 \le |s|<2 \\ 0 & \text{otherwise} \end{matrix} \right. $$와 같이 쓸 수 있다. 계수를 결정하기 위해서는 8개의 조건이 필요한다. 우선 \(K(s)\)가 우함수이므로 원점에서 미분값이 제대로 정의되려면 $$C_1=0$$도 만족해야 한다, 그리고 각 node에서 연속성을 요구하면

\begin{align} s=1^\pm:~~~& A_1+B_1 +D_1 = A_2 +B_2+C_2 +D_2 \\ s=2^\pm:~~~& 8A_2 +4B_2 +2C_2 +D_2 =0 \end{align}을 만족해야 한다. 곡선이 부드럽게 연결되도록 요구하면

\begin{align} s = 1^\pm:~~~& 3A_1 + 2B_1 = 3A_2 + 2B_2+C_2 \end{align} 여기에 더하여 2차 도함수 역시 \(s=1\)에서 연속이라고 요구하면

\begin{align} 6A_1 + 2B_1 = 6A_2 + 2B_2\end{align}이 추가된다.

또 하나의 중요한 조건으로는 샘플링된 데이터가 모두 같은 경우 보간함수도 상수함수가 되는 것이 타당하므로

$$ g(x) = \sum_k K \left( \frac{x-x_k}{h}\right) = 1~~~\text{if} ~~\forall f_k = 1$$

을 만족시켜야 한다. $x_j <x<x_{j+1}$일 때 $x  =  x_j + sh,  ~(0< s<1)$이라고 쓰면, kernel의 support가 \([-2, 2]\)이므로 

$$K(s+1) + K(s) +  K(s-1) + K(s-2)=1$$을 만족해야 한다. 이를 kernel 식에 대입해 정리하면 다음과 같은 항등식을 얻는다.

\begin{gather} -1 + A_1 + 9 A_2 + B_1 + 5 B_2 + 3 C_2 + 2 D_1 +2 D_2 \\ +(-3 A_1 - 9 A_2 - 2 B_1 - 2 B_2) s + (3 A_1 + 9 A_2 + 2 B_1 + 2 B_2) s^2 =0\end{gather}

항등식으로부터 세 개의 조건을 얻지만, 이 중 하나는 앞에서 얻은 조건과 독립적이지 않다. 따라서 최종적으로 독립적인 조건은 모두 여섯 개가 되고, 두 개의 자유 매개변수가 남는다.  보통 이 두 계수는 $D_1=1-B/3$,  $D_2=4C + 4B/3$처럼 매개화한다. 이 경우 kernel은

$$ K(s)  = \frac{1}{6} \left\{ \begin{matrix}   (12-9B-6C)|s|^3 +(-18+12B+6C) |s|^2 + (6-2B) & |s|< 1\\ (-B-6C)|s|^3 +(6B+30C)|s|^2 + (-12B-48C) |s| + (8B+24C) & 1\le |s|<2\\0 & \text{otherwise} \end{matrix} \right.$$

따라서 cubic spline kernel은 두 개의 파라미터 \((B,C)\)에 의해서 정해진다. 또한 kernel의 적분은 (B,C)(B,C)(B,C)와 무관하게 항상 1이므로 정규화(normalization)되어 있으며, 앞에서 유도한 partition of unity 조건에 의해 재구성 시 가중치의 합 역시 항상 1이 된다.

$$\int_{-\infty}^\infty K(s)ds =1$$

이 중에는 이미지의 resampling에서 많이 사용되는 커널도 있는데, 잘 알려진 경우를 보면

$$ \begin{matrix} (B,C)=(0,1) & \text{Cardinal spline} \\ (B,C)=(0,1/2) & \text{Catmull-Rom spline } \\ (B,C)=(0,3/4) & \text{used in photoshop} \\ (B,C)=(1/3,1/3) & \text{ Mitchell-Netravali spline}  \\ (B,C)=(1,0) & \text{B-spline}\end{matrix}$$

특히 $B=0$인 경우는 $K(0)=1$이고, $K(\pm1)=K(\pm2)=0$이므로 $$K(i-j)=\delta_{ij}$$를 만족하는 interpolation kernel에 해당한다. 그중에서도 $(B,C)=(0, \frac{1}{2})$인 Catmull-Rom spline은 가장 널리 사용되는 cubic interpolation kernel 중 하나이다. 이 kernel은 충분히 매끄러운 원래 함수에 대해 3차 정확도를 가지므로 재구성 오차가 $O(h^4)$이며, 계산량과 정확도 사이의 균형이 뛰어나 이미지의 확대, 축소, 회전과 같은 resampling에서 널리 사용된다. 반면 $(B,C)=(1,0)$인 B-spline은 매우 부드러운 결과를 제공하지만 interpolation kernel이 아니므로 원래 샘플값을 정확히 재현하지는 않는다. 따라서 보다 매끄러운 영상을 얻을 수 있는 대신, 영상이 다소 흐려지는 특성을 갖는다.

// Mitchell Netravali Reconstruction Filter
// B = 0    C = 0   - Hermite B-Spline interpolator 
// B = 0,   C = 1/2 - Catmull-Rom spline
// B = 1/3, C = 1/3 - Mitchell Netravali spline
// B = 1,   C = 0   - cubic B-spline
double MitchellNetravali(double x, double B, double C) {
    x = fabs(x);
    if (x >= 2) return 0;
    double xx = x*x;
    if (x >= 1) return ((-B - 6*C)*xx*x 
                + (6*B + 30*C)*xx + (-12*B - 48*C)*x 
                + (8*B + 24*C))/6;
    if (x < 1) return ((12 - 9*B - 6*C)*xx*x +
        (-18 + 12*B + 6*C) * xx + (6 - 2*B))/6;
}