Merge branch 'release-0.12.x'

Author: gregjohnson

CHANGELOG.md

@@ -1,6 +1,11 @@
 Version History
 ---------------
 
+### Open VKL 0.12.1
+
+-   Fixed bug in VDB volume interval iterator implementation which could lead to
+    missed intervals or incorrect value ranges in returned intervals
+
 ### Open VKL 0.12.0
 
 -   Added support for temporally varying volumes with associated API changes for

CMakeLists.txt

@@ -19,7 +19,7 @@ set(CMAKE_POSITION_INDEPENDENT_CODE ON)
 
 ## Establish project ##
 
-project(openvkl VERSION 0.12.0 LANGUAGES C CXX)
+project(openvkl VERSION 0.12.1 LANGUAGES C CXX)
 
 ## Add openvkl specific macros ##
 

README.md

@@ -1,6 +1,6 @@
 # Intel® Open Volume Kernel Library
 
-This is release v0.12.0 of Intel® Open VKL. For changes and new features
+This is release v0.12.1 of Intel® Open VKL. For changes and new features
 see the [changelog](CHANGELOG.md). Visit http://www.openvkl.org for more
 information.
 
@@ -32,6 +32,12 @@ example renderers to demonstrate how to best use the Open VKL API.
 
 ## Version History
 
+### Open VKL 0.12.1
+
+  - Fixed bug in VDB volume interval iterator implementation which could
+    lead to missed intervals or incorrect value ranges in returned
+    intervals
+
 ### Open VKL 0.12.0
 
   - Added support for temporally varying volumes with associated API

gitlab/.gitlab-ci.yml

@@ -4,7 +4,8 @@
 variables:
   GIT_DEPTH: "15"
   KW_PROJECT_NAME: openvkl
-  OPENVKL_RELEASE_PACKAGE_VERSION: "0.12.0"
+  OPENVKL_RELEASE_PACKAGE_VERSION: "0.12.1"
+  MACOSX_DEPLOYMENT_TARGET: "10.13"
 
 stages:
   - build
@@ -164,6 +165,7 @@ build-arch-devel-rkdeps:
   image: $DOCKER_REGISTRY/ospray/docker-images:arch
   script:
     - gitlab/build.sh -D RKCOMMON_VERSION=devel -D RKCOMMON_HASH="" -D EMBREE_VERSION=devel -D EMBREE_HASH=""
+  allow_failure: true
 
 build-arch-TBB2021:
   <<: *build_job_docker
@@ -298,6 +300,7 @@ test-functional-arch-devel-rkdeps:
   image: $DOCKER_REGISTRY/ospray/docker-images:arch
   needs:
     - build-arch-devel-rkdeps
+  allow_failure: true
 
 test-functional-arch-TBB2021:
   <<: *test_functional_job_docker

openvkl/drivers/ispc/volume/vdb/VdbSampler.ispc

@@ -308,85 +308,152 @@ inline uniform float VdbSampler_sample(const VdbSampler *uniform sampler,
 // Value range computation.
 // ---------------------------------------------------------------------------
 
-inline void extendValueRangeFilterFloat(const void *uniform _grid,
-                                        const vec3ui *uniform offset,
-                                        uniform uint32 level,
-                                        uniform uint32 attributeIndex,
-                                        uniform box1f *uniform range)
-{
-  uniform VdbSampler sampler;
-  memset(&sampler, 0, sizeof(uniform VdbSampler));
-  sampler.grid             = (const VdbGrid *uniform)_grid;
-  sampler.maxSamplingDepth = VKL_VDB_NUM_LEVELS - 1;
-  // We set none of the other members as the internal code does not require
-  // them.
-
-  float vminFilter            = pos_inf;
-  float vmaxFilter            = neg_inf;
-  const uniform uint32 res    = vklVdbLevelRes(level);
-  const VdbGrid *uniform grid = (const VdbGrid *uniform)_grid;
-  const uniform vec3ui os     = *offset;
-  const uniform int radius    = 1;  // TRILINEAR INTERPOLATION!
-  // We compute min/max over the voxels (one filter radius) after upper bound.
-  // These are voxels that can be interpolated into our node, so we must include
-  // them for conservative value ranges.
-  foreach (z = -radius... res + radius, y = -radius... res + radius) {
-    const int32 x      = res;
-    const float sample = VdbSampler_sample(
-        &sampler, make_vec3i(x + os.x, y + os.y, z + os.z), attributeIndex);
-    vminFilter = min(vminFilter, sample);
-    vmaxFilter = max(vmaxFilter, sample);
-  }
-  // We only have to go to x = res+radius as we added those voxels above.
-  foreach (z = -radius... res + radius, x = -radius... res) {
-    const int32 y      = res;
-    const float sample = VdbSampler_sample(
-        &sampler, make_vec3i(x + os.x, y + os.y, z + os.z), attributeIndex);
-    vminFilter = min(vminFilter, sample);
-    vmaxFilter = max(vmaxFilter, sample);
-  }
-  // We only have to go to x = res+radius and y = res+radius as we added those
-  // voxels above.
-  foreach (y = -radius... res, x = -radius... res) {
-    const int32 z      = res;
-    const float sample = VdbSampler_sample(
-        &sampler, make_vec3i(x + os.x, y + os.y, z + os.z), attributeIndex);
-    vminFilter = min(vminFilter, sample);
-    vmaxFilter = max(vmaxFilter, sample);
-  }
-
-  range->lower = min(reduce_min(vminFilter), range->lower);
-  range->upper = max(reduce_max(vmaxFilter), range->upper);
-}
-
 /*
- * Compute the value range on the given tile
+ * Compute the value range for the given node and index range.
  */
-export void EXPORT_UNIQUE(VdbSampler_valueRangeTileFloat,
-                          const void *uniform _grid,
-                          const vec3ui *uniform offset,
-                          uniform uint32 level,
-                          uniform uint32 attributeIndex,
-                          uniform box1f *uniform range)
+inline uniform box1f
+VdbSampler_computeNodeValueRange(const VdbGrid *uniform grid,
+                                 uniform uint64 voxel,
+                                 const uniform vec2ui &rangeX,
+                                 const uniform vec2ui &rangeY,
+                                 const uniform vec2ui &rangeZ,
+                                 uniform uint32 attributeIndex)
 {
-  range->lower = pos_inf;
-  range->upper = neg_inf;
-  extendValueRangeFilterFloat(_grid, offset, level, attributeIndex, range);
+  const uniform VKLFormat leafFormat = vklVdbVoxelIsLeafPtr(voxel)
+                                           ? vklVdbVoxelLeafGetFormat(voxel)
+                                           : VKL_FORMAT_INVALID;
+  const void *uniform leafPtr        = vklVdbVoxelLeafGetPtr(voxel);
+  const uniform float background     = 0.f;
+  uniform box1f valueRange           = make_box1f(background, background);
+
+  switch (leafFormat) {
+  case VKL_FORMAT_TILE:
+    valueRange.lower = VdbSampler_sampleTileFloatLeaf(
+        ((const Data1D *uniform)leafPtr), attributeIndex);
+    valueRange.upper = valueRange.lower;
+    break;
+  case VKL_FORMAT_CONSTANT_ZYX:
+    valueRange = make_box1f(pos_inf, neg_inf);
+    if (grid->allLeavesCompact) {
+      for (uniform unsigned int x = rangeX.x; x < rangeX.y; ++x) {
+        for (uniform unsigned int y = rangeY.x; y < rangeY.y; ++y) {
+          for (uniform unsigned int z = rangeZ.x; z < rangeZ.y; ++z) {
+            const uniform float sample =
+                VdbSampler_sampleConstantFloatLeafCompact(
+                    ((const uniform Data1D *uniform)leafPtr),
+                    make_vec3ui(x, y, z),
+                    attributeIndex);
+            extend(valueRange, sample);
+          }
+        }
+      }
+    } else {
+      for (uniform unsigned int x = rangeX.x; x < rangeX.y; ++x) {
+        for (uniform unsigned int y = rangeY.x; y < rangeY.y; ++y) {
+          for (uniform unsigned int z = rangeZ.x; z < rangeZ.y; ++z) {
+            const uniform float sample =
+                VdbSampler_sampleConstantFloatLeafStrided(
+                    ((const uniform Data1D *uniform)leafPtr),
+                    make_vec3ui(x, y, z),
+                    attributeIndex);
+            extend(valueRange, sample);
+          }
+        }
+      }
+    }
+    break;
+  default: 
+    break;
+  }
+
+  return valueRange;
 }
 
 /*
  * Compute the value range on the given constant float leaf.
  */
-export void EXPORT_UNIQUE(VdbSampler_valueRangeConstantFloat,
+export void EXPORT_UNIQUE(VdbSampler_computeValueRange,
                           const void *uniform _grid,
-                          const vec3ui *uniform offset,
+                          const vec3ui *uniform centerNodeOffset,
                           uniform uint32 level,
                           uniform uint32 attributeIndex,
                           uniform box1f *uniform range)
 {
+  const VdbGrid *uniform grid = (const VdbGrid *uniform)_grid;
+
+  // We need to consider voxels outside this node in the value range because of
+  // reconstruction filter support.
+  //
+  // In particular, we support filters that use up to a 4x4x4 neighborhood.
+  //
+  // In the following 1D illustration, we show voxels in the current leaf node
+  // marked as a '#'. We also show voxels from neighboring nodes, and
+  // mark those that could be used for interpolation with a '?'.
+  //
+  // | | | | | | | |?|#|#|#|#|#|#|#|#|?|?| | | | | | |
+  //
+  // We therefore iterate over a 3x3x3 neighborhood of nodes and compute the
+  // value range over the following voxel range (in each dimension):
+  //
+  // "previous node" : [VKL_VDB_RES_LEAF - 1, VKL_VDB_RES_LEAF[
+  // "central node"  : [0, VKL_VDB_RES_LEAF[
+  // "next node"     : [0, 2[
+  //
+  // Note that the upper limit is exclusive. In 3D, each dimension follows
+  // the same pattern.
+  //
+  // The following array encodes the iteration range for the previous,
+  // central, and next leaf node:
+  const uniform vec2ui voxelRange[] = {
+    make_vec2ui(VKL_VDB_RES_LEAF - 1, VKL_VDB_RES_LEAF),
+    make_vec2ui(0, VKL_VDB_RES_LEAF),
+    make_vec2ui(0, 2)
+  };
+
   range->lower = pos_inf;
   range->upper = neg_inf;
-  extendValueRangeFilterFloat(_grid, offset, level, attributeIndex, range);
+
+  uniform VdbSampler sampler;
+  memset(&sampler, 0, sizeof(uniform VdbSampler));
+  sampler.grid             = (const VdbGrid *uniform)_grid;
+  sampler.maxSamplingDepth = VKL_VDB_NUM_LEVELS - 1;
+
+  // Iterate over a 3x3x3 neighborhood of leaf nodes, and compute the value
+  // range for each of those nodes.
+  for (uniform int ix = 0; ix < 3; ++ix) {
+    for (uniform int iy = 0; iy < 3; ++iy) {
+      for (uniform int iz = 0; iz < 3; ++iz) {
+        // This is the offset of the current node relative to the central node.
+        const uniform vec3i nodeOffset =
+            make_vec3i((ix - 1) * VKL_VDB_RES_LEAF,
+                       (iy - 1) * VKL_VDB_RES_LEAF,
+                       (iz - 1) * VKL_VDB_RES_LEAF);
+
+        // *centerNodeOffset is relative to the rootOrigin, but 
+        // VdbSampler_traverse expects an index space coordinate. 
+        const uniform vec3i nodeOrigin = make_vec3i(
+          grid->rootOrigin.x + centerNodeOffset->x + nodeOffset.x,
+          grid->rootOrigin.y + centerNodeOffset->y + nodeOffset.y,
+          grid->rootOrigin.z + centerNodeOffset->z + nodeOffset.z
+        );
+
+        uniform uint64 voxel;
+        uniform vec3ui domainOffset;
+        VdbSampler_traverse(&sampler, nodeOrigin, voxel, domainOffset);
+
+        const uniform range1f nodeValueRange =
+            VdbSampler_computeNodeValueRange(grid,
+                                             voxel,
+                                             voxelRange[ix],
+                                             voxelRange[iy],
+                                             voxelRange[iz],
+                                             attributeIndex);
+
+        range->lower = min(range->lower, nodeValueRange.lower);
+        range->upper = max(range->upper, nodeValueRange.upper);
+      }
+    }
+  }
 }
 
 // ---------------------------------------------------------------------------

openvkl/drivers/ispc/volume/vdb/VdbVolume.cpp

@@ -251,32 +251,12 @@ namespace openvkl {
     {
       range1f range;
 
-      switch (format) {
-      case VKL_FORMAT_TILE: {
-        CALL_ISPC(VdbSampler_valueRangeTileFloat,
-                  grid,
-                  reinterpret_cast<const ispc::vec3ui *>(&offset),
-                  level,
-                  attributeIndex,
-                  reinterpret_cast<ispc::box1f *>(&range));
-        break;
-      }
-
-      case VKL_FORMAT_CONSTANT_ZYX: {
-        range1f leafRange;
-        CALL_ISPC(VdbSampler_valueRangeConstantFloat,
-                  grid,
-                  reinterpret_cast<const ispc::vec3ui *>(&offset),
-                  level,
-                  attributeIndex,
-                  reinterpret_cast<ispc::box1f *>(&range));
-        break;
-      }
-
-      default:
-        runtimeError(
-            "Only VKL_FORMAT_TILE and VKL_FORMAT_CONSTANT_ZYX are supported.");
-      }
+      CALL_ISPC(VdbSampler_computeValueRange,
+                grid,
+                reinterpret_cast<const ispc::vec3ui *>(&offset),
+                level,
+                attributeIndex,
+                reinterpret_cast<ispc::box1f *>(&range));
 
       return range;
     }
@@ -373,9 +353,10 @@ namespace openvkl {
         const auto format    = static_cast<VKLFormat>(leafFormat[idx]);
         const vec3ui &offset = leafOffsets[idx];
 
-        for (unsigned int j = 0; j < grid->numAttributes; j++)
+        for (unsigned int j = 0; j < grid->numAttributes; j++) {
           valueRanges[idx][j] =
               computeValueRangeFloat(grid, format, leafLevel[idx], offset, j);
+        }
       });
 
       for (size_t idx = 0; idx < numLeaves; ++idx) {
@@ -388,20 +369,22 @@ namespace openvkl {
           assert(nodeIndex < level.numNodes);
 
           const uint64_t voxelIndex = offsetToLinearVoxelIndex(offset, l);
-          // NOTE: If this is every greater than 2^32-1 then we will have to
+          // NOTE: If this is ever greater than 2^32-1 then we will have to
           // use 64 bit addressing.
           const uint64_t v = nodeIndex * vklVdbLevelNumVoxels(l) + voxelIndex;
           assert(v < ((uint64_t)1) << 32);
 
-          for (unsigned int j = 0; j < grid->numAttributes; j++)
+          for (unsigned int j = 0; j < grid->numAttributes; j++) {
             level.valueRange[v * grid->numAttributes + j].extend(
                 valueRanges[idx][j]);
+          }
 
           uint64_t &voxel = level.voxels[v];
           assert(!vklVdbVoxelIsEmpty(voxel));
 
-          if (vklVdbVoxelIsLeafPtr(voxel))
+          if (vklVdbVoxelIsLeafPtr(voxel)) {
             break;
+          }
 
           nodeIndex = vklVdbVoxelChildGetIndex(voxel);
           assert(nodeIndex < grid->levels[l + 1].numNodes);