From 071353eddc64ecbb71939cfd2b0f6ee5a0c07125 Mon Sep 17 00:00:00 2001
From: Thomas Krijnen <t.krijnen@gmail.com>
Date: Wed, 6 Aug 2025 05:46:07 +0200
Subject: [PATCH] Vertex normal smoothing in CGAL --cgal-smooth-angle-degrees
 #6997

---
 src/ifcgeom/ConversionSettings.h              |  8 ++-
 .../kernels/cgal/CgalConversionResult.cpp     | 56 ++++++++++++++++---
 2 files changed, 55 insertions(+), 9 deletions(-)

diff --git a/src/ifcgeom/ConversionSettings.h b/src/ifcgeom/ConversionSettings.h
index a9e24d5c4c..2589820ab6 100644
--- a/src/ifcgeom/ConversionSettings.h
+++ b/src/ifcgeom/ConversionSettings.h
@@ -364,6 +364,12 @@ namespace ifcopenshell {
 				static constexpr int defaultvalue = 16;
 			};
 
+			struct CgalSmoothAngleDegrees : public SettingBase<CgalSmoothAngleDegrees, double> {
+				static constexpr const char* const name = "cgal-smooth-angle-degrees";
+				static constexpr const char* const description = "Angle in degrees under which adjacent facets will have averaged vertex normals in CGAL output. NB irrespective of original IFC geometry types. Defaults to -1 to disable smoothing.";
+				static constexpr int defaultvalue = -1.;
+			};
+
 			struct KeepBoundingBoxes : public SettingBase<KeepBoundingBoxes, bool> {
 				static constexpr const char* const name = "keep-bounding-boxes";
 				static constexpr const char* const description =
@@ -640,7 +646,7 @@ namespace ifcopenshell {
 		};
 
 		class IFC_GEOM_API Settings : public SettingsContainer<
-                                          std::tuple<MesherLinearDeflection, MesherAngularDeflection, ReorientShells, LengthUnit, PlaneUnit, Precision, OutputDimensionality, LayersetFirst, DisableBooleanResult, NoWireIntersectionCheck, NoWireIntersectionTolerance, PrecisionFactor, DebugBooleanOperations, BooleanAttempt2d, SurfaceColour, WeldVertices, UseWorldCoords, UnifyShapes, UseMaterialNames, ConvertBackUnits, ContextIds, ContextTypes, ContextIdentifiers, IteratorOutput, DisableOpeningSubtractions, ApplyDefaultMaterials, DontEmitNormals, GenerateUvs, ApplyLayerSets, UseElementHierarchy, ValidateQuantities, EdgeArrows, BuildingLocalPlacement, SiteLocalPlacement, ForceSpaceTransparency, CircleSegments, KeepBoundingBoxes, ComputeCurvature, FunctionStepType, FunctionStepParam, NoParallelMapping, PermissiveShapeReuse, ModelOffset, ModelRotation, TriangulationType, CgalEmitOriginalEdges, OcctNoCleanTriangulation, CacheShapes, DeferProcessingFirstElement, MaxOffset, MaxOffsetDeviation, ApplyOffset>
+                                          std::tuple<MesherLinearDeflection, MesherAngularDeflection, ReorientShells, LengthUnit, PlaneUnit, Precision, OutputDimensionality, LayersetFirst, DisableBooleanResult, NoWireIntersectionCheck, NoWireIntersectionTolerance, PrecisionFactor, DebugBooleanOperations, BooleanAttempt2d, SurfaceColour, WeldVertices, UseWorldCoords, UnifyShapes, UseMaterialNames, ConvertBackUnits, ContextIds, ContextTypes, ContextIdentifiers, IteratorOutput, DisableOpeningSubtractions, ApplyDefaultMaterials, DontEmitNormals, GenerateUvs, ApplyLayerSets, UseElementHierarchy, ValidateQuantities, EdgeArrows, BuildingLocalPlacement, SiteLocalPlacement, ForceSpaceTransparency, CircleSegments, CgalSmoothAngleDegrees, KeepBoundingBoxes, ComputeCurvature, FunctionStepType, FunctionStepParam, NoParallelMapping, PermissiveShapeReuse, ModelOffset, ModelRotation, TriangulationType, CgalEmitOriginalEdges, OcctNoCleanTriangulation, CacheShapes, DeferProcessingFirstElement, MaxOffset, MaxOffsetDeviation, ApplyOffset>
 		>
 		{};
 }
diff --git a/src/ifcgeom/kernels/cgal/CgalConversionResult.cpp b/src/ifcgeom/kernels/cgal/CgalConversionResult.cpp
index eaa6dad6f9..ef2842107e 100644
--- a/src/ifcgeom/kernels/cgal/CgalConversionResult.cpp
+++ b/src/ifcgeom/kernels/cgal/CgalConversionResult.cpp
@@ -217,6 +217,14 @@ void ifcopenshell::geometry::CgalShape::Triangulate(ifcopenshell::geometry::Sett
 		}
 	}
 
+	boost::optional<double> smooth_treshold;
+	{
+		auto setting_value = settings.get<ifcopenshell::geometry::settings::CgalSmoothAngleDegrees>().get();
+		if (setting_value > 0.) {
+			smooth_treshold = std::cos(setting_value * boost::math::constants::pi<double>() / 180.0);
+		}
+	}
+
 	if (!all_triangles) {
 		if (!shape_to_use->is_valid()) {
 			Logger::Message(Logger::LOG_ERROR, "Invalid Polyhedron_3 in object (before triangulation)");
@@ -261,8 +269,8 @@ void ifcopenshell::geometry::CgalShape::Triangulate(ifcopenshell::geometry::Sett
 	// boost::associative_property_map<std::map<cgal_vertex_descriptor_t, Kernel_::Vector_3>> vertex_normals_map(vertex_normals);
 	
 	// Triangulate the shape and compute the normals
-	std::map<cgal_face_descriptor_t, Kernel_::Vector_3> face_normals;
-	boost::associative_property_map<std::map<cgal_face_descriptor_t, Kernel_::Vector_3>> face_normals_map(face_normals);
+	std::map<Facet_const_handle, Kernel_::Vector_3> face_normals;
+	boost::associative_property_map<std::map<Facet_const_handle, Kernel_::Vector_3>> face_normals_map(face_normals);
 
 	//  CGAL::Polygon_mesh_processing::compute_normals(s, vertex_normals_map, face_normals_map);
 	try {
@@ -286,17 +294,49 @@ void ifcopenshell::geometry::CgalShape::Triangulate(ifcopenshell::geometry::Sett
 			continue;
 		}
 		CGAL::Polyhedron_3<Kernel_>::Halfedge_around_facet_const_circulator current_halfedge = face->facet_begin();
+
+		const Kernel_::Vector_3 facet_normal = face_normals_map[face];
+
 		int vertexidx[3];
 		bool is_face_boundary[3];
 		int i = 0;
 		do {
+			auto v = current_halfedge->vertex();
+
+			auto vertex_norm = facet_normal;
+
+			if (smooth_treshold) {
+				Kernel_::Vector_3 normal_accum(0, 0, 0);
+				{
+					// circulator around the vertex
+					auto vh_begin = v->vertex_begin();
+					if (vh_begin != nullptr) {
+						auto vh = vh_begin;
+						do {
+							if (!vh->is_border()) {
+								Facet_const_handle adj_f = vh->facet();
+								const auto fn2 = face_normals_map[adj_f];
+								if ((fn2 * facet_normal) >= *smooth_treshold) {
+									normal_accum = normal_accum + fn2;
+								}
+								++vh;
+							}
+						} while (vh != vh_begin);
+					}
+				}
+				const double len = std::sqrt(CGAL::to_double(normal_accum.squared_length()));
+				if (len > 0) {
+					vertex_norm = normal_accum / len;
+				}
+			}
+
 			postion_normal pn = {
-				current_halfedge->vertex()->point().cartesian(0),
-				current_halfedge->vertex()->point().cartesian(1),
-				current_halfedge->vertex()->point().cartesian(2),
-				face_normals_map[face].cartesian(0),
-				face_normals_map[face].cartesian(1),
-				face_normals_map[face].cartesian(2)
+				v->point().cartesian(0),
+				v->point().cartesian(1),
+				v->point().cartesian(2),
+				vertex_norm.cartesian(0),
+				vertex_norm.cartesian(1),
+				vertex_norm.cartesian(2)
 			};
 
 			// @todo normalzie based on largest component?