Fix T99004: scaling volume down results in crash

OpenVDB crashes when the determinant of the grid transformation is
too small. The solution is too detect when the determinant is too small
and to replace the grid with an empty one. If possible the translation
and rotation of the grid remains unchanged.

Differential Revision: https://developer.blender.org/D15806

6 files changed, 114 insertions, 30 deletions

diff --git a/source/blender/blenkernel/BKE_volume.h b/source/blender/blenkernel/BKE_volume.h
index bc578ef8b28..4215cad5858 100644
--- a/source/blender/blenkernel/BKE_volume.h
+++ b/source/blender/blenkernel/BKE_volume.h
@@ -114,6 +114,7 @@ int BKE_volume_grid_channels(const struct VolumeGrid *grid);
  * Transformation from index space to object space.
  */
 void BKE_volume_grid_transform_matrix(const struct VolumeGrid *grid, float mat[4][4]);
+void BKE_volume_grid_transform_matrix_set(struct VolumeGrid *volume_grid, const float mat[4][4]);
 
 /* Volume Editing
  *
@@ -133,6 +134,11 @@ struct VolumeGrid *BKE_volume_grid_add(struct Volume *volume,
                                        VolumeGridType type);
 void BKE_volume_grid_remove(struct Volume *volume, struct VolumeGrid *grid);
 
+/**
+ * Openvdb crashes when the determinant of the transform matrix becomes too small.
+ */
+bool BKE_volume_grid_determinant_valid(double determinant);
+
 /* Simplify */
 int BKE_volume_simplify_level(const struct Depsgraph *depsgraph);
 float BKE_volume_simplify_factor(const struct Depsgraph *depsgraph);
@@ -186,6 +192,9 @@ openvdb::GridBase::Ptr BKE_volume_grid_openvdb_for_write(const struct Volume *vo
                                                          struct VolumeGrid *grid,
                                                          bool clear);
 
+void BKE_volume_grid_clear_tree(Volume &volume, VolumeGrid &volume_grid);
+void BKE_volume_grid_clear_tree(openvdb::GridBase &grid);
+
 VolumeGridType BKE_volume_grid_type_openvdb(const openvdb::GridBase &grid);
 
 template<typename OpType>
diff --git a/source/blender/blenkernel/intern/volume.cc b/source/blender/blenkernel/intern/volume.cc
index 2f5bf0de58f..6ba396259aa 100644
--- a/source/blender/blenkernel/intern/volume.cc
+++ b/source/blender/blenkernel/intern/volume.cc
@@ -1361,6 +1361,26 @@ const char *BKE_volume_grid_name(const VolumeGrid *volume_grid)
 }
 
 #ifdef WITH_OPENVDB
+struct ClearGridOp {
+  openvdb::GridBase &grid;
+
+  template<typename Grid> void operator()()
+  {
+    static_cast<Grid &>(grid).clear();
+  }
+};
+void BKE_volume_grid_clear_tree(openvdb::GridBase &grid)
+{
+  const VolumeGridType grid_type = BKE_volume_grid_type_openvdb(grid);
+  ClearGridOp op{grid};
+  BKE_volume_grid_type_operation(grid_type, op);
+}
+void BKE_volume_grid_clear_tree(Volume &volume, VolumeGrid &volume_grid)
+{
+  openvdb::GridBase::Ptr grid = BKE_volume_grid_openvdb_for_write(&volume, &volume_grid, false);
+  BKE_volume_grid_clear_tree(*grid);
+}
+
 VolumeGridType BKE_volume_grid_type_openvdb(const openvdb::GridBase &grid)
 {
   if (grid.isType<openvdb::FloatGrid>()) {
@@ -1451,6 +1471,23 @@ void BKE_volume_grid_transform_matrix(const VolumeGrid *volume_grid, float mat[4
 #endif
 }
 
+void BKE_volume_grid_transform_matrix_set(struct VolumeGrid *volume_grid, const float mat[4][4])
+{
+#ifdef WITH_OPENVDB
+  openvdb::math::Mat4f mat_openvdb;
+  for (int col = 0; col < 4; col++) {
+    for (int row = 0; row < 4; row++) {
+      mat_openvdb(col, row) = mat[col][row];
+    }
+  }
+  openvdb::GridBase::Ptr grid = volume_grid->grid();
+  grid->setTransform(std::make_shared<openvdb::math::Transform>(
+      std::make_shared<openvdb::math::AffineMap>(mat_openvdb)));
+#else
+  UNUSED_VARS(grid, mat);
+#endif
+}
+
 /* Grid Tree and Voxels */
 
 /* Volume Editing */
@@ -1547,6 +1584,16 @@ void BKE_volume_grid_remove(Volume *volume, VolumeGrid *grid)
 #endif
 }
 
+bool BKE_volume_grid_determinant_valid(const double determinant)
+{
+#ifdef WITH_OPENVDB
+  /* Limit taken from openvdb/math/Maps.h. */
+  return std::abs(determinant) >= 3.0 * openvdb::math::Tolerance<double>::value();
+#else
+  return true;
+#endif
+}
+
 int BKE_volume_simplify_level(const Depsgraph *depsgraph)
 {
   if (DEG_get_mode(depsgraph) != DAG_EVAL_RENDER) {
diff --git a/source/blender/nodes/geometry/node_geometry_util.hh b/source/blender/nodes/geometry/node_geometry_util.hh
index efb7efaf1cc..a4af608a40e 100644
--- a/source/blender/nodes/geometry/node_geometry_util.hh
+++ b/source/blender/nodes/geometry/node_geometry_util.hh
@@ -34,7 +34,8 @@ void transform_mesh(Mesh &mesh,
                     const float3 rotation,
                     const float3 scale);
 
-void transform_geometry_set(GeometrySet &geometry,
+void transform_geometry_set(GeoNodeExecParams &params,
+                            GeometrySet &geometry,
                             const float4x4 &transform,
                             const Depsgraph &depsgraph);
 
diff --git a/source/blender/nodes/geometry/nodes/node_geo_object_info.cc b/source/blender/nodes/geometry/nodes/node_geo_object_info.cc
index 0b2159364f1..691988249df 100644
--- a/source/blender/nodes/geometry/nodes/node_geo_object_info.cc
+++ b/source/blender/nodes/geometry/nodes/node_geo_object_info.cc
@@ -80,7 +80,7 @@ static void node_geo_exec(GeoNodeExecParams params)
     else {
       geometry_set = bke::object_get_evaluated_geometry_set(*object);
       if (transform_space_relative) {
-        transform_geometry_set(geometry_set, transform, *params.depsgraph());
+        transform_geometry_set(params, geometry_set, transform, *params.depsgraph());
       }
     }
 
diff --git a/source/blender/nodes/geometry/nodes/node_geo_transform.cc b/source/blender/nodes/geometry/nodes/node_geo_transform.cc
index 945d5fbdcac..adf55abbbec 100644
--- a/source/blender/nodes/geometry/nodes/node_geo_transform.cc
+++ b/source/blender/nodes/geometry/nodes/node_geo_transform.cc
@@ -83,44 +83,61 @@ static void transform_instances(InstancesComponent &instances, const float4x4 &t
   }
 }
 
-static void transform_volume(Volume &volume, const float4x4 &transform, const Depsgraph &depsgraph)
+static void transform_volume(GeoNodeExecParams &params,
+                             Volume &volume,
+                             const float4x4 &transform,
+                             const Depsgraph &depsgraph)
 {
 #ifdef WITH_OPENVDB
-  /* Scaling an axis to zero is not supported for volumes. */
-  const float3 translation = transform.translation();
-  const float3 rotation = transform.to_euler();
-  const float3 scale = transform.scale();
-  const float3 limited_scale = {
-      (scale.x == 0.0f) ? FLT_EPSILON : scale.x,
-      (scale.y == 0.0f) ? FLT_EPSILON : scale.y,
-      (scale.z == 0.0f) ? FLT_EPSILON : scale.z,
-  };
-  const float4x4 scale_limited_transform = float4x4::from_loc_eul_scale(
-      translation, rotation, limited_scale);
-
   const Main *bmain = DEG_get_bmain(&depsgraph);
   BKE_volume_load(&volume, bmain);
 
   openvdb::Mat4s vdb_matrix;
-  memcpy(vdb_matrix.asPointer(), &scale_limited_transform, sizeof(float[4][4]));
+  memcpy(vdb_matrix.asPointer(), &transform, sizeof(float[4][4]));
   openvdb::Mat4d vdb_matrix_d{vdb_matrix};
 
+  bool found_too_small_scale = false;
   const int grids_num = BKE_volume_num_grids(&volume);
   for (const int i : IndexRange(grids_num)) {
     VolumeGrid *volume_grid = BKE_volume_grid_get_for_write(&volume, i);
-
-    openvdb::GridBase::Ptr grid = BKE_volume_grid_openvdb_for_write(&volume, volume_grid, false);
-    openvdb::math::Transform &grid_transform = grid->transform();
-    grid_transform.postMult(vdb_matrix_d);
+    float4x4 grid_matrix;
+    BKE_volume_grid_transform_matrix(volume_grid, grid_matrix.values);
+    mul_m4_m4_pre(grid_matrix.values, transform.values);
+    const float determinant = determinant_m4(grid_matrix.values);
+    if (!BKE_volume_grid_determinant_valid(determinant)) {
+      found_too_small_scale = true;
+      /* Clear the tree because it is too small. */
+      BKE_volume_grid_clear_tree(volume, *volume_grid);
+      if (determinant == 0) {
+        /* Reset rotation and scale. */
+        copy_v3_fl3(grid_matrix.values[0], 1, 0, 0);
+        copy_v3_fl3(grid_matrix.values[1], 0, 1, 0);
+        copy_v3_fl3(grid_matrix.values[2], 0, 0, 1);
+      }
+      else {
+        /* Keep rotation but reset scale. */
+        normalize_v3(grid_matrix.values[0]);
+        normalize_v3(grid_matrix.values[1]);
+        normalize_v3(grid_matrix.values[2]);
+      }
+    }
+    BKE_volume_grid_transform_matrix_set(volume_grid, grid_matrix.values);
+  }
+  if (found_too_small_scale) {
+    params.error_message_add(NodeWarningType::Warning,
+                             TIP_("Volume scale is lower than permitted by OpenVDB"));
   }
 #else
   UNUSED_VARS(volume, transform, depsgraph);
 #endif
 }
 
-static void translate_volume(Volume &volume, const float3 translation, const Depsgraph &depsgraph)
+static void translate_volume(GeoNodeExecParams &params,
+                             Volume &volume,
+                             const float3 translation,
+                             const Depsgraph &depsgraph)
 {
-  transform_volume(volume, float4x4::from_location(translation), depsgraph);
+  transform_volume(params, volume, float4x4::from_location(translation), depsgraph);
 }
 
 static void transform_curve_edit_hints(bke::CurvesEditHints &edit_hints, const float4x4 &transform)
@@ -151,7 +168,8 @@ static void translate_curve_edit_hints(bke::CurvesEditHints &edit_hints, const f
   }
 }
 
-static void translate_geometry_set(GeometrySet &geometry,
+static void translate_geometry_set(GeoNodeExecParams &params,
+                                   GeometrySet &geometry,
                                    const float3 translation,
                                    const Depsgraph &depsgraph)
 {
@@ -165,7 +183,7 @@ static void translate_geometry_set(GeometrySet &geometry,
     translate_pointcloud(*pointcloud, translation);
   }
   if (Volume *volume = geometry.get_volume_for_write()) {
-    translate_volume(*volume, translation, depsgraph);
+    translate_volume(params, *volume, translation, depsgraph);
   }
   if (geometry.has_instances()) {
     translate_instances(geometry.get_component_for_write<InstancesComponent>(), translation);
@@ -175,7 +193,8 @@ static void translate_geometry_set(GeometrySet &geometry,
   }
 }
 
-void transform_geometry_set(GeometrySet &geometry,
+void transform_geometry_set(GeoNodeExecParams &params,
+                            GeometrySet &geometry,
                             const float4x4 &transform,
                             const Depsgraph &depsgraph)
 {
@@ -189,7 +208,7 @@ void transform_geometry_set(GeometrySet &geometry,
     transform_pointcloud(*pointcloud, transform);
   }
   if (Volume *volume = geometry.get_volume_for_write()) {
-    transform_volume(*volume, transform, depsgraph);
+    transform_volume(params, *volume, transform, depsgraph);
   }
   if (geometry.has_instances()) {
     transform_instances(geometry.get_component_for_write<InstancesComponent>(), transform);
@@ -230,10 +249,11 @@ static void node_geo_exec(GeoNodeExecParams params)
 
   /* Use only translation if rotation and scale don't apply. */
   if (use_translate(rotation, scale)) {
-    translate_geometry_set(geometry_set, translation, *params.depsgraph());
+    translate_geometry_set(params, geometry_set, translation, *params.depsgraph());
   }
   else {
-    transform_geometry_set(geometry_set,
+    transform_geometry_set(params,
+                           geometry_set,
                            float4x4::from_loc_eul_scale(translation, rotation, scale),
                            *params.depsgraph());
   }
diff --git a/source/blender/nodes/geometry/nodes/node_geo_volume_cube.cc b/source/blender/nodes/geometry/nodes/node_geo_volume_cube.cc
index d7e9e38ea0d..e964bf03ed2 100644
--- a/source/blender/nodes/geometry/nodes/node_geo_volume_cube.cc
+++ b/source/blender/nodes/geometry/nodes/node_geo_volume_cube.cc
@@ -137,6 +137,14 @@ static void node_geo_exec(GeoNodeExecParams params)
     return;
   }
 
+  const double3 scale_fac = double3(bounds_max - bounds_min) / double3(resolution - 1);
+  if (!BKE_volume_grid_determinant_valid(scale_fac.x * scale_fac.y * scale_fac.z)) {
+    params.error_message_add(NodeWarningType::Warning,
+                             TIP_("Volume scale is lower than permitted by OpenVDB"));
+    params.set_default_remaining_outputs();
+    return;
+  }
+
   Field<float> input_field = params.extract_input<Field<float>>("Density");
 
   /* Evaluate input field on a 3D grid. */
@@ -157,8 +165,7 @@ static void node_geo_exec(GeoNodeExecParams params)
   openvdb::tools::copyFromDense(dense_grid, *grid, 0.0f);
 
   grid->transform().preTranslate(openvdb::math::Vec3<float>(-0.5f));
-  const float3 scale_fac = (bounds_max - bounds_min) / float3(resolution - 1);
-  grid->transform().postScale(openvdb::math::Vec3<float>(scale_fac.x, scale_fac.y, scale_fac.z));
+  grid->transform().postScale(openvdb::math::Vec3<double>(scale_fac.x, scale_fac.y, scale_fac.z));
   grid->transform().postTranslate(
       openvdb::math::Vec3<float>(bounds_min.x, bounds_min.y, bounds_min.z));