image processing - Bicubic Interpolation

I'm trying to do bicubic interpolation on an 8*8 matrix(image) shown below.

A = 10   20  30  40  50  60  70  80
    90  100 110 120 130 140 150 160
    170 180 190 200 210 220 230 240
    250  10  20  30  40  50  60  70
     80  90 100 110 120 130 140 150
    160 170 180 190 200 210 220 230
    240 250  10  20  30  40  50  60

     70  80  90 100 110 120 130 140

I used the following matlab code for it and I got an output

interp2(A,'bicubic')

But even after googling I'm not getting how to do it manually. Can anyone help me to get the idea using this example(only for few pixels)?

Answer

Here pretty good explanation.

I'll start by consider the 1D case of cubic interpolation (because its easier to explain) and then go onto the 2D case.

The basic idea of cubic interpolation is to estimate the values between any two points by as a cubic function, $f(x)$, with first derivative, $f'(x)$

$$f(x) = ax^3 + bx^2 + cx + d$$ $$f'(x) = 3ax^2 + 2bx + c$$

so that, when spliced to the adjacent cubic splines on the left and on the right, the interpolated function (the "zeroth derivative") and its first derivative are continuous. That means, at the left boundary of some given spline, the zeroth and first derivatives are equal to the those of the right boundary of the adjacent spline immediately to the left.

Without loss of generality, we can assume that we are interpolating between $x=0$ and $x=1$. On the left and right boundaries, we see that

$$f(0) = d$$ $$f(1) = a+b+c+d$$ $$f'(0) = c$$ $$f'(1) = 3a+2b+c$$

Solving for the polynomial coefficients

$$a = 2f(0) - 2f(1) + f'(0) + f'(1)$$ $$b = -3f(0) + 3f(1) - 2f'(0) - f'(1)$$ $$c = f'(0)$$ $$d = f(0)$$

Now consider the four points around the region being estimated, which are $(-1,y_{-1})$, $(0,y_0)$, $(1,y_1)$, $(2,y_2)$, and impose boundary conditions that

$$f(0) = y_0$$ $$f(1) = y_1$$

and define the two boundary derivatives as

$$f'(0) = \frac{y_1-y_{-1}}{2}$$ $$f'(1) = \frac{y_2-y_0}{2}$$

which is consistent to how they will be defined for the adjacent splines.

Then

$$a = -\frac{y_{-1}}{2} + \frac{3y_0}{2} - \frac{3y_1}{2} + \frac{y_2}{2}$$ $$b = y_{-1} - \frac{5y_0}{2} + 2y_1 - \frac{y_2}{2}$$ $$c = -\frac{y_{-1}}{2} + \frac{y_1}{2}$$ $$d = y_0$$

once the coefficients are known you can interpolate $f(x)$ for any point between the central 2 points, $(0,y_0)$ and $(1,y_1)$.

For bicubic interpolation the principle is pretty much the same but you estimate a surface using 16 points (4x4 grid) rather than just a curve. One simple way to do this is to first interpolate the columns and then interpolate the resulting rows.