Correctness fixes. The port now produces equivalent results to the original CUDA implementation, and further changes should be checked against regressions.

Author: sfn

src/aliceVision/depthMap_sycl/NormalMapEstimator.cpp

@@ -42,6 +42,9 @@ void NormalMapEstimator::compute(sycl::queue& queue, const std::vector<int>& cam
             // get R camera parameters
             CameraParams camParams = getCameraParameters(rc, 1 /*downscale*/, _mp);
 
+            // init prerequisite
+            sycl::event event{};
+
             // read input depth map
             image::Image<float> in_depthMap;
             mvsUtils::readMap(rc, _mp, mvsUtils::EFileType::depthMapFiltered, in_depthMap);
@@ -65,22 +68,23 @@ void NormalMapEstimator::compute(sycl::queue& queue, const std::vector<int>& cam
             // copy input depth map into depth/sim map in device memory
             // note: we don't need similarity for normal map computation
             //       we use depth/sim map in order to avoid code duplication
-            SyclDeviceMemoryPitched<sycl::float2, 2> in_depthSimMap_dmp(out_normalMap_dmp.getSize(), allocSuccess, queue);
+            SyclDeviceMemoryPitched<sycl::float2, 2> in_depthMap_dmp(out_normalMap_dmp.getSize(), allocSuccess, queue);
 
             if (!allocSuccess) ALICEVISION_THROW_ERROR("Not enough device memory to compute normal map!")
 
             {
-                SyclHostMemoryHeap<sycl::float2, 2> in_depthSimMap_hmh(in_depthSimMap_dmp.getSize(), queue);
+                SyclHostMemoryHeap<sycl::float2, 2> in_depthMap_hmh(in_depthMap_dmp.getSize(), queue);
 
-                for (int x = 0; x < width; ++x)
-                    for (int y = 0; y < height; ++y)
-                        in_depthSimMap_hmh(size_t(x), size_t(y)) = sycl::float2(in_depthMap(y, x), 1.f);
+                for (int y = 0; y < height; ++y)
+                    for (int x = 0; x < width; ++x)
+                        in_depthMap_hmh(size_t(x), size_t(y)).x() = in_depthMap(y, x);
 
-                in_depthSimMap_dmp.copyFrom(in_depthSimMap_hmh, sycl::event()).wait();
+                event = in_depthMap_dmp.copyFrom(in_depthMap_hmh, event);
             }
 
             // compute normal map synchronosly
-            sycl_depthSimMapComputeNormal(out_normalMap_dmp, in_depthSimMap_dmp, camParams, 1 /*step*/, roi, queue, sycl::event()).wait();
+            event = sycl_depthMapComputeNormal(out_normalMap_dmp, in_depthMap_dmp, camParams, 1 /*step*/, roi, queue, event);
+            event.wait();
 
             // write output normal map
             writeNormalMapFiltered(rc, _mp, tileParams, roi, out_normalMap_dmp, queue);

src/aliceVision/depthMap_sycl/Refine.cpp

@@ -30,7 +30,6 @@ namespace depthMap {
     _sgmNormalMap_dmp(queue),
     _normalMap_dmp(queue),
     _volumeRefineSim_dmp(queue),
-    _optTmpDepthMap_dmp(queue),
     _optImgVariance_dmp(queue),
     _queue(queue)
 {
@@ -64,42 +63,7 @@ namespace depthMap {
 
     // allocate depth/sim map optimization buffers
     if (_refineParams.useColorOptimization)
-    {
-        allocationSuccess &= _optTmpDepthMap_dmp.allocate(depthSimMapDim);
         allocationSuccess &= _optImgVariance_dmp.allocate(depthSimMapDim);
-    }
-}
-
-double Refine::getDeviceMemoryConsumption() const
-{
-    size_t bytes = 0;
-
-    bytes += _sgmDepthPixSizeMap_dmp.getBytesPadded();
-    bytes += _refinedDepthSimMap_dmp.getBytesPadded();
-    bytes += _optimizedDepthSimMap_dmp.getBytesPadded();
-    bytes += _sgmNormalMap_dmp.getBytesPadded();
-    bytes += _normalMap_dmp.getBytesPadded();
-    bytes += _volumeRefineSim_dmp.getBytesPadded();
-    bytes += _optTmpDepthMap_dmp.getBytesPadded();
-    bytes += _optImgVariance_dmp.getBytesPadded();
-
-    return (double(bytes) / (1024.0 * 1024.0));
-}
-
-double Refine::getDeviceMemoryConsumptionUnpadded() const
-{
-    size_t bytes = 0;
-
-    bytes += _sgmDepthPixSizeMap_dmp.getBytesUnpadded();
-    bytes += _refinedDepthSimMap_dmp.getBytesUnpadded();
-    bytes += _optimizedDepthSimMap_dmp.getBytesUnpadded();
-    bytes += _sgmNormalMap_dmp.getBytesUnpadded();
-    bytes += _normalMap_dmp.getBytesUnpadded();
-    bytes += _volumeRefineSim_dmp.getBytesUnpadded();
-    bytes += _optTmpDepthMap_dmp.getBytesUnpadded();
-    bytes += _optImgVariance_dmp.getBytesUnpadded();
-
-    return (double(bytes) / (1024.0 * 1024.0));
 }
 
 sycl::event Refine::refineRc(const Tile& tile,
@@ -206,7 +170,7 @@ sycl::event Refine::refineAndFuseDepthSimMap(const Tile& tile, sycl::event prere
 
     // initialize the similarity volume at 0
     // each tc filtered and inverted similarity value will be summed in this volume
-    prerequisite = sycl_volumeInitialize(_volumeRefineSim_dmp, TSimRefine(0.f), _queue, prerequisite);
+    prerequisite = sycl_volumeInitialize(_volumeRefineSim_dmp, TSimRefine(0.f), prerequisite);
 
     // get device cache instance
     DeviceCache& deviceCache = DeviceCache::getInstance();
@@ -234,37 +198,20 @@ sycl::event Refine::refineAndFuseDepthSimMap(const Tile& tile, sycl::event prere
                                    << "\t- tc: " << tc << " (" << (tci + 1) << "/" << tile.refineTCams.size() << ")" << std::endl
                                    << "\t- tile range x: [" << downscaledRoi.x.begin << " - " << downscaledRoi.x.end << "]" << std::endl
                                    << "\t- tile range y: [" << downscaledRoi.y.begin << " - " << downscaledRoi.y.end << "]" << std::endl);
-        if (_refineParams.useSgmNormalMap)
-        {
-            prerequisite = sycl_volumeRefineSimilarityNormal(
-                                    _volumeRefineSim_dmp,
-                                    _sgmDepthPixSizeMap_dmp,
-                                    _sgmNormalMap_dmp,
-                                    rcParams,
-                                    tcParams,
-                                    rcDeviceMipmapImage,
-                                    tcDeviceMipmapImage,
-                                    _refineParams,
-                                    depthRange,
-                                    downscaledRoi,
-                                    _queue,
-                                    prerequisite);
-        }
-        else
-        {
-            prerequisite = sycl_volumeRefineSimilarityNoNormal(
-                                    _volumeRefineSim_dmp,
-                                    _sgmDepthPixSizeMap_dmp,
-                                    rcParams,
-                                    tcParams,
-                                    rcDeviceMipmapImage,
-                                    tcDeviceMipmapImage,
-                                    _refineParams,
-                                    depthRange,
-                                    downscaledRoi,
-                                    _queue,
-                                    prerequisite);
-        }
+
+        prerequisite = sycl_volumeRefineSimilarity(
+                            _volumeRefineSim_dmp,
+                            _sgmDepthPixSizeMap_dmp,
+                            _refineParams.useSgmNormalMap ? &_sgmNormalMap_dmp : nullptr,
+                            rcParams,
+                            tcParams,
+                            rcDeviceMipmapImage,
+                            tcDeviceMipmapImage,
+                            _refineParams,
+                            depthRange,
+                            downscaledRoi,
+                            _queue,
+                            prerequisite);
     }
 
     // export intermediate volume information (if requested by user)
@@ -302,7 +249,6 @@ sycl::event Refine::optimizeDepthSimMap(const Tile& tile, sycl::event prerequisi
     sycl::event gradientDescent = sycl_depthSimMapOptimizeGradientDescent(
                                             _optimizedDepthSimMap_dmp,  // output depth/sim map optimized
                                             _optImgVariance_dmp,        // image variance buffer pre-allocate
-                                            _optTmpDepthMap_dmp,        // temporary depth map buffer pre-allocate
                                             _sgmDepthPixSizeMap_dmp,    // input SGM upscaled depth/pixSize map
                                             _refinedDepthSimMap_dmp,    // input refined and fused depth/sim map
                                             camParams,
@@ -328,7 +274,7 @@ void Refine::computeAndWriteNormalMap(const Tile& tile, const SyclDeviceMemoryPi
     ALICEVISION_LOG_INFO(tile << "Refine compute normal map of view id: " << _mp.getViewId(tile.rc) << ", rc: " << tile.rc << " (" << (tile.rc + 1)
                               << " / " << _mp.ncams << ").");
 
-    sycl_depthSimMapComputeNormal(_normalMap_dmp, in_depthSimMap_dmp, camParams, _refineParams.stepXY, downscaledRoi, _queue, prerequisite).wait();
+    sycl_depthMapComputeNormal(_normalMap_dmp, in_depthSimMap_dmp, camParams, _refineParams.stepXY, downscaledRoi, _queue, prerequisite).wait();
 
     writeNormalMap(tile.rc, _mp, _tileParams, tile.roi, _normalMap_dmp, _queue, _refineParams.scale, _refineParams.stepXY, name);
 }

src/aliceVision/depthMap_sycl/Refine.hpp

@@ -46,18 +46,6 @@ class Refine
     // final depth/similarity map getter
     inline const SyclDeviceMemoryPitched<sycl::float2, 2>& getDeviceDepthSimMap() const { return _optimizedDepthSimMap_dmp; }
 
-    /**
-     * @brief Get memory consumpyion in device memory.
-     * @return device memory consumpyion (in MB)
-     */
-    double getDeviceMemoryConsumption() const;
-
-    /**
-     * @brief Get unpadded memory consumpyion in device memory.
-     * @return unpadded device memory consumpyion (in MB)
-     */
-    double getDeviceMemoryConsumptionUnpadded() const;
-
     /**
      * @brief Refine for a single R camera the Semi-Global Matching depth/sim map.
      * @param[in] tile The given tile for Refine computation
@@ -113,7 +101,6 @@ class Refine
     SyclDeviceMemoryPitched<sycl::float3, 2> _sgmNormalMap_dmp;          //< rc upscaled SGM normal map (for experimentation purposes)
     SyclDeviceMemoryPitched<sycl::float3, 2> _normalMap_dmp;             //< rc normal map (for debug / intermediate results purposes)
     SyclDeviceMemoryPitched<TSimRefine, 3> _volumeRefineSim_dmp;   //< rc refine similarity volume
-    SyclDeviceMemoryPitched<float, 2> _optTmpDepthMap_dmp;         //< for color optimization: temporary depth map buffer
     SyclDeviceMemoryPitched<float, 2> _optImgVariance_dmp;         //< for color optimization: image variance buffer
     sycl::queue _queue;                                            //< queue for device execution
 };

src/aliceVision/depthMap_sycl/Sgm.cpp

@@ -90,40 +90,6 @@ Sgm::Sgm(const mvsUtils::MultiViewParams& mp,
     }
 }
 
-double Sgm::getDeviceMemoryConsumption() const
-{
-    size_t bytes = 0;
-
-    bytes += _depths_dmp.getBytesPadded();
-    bytes += _depthThicknessMap_dmp.getBytesPadded();
-    bytes += _depthSimMap_dmp.getBytesPadded();
-    bytes += _normalMap_dmp.getBytesPadded();
-    bytes += _volumeBestSim_dmp.getBytesPadded();
-    bytes += _volumeSecBestSim_dmp.getBytesPadded();
-    bytes += _volumeSliceAccA_dmp.getBytesPadded();
-    bytes += _volumeSliceAccB_dmp.getBytesPadded();
-    bytes += _volumeAxisAcc_dmp.getBytesPadded();
-
-    return (double(bytes) / (1024.0 * 1024.0));
-}
-
-double Sgm::getDeviceMemoryConsumptionUnpadded() const
-{
-    size_t bytes = 0;
-
-    bytes += _depths_dmp.getBytesUnpadded();
-    bytes += _depthThicknessMap_dmp.getBytesUnpadded();
-    bytes += _depthSimMap_dmp.getBytesUnpadded();
-    bytes += _normalMap_dmp.getBytesUnpadded();
-    bytes += _volumeBestSim_dmp.getBytesUnpadded();
-    bytes += _volumeSecBestSim_dmp.getBytesUnpadded();
-    bytes += _volumeSliceAccA_dmp.getBytesUnpadded();
-    bytes += _volumeSliceAccB_dmp.getBytesUnpadded();
-    bytes += _volumeAxisAcc_dmp.getBytesUnpadded();
-
-    return (double(bytes) / (1024.0 * 1024.0));
-}
-
 sycl::event Sgm::sgmRc(const Tile& tile, const SgmDepthList& tileDepthList, sycl::event prerequisite)
 {
     const IndexT viewId = _mp.getViewId(tile.rc);
@@ -184,7 +150,7 @@ sycl::event Sgm::sgmRc(const Tile& tile, const SgmDepthList& tileDepthList, sycl
 
         ALICEVISION_LOG_INFO(tile << "SGM compute normal map of view id: " << viewId << ", rc: " << tile.rc << " (" << (tile.rc + 1) << " / "
                                   << _mp.ncams << ").");
-        finished = sycl_depthSimMapComputeNormal(_normalMap_dmp, _depthSimMap_dmp, getCameraParameters(tile.rc, _sgmParams.scale, _mp), _sgmParams.stepXY, downscaledRoi, _queue, retrieveDepth);
+        finished = sycl_depthMapComputeNormal(_normalMap_dmp, _depthSimMap_dmp, getCameraParameters(tile.rc, _sgmParams.scale, _mp), _sgmParams.stepXY, downscaledRoi, _queue, retrieveDepth);
 
         // export intermediate normal map (if requested by user)
         if (_sgmParams.exportIntermediateNormalMaps)
@@ -222,8 +188,8 @@ sycl::event Sgm::computeSimilarityVolumes(const Tile& tile, const SgmDepthList&
     const ROI downscaledRoi = downscaleROI(tile.roi, _sgmParams.scale * _sgmParams.stepXY);
 
     // initialize the two similarity volumes at 255
-    prerequisite = sycl_volumeInitialize<TSim>(_volumeBestSim_dmp, 255., _queue, prerequisite);
-    prerequisite = sycl_volumeInitialize<TSim>(_volumeSecBestSim_dmp, 255., _queue, prerequisite);
+    prerequisite = sycl_volumeInitialize<TSim>(_volumeBestSim_dmp, 255., prerequisite);
+    prerequisite = sycl_volumeInitialize<TSim>(_volumeSecBestSim_dmp, 255., prerequisite);
 
     // get device cache instance
     DeviceCache& deviceCache = DeviceCache::getInstance();

src/aliceVision/depthMap_sycl/Sgm.hpp

@@ -62,18 +62,6 @@ class Sgm
     // final normal map getter (optional: could be empty)
     inline const SyclDeviceMemoryPitched<sycl::float3, 2>& getDeviceNormalMap() const { return _normalMap_dmp; }
 
-    /**
-     * @brief Get memory consumpyion in device memory.
-     * @return device memory consumpyion (in MB)
-     */
-    double getDeviceMemoryConsumption() const;
-
-    /**
-     * @brief Get unpadded memory consumpyion in device memory.
-     * @return unpadded device memory consumpyion (in MB)
-     */
-    double getDeviceMemoryConsumptionUnpadded() const;
-
     /**
      * @brief Compute for a single R camera the Semi-Global Matching.
      * @param[in] tile The given tile for SGM computation

src/aliceVision/depthMap_sycl/sycl/Patch.hpp

@@ -25,14 +25,14 @@ struct Patch
     float d;  //< pixel size
 };
 
-inline void rotPointAroundVect(sycl::float3& out, const sycl::float3& X, const sycl::float3& vect, const int& angle)
+inline void rotPointAroundVect(sycl::float3& point, const sycl::float3& vect, const int& angle)
 {
     double ux, uy, uz, vx, vy, vz, wx, wy, wz, sa, ca, x, y, z, u, v, w;
 
-    double sizeX = sycl::length(X);
-    x = X.x() / sizeX;
-    y = X.y() / sizeX;
-    z = X.z() / sizeX;
+    const double size = sycl::length(point);
+    x = point.x() / size;
+    y = point.y() / size;
+    z = point.z() / size;
     u = vect.x();
     v = vect.y();
     w = vect.z();
@@ -53,23 +53,21 @@ inline void rotPointAroundVect(sycl::float3& out, const sycl::float3& X, const s
     y = v * (ux + vy + wz) + (y * (u * u + w * w) - v * (ux + wz)) * ca + (wx - uz) * sa;
     z = w * (ux + vy + wz) + (z * (u * u + v * v) - w * (ux + vy)) * ca + (-vx + uy) * sa;
 
-    u = sycl::sqrt(x * x + y * y + z * z);
-
-    out[0] = x;
-    out[1] = y;
-    out[2] = z;
-    out *= sizeX/u;
+    point[0] = x;
+    point[1] = y;
+    point[2] = z;
+    point = sycl::normalize(point) * size;
 }
 
 inline static constexpr void rotatePatch(Patch& ptch, const int& rx, const int& ry)
 {
     // rotate patch around x axis by angle rx
-    rotPointAroundVect(ptch.n, ptch.n, ptch.x, rx);
-    rotPointAroundVect(ptch.y, ptch.y, ptch.x, rx);
+    rotPointAroundVect(ptch.n, ptch.x, rx);
+    rotPointAroundVect(ptch.y, ptch.x, rx);
 
     // rotate new patch around y axis by angle ry
-    rotPointAroundVect(ptch.n, ptch.n, ptch.y, ry);
-    rotPointAroundVect(ptch.x, ptch.x, ptch.y, ry);
+    rotPointAroundVect(ptch.n, ptch.y, ry);
+    rotPointAroundVect(ptch.x, ptch.y, ry);
 }
 
 inline void movePatch(Patch& ptch, const int& pt)

src/aliceVision/depthMap_sycl/sycl/color.hpp

@@ -100,11 +100,11 @@ inline sycl::float3 rgb2hsl(const sycl::float3& c)
 inline sycl::float3 xyz2lab(const sycl::float3& c)
 {
     // assuming whitepoint D65, XYZ=(0.95047, 1.00000, 1.08883)
-    sycl::float3 r = sycl::float3(c.x() / 0.95047f, c.y(), c.z() / 1.08883f);
+    const sycl::float3 r = c * sycl::float3(1.f / 0.95047f, 1.f, 1.f / 1.08883f);
 
-    sycl::float3 f = sycl::float3((r.x() > 216.0f / 24389.0f ? sycl::cbrt(r.x()) : (24389.0f / 27.0f * r.x() + 16.0f) / 116.0f),
-                                  (r.y() > 216.0f / 24389.0f ? sycl::cbrt(r.y()) : (24389.0f / 27.0f * r.y() + 16.0f) / 116.0f),
-                                  (r.z() > 216.0f / 24389.0f ? sycl::cbrt(r.z()) : (24389.0f / 27.0f * r.z() + 16.0f) / 116.0f));
+    sycl::float3 f = sycl::float3(r.x() > 216.0f / 24389.0f ? sycl::cbrt(r.x()) : (24389.0f / 27.0f * r.x() + 16.0f) / 116.0f,
+                                  r.y() > 216.0f / 24389.0f ? sycl::cbrt(r.y()) : (24389.0f / 27.0f * r.y() + 16.0f) / 116.0f,
+                                  r.z() > 216.0f / 24389.0f ? sycl::cbrt(r.z()) : (24389.0f / 27.0f * r.z() + 16.0f) / 116.0f);
 
     sycl::float3 out = sycl::float3(116.0f * f.y() - 16.0f, 500.0f * (f.x() - f.y()), 200.0f * (f.y() - f.z()));
 
@@ -122,7 +122,7 @@ inline sycl::float3 xyz2lab(const sycl::float3& c)
 inline float rgb2gray(const sycl::uchar4& c_char)
 {
     const sycl::float4 c = c_char.convert<float>();
-    return sycl::dot(sycl::float3(c.x(), c.y(), c.z()), sycl::float3(0.2989f, 0.5870f, 0.1140f));
+    return c.x() * 0.2989f + c.y() * 0.5870f + c.z() * 0.1140f;
 }
 
 /**