// Copyright 2021 The Manifold Authors.
//
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#pragma once
#include <map>
#include "collider.h"
#include "manifold/common.h"
#include "manifold/manifold.h"
#include "shared.h"
#include "vec.h"
namespace manifold {
/** @ingroup Private */
struct Manifold::Impl {
struct Relation {
int originalID = -1;
mat3x4 transform = la::identity;
bool backSide = false;
};
struct MeshRelationD {
/// The originalID of this Manifold if it is an original; -1 otherwise.
int originalID = -1;
std::map<int, Relation> meshIDtransform;
Vec<TriRef> triRef;
};
struct BaryIndices {
int tri, start4, end4;
};
Box bBox_;
double epsilon_ = -1;
double tolerance_ = -1;
int numProp_ = 0;
Error status_ = Error::NoError;
Vec<vec3> vertPos_;
Vec<Halfedge> halfedge_;
Vec<double> properties_;
// Note that vertNormal_ is not precise due to the use of an approximated acos
// function
Vec<vec3> vertNormal_;
Vec<vec3> faceNormal_;
Vec<vec4> halfedgeTangent_;
MeshRelationD meshRelation_;
Collider collider_;
static std::atomic<uint32_t> meshIDCounter_;
static uint32_t ReserveIDs(uint32_t);
Impl() {}
enum class Shape { Tetrahedron, Cube, Octahedron };
Impl(Shape, const mat3x4 = la::identity);
template <typename Precision, typename I>
Impl(const MeshGLP<Precision, I>& meshGL) {
const uint32_t numVert = meshGL.NumVert();
const uint32_t numTri = meshGL.NumTri();
if (numVert == 0 && numTri == 0) {
MakeEmpty(Error::NoError);
return;
}
if (numVert < 4 || numTri < 4) {
MakeEmpty(Error::NotManifold);
return;
}
if (meshGL.numProp < 3) {
MakeEmpty(Error::MissingPositionProperties);
return;
}
if (meshGL.mergeFromVert.size() != meshGL.mergeToVert.size()) {
MakeEmpty(Error::MergeVectorsDifferentLengths);
return;
}
if (!meshGL.runTransform.empty() &&
12 * meshGL.runOriginalID.size() != meshGL.runTransform.size()) {
MakeEmpty(Error::TransformWrongLength);
return;
}
if (!meshGL.runOriginalID.empty() && !meshGL.runIndex.empty() &&
meshGL.runOriginalID.size() + 1 != meshGL.runIndex.size() &&
meshGL.runOriginalID.size() != meshGL.runIndex.size()) {
MakeEmpty(Error::RunIndexWrongLength);
return;
}
if (!meshGL.faceID.empty() && meshGL.faceID.size() != meshGL.NumTri()) {
MakeEmpty(Error::FaceIDWrongLength);
return;
}
if (!manifold::all_of(meshGL.vertProperties.begin(),
meshGL.vertProperties.end(),
[](Precision x) { return std::isfinite(x); })) {
MakeEmpty(Error::NonFiniteVertex);
return;
}
if (!manifold::all_of(meshGL.runTransform.begin(),
meshGL.runTransform.end(),
[](Precision x) { return std::isfinite(x); })) {
MakeEmpty(Error::InvalidConstruction);
return;
}
if (!manifold::all_of(meshGL.halfedgeTangent.begin(),
meshGL.halfedgeTangent.end(),
[](Precision x) { return std::isfinite(x); })) {
MakeEmpty(Error::InvalidConstruction);
return;
}
std::vector<int> prop2vert;
if (!meshGL.mergeFromVert.empty()) {
prop2vert.resize(numVert);
std::iota(prop2vert.begin(), prop2vert.end(), 0);
for (size_t i = 0; i < meshGL.mergeFromVert.size(); ++i) {
const uint32_t from = meshGL.mergeFromVert[i];
const uint32_t to = meshGL.mergeToVert[i];
if (from >= numVert || to >= numVert) {
MakeEmpty(Error::MergeIndexOutOfBounds);
return;
}
prop2vert[from] = to;
}
}
const auto numProp = meshGL.numProp - 3;
numProp_ = numProp;
properties_.resize_nofill(meshGL.NumVert() * numProp);
tolerance_ = meshGL.tolerance;
// This will have unreferenced duplicate positions that will be removed by
// Impl::RemoveUnreferencedVerts().
vertPos_.resize_nofill(meshGL.NumVert());
for (size_t i = 0; i < meshGL.NumVert(); ++i) {
for (const int j : {0, 1, 2})
vertPos_[i][j] = meshGL.vertProperties[meshGL.numProp * i + j];
for (size_t j = 0; j < numProp; ++j)
properties_[i * numProp + j] =
meshGL.vertProperties[meshGL.numProp * i + 3 + j];
}
halfedgeTangent_.resize_nofill(meshGL.halfedgeTangent.size() / 4);
for (size_t i = 0; i < halfedgeTangent_.size(); ++i) {
for (const int j : {0, 1, 2, 3})
halfedgeTangent_[i][j] = meshGL.halfedgeTangent[4 * i + j];
}
Vec<TriRef> triRef;
triRef.resize_nofill(meshGL.NumTri());
auto runIndex = meshGL.runIndex;
const auto runEnd = meshGL.triVerts.size();
if (runIndex.empty()) {
runIndex = {0, static_cast<I>(runEnd)};
} else if (runIndex.size() == meshGL.runOriginalID.size()) {
runIndex.push_back(runEnd);
} else if (runIndex.size() == 1) {
runIndex.push_back(runEnd);
}
const auto startID =
Impl::ReserveIDs(std::max(1_uz, meshGL.runOriginalID.size()));
auto runOriginalID = meshGL.runOriginalID;
if (runOriginalID.empty()) {
runOriginalID.push_back(startID);
}
for (size_t i = 0; i < runOriginalID.size(); ++i) {
const int meshID = startID + i;
const int originalID = runOriginalID[i];
for (size_t tri = runIndex[i] / 3; tri < runIndex[i + 1] / 3; ++tri) {
TriRef& ref = triRef[tri];
ref.meshID = meshID;
ref.originalID = originalID;
ref.faceID = meshGL.faceID.empty() ? -1 : meshGL.faceID[tri];
ref.coplanarID = tri;
}
if (meshGL.runTransform.empty()) {
meshRelation_.meshIDtransform[meshID] = {originalID};
} else {
const Precision* m = meshGL.runTransform.data() + 12 * i;
meshRelation_.meshIDtransform[meshID] = {originalID,
{{m[0], m[1], m[2]},
{m[3], m[4], m[5]},
{m[6], m[7], m[8]},
{m[9], m[10], m[11]}}};
}
}
Vec<ivec3> triProp;
triProp.reserve(numTri);
Vec<ivec3> triVert;
const bool needsPropMap = numProp > 0 && !prop2vert.empty();
if (needsPropMap) triVert.reserve(numTri);
if (triRef.size() > 0) meshRelation_.triRef.reserve(numTri);
for (size_t i = 0; i < numTri; ++i) {
ivec3 triP, triV;
for (const size_t j : {0, 1, 2}) {
uint32_t vert = (uint32_t)meshGL.triVerts[3 * i + j];
if (vert >= numVert) {
MakeEmpty(Error::VertexOutOfBounds);
return;
}
triP[j] = vert;
triV[j] = prop2vert.empty() ? vert : prop2vert[vert];
}
if (triV[0] != triV[1] && triV[1] != triV[2] && triV[2] != triV[0]) {
if (needsPropMap) {
triProp.push_back(triP);
triVert.push_back(triV);
} else {
triProp.push_back(triV);
}
if (triRef.size() > 0) {
meshRelation_.triRef.push_back(triRef[i]);
}
}
}
CreateHalfedges(triProp, triVert);
if (!IsManifold()) {
MakeEmpty(Error::NotManifold);
return;
}
CalculateBBox();
SetEpsilon(-1, std::is_same<Precision, float>::value);
// we need to split pinched verts before calculating vertex normals, because
// the algorithm doesn't work with pinched verts
CleanupTopology();
CalculateNormals();
DedupePropVerts();
MarkCoplanar();
RemoveDegenerates();
RemoveUnreferencedVerts();
Finish();
if (!IsFinite()) {
MakeEmpty(Error::NonFiniteVertex);
return;
}
// A Manifold created from an input mesh is never an original - the input is
// the original.
meshRelation_.originalID = -1;
}
template <typename F>
inline void ForVert(int halfedge, F func) {
int current = halfedge;
do {
current = NextHalfedge(halfedge_[current].pairedHalfedge);
func(current);
} while (current != halfedge);
}
template <typename T>
void ForVert(
int halfedge, std::function<T(int halfedge)> transform,
std::function<void(int halfedge, const T& here, T& next)> binaryOp) {
T here = transform(halfedge);
int current = halfedge;
do {
const int nextHalfedge = NextHalfedge(halfedge_[current].pairedHalfedge);
T next = transform(nextHalfedge);
binaryOp(current, here, next);
here = next;
current = nextHalfedge;
} while (current != halfedge);
}
void MarkCoplanar();
void DedupePropVerts();
void RemoveUnreferencedVerts();
void InitializeOriginal(bool keepFaceID = false);
void CreateHalfedges(const Vec<ivec3>& triProp,
const Vec<ivec3>& triVert = {});
void CalculateNormals();
void IncrementMeshIDs();
void Update();
void MakeEmpty(Error status);
void Warp(std::function<void(vec3&)> warpFunc);
void WarpBatch(std::function<void(VecView<vec3>)> warpFunc);
Impl Transform(const mat3x4& transform) const;
bool IsEmpty() const { return NumTri() == 0; }
size_t NumVert() const { return vertPos_.size(); }
size_t NumEdge() const { return halfedge_.size() / 2; }
size_t NumTri() const { return halfedge_.size() / 3; }
size_t NumProp() const { return numProp_; }
size_t NumPropVert() const {
return NumProp() == 0 ? NumVert() : properties_.size() / NumProp();
}
// properties.cpp
enum class Property { Volume, SurfaceArea };
double GetProperty(Property prop) const;
void CalculateCurvature(int gaussianIdx, int meanIdx);
void CalculateBBox();
bool IsFinite() const;
bool IsIndexInBounds(VecView<const ivec3> triVerts) const;
void SetEpsilon(double minEpsilon = -1, bool useSingle = false);
bool IsManifold() const;
bool Is2Manifold() const;
bool IsSelfIntersecting() const;
bool MatchesTriNormals() const;
int NumDegenerateTris() const;
double MinGap(const Impl& other, double searchLength) const;
// sort.cpp
void Finish();
void SortVerts();
void ReindexVerts(const Vec<int>& vertNew2Old, size_t numOldVert);
void CompactProps();
void GetFaceBoxMorton(Vec<Box>& faceBox, Vec<uint32_t>& faceMorton) const;
void SortFaces(Vec<Box>& faceBox, Vec<uint32_t>& faceMorton);
void GatherFaces(const Vec<int>& faceNew2Old);
void GatherFaces(const Impl& old, const Vec<int>& faceNew2Old);
// face_op.cpp
void Face2Tri(const Vec<int>& faceEdge, const Vec<TriRef>& halfedgeRef,
bool allowConvex = false);
Polygons Slice(double height) const;
Polygons Project() const;
// edge_op.cpp
void CleanupTopology();
void SimplifyTopology(int firstNewVert = 0);
void RemoveDegenerates(int firstNewVert = 0);
void CollapseShortEdges(int firstNewVert = 0);
void CollapseColinearEdges(int firstNewVert = 0);
void SwapDegenerates(int firstNewVert = 0);
void DedupeEdge(int edge);
bool CollapseEdge(int edge, std::vector<int>& edges);
void RecursiveEdgeSwap(int edge, int& tag, std::vector<int>& visited,
std::vector<int>& edgeSwapStack,
std::vector<int>& edges);
void RemoveIfFolded(int edge);
void PairUp(int edge0, int edge1);
void UpdateVert(int vert, int startEdge, int endEdge);
void FormLoop(int current, int end);
void CollapseTri(const ivec3& triEdge);
void SplitPinchedVerts();
void DedupeEdges();
// subdivision.cpp
int GetNeighbor(int tri) const;
ivec4 GetHalfedges(int tri) const;
BaryIndices GetIndices(int halfedge) const;
void FillRetainedVerts(Vec<Barycentric>& vertBary) const;
Vec<Barycentric> Subdivide(std::function<int(vec3, vec4, vec4)>,
bool = false);
// smoothing.cpp
bool IsInsideQuad(int halfedge) const;
bool IsMarkedInsideQuad(int halfedge) const;
vec3 GetNormal(int halfedge, int normalIdx) const;
vec4 TangentFromNormal(const vec3& normal, int halfedge) const;
std::vector<Smoothness> UpdateSharpenedEdges(
const std::vector<Smoothness>&) const;
Vec<bool> FlatFaces() const;
Vec<int> VertFlatFace(const Vec<bool>&) const;
Vec<int> VertHalfedge() const;
std::vector<Smoothness> SharpenEdges(double minSharpAngle,
double minSmoothness) const;
void SharpenTangent(int halfedge, double smoothness);
void SetNormals(int normalIdx, double minSharpAngle);
void LinearizeFlatTangents();
void DistributeTangents(const Vec<bool>& fixedHalfedges);
void CreateTangents(int normalIdx);
void CreateTangents(std::vector<Smoothness>);
void Refine(std::function<int(vec3, vec4, vec4)>, bool = false);
// quickhull.cpp
void Hull(VecView<vec3> vertPos);
};
extern std::mutex dump_lock;
std::ostream& operator<<(std::ostream& stream, const Manifold::Impl& impl);
} // namespace manifold