Artifact in Bicubic interpolation

Question

I wrote an algorithm to do bicubic interpolation of an image. I used the method desribed in the wikipedia page. On simple images, the result looks good, but on more complex ones, I got strange artifacts in non-smooth zones.

What could be a probable source of those artifacts, and what can I try to get rid of them ?

My code is on github, but I checked it several times, so I don't think the problem comes from it but rather from the simple implementation suggested in wikipedia

I'm using this particular equation (3rd-order Hermite polynomial interpolation) to do the interpolation (in one dimension):

$$x(n+t)=\frac{1}{2} \begin{bmatrix} 1 & t & t^2 & t^3 \\ \end{bmatrix} \begin{bmatrix} 0 & 2 & 0 & 0 \\ -1 & 0 & 1 & 0 \\ 2 & -5 & 4 & -1 \\ -1 & 3 & -3 & 1 \\ \end{bmatrix} \begin{bmatrix} x[n-1] \\ x[n] \\ x[n+1] \\ x[n+2] \\ \end{bmatrix}$$

where $x[n-1]$, ... $x[n+2]$, .. are the values of the surrounding pixels (in one dimension), and $t$ is the fractional coordinate (where $0 \le t < 1$) of the interpolated pixel between the adjacent pixels $x[n]$, $x[n+1]$. $n$ is an integer and $x(n)=x[n]$ for any integer $n$.

The expected output:

My result, look closely at the bird head,a lot of white pixels have appeared.

---EDIT--- To show that my code does produce coherent results, here is the example of my scaling up on simple image with 16 red pixels.

And here is the ouput of the library (which also use 3rd-order Hermite polynomial interpolation):

And here is my results:

In this case the difference are IMO totally acceptable, but it shows that my interpolation does something legit.

Answer 1

What you're seeing is acutance, as described in the wikipedia article on bicubic interpolation.

However, for downsampling an image, cubic interpolation is perhaps not the ideal choice, but makes more sense for upsampling.

Note that if you're downsampling by more than a facor of 4, there will be many pixels in the original image that are simply never used in the downsampled version.

The easiest way to get around this is to downsample by iteratively taking the mean value of 4 neighbouring pixels, until you're within a factor of 2 of the desired image size. Then you can do the last part by cubic interpolation.

If you're want high quality resampling, try taking a look at the ImageMagik documentation for inspiration. In particular they recommend using windowed filters for downsampling, such as Lanczos or Sinc.

Answer 2

I'm not sure what the problem is, but my first guess is that somehow the 16 pixels that the code is pulling from the image is getting interpreted in the wrong order.

In other words, perhaps the block of code currently

val b_1 = q(x1, x2, x3)(extracted(0), extracted(1), extracted(2), extracted(3))
val b0  = q(x1, x2, x3)(extracted(4), extracted(5), extracted(6), extracted(7))
val b1  = q(x1, x2, x3)(extracted(8), extracted(9), extracted(10), extracted(11))
val b2  = q(x1, x2, x3)(extracted(12), extracted(13), extracted(14), extracted(15))

needs to be something else, perhaps something like

val b_1 = q(x1, x2, x3)(extracted(0), extracted(4), extracted(8), extracted(12))
val b0  = q(x1, x2, x3)(extracted(1), extracted(5), extracted(9), extracted(13))
val b1  = q(x1, x2, x3)(extracted(2), extracted(6), extracted(10), extracted(14))
val b2  = q(x1, x2, x3)(extracted(3), extracted(7), extracted(11), extracted(15))

?

That sort of shuffling pixels around would give results similar to what you are seeing -- the output looks fine in "smooth" regions, but the results are more obviously wrong in non-smooth regions.

This "nonscaling" test (run the "interpolation" to generate an image of exactly the same size) is a great first test -- it should give pixel-for-pixel exactly the same output image as the input image. (Tests that expand or shrink the image, it's not as clear exactly what the output image ought to be).

I would also consider running a nonscaling test (or perhaps even a scaling test) on a very simple image, small enough that a human can single-step through every pixel in the image -- perhaps something like

0 0 0
0 5 6
0 8 9

Once you get this "nonscaling" test working, I agree with Robert that to really test the interpolation, you need to also run some tests that expand or shrink the image.