updated test

src/manifold/include/manifold.h

@@ -100,7 +100,7 @@ class Manifold {
    */
   ///@{
   Mesh GetMesh() const;
-  MeshGL GetMeshGL() const;
+  MeshGL GetMeshGL(glm::ivec3 normalIdx = glm::ivec3(0)) const;
   bool IsEmpty() const;
   enum class Error {
     NO_ERROR,

src/manifold/src/boolean_result.cpp

@@ -412,7 +412,7 @@ struct MapTriRef {
 };
 
 VecDH<TriRef> UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP,
-                              const Manifold::Impl &inQ,
+                              const Manifold::Impl &inQ, bool invertQ,
                               ExecutionPolicy policy) {
   VecDH<TriRef> refPQ = outR.meshRelation_.triRef;
   const int offsetQ = Manifold::Impl::meshIDCounter_;
@@ -425,6 +425,8 @@ VecDH<TriRef> UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP,
   }
   for (const auto &pair : inQ.meshRelation_.meshIDtransform) {
     outR.meshRelation_.meshIDtransform[pair.first + offsetQ] = pair.second;
+    outR.meshRelation_.meshIDtransform[pair.first + offsetQ].backSide ^=
+        invertQ;
   }
   return refPQ;
 }
@@ -692,7 +694,7 @@ Manifold::Impl Boolean3::Result(Manifold::OpType op) const {
   if (ManifoldParams().intermediateChecks)
     ASSERT(outR.IsManifold(), logicErr, "triangulated mesh is not manifold!");
 
-  VecDH<TriRef> refPQ = UpdateReference(outR, inP_, inQ_, policy_);
+  VecDH<TriRef> refPQ = UpdateReference(outR, inP_, inQ_, invertQ, policy_);
 
   outR.SimplifyTopology();
 

src/manifold/src/impl.cpp

@@ -491,12 +491,12 @@ Manifold::Impl::Impl(const MeshGL& meshGL,
       }
 
       if (meshGL.transform.empty()) {
-        relation.meshIDtransform[meshID] = glm::mat4x3(1);
+        relation.meshIDtransform[meshID] = {};
       } else {
         const float* m = meshGL.transform.data() + 12 * i;
-        relation.meshIDtransform[meshID] = {m[0], m[1], m[2],  m[3],
-                                            m[4], m[5], m[6],  m[7],
-                                            m[8], m[9], m[10], m[11]};
+        relation.meshIDtransform[meshID] = {{m[0], m[1], m[2], m[3], m[4], m[5],
+                                             m[6], m[7], m[8], m[9], m[10],
+                                             m[11]}};
       }
     }
   }
@@ -630,7 +630,7 @@ void Manifold::Impl::InitializeOriginal() {
              zip(meshRelation_.triRef.begin(), countAt(0)), NumTri(),
              InitializeTriRef({meshID, halfedge_.cptrD()}));
   meshRelation_.meshIDtransform.clear();
-  meshRelation_.meshIDtransform[meshID] = glm::mat4x3(1);
+  meshRelation_.meshIDtransform[meshID] = {};
 }
 
 void Manifold::Impl::CreateFaces(const std::vector<float>& propertyTolerance) {
@@ -736,7 +736,7 @@ Manifold::Impl Manifold::Impl::Transform(const glm::mat4x3& transform_) const {
   result.halfedgeTangent_.resize(halfedgeTangent_.size());
 
   for (auto& m : result.meshRelation_.meshIDtransform) {
-    m.second = transform_ * glm::mat4(m.second);
+    m.second.transform = transform_ * glm::mat4(m.second.transform);
   }
 
   result.vertPos_.resize(NumVert());
@@ -745,8 +745,7 @@ Manifold::Impl Manifold::Impl::Transform(const glm::mat4x3& transform_) const {
   transform(policy, vertPos_.begin(), vertPos_.end(), result.vertPos_.begin(),
             Transform4x3({transform_}));
 
-  glm::mat3 normalTransform =
-      glm::inverse(glm::transpose(glm::mat3(transform_)));
+  glm::mat3 normalTransform = NormalTransform(transform_);
   transform(policy, faceNormal_.begin(), faceNormal_.end(),
             result.faceNormal_.begin(), TransformNormals({normalTransform}));
   transform(policy, vertNormal_.begin(), vertNormal_.end(),
@@ -845,7 +844,7 @@ void Manifold::Impl::IncrementMeshIDs() {
   for_each_n(policy, meshRelation_.triRef.begin(), numTri,
              UpdateMeshID({meshIDold2new.D(), meshIDstart}));
   // Update keys of the transform map
-  std::map<int, glm::mat4x3> oldTransforms;
+  std::map<int, Relation> oldTransforms;
   std::swap(meshRelation_.meshIDtransform, oldTransforms);
   const int tableSize = meshIDold2new.Size();
   for (int i = 0; i < tableSize; ++i) {

src/manifold/src/impl.h

@@ -27,6 +27,10 @@ namespace manifold {
 
 /** @ingroup Private */
 struct Manifold::Impl {
+  struct Relation {
+    glm::mat4x3 transform = glm::mat4x3(1);
+    bool backSide = false;
+  };
   struct MeshRelationD {
     /// The originalID of this Manifold if it is an original; -1 otherwise.
     int originalID = -1;
@@ -34,7 +38,7 @@ struct Manifold::Impl {
     VecDH<TriRef> triRef;
     VecDH<glm::ivec3> triProperties;
     VecDH<float> properties;
-    std::map<int, glm::mat4x3> meshIDtransform;
+    std::map<int, Relation> meshIDtransform;
   };
 
   Box bBox_;

src/manifold/src/manifold.cpp

@@ -152,14 +152,24 @@ Mesh Manifold::GetMesh() const {
  * to easily push into a renderer, including all interleaved vertex properties
  * that may have been input. It also includes relations to all the input meshes
  * that form a part of this result and the transforms applied to each.
+ *
+ * @param normalIdx If the original MeshGL inputs that formed this manifold had
+ * properties corresponding to normal vectors, you can specify which property
+ * channels these are (x, y, z), which will cause this output MeshGL to
+ * automatically update these normals according to the applied transforms and
+ * front/back side. Each channel must be >= 3 and < numProp, and all original
+ * MeshGLs must use the same channels for their normals.
  */
-MeshGL Manifold::GetMeshGL() const {
+MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const {
   const Impl& impl = *GetCsgLeafNode().GetImpl();
 
   const int numProp = NumProp();
   const int numVert = NumPropVert();
   const int numTri = NumTri();
 
+  const bool updateNormals =
+      glm::all(glm::greaterThan(normalIdx, glm::ivec3(2)));
+
   MeshGL out;
   out.numProp = 3 + numProp;
   out.triVerts.resize(3 * numTri);
@@ -187,6 +197,8 @@ MeshGL Manifold::GetMeshGL() const {
                  : triRef[a].originalID < triRef[b].originalID;
     });
   }
+
+  std::vector<glm::mat3> runNormalTransform;
   int lastID = -1;
   for (int tri = 0; tri < numTri; ++tri) {
     const int oldTri = triNew2Old[tri];
@@ -195,11 +207,15 @@ MeshGL Manifold::GetMeshGL() const {
     if (meshID != lastID) {
       out.runIndex.push_back(3 * tri);
       out.originalID.push_back(ref.originalID);
+      const Impl::Relation& m = impl.meshRelation_.meshIDtransform.at(meshID);
+      if (updateNormals) {
+        runNormalTransform.push_back(NormalTransform(m.transform) *
+                                     (m.backSide ? -1.0f : 1.0f));
+      }
       if (impl.meshRelation_.originalID < 0) {
-        const glm::mat4x3& m = impl.meshRelation_.meshIDtransform.at(meshID);
         for (const int col : {0, 1, 2, 3}) {
           for (const int row : {0, 1, 2}) {
-            out.transform.push_back(m[col][row]);
+            out.transform.push_back(m.transform[col][row]);
           }
         }
       }
@@ -226,35 +242,49 @@ MeshGL Manifold::GetMeshGL() const {
   // Duplicate verts with different props
   std::vector<int> vert2idx(impl.NumVert(), -1);
   std::map<std::pair<int, int>, int> vertPropPair;
-  for (int tri = 0; tri < numTri; ++tri) {
-    const glm::ivec3 triProp =
-        impl.meshRelation_.triProperties[triNew2Old[tri]];
-    for (const int i : {0, 1, 2}) {
-      const int prop = triProp[i];
-      const int vert = out.triVerts[3 * tri + i];
+  for (int run = 0; run < out.originalID.size(); ++run) {
+    for (int tri = out.runIndex[run] / 3; tri < out.runIndex[run + 1] / 3;
+         ++tri) {
+      const glm::ivec3 triProp =
+          impl.meshRelation_.triProperties[triNew2Old[tri]];
+      for (const int i : {0, 1, 2}) {
+        const int prop = triProp[i];
+        const int vert = out.triVerts[3 * tri + i];
 
-      const auto it = vertPropPair.find({vert, prop});
-      if (it != vertPropPair.end()) {
-        out.triVerts[3 * tri + i] = it->second;
-        continue;
-      }
-      const int idx = out.vertProperties.size() / out.numProp;
-      vertPropPair.insert({{vert, prop}, idx});
-      out.triVerts[3 * tri + i] = idx;
+        const auto it = vertPropPair.find({vert, prop});
+        if (it != vertPropPair.end()) {
+          out.triVerts[3 * tri + i] = it->second;
+          continue;
+        }
+        const int idx = out.vertProperties.size() / out.numProp;
+        vertPropPair.insert({{vert, prop}, idx});
+        out.triVerts[3 * tri + i] = idx;
 
-      for (int p : {0, 1, 2}) {
-        out.vertProperties.push_back(impl.vertPos_[vert][p]);
-      }
-      for (int p = 0; p < numProp; ++p) {
-        out.vertProperties.push_back(
-            impl.meshRelation_.properties[prop * numProp + p]);
-      }
+        for (int p : {0, 1, 2}) {
+          out.vertProperties.push_back(impl.vertPos_[vert][p]);
+        }
+        for (int p = 0; p < numProp; ++p) {
+          out.vertProperties.push_back(
+              impl.meshRelation_.properties[prop * numProp + p]);
+        }
+        if (updateNormals) {
+          glm::vec3 normal;
+          const int start = out.vertProperties.size() - out.numProp;
+          for (int i : {0, 1, 2}) {
+            normal[i] = out.vertProperties[start + normalIdx[i]];
+          }
+          normal = glm::normalize(runNormalTransform[run] * normal);
+          for (int i : {0, 1, 2}) {
+            out.vertProperties[start + normalIdx[i]] = normal[i];
+          }
+        }
 
-      if (vert2idx[vert] == -1) {
-        vert2idx[vert] = idx;
-      } else {
-        out.mergeFromVert.push_back(idx);
-        out.mergeToVert.push_back(vert2idx[vert]);
+        if (vert2idx[vert] == -1) {
+          vert2idx[vert] = idx;
+        } else {
+          out.mergeFromVert.push_back(idx);
+          out.mergeToVert.push_back(vert2idx[vert]);
+        }
       }
     }
   }

src/manifold/src/shared.h

@@ -45,6 +45,11 @@ __host__ __device__ inline glm::vec3 UVW(int vert,
   return uvw;
 }
 
+__host__ __device__ inline glm::mat3 NormalTransform(
+    const glm::mat4x3& transform) {
+  return glm::inverse(glm::transpose(glm::mat3(transform)));
+}
+
 /**
  * By using the closest axis-aligned projection to the normal instead of a
  * projection along the normal, we avoid introducing any rounding error.

test/meshIO/src/meshIO.cpp

@@ -190,7 +190,7 @@ void ExportMesh(const std::string& filename, const MeshGL& mesh,
     bool validChannels = true;
     for (int i : {0, 1, 2}) {
       int c = options.mat.normalChannels[i];
-      validChannels &= c >= 3 && c < mesh.numProp;
+      validChannels &= c >= 3 && c < (int)mesh.numProp;
     }
     ASSERT(validChannels, userErr,
            "When faceted is false, valid normalChannels must be supplied.");
@@ -201,7 +201,7 @@ void ExportMesh(const std::string& filename, const MeshGL& mesh,
   bool hasColor = false;
   for (int i : {0, 1, 2, 3}) {
     int c = options.mat.colorChannels[i];
-    validChannels &= c < mesh.numProp;
+    validChannels &= c < (int)mesh.numProp;
     hasColor |= c >= 0;
   }
   ASSERT(validChannels, userErr, "Invalid colorChannels.");

test/mesh_test.cpp

@@ -79,6 +79,37 @@ MeshGL TetGL() {
   return tet;
 }
 
+MeshGL CubeGL() {
+  MeshGL cube;
+  cube.numProp = 6;
+  glm::mat3 t[3] = {
+      glm::mat3(1), {1, 0, 0, 0, 0, -1, 0, 1, 0}, {0, 1, 0, -1, 0, 0, 0, 0, 1}};
+  glm::vec3 n(0, 1, 0);
+  glm::mat4x3 p(1, 1, 1,    //
+                1, 1, -1,   //
+                -1, 1, -1,  //
+                -1, 1, 1);
+  for (const int dir : {0, 1, 2}) {
+    for (const float side : {1, -1}) {
+      const int start = cube.vertProperties.size();
+      cube.triVerts.push_back(start);
+      cube.triVerts.push_back(start + 1);
+      cube.triVerts.push_back(start + 2);
+      cube.triVerts.push_back(start);
+      cube.triVerts.push_back(start + 2);
+      cube.triVerts.push_back(start + 3);
+      const glm::vec3 norm = side * t[dir] * n;
+      for (const int i : {0, 1, 2, 3}) {
+        const glm::vec3 v = side * t[dir] * p[i];
+        const glm::vec3 ni = side * t[dir] * n;
+        for (const int j : {0, 1, 2}) cube.vertProperties.push_back(v[j]);
+        for (const int j : {0, 1, 2}) cube.vertProperties.push_back(ni[j]);
+      }
+    }
+  }
+  return Weld(cube);
+}
+
 MeshGL WithIndexColors(const Mesh& in) {
   MeshGL inGL(in);
   inGL.originalID = {Manifold::ReserveIDs(1)};
@@ -119,6 +150,25 @@ MeshGL WithPositionColors(const Manifold& in) {
   return inGL;
 }
 
+MeshGL WithNormals(const Manifold& in) {
+  const Mesh mesh = in.GetMesh();
+  MeshGL out;
+  out.originalID = {Manifold::ReserveIDs(1)};
+  out.numProp = 6;
+  out.vertProperties.resize(out.numProp * mesh.vertPos.size());
+  for (int i = 0; i < mesh.vertPos.size(); ++i) {
+    for (int j : {0, 1, 2}) {
+      out.vertProperties[6 * i + j] = mesh.vertPos[i][j];
+      out.vertProperties[6 * i + 3 + j] = mesh.vertNormal[i][j];
+    }
+  }
+  out.triVerts.resize(3 * mesh.triVerts.size());
+  for (int i = 0; i < mesh.triVerts.size(); ++i) {
+    for (int j : {0, 1, 2}) out.triVerts[3 * i + j] = mesh.triVerts[i][j];
+  }
+  return out;
+}
+
 void Identical(const Mesh& mesh1, const Mesh& mesh2) {
   ASSERT_EQ(mesh1.vertPos.size(), mesh2.vertPos.size());
   for (int i = 0; i < mesh1.vertPos.size(); ++i)
@@ -139,9 +189,12 @@ int NumUnique(const std::vector<T>& in) {
   return unique.size();
 }
 
-void RelatedGL(const Manifold& out, const std::vector<MeshGL>& originals) {
+void RelatedGL(const Manifold& out, const std::vector<MeshGL>& originals,
+               bool checkNormals = false) {
   ASSERT_FALSE(out.IsEmpty());
-  MeshGL output = out.GetMeshGL();
+  const glm::ivec3 normalIdx =
+      checkNormals ? glm::ivec3(3, 4, 5) : glm::ivec3(0);
+  MeshGL output = out.GetMeshGL(normalIdx);
   for (int run = 0; run < output.originalID.size(); ++run) {
     const float* m = output.transform.data() + 12 * run;
     const glm::mat4x3 transform =
@@ -169,45 +222,56 @@ void RelatedGL(const Manifold& out, const std::vector<MeshGL>& originals) {
       inTriangle *= inMesh.numProp;
 
       glm::mat3 inTriPos;
+      glm::mat3 outTriPos;
       for (int j : {0, 1, 2}) {
+        const int vert = output.triVerts[3 * tri + j];
         glm::vec4 pos;
-        for (int k : {0, 1, 2})
+        for (int k : {0, 1, 2}) {
           pos[k] = inMesh.vertProperties[inTriangle[j] + k];
+          outTriPos[j][k] = output.vertProperties[vert * output.numProp + k];
+        }
         pos[3] = 1;
         inTriPos[j] = transform * pos;
       }
-      glm::vec3 normal =
+      glm::vec3 outNormal =
+          glm::cross(outTriPos[1] - outTriPos[0], outTriPos[2] - outTriPos[0]);
+      glm::vec3 inNormal =
           glm::cross(inTriPos[1] - inTriPos[0], inTriPos[2] - inTriPos[0]);
-      const float area = glm::length(normal);
+      const float area = glm::length(inNormal);
       if (area == 0) continue;
-      normal /= area;
+      inNormal /= area;
 
       for (int j : {0, 1, 2}) {
         const int vert = output.triVerts[3 * tri + j];
-        glm::vec3 outPos;
-        for (int k : {0, 1, 2})
-          outPos[k] = output.vertProperties[vert * output.numProp + k];
-
         glm::vec3 edges[3];
-        for (int k : {0, 1, 2}) edges[k] = inTriPos[k] - outPos;
+        for (int k : {0, 1, 2}) edges[k] = inTriPos[k] - outTriPos[j];
         const float volume = glm::dot(edges[0], glm::cross(edges[1], edges[2]));
         ASSERT_LE(volume, area * 100 * out.Precision());
 
-        for (int p = 3; p < output.numProp; ++p) {
-          const float propOut =
-              output.vertProperties[vert * output.numProp + p];
+        if (checkNormals) {
+          glm::vec3 normal;
+          for (int k : {0, 1, 2})
+            normal[k] =
+                output.vertProperties[vert * output.numProp + normalIdx[k]];
+          ASSERT_NEAR(glm::length(normal), 1, 0.0001);
+          ASSERT_GT(glm::dot(normal, outNormal), 0);
+        } else {
+          for (int p = 3; p < output.numProp; ++p) {
+            const float propOut =
+                output.vertProperties[vert * output.numProp + p];
 
-          glm::vec3 inProp = {inMesh.vertProperties[inTriangle[0] + p],
-                              inMesh.vertProperties[inTriangle[1] + p],
-                              inMesh.vertProperties[inTriangle[2] + p]};
-          glm::vec3 edgesP[3];
-          for (int k : {0, 1, 2}) {
-            edgesP[k] = edges[k] + normal * inProp[k] - normal * propOut;
+            glm::vec3 inProp = {inMesh.vertProperties[inTriangle[0] + p],
+                                inMesh.vertProperties[inTriangle[1] + p],
+                                inMesh.vertProperties[inTriangle[2] + p]};
+            glm::vec3 edgesP[3];
+            for (int k : {0, 1, 2}) {
+              edgesP[k] = edges[k] + inNormal * inProp[k] - inNormal * propOut;
+            }
+            const float volumeP =
+                glm::dot(edgesP[0], glm::cross(edgesP[1], edgesP[2]));
+
+            ASSERT_LE(volumeP, area * 100 * out.Precision());
           }
-          const float volumeP =
-              glm::dot(edgesP[0], glm::cross(edgesP[1], edgesP[2]));
-
-          ASSERT_LE(volumeP, area * 100 * out.Precision());
         }
       }
     }
@@ -777,6 +841,28 @@ TEST(Boolean, MeshGLRoundTrip) {
   ASSERT_EQ(outGL.originalID.size(), 2);
 }
 
+TEST(Boolean, Normals) {
+  const MeshGL cubeGL = CubeGL();
+  const Manifold cube(cubeGL);
+  const MeshGL sphereGL = WithNormals(Manifold::Sphere(60));
+  const Manifold sphere(sphereGL);
+
+  Manifold result =
+      cube - (sphere.Rotate(180) -
+              sphere.Scale(glm::vec3(0.5)).Rotate(90).Translate({40, 40, 40}));
+
+#ifdef MANIFOLD_EXPORT
+  ExportOptions opt;
+  opt.faceted = false;
+  opt.mat.roughness = 0;
+  opt.mat.normalChannels = {3, 4, 5};
+  if (options.exportModels)
+    ExportMesh("normals.glb", result.GetMeshGL({3, 4, 5}), opt);
+#endif
+
+  RelatedGL(result, {cubeGL, sphereGL}, true);
+}
+
 TEST(Boolean, Mirrored) {
   Manifold cube = Manifold::Cube(glm::vec3(1)).Scale({1, -1, 1});
   EXPECT_TRUE(cube.IsManifold());