Nearest neighbor interpolation does not give expected results

[gerhardneuhold]

System information (version)

• OpenCV => 2.4.9
• Operating System / Platform => Ubuntu 16.04
• Compiler => gcc 5.4.0

Detailed description

Nearest neighbor interpolation using cv2.resize does not give expected results.

Steps to reproduce

import cv2
import numpy as np

a = np.array([[0, 1, 2, 3, 4]], dtype=np.uint8)
print 'nearest', cv2.resize(a, dsize=(3, 1), interpolation=cv2.INTER_NEAREST)

gives

nearest [[0 1 3]]

expected

nearest [[0 2 4]]

[sergiud]

How did you come up with expected result?

The downsampling factor is 3/5 (in horizontal direction) meaning destination columns 0, 1, 2 are mapped to source columns with indices 0, 1 * 5/3 and 2 * 5/3 (i.e., 0, 1.67 and 3.33). After rounding half up (probably what you expect), the values of the corresponding columns are 0, 2, 3.

Either way, it seems OpenCV rounds towards zero instead of rounding half away from zero (not verified).

[gerhardneuhold]

You can compute the coordinates as following:

np.floor((np.arange(3) + 0.5) / 3. * 5.)

to get the expected result. PIL works accordingly:

from PIL import Image
print np.asarray(Image.fromarray(a).resize((3,1), Image.NEAREST))

gives:

[[0 2 4]]

[sergiud]

Why do you shift by 0.5 before transforming the coordinates? Do you assume the center of the pixel to be at .5 fraction of the pixel coordinate? I don't know much about conventions used by PIL.

np.round((np.arange(3)) / 3. * 5.) produces [[0, 2, 3]].

[gerhardneuhold]

Exactly, treating pixels as squares with center at .5.

[gerhardneuhold]

@sergiud The issue is that opencv's resize function and your coordinate example lead to a slight shift of the pixels to the right/bottom of the image in addition to scaling, hence not only performs scaling but also a translation of the pixels.

@alalek What is the reason for the question (invalid) label?

[sergiud]

@gerhardneuhold I agree.

[alalek]

There is no exact formula in documentation (and probably will not be there for general case - for performance reason and rounding issues). Existed resize tests are passed. So this issue looks like the question about the used resize formula in OpenCV.

[homm]

Just copy here from this thread:

Results of reading and resizing can be different in cv2 and Pilllow. This creates problems when you want to reuse a model (neural network) trained using cv2 with Pillow.

import cv2
from PIL import Image

im = cv2.imread('_in.png', cv2.IMREAD_COLOR)
cv2.resize(im, (3, 3), interpolation=cv2.INTER_NEAREST)
cv2.imwrite('_out.cv2.png', im)

im = Image.open('_in.png')
im = im.resize((3, 3), Image.NEAREST)
im.save('_out.pil.png')

Please, look at the sample:

This image with the uniform gradient (from 100% white to 100% black) allows us to find out which pixels are used by each library. This is the same image after resizing to (3, 3). Left is CV2, right is Pillow:

OpenCV uses the topmost left white pixel from the source image, but the bottommost right pixel on the result is too bright. Here are maps of pixels on the source image:

It's clear to me there are some problems with rounding in OpenCV there.

[Kirill888]

@sergiud, @alalek I just tried this on my machine (version 3.3.0 on a mac via python bindings) and the behaviour is as described by @gerhardneuhold. He is quite right that your earlier example from July introduces sub-pixel translation, and if this is the logic that is implemented in the code than it needs to be fixed.

When resizing the image it's not helpful to think in terms of pixel centers, you have to think in terms of image edges, and pixel spans. Pixel 0 of the destination image spans from the left edge of the source image pixel 0 and into pixel 1, middle pixel of the destination image spans spans from pixel 1 all the way into 3, see diagram below:

         0       1       2       3       4
  5: |---+---|---+---|---+---|---+---|---+---|
phy: |~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~|
  3: |------+-----|------+------|------+-----|
            0'           1'             2'

So if we assume pixel centers are at 0,0, then expected coords should be 1/3, 2, 3+2/3. Remember that source and destination images should span the same "physical region".

What current code seems to be doing is this, instead

            0       1       2       3       4
     5: |---+---|---+---|---+---|---+---|---+---|
   phy: |~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~|
  3: |------+-----|------+------|------+-----|
            0'           1'             2'

So you are sampling from a slightly different "physical region", introducing sub-pixel translation.

I don't think those differences come from rounding behaviour differences or anything like that, it comes from subpixel translation introduced by resize method. The error is in your coordinate computation math. And yes resize should return 0,2,4 in the example above and not 0,1,3.

[crackwitz]

using warpAffine with custom matrix respecting pixel coordinate flavors would be a workaround. notice that I use INTER_LINEAR so you can see roughly where sampling is done. INTER_LINEAR appears to quantize a little; the values should be exact multiples of thirds.

The convention whether pixels are points or whether pixels are areas is a common debate in computer graphics. It's the same issue in OpenGL, Direct3D, ...

import numpy as np
import cv2 as cv

src = np.float32([[0, 1, 2, 3, 4]])
(srows, scols) = src.shape

(dcols, drows) = dsize = (3,1)

# some composition functions

def translate(tvec):
	H = np.eye(3)
	H[0,2] = tvec[0]
	H[1,2] = tvec[1]
	return H

def scale(svec):
	H = np.eye(3)
	H[0,0] = svec[0]
	H[1,1] = svec[1]
	return H

# don't need rotation, diy
# or use getRotationMatrix2D(center, angle, scale)

def project(H):
	# it's a projective space,
	# where (x,y) is represented by all (x,y,1)*w coordinates, i.e. a ray through (x,y,1)
	# these are homogeneous coordinates
	# project onto w=1
	M = H / H[2,2]
	
	# check that this isn't a perspective transform, affine only
	assert np.allclose(M[2,0:2], 0)
	# if you want perspective transforms, return the full 3x3 matrix and use warpPerspective

	return M[0:2, 0:3]

#                 0.0         1.0         2.0         3.0         4.0
#     5:     |-----+-----|-----+-----|-----+-----|-----+-----|-----+-----|
#   phy:     |~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~|
#     3: |---------+---------|---------+---------|---------+---------|
#                  0'                  1'                  2'

# define dst->src transformation
# for sampling dest point value in src space
M = project(
	scale([scols/dcols, srows/drows])
)

dst1 = cv.warpAffine(src, M, dsize=dsize, flags=cv.WARP_INVERSE_MAP | cv.INTER_LINEAR)
print(dst1)
# => [[0.      1.65625 3.34375]]
dst1 = cv.warpAffine(src, M, dsize=dsize, flags=cv.WARP_INVERSE_MAP | cv.INTER_NEAREST)
print(dst1)
# => [[0. 2. 3.]]


#             0.0         1.0         2.0         3.0         4.0
#     5: |-----+-----|-----+-----|-----+-----|-----+-----|-----+-----|
#   phy: |~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~|
#     3: |---------+---------|---------+---------|---------+---------|
#                  0'                  1'                  2'

# define dst->src transformation
# for sampling dest point value in src space
M = project(
	translate([-0.5, -0.5]) @
	scale([scols/dcols, srows/drows]) @
	translate([+0.5, +0.5])
)

dst2 = cv.warpAffine(src, M, dsize=dsize, flags=cv.WARP_INVERSE_MAP | cv.INTER_LINEAR)
print(dst2)
# => [[0.34375 2.      3.65625]]
dst2 = cv.warpAffine(src, M, dsize=dsize, flags=cv.WARP_INVERSE_MAP | cv.INTER_NEAREST)
print(dst2)
# => [[0. 2. 4.]]

[angelo-peronio]

Possibly fixed by cv::InterpolationFlags::INTER_LINEAR_EXACT cv::InterpolationFlags::INTER_NEAREST_EXACT in #18053 .

[gerhardneuhold]

Possibly fixed by cv::InterpolationFlags::INTER_LINEAR_EXACT in #18053 .

Great - just for completeness, the related interpolation mode for this issue here is cv2.INTER_NEAREST_EXACT (instead of INTER_LINEAR_EXACT).

[crackwitz]

User error. You aren't passing the INTER value in the correct position. Please direct your usage questions to the forum or Stack Overflow.

[KongCang]

User error. You aren't passing the INTER value in the correct position. Please direct your usage questions to the forum or Stack Overflow.

Thanks for your reply! it's really my mistake, i forgot to add fx,fy parms this leading INTER value in wrong position.

[dkurt]

Can this issue be closed? The story about resize compatibility across the libraries is pretty old. However I believe that we have finished with the described problem by this code:

import cv2
import numpy as np

a = np.array([[0, 1, 2, 3, 4]], dtype=np.uint8)
print('nearest', cv2.resize(a, dsize=(3, 1), interpolation=cv2.INTER_NEAREST_EXACT))

INTER_NEAREST_EXACT | Bit exact nearest neighbor interpolation. This will produce same results as the nearest neighbor method in PIL, scikit-image or Matlab.

https://docs.opencv.org/4.x/da/d54/group__imgproc__transform.html#gga5bb5a1fea74ea38e1a5445ca803ff121a9c081e5a4d3685625a5e8a209d8cdab5

[homm]

The story about resize compatibility across the libraries is pretty old.

Basically this is not a story about compatibility, this is story of wrong behavior.

Bit exact nearest neighbor interpolation.

What does mean by "bit exact" if neighbor interpolation doesn't isn't interpolation actually and don't operates pixels' values? It should be "Correct nearest neighbor interpolation".

This will produce same results as the nearest neighbor method in PIL, scikit-image or Matlab.

It will not produce the same results as PIL since PIL originally had the same issue and it was fixed in Pillow (python-pillow/Pillow#2022) which is different library.

[dkurt]

Oh, despite #9096 (comment), it worked only for odd dimensions at input:

import cv2
import numpy as np
from skimage.transform import resize

for l in range(5, 9):
    for t in (3, 4):
        a = np.arange(l, dtype=np.uint8).reshape(1, l)
        outCV = cv2.resize(a, dsize=(t, 1), interpolation=cv2.INTER_NEAREST_EXACT)
        outScikit = resize(a, output_shape=(1, t), order=0, preserve_range=True, anti_aliasing=False)
        print(a, outCV, outScikit)

[[0 1 2 3 4]] [[0 2 4]] [[0 2 4]]
[[0 1 2 3 4]] [[0 1 3 4]] [[0 1 3 4]]
[[0 1 2 3 4 5]] [[0 2 4]] [[1 3 5]]
[[0 1 2 3 4 5]] [[0 2 3 5]] [[0 2 3 5]]
[[0 1 2 3 4 5 6]] [[1 3 5]] [[1 3 5]]
[[0 1 2 3 4 5 6]] [[0 2 4 6]] [[0 2 4 6]]
[[0 1 2 3 4 5 6 7]] [[1 3 6]] [[1 4 6]]
[[0 1 2 3 4 5 6 7]] [[0 2 4 6]] [[1 3 5 7]]

[ppwwyyxx]

[[0 1 2 3 4 5]] [[0 2 4]] [[1 3 5]]

Both results are correct. The sampling points are in the exact middle of two pixels, so the "nearest" pixel is ambiguous and different rounding conventions can produce different results.

[dkurt]

@ppwwyyxx, there is a statement that

INTER_NEAREST_EXACT will produce same results as the nearest neighbor method in PIL, scikit-image or Matlab.

So the goal of #23634 is to resolve #22204. However, the results are not the same for all the scales. For example:

Input	scale	OpenCV 4.x	#23634	Scikit-Image 0.20.0	PIL 9.4.0
[0 1 2 3 4 5 6 7 8 9]	0.3	[1 4 8]	[1 4 8]	[1 5 8]	[1 5 8]
[0 1 2 3 4 5]	5/6	[0 1 2 4 5]	[0 1 2 4 5]	[0 1 3 4 5]	[0 1 3 4 5]
[0 1 2 3 4 5 6 7]	0.75	[0 1 3 4 5 7]	[0 1 3 4 5 7]	[0 2 3 4 6 7]	[0 2 3 4 5 7]
[0 1 2 3 4 5]	0.5	[0 2 4]	[1 3 5]	[1 3 5]	[1 3 5]

But I believe that at least the case with x2 nearest downsampling from even sizes should be deterministic. So it's not about correctness but portability.

[scott-vsi]

FYI using project, translate, scale from this comment and INTER_NEAREST (because INTER_NEAREST_EXACT isn't recognized by warpAffine it seems), cv.warpAffine matches skimage for these tests:

M = project(
	translate([-0.5, -0.5]) @
	scale(s,1) @
	translate([+0.5, +0.5])
)

Input	scale	cv2.warpAffine
[0 1 2 3 4 5 6 7 8 9]	0.3	[1 5 8]
[0 1 2 3 4 5]	5/6	[0 1 3 4 5]
[0 1 2 3 4 5 6 7]	0.75	[0 2 3 4 6 7]
[0 1 2 3 4 5]	0.5	[1 3 5]

However, for a slightly different scale, it does not:

Input	scale	cv2.warpAffine	skimage (v0.20.0)
[0 1 2 3 4 5 6 7 8]	0.4	[1 3 6 8]	[1 3 5 7]
[0 1 2 3 4 5 6 7]	0.4	[1 3 6]	[1 4 6]

Is that expected?

It seems to usually be the same as skimage when the scale divides the shape evenly (i.e., N*s is a whole number). However, these are the same, even though scale does not divide the input length evenly (7*0.3 = 2.1, which is not a whole number):

Input	scale	cv2.warpAffine	skimage (v0.20.0)
[0 1 2 3 4 5 6]	0.3	[1 5]	[1 5]

[crackwitz]

You used integer image inputs? Then rounding may accidentally make results the same.

I would advise, in order of insight, highest to lowest:

1. investigating the source code of each implementation for its behavior
2. blackbox-testing each implementation to figure out its behavior
3. Comparing, if at all

The goal should be to behave in sensible ways, regardless of what some other library happens to do. They are not the benchmark.

[scott-vsi]

Yes. I am using integer images.

Since this thread is about making sure NN interpolation works correctly, I was simply commenting on a difference in behavior in an area that I hadn't seen mentioned in this thread.