From d51896639de07ec381ee068ef4645f23217d2c8d Mon Sep 17 00:00:00 2001 From: Povilas Kanapickas Date: Sun, 2 Oct 2022 02:56:11 +0300 Subject: [PATCH 1/2] [numeric] Introduce divideRoundUp() for integer division with rounding --- src/aliceVision/numeric/numeric.hpp | 14 ++++++ src/aliceVision/numeric/numeric_test.cpp | 61 ++++++++++++++++++++++++ 2 files changed, 75 insertions(+) diff --git a/src/aliceVision/numeric/numeric.hpp b/src/aliceVision/numeric/numeric.hpp index 00d13a5ed4..5c466f1d39 100644 --- a/src/aliceVision/numeric/numeric.hpp +++ b/src/aliceVision/numeric/numeric.hpp @@ -570,6 +570,20 @@ void SplitRange(const T range_start, const T range_end, const int nb_split, } } +template +constexpr T divideRoundUp(T x, T y) +{ + static_assert(std::is_integral::value, "divideRoundUp only works with integer arguments"); + const auto xPos = x >= 0; + const auto yPos = y >= 0; + if (xPos == yPos) { + return x / y + T((x % y) != 0); + } else { + // negative result, rounds towards zero anyways + return x / y; + } +} + /** * This function initializes the global state of random number generators that e.g. our tests * depend on. This makes it possible to have exactly reproducible program runtime behavior diff --git a/src/aliceVision/numeric/numeric_test.cpp b/src/aliceVision/numeric/numeric_test.cpp index b99a9d3d9b..b7a367160a 100644 --- a/src/aliceVision/numeric/numeric_test.cpp +++ b/src/aliceVision/numeric/numeric_test.cpp @@ -113,3 +113,64 @@ BOOST_AUTO_TEST_CASE(Numeric_MeanAndVarianceAlongRows) { BOOST_CHECK_SMALL(0.25-variance(0), 1e-8); BOOST_CHECK_SMALL(1.25-variance(1), 1e-8); } + +BOOST_AUTO_TEST_CASE(Numeric_divideRoundUp) +{ + BOOST_CHECK_EQUAL(divideRoundUp(0, 1), 0); + BOOST_CHECK_EQUAL(divideRoundUp(1, 1), 1); + BOOST_CHECK_EQUAL(divideRoundUp(2, 1), 2); + BOOST_CHECK_EQUAL(divideRoundUp(0, 2), 0); + BOOST_CHECK_EQUAL(divideRoundUp(1, 2), 1); + BOOST_CHECK_EQUAL(divideRoundUp(2, 2), 1); + BOOST_CHECK_EQUAL(divideRoundUp(3, 2), 2); + BOOST_CHECK_EQUAL(divideRoundUp(4, 2), 2); + BOOST_CHECK_EQUAL(divideRoundUp(999, 1000), 1); + BOOST_CHECK_EQUAL(divideRoundUp(1000, 1000), 1); + BOOST_CHECK_EQUAL(divideRoundUp(1001, 1000), 2); + BOOST_CHECK_EQUAL(divideRoundUp(1000999, 1000), 1001); + BOOST_CHECK_EQUAL(divideRoundUp(1001000, 1000), 1001); + BOOST_CHECK_EQUAL(divideRoundUp(1001001, 1000), 1002); + + BOOST_CHECK_EQUAL(divideRoundUp(-1, 1), -1); + BOOST_CHECK_EQUAL(divideRoundUp(-2, 1), -2); + BOOST_CHECK_EQUAL(divideRoundUp(-0, 2), 0); + BOOST_CHECK_EQUAL(divideRoundUp(-1, 2), 0); + BOOST_CHECK_EQUAL(divideRoundUp(-2, 2), -1); + BOOST_CHECK_EQUAL(divideRoundUp(-3, 2), -1); + BOOST_CHECK_EQUAL(divideRoundUp(-4, 2), -2); + BOOST_CHECK_EQUAL(divideRoundUp(-999, 1000), 0); + BOOST_CHECK_EQUAL(divideRoundUp(-1000, 1000), -1); + BOOST_CHECK_EQUAL(divideRoundUp(-1001, 1000), -1); + BOOST_CHECK_EQUAL(divideRoundUp(-1000999, 1000), -1000); + BOOST_CHECK_EQUAL(divideRoundUp(-1001000, 1000), -1001); + BOOST_CHECK_EQUAL(divideRoundUp(-1001001, 1000), -1001); + + BOOST_CHECK_EQUAL(divideRoundUp(0, -1), 0); + BOOST_CHECK_EQUAL(divideRoundUp(1, -1), -1); + BOOST_CHECK_EQUAL(divideRoundUp(2, -1), -2); + BOOST_CHECK_EQUAL(divideRoundUp(0, -2), 0); + BOOST_CHECK_EQUAL(divideRoundUp(1, -2), 0); + BOOST_CHECK_EQUAL(divideRoundUp(2, -2), -1); + BOOST_CHECK_EQUAL(divideRoundUp(3, -2), -1); + BOOST_CHECK_EQUAL(divideRoundUp(4, -2), -2); + BOOST_CHECK_EQUAL(divideRoundUp(999, -1000), 0); + BOOST_CHECK_EQUAL(divideRoundUp(1000, -1000), -1); + BOOST_CHECK_EQUAL(divideRoundUp(1001, -1000), -1); + BOOST_CHECK_EQUAL(divideRoundUp(1000999, -1000), -1000); + BOOST_CHECK_EQUAL(divideRoundUp(1001000, -1000), -1001); + BOOST_CHECK_EQUAL(divideRoundUp(1001001, -1000), -1001); + + BOOST_CHECK_EQUAL(divideRoundUp(-1, -1), 1); + BOOST_CHECK_EQUAL(divideRoundUp(-2, -1), 2); + BOOST_CHECK_EQUAL(divideRoundUp(0, -2), 0); + BOOST_CHECK_EQUAL(divideRoundUp(-1, -2), 1); + BOOST_CHECK_EQUAL(divideRoundUp(-2, -2), 1); + BOOST_CHECK_EQUAL(divideRoundUp(-3, -2), 2); + BOOST_CHECK_EQUAL(divideRoundUp(-4, -2), 2); + BOOST_CHECK_EQUAL(divideRoundUp(-999, -1000), 1); + BOOST_CHECK_EQUAL(divideRoundUp(-1000, -1000), 1); + BOOST_CHECK_EQUAL(divideRoundUp(-1001, -1000), 2); + BOOST_CHECK_EQUAL(divideRoundUp(-1000999, -1000), 1001); + BOOST_CHECK_EQUAL(divideRoundUp(-1001000, -1000), 1001); + BOOST_CHECK_EQUAL(divideRoundUp(-1001001, -1000), 1002); +} From 2caea252cb93de9d4fb08cc52aae2999b4a678dc Mon Sep 17 00:00:00 2001 From: Povilas Kanapickas Date: Sun, 2 Oct 2022 02:56:12 +0300 Subject: [PATCH 2/2] [fuseCut] Fix incorrect rounding of integer division result The intention here likely was to get floating-point division result and round it upwards. However, the division was performed with integer arguments, fractional part was lost and thus std::ceil did nothing. --- src/aliceVision/fuseCut/DelaunayGraphCut.cpp | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/src/aliceVision/fuseCut/DelaunayGraphCut.cpp b/src/aliceVision/fuseCut/DelaunayGraphCut.cpp index 12562c8c48..0627cf5efe 100644 --- a/src/aliceVision/fuseCut/DelaunayGraphCut.cpp +++ b/src/aliceVision/fuseCut/DelaunayGraphCut.cpp @@ -958,8 +958,8 @@ void DelaunayGraphCut::addMaskHelperPoints(const Point3d voxel[8], const StaticV } } - int syMax = std::ceil(height / step); - int sxMax = std::ceil(width / step); + int syMax = divideRoundUp(height, step); + int sxMax = divideRoundUp(width, step); for(int sy = 0; sy < syMax; ++sy) { @@ -1063,7 +1063,8 @@ void DelaunayGraphCut::fuseFromDepthMaps(const StaticVector& cams, const Po { const auto& imgParams = _mp.getImageParams(i); startIndex[i] = realMaxVertices; - realMaxVertices += std::ceil(imgParams.width / step) * std::ceil(imgParams.height / step); + realMaxVertices += divideRoundUp(imgParams.width, step) * + divideRoundUp(imgParams.height, step); } std::vector verticesCoordsPrepare(realMaxVertices); std::vector pixSizePrepare(realMaxVertices); @@ -1135,8 +1136,8 @@ void DelaunayGraphCut::fuseFromDepthMaps(const StaticVector& cams, const Po } } - int syMax = std::ceil(height/step); - int sxMax = std::ceil(width/step); + int syMax = divideRoundUp(height, step); + int sxMax = divideRoundUp(width, step); #pragma omp parallel for for(int sy = 0; sy < syMax; ++sy) {