From abe5adf7ecdd739effa167eda7e4233604c6bb1a Mon Sep 17 00:00:00 2001
From: ghurstunither <62885595+ghurstunither@users.noreply.github.com>
Date: Tue, 16 Apr 2024 11:46:14 -0400
Subject: [PATCH] Update LevelSetFilter.h

Adds fillet filter, which rounds off concave edges to create a smoother transition between surfaces.

Signed-off-by: ghurstunither <62885595+ghurstunither@users.noreply.github.com>
---
 openvdb/openvdb/tools/LevelSetFilter.h | 77 ++++++++++++++++++++++++++
 1 file changed, 77 insertions(+)

diff --git a/openvdb/openvdb/tools/LevelSetFilter.h b/openvdb/openvdb/tools/LevelSetFilter.h
index a31b7344be..08d96e60dc 100644
--- a/openvdb/openvdb/tools/LevelSetFilter.h
+++ b/openvdb/openvdb/tools/LevelSetFilter.h
@@ -95,6 +95,19 @@ class LevelSetFilter : public LevelSetTracker<GridT, InterruptT>
         Filter f(this, mask); f.meanCurvature();
     }
 
+    /// @brief One iteration of filleting on the level set.
+    /// @param mask Optional alpha mask.
+    ///
+    /// @note This filter rounds off concave edges to create a smoother transition between surfaces
+    /// Fillets a level set by
+    ///   offsetting at locations with a negative principal curvature, proportional to its magnitude
+    ///   leaves locations with non-negative principal curvatures untouched
+    /// This filter converges to the convex hull if iterated enough times
+    void fillet(const MaskType* mask = nullptr)
+    {
+        Filter f(this, mask); f.fillet();
+    }
+
     /// @brief One iteration of Laplacian flow of the level set.
     /// @param mask Optional alpha mask.
     void laplacian(const MaskType* mask = nullptr)
@@ -168,6 +181,7 @@ class LevelSetFilter : public LevelSetTracker<GridT, InterruptT>
         void gaussian(int width);
         void laplacian();
         void meanCurvature();
+        void fillet();
         void offset(ValueType value);
         void operator()(const LeafRange& r) const
         {
@@ -210,6 +224,7 @@ class LevelSetFilter : public LevelSetTracker<GridT, InterruptT>
 
         void medianImpl(const LeafRange&, int);
         void meanCurvatureImpl(const LeafRange&);
+        void filletImpl(const LeafRange&);
         void laplacianImpl(const LeafRange&);
         void offsetImpl(const LeafRange&, ValueType);
 
@@ -301,6 +316,22 @@ LevelSetFilter<GridT, MaskT, InterruptT>::Filter::meanCurvature()
     mParent->endInterrupter();
 }
 
+template<typename GridT, typename MaskT, typename InterruptT>
+inline void
+LevelSetFilter<GridT, MaskT, InterruptT>::Filter::fillet()
+{
+    mParent->startInterrupter("Filleting level set");
+
+    mParent->leafs().rebuildAuxBuffers(1, mParent->getGrainSize()==0);
+
+    mTask = std::bind(&Filter::filletImpl, std::placeholders::_1, std::placeholders::_2);
+    this->cook(true);
+
+    mParent->track();
+
+    mParent->endInterrupter();
+}
+
 template<typename GridT, typename MaskT, typename InterruptT>
 inline void
 LevelSetFilter<GridT, MaskT, InterruptT>::Filter::laplacian()
@@ -378,6 +409,52 @@ LevelSetFilter<GridT, MaskT, InterruptT>::Filter::meanCurvatureImpl(const LeafRa
     }
 }
 
+/// Fillets a level set by
+///   offsetting at locations with a negative principal curvature, proportional to its magnitude
+///   leaves locations with non-negative principal curvatures untouched
+/// This filter converges to the convex hull if iterated enough times
+template<typename GridT, typename MaskT, typename InterruptT>
+inline void
+LevelSetFilter<GridT, MaskT, InterruptT>::Filter::filletImpl(const LeafRange& range)
+{
+    mParent->checkInterrupter();
+
+    const ValueType dx = mParent->voxelSize(), dt = math::Pow2(dx) / ValueType(3);
+    math::CurvatureStencil<GridType> stencil(mParent->grid(), dx);
+
+    if (mMask) {
+        typename AlphaMaskT::FloatType a, b;
+        AlphaMaskT alpha(mParent->grid(), *mMask, mParent->minMask(),
+                         mParent->maxMask(), mParent->isMaskInverted());
+        for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
+            ValueType* buffer = leafIter.buffer(1).data();
+            for (VoxelCIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
+                if (alpha(iter.getCoord(), a, b)) {
+                    stencil.moveTo(iter);
+
+                    const ValueType kappa = stencil.principalCurvatures().first;
+
+                    const ValueType phi0 = *iter,
+                                    phi1 = phi0 + math::Min(ValueType(0), dt*kappa);
+                    buffer[iter.pos()] = b * phi0 + a * phi1;
+                }
+            }
+        }
+    } else {
+        for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
+            ValueType* buffer = leafIter.buffer(1).data();
+            for (VoxelCIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
+                stencil.moveTo(iter);
+
+                const ValueType kappa = stencil.principalCurvatures().first;
+
+                if (math::isNegative(kappa))
+                    buffer[iter.pos()] = *iter + dt*kappa;
+            }
+        }
+    }
+}
+
 /// Performs Laplacian diffusion. Note if the grids contains a true
 /// signed distance field (e.g. a solution to the Eikonal equation)
 /// Laplacian diffusions (e.g. geometric heat equation) is actually