songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Harmonize MC/TMC vertex point computations

This commit is contained in:
Mael Rouxel-Labbé 2025-03-26 22:29:01 +01:00
parent 08109aa12b
commit e843cefce0
2 changed files with 7 additions and 10 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
Isosurfacing_3/include/CGAL/Isosurfacing_3/internal/marching_cubes_functors.h
View File

@ -78,20 +78,14 @@ typename GeomTraits::Point_3 vertex_interpolation(const typename GeomTraits::Poi
typename GeomTraits::Compute_z_3 z_coord = gt.compute_z_3_object(); typename GeomTraits::Compute_z_3 z_coord = gt.compute_z_3_object();
typename GeomTraits::Construct_point_3 point = gt.construct_point_3_object(); typename GeomTraits::Construct_point_3 point = gt.construct_point_3_object();
FT mu = FT(0);
// @todo, technically we should be using the edge intersection oracle here, but there is a nuance // @todo, technically we should be using the edge intersection oracle here, but there is a nuance
// between MC and DC on the handling of edges that have val0 = val1 = isovalue: in MC we assume // between MC and DC on the handling of edges that have val0 = val1 = isovalue: in MC we assume
// the isosurface is in the middle, in DC we assume the isosurface is not intersecting the edge. // the isosurface is in the middle, in DC we assume the isosurface is not intersecting the edge.
// In the oracle, we follow DC right now. Could put a Boolean parameter, but it's ugly. // In the oracle, we follow DC right now. Could put a Boolean parameter, but it's ugly.
// don't divide by 0 const FT den = d1 - d0;
if(abs(d1 - d0) < 0.000001) // @fixme hardcoded bound FT mu = is_zero(den) ? FT(1) / FT(2) : (isovalue - d0) / den;
mu = FT(0.5); // if both points have the same value, assume isolevel is in the middle mu = std::clamp<FT>(mu, FT(0), FT(1));
else
mu = (isovalue - d0) / (d1 - d0);
CGAL_assertion(mu >= FT(0.0) || mu <= FT(1.0));
// linear interpolation // linear interpolation
return point((FT(1) - mu) * x_coord(p0) + mu * x_coord(p1), return point((FT(1) - mu) * x_coord(p0) + mu * x_coord(p1),

5
Isosurfacing_3/include/CGAL/Isosurfacing_3/internal/topologically_correct_marching_cubes_functors.h
View File

@ -400,7 +400,10 @@ private:
get_edge_vertex(eg, v0, v1, l_edges_); get_edge_vertex(eg, v0, v1, l_edges_);
// @todo use the domain's interpolation scheme? // @todo use the domain's interpolation scheme?
FT l = (i0 - values[v0]) / (values[v1] - values[v0]); const FT den = values[v1] - values[v0];
FT l = is_zero(den) ? FT(1) / FT(2) : (i0 - values[v0]) / den;
l = std::clamp<FT>(l, FT(0), FT(1));
ecoord[eg] = l; ecoord[eg] = l;
const FT px = (FT(1) - l) * x_coord(corners[v0]) + l * x_coord(corners[v1]); const FT px = (FT(1) - l) * x_coord(corners[v0]) + l * x_coord(corners[v1]);