Merge pull request opencv#18053 from Yosshi999:bit-exact-resizeNN

Bit-exact Nearest Neighbor Resizing

* bit exact resizeNN

* change the value of method enum

* add bitexact-nn to ResizeExactTest

* test to compare with non-exact version

* add perf for bit-exact resizenn

* use cvFloor-equivalent

* 1/3 scaling is not stable for floating calculation

* stricter test

* bugfix: broken data in case of 6 or 12bytes elements

* bugfix: broken data in default pix_size

* stricter threshold

* use raw() for floor

* use double instead of int

* follow code reviews

* fewer cases in perf test

* center pixel convention

Author: Yosshi999

modules/imgproc/include/opencv2/imgproc.hpp

@@ -252,6 +252,9 @@ enum InterpolationFlags{
     INTER_LANCZOS4       = 4,
     /** Bit exact bilinear interpolation */
     INTER_LINEAR_EXACT = 5,
+    /** Bit exact nearest neighbor interpolation. This will produce same results as
+    the nearest neighbor method in PIL, scikit-image or Matlab. */
+    INTER_NEAREST_EXACT  = 6,
     /** mask for interpolation codes */
     INTER_MAX            = 7,
     /** flag, fills all of the destination image pixels. If some of them correspond to outliers in the

modules/imgproc/perf/perf_resize.cpp

@@ -254,4 +254,30 @@ PERF_TEST_P(MatInfo_Size_Scale_NN, ResizeNN,
     SANITY_CHECK_NOTHING();
 }
 
+PERF_TEST_P(MatInfo_Size_Scale_NN, ResizeNNExact,
+    testing::Combine(
+        testing::Values(CV_8UC1, CV_8UC3, CV_8UC4),
+        testing::Values(sz720p, sz1080p),
+        testing::Values(0.25, 0.5, 2.0)
+    )
+)
+{
+    int matType = get<0>(GetParam());
+    Size from = get<1>(GetParam());
+    double scale = get<2>(GetParam());
+
+    cv::Mat src(from, matType);
+
+    Size to(cvRound(from.width * scale), cvRound(from.height * scale));
+    cv::Mat dst(to, matType);
+
+    declare.in(src, WARMUP_RNG).out(dst);
+    declare.time(100);
+
+    TEST_CYCLE() resize(src, dst, dst.size(), 0, 0, INTER_NEAREST_EXACT);
+
+    EXPECT_GT(countNonZero(dst.reshape(1)), 0);
+    SANITY_CHECK_NOTHING();
+}
+
 } // namespace

modules/imgproc/src/fixedpoint.inl.hpp

@@ -157,6 +157,7 @@ class ufixedpoint64
     CV_ALWAYS_INLINE bool isZero() { return val == 0; }
     static CV_ALWAYS_INLINE ufixedpoint64 zero() { return ufixedpoint64(); }
     static CV_ALWAYS_INLINE ufixedpoint64 one() { return ufixedpoint64((uint64_t)(1ULL << fixedShift)); }
+    CV_ALWAYS_INLINE uint32_t cvFloor() const { return cv::saturate_cast<uint32_t>(val >> fixedShift); }
     friend class ufixedpoint32;
 };
 

modules/imgproc/src/resize.cpp

@@ -51,6 +51,7 @@
 #include "opencl_kernels_imgproc.hpp"
 #include "hal_replacement.hpp"
 #include "opencv2/core/hal/intrin.hpp"
+#include "opencv2/core/utils/buffer_area.private.hpp"
 
 #include "opencv2/core/openvx/ovx_defs.hpp"
 #include "resize.hpp"
@@ -1104,6 +1105,121 @@ resizeNN( const Mat& src, Mat& dst, double fx, double fy )
     }
 }
 
+class resizeNN_bitexactInvoker : public ParallelLoopBody
+{
+public:
+    resizeNN_bitexactInvoker(const Mat& _src, Mat& _dst, int* _x_ofse, int _ify, int _ify0)
+        : src(_src), dst(_dst), x_ofse(_x_ofse), ify(_ify), ify0(_ify0) {}
+
+    virtual void operator() (const Range& range) const CV_OVERRIDE
+    {
+        Size ssize = src.size(), dsize = dst.size();
+        int pix_size = (int)src.elemSize();
+        for( int y = range.start; y < range.end; y++ )
+        {
+            uchar* D = dst.ptr(y);
+            int _sy = (ify * y + ify0) >> 16;
+            int sy = std::min(_sy, ssize.height-1);
+            const uchar* S = src.ptr(sy);
+
+            int x = 0;
+            switch( pix_size )
+            {
+            case 1:
+#if CV_SIMD
+                for( ; x <= dsize.width - v_uint8::nlanes; x += v_uint8::nlanes )
+                    v_store(D + x, vx_lut(S, x_ofse + x));
+#endif
+                for( ; x < dsize.width; x++ )
+                    D[x] = S[x_ofse[x]];
+                break;
+            case 2:
+#if CV_SIMD
+                for( ; x <= dsize.width - v_uint16::nlanes; x += v_uint16::nlanes )
+                    v_store((ushort*)D + x, vx_lut((ushort*)S, x_ofse + x));
+#endif
+                for( ; x < dsize.width; x++ )
+                    *((ushort*)D + x) = *((ushort*)S + x_ofse[x]);
+                break;
+            case 3:
+                for( ; x < dsize.width; x++, D += 3 )
+                {
+                    const uchar* _tS = S + x_ofse[x] * 3;
+                    D[0] = _tS[0]; D[1] = _tS[1]; D[2] = _tS[2];
+                }
+                break;
+            case 4:
+#if CV_SIMD
+                for( ; x <= dsize.width - v_uint32::nlanes; x += v_uint32::nlanes )
+                    v_store((uint32_t*)D + x, vx_lut((uint32_t*)S, x_ofse + x));
+#endif
+                for( ; x < dsize.width; x++ )
+                    *((uint32_t*)D + x) = *((uint32_t*)S + x_ofse[x]);
+                break;
+            case 6:
+                for( ; x < dsize.width; x++, D += 6 )
+                {
+                    const ushort* _tS = (const ushort*)(S + x_ofse[x]*6);
+                    ushort* _tD = (ushort*)D;
+                    _tD[0] = _tS[0]; _tD[1] = _tS[1]; _tD[2] = _tS[2];
+                }
+                break;
+            case 8:
+#if CV_SIMD
+                for( ; x <= dsize.width - v_uint64::nlanes; x += v_uint64::nlanes )
+                    v_store((uint64_t*)D + x, vx_lut((uint64_t*)S, x_ofse + x));
+#endif
+                for( ; x < dsize.width; x++ )
+                    *((uint64_t*)D + x) = *((uint64_t*)S + x_ofse[x]);
+                break;
+            case 12:
+                for( ; x < dsize.width; x++, D += 12 )
+                {
+                    const int* _tS = (const int*)(S + x_ofse[x]*12);
+                    int* _tD = (int*)D;
+                    _tD[0] = _tS[0]; _tD[1] = _tS[1]; _tD[2] = _tS[2];
+                }
+                break;
+            default:
+                for( x = 0; x < dsize.width; x++, D += pix_size )
+                {
+                    const uchar* _tS = S + x_ofse[x] * pix_size;
+                    for (int k = 0; k < pix_size; k++)
+                        D[k] = _tS[k];
+                }
+            }
+        }
+    }
+private:
+    const Mat& src;
+    Mat& dst;
+    int* x_ofse;
+    const int ify;
+    const int ify0;
+};
+
+static void resizeNN_bitexact( const Mat& src, Mat& dst, double /*fx*/, double /*fy*/ )
+{
+    Size ssize = src.size(), dsize = dst.size();
+    int ifx = ((ssize.width << 16) + dsize.width / 2) / dsize.width; // 16bit fixed-point arithmetic
+    int ifx0 = ifx / 2 - 1;                                     // This method uses center pixel coordinate as Pillow and scikit-images do.
+    int ify = ((ssize.height << 16) + dsize.height / 2) / dsize.height;
+    int ify0 = ify / 2 - 1;
+
+    cv::utils::BufferArea area;
+    int* x_ofse = 0;
+    area.allocate(x_ofse, dsize.width, CV_SIMD_WIDTH);
+    area.commit();
+
+    for( int x = 0; x < dsize.width; x++ )
+    {
+        int sx = (ifx * x + ifx0) >> 16;
+        x_ofse[x] = std::min(sx, ssize.width-1);    // offset in element (not byte)
+    }
+    Range range(0, dsize.height);
+    resizeNN_bitexactInvoker invoker(src, dst, x_ofse, ify, ify0);
+    parallel_for_(range, invoker, dst.total()/(double)(1<<16));
+}
 
 struct VResizeNoVec
 {
@@ -3723,6 +3839,12 @@ void resize(int src_type,
         return;
     }
 
+    if( interpolation == INTER_NEAREST_EXACT )
+    {
+        resizeNN_bitexact( src, dst, inv_scale_x, inv_scale_y );
+        return;
+    }
+
     int k, sx, sy, dx, dy;
 
 

modules/imgproc/test/test_imgwarp.cpp

@@ -346,14 +346,24 @@ class CV_ResizeExactTest : public CV_ResizeTest
 
 CV_ResizeExactTest::CV_ResizeExactTest() : CV_ResizeTest()
 {
-    max_interpolation = 1;
+    max_interpolation = 2;
 }
 
 
 void CV_ResizeExactTest::get_test_array_types_and_sizes(int test_case_idx, vector<vector<Size> >& sizes, vector<vector<int> >& types)
 {
     CV_ResizeTest::get_test_array_types_and_sizes(test_case_idx, sizes, types);
-    interpolation = INTER_LINEAR_EXACT;
+    switch (interpolation)
+    {
+    case 0:
+        interpolation = INTER_LINEAR_EXACT;
+        break;
+    case 1:
+        interpolation = INTER_NEAREST_EXACT;
+        break;
+    default:
+        CV_Assert(interpolation < max_interpolation);
+    }
     if (CV_MAT_DEPTH(types[INPUT][0]) == CV_32F ||
         CV_MAT_DEPTH(types[INPUT][0]) == CV_64F)
         types[INPUT][0] = types[INPUT_OUTPUT][0] = types[REF_INPUT_OUTPUT][0] = CV_MAKETYPE(CV_8U, CV_MAT_CN(types[INPUT][0]));

modules/imgproc/test/test_resize_bitexact.cpp

@@ -152,4 +152,89 @@ TEST(Resize_Bitexact, Linear8U)
         }
 }
 
+PARAM_TEST_CASE(Resize_Bitexact, int)
+{
+public:
+    int depth;
+
+    virtual void SetUp()
+    {
+        depth = GET_PARAM(0);
+    }
+
+    double CountDiff(const Mat& src)
+    {
+        Mat dstExact; cv::resize(src, dstExact, Size(), 2, 1, INTER_NEAREST_EXACT);
+        Mat dstNonExact; cv::resize(src, dstNonExact, Size(), 2, 1, INTER_NEAREST);
+
+        return cv::norm(dstExact, dstNonExact, NORM_INF);
+    }
+};
+
+TEST_P(Resize_Bitexact, Nearest8U_vsNonExact)
+{
+    Mat mat_color, mat_gray;
+    Mat src_color = imread(cvtest::findDataFile("shared/lena.png"));
+    Mat src_gray; cv::cvtColor(src_color, src_gray, COLOR_BGR2GRAY);
+    src_color.convertTo(mat_color, depth);
+    src_gray.convertTo(mat_gray, depth);
+
+    EXPECT_EQ(CountDiff(mat_color), 0) << "color, type: " << depth;
+    EXPECT_EQ(CountDiff(mat_gray), 0) << "gray, type: " << depth;
+}
+
+// Now INTER_NEAREST's convention and INTER_NEAREST_EXACT's one are different.
+INSTANTIATE_TEST_CASE_P(DISABLED_Imgproc, Resize_Bitexact,
+    testing::Values(CV_8U, CV_16U, CV_32F, CV_64F)
+);
+
+TEST(Resize_Bitexact, Nearest8U)
+{
+    Mat src[6], dst[6];
+
+    // 2x decimation
+    src[0] = (Mat_<uint8_t>(1, 6) << 0, 1, 2, 3, 4, 5);
+    dst[0] = (Mat_<uint8_t>(1, 3) << 0, 2, 4);
+
+    // decimation odd to 1
+    src[1] = (Mat_<uint8_t>(1, 5) << 0, 1, 2, 3, 4);
+    dst[1] = (Mat_<uint8_t>(1, 1) << 2);
+
+    // decimation n*2-1 to n
+    src[2] = (Mat_<uint8_t>(1, 5) << 0, 1, 2, 3, 4);
+    dst[2] = (Mat_<uint8_t>(1, 3) << 0, 2, 4);
+
+    // decimation n*2+1 to n
+    src[3] = (Mat_<uint8_t>(1, 5) << 0, 1, 2, 3, 4);
+    dst[3] = (Mat_<uint8_t>(1, 2) << 1, 3);
+
+    // zoom
+    src[4] = (Mat_<uint8_t>(3, 5) <<
+        0, 1, 2, 3, 4,
+        5, 6, 7, 8, 9,
+        10, 11, 12, 13, 14);
+    dst[4] = (Mat_<uint8_t>(5, 7) <<
+        0, 1, 1, 2, 3, 3, 4,
+        0, 1, 1, 2, 3, 3, 4,
+        5, 6, 6, 7, 8, 8, 9,
+        10, 11, 11, 12, 13, 13, 14,
+        10, 11, 11, 12, 13, 13, 14);
+
+    src[5] = (Mat_<uint8_t>(2, 3) <<
+        0, 1, 2,
+        3, 4, 5);
+    dst[5] = (Mat_<uint8_t>(4, 6) <<
+        0, 0, 1, 1, 2, 2,
+        0, 0, 1, 1, 2, 2,
+        3, 3, 4, 4, 5, 5,
+        3, 3, 4, 4, 5, 5);
+
+    for (int i = 0; i < 6; i++)
+    {
+        Mat calc;
+        resize(src[i], calc, dst[i].size(), 0, 0, INTER_NEAREST_EXACT);
+        EXPECT_EQ(cvtest::norm(calc, dst[i], cv::NORM_L1), 0);
+    }
+}
+
 }} // namespace