3D-printing/OrcaSlicer
1
0
Fork 0
mirror of https://github.com/SoftFever/OrcaSlicer.git synced 2025-07-08 07:27:41 -06:00
Code Issues Projects Releases Packages Wiki Activity Actions 3

ENH: revert boolean tool to mcut

Browse source 
Change-Id: I1aca763869e107a996519cb74e025043407005f2
(cherry picked from commit f7865828cf4b7b3ab8987bf0fc3d45efad2f08fe)
This commit is contained in:
Arthur 2023-05-31 15:05:28 +08:00 • committed by Lane.Wei
parent b55a5b7556
commit c7ed4e7e14
38 changed files with 39957 additions and 5 deletions
Download patch file Download diff file Expand all files Collapse all files

1
CMakeLists.txt
View file

@ -646,7 +646,6 @@ function(bambustudio_copy_dlls target config postfix output_dlls)
${_out_dir}/TKXSBase.dll
${_out_dir}/freetype.dll
PARENT_SCOPE
)

1
src/CMakeLists.txt
View file

@ -16,6 +16,7 @@ add_subdirectory(Shiny)
add_subdirectory(semver)
add_subdirectory(libigl)
add_subdirectory(hints)
add_subdirectory(mcut)
# Adding libnest2d project for bin packing...
add_subdirectory(libnest2d)

3
src/libslic3r/CMakeLists.txt
View file

@ -425,7 +425,7 @@ if (_opts)
target_compile_options(libslic3r_cgal PRIVATE "${_opts_bad}")
endif()
target_link_libraries(libslic3r_cgal PRIVATE ${_cgal_tgt} libigl)
target_link_libraries(libslic3r_cgal PRIVATE ${_cgal_tgt} libigl mcut)
if (MSVC AND "${CMAKE_SIZEOF_VOID_P}" STREQUAL "4") # 32 bit MSVC workaround
target_compile_definitions(libslic3r_cgal PRIVATE CGAL_DO_NOT_USE_MPZF)
@ -492,6 +492,7 @@ target_link_libraries(libslic3r
ZLIB::ZLIB
${OCCT_LIBS}
Clipper2
mcut
)
if(NOT WIN32)

174
src/libslic3r/CSGMesh/PerformCSGMeshBooleans.hpp
View file

@ -38,6 +38,32 @@ MeshBoolean::cgal::CGALMeshPtr get_cgalmesh(const CSGPartT &csgpart)
return ret;
}
// This method can be overriden when a specific CSGPart type supports caching
// of the voxel grid
template<class CSGPartT>
MeshBoolean::mcut::McutMeshPtr get_mcutmesh(const CSGPartT& csgpart)
{
const indexed_triangle_set* its = csg::get_mesh(csgpart);
indexed_triangle_set dummy;
if (!its)
its = &dummy;
MeshBoolean::mcut::McutMeshPtr ret;
indexed_triangle_set m = *its;
its_transform(m, get_transform(csgpart), true);
try {
ret = MeshBoolean::mcut::triangle_mesh_to_mcut(m);
}
catch (...) {
// errors are ignored, simply return null
ret = nullptr;
}
return ret;
}
namespace detail_cgal {
@ -83,6 +109,50 @@ std::vector<CGALMeshPtr> get_cgalptrs(Ex policy, const Range<It> &csgrange)
} // namespace detail
namespace detail_mcut {
using MeshBoolean::mcut::McutMeshPtr;
inline void perform_csg(CSGType op, McutMeshPtr& dst, McutMeshPtr& src)
{
if (!dst && op == CSGType::Union && src) {
dst = std::move(src);
return;
}
if (!dst || !src)
return;
switch (op) {
case CSGType::Union:
MeshBoolean::mcut::do_boolean(*dst, *src,"UNION");
break;
case CSGType::Difference:
MeshBoolean::mcut::do_boolean(*dst, *src,"A_NOT_B");
break;
case CSGType::Intersection:
MeshBoolean::mcut::do_boolean(*dst, *src,"INTERSECTION");
break;
}
}
template<class Ex, class It>
std::vector<McutMeshPtr> get_mcutptrs(Ex policy, const Range<It>& csgrange)
{
std::vector<McutMeshPtr> ret(csgrange.size());
execution::for_each(policy, size_t(0), csgrange.size(),
[&csgrange, &ret](size_t i) {
auto it = csgrange.begin();
std::advance(it, i);
auto& csgpart = *it;
ret[i] = get_mcutmesh(csgpart);
});
return ret;
}
} // namespace mcut_detail
// Process the sequence of CSG parts with CGAL.
template<class It>
void perform_csgmesh_booleans_cgal(MeshBoolean::cgal::CGALMeshPtr &cgalm,
@ -133,6 +203,58 @@ void perform_csgmesh_booleans_cgal(MeshBoolean::cgal::CGALMeshPtr &cgalm,
cgalm = std::move(opstack.top().cgalptr);
}
// Process the sequence of CSG parts with mcut.
template<class It>
void perform_csgmesh_booleans_mcut(MeshBoolean::mcut::McutMeshPtr& mcutm,
const Range<It>& csgrange)
{
using MeshBoolean::mcut::McutMesh;
using MeshBoolean::mcut::McutMeshPtr;
using namespace detail_mcut;
struct Frame {
CSGType op; McutMeshPtr mcutptr;
explicit Frame(CSGType csgop = CSGType::Union)
: op{ csgop }
, mcutptr{ MeshBoolean::mcut::triangle_mesh_to_mcut(indexed_triangle_set{}) }
{}
};
std::stack opstack{ std::vector<Frame>{} };
opstack.push(Frame{});
std::vector<McutMeshPtr> McutMeshes = get_mcutptrs(ex_tbb, csgrange);
size_t csgidx = 0;
for (auto& csgpart : csgrange) {
auto op = get_operation(csgpart);
McutMeshPtr& mcutptr = McutMeshes[csgidx++];
if (get_stack_operation(csgpart) == CSGStackOp::Push) {
opstack.push(Frame{ op });
op = CSGType::Union;
}
Frame* top = &opstack.top();
perform_csg(get_operation(csgpart), top->mcutptr, mcutptr);
if (get_stack_operation(csgpart) == CSGStackOp::Pop) {
McutMeshPtr src = std::move(top->mcutptr);
auto popop = opstack.top().op;
opstack.pop();
McutMeshPtr& dst = opstack.top().mcutptr;
perform_csg(popop, dst, src);
}
}
mcutm = std::move(opstack.top().mcutptr);
}
template<class It, class Visitor>
It check_csgmesh_booleans(const Range<It> &csgrange, Visitor &&vfn)
{
@ -184,18 +306,66 @@ It check_csgmesh_booleans(const Range<It> &csgrange, Visitor &&vfn)
}
template<class It>
It check_csgmesh_booleans(const Range<It> &csgrange)
It check_csgmesh_booleans(const Range<It> &csgrange, bool use_mcut=false)
{
if(!use_mcut)
return check_csgmesh_booleans(csgrange, [](auto &) {});
else {
using namespace detail_mcut;
std::vector<McutMeshPtr> McutMeshes(csgrange.size());
auto check_part = [&csgrange, &McutMeshes](size_t i) {
auto it = csgrange.begin();
std::advance(it, i);
auto& csgpart = *it;
auto m = get_mcutmesh(csgpart);
// mesh can be nullptr if this is a stack push or pull
if (!get_mesh(csgpart) && get_stack_operation(csgpart) != CSGStackOp::Continue) {
McutMeshes[i] = MeshBoolean::mcut::triangle_mesh_to_mcut(indexed_triangle_set{});
return;
}
try {
if (!m || MeshBoolean::mcut::empty(*m))
return;
}
catch (...) { return; }
McutMeshes[i] = std::move(m);
};
execution::for_each(ex_tbb, size_t(0), csgrange.size(), check_part);
It ret = csgrange.end();
for (size_t i = 0; i < csgrange.size(); ++i) {
if (!McutMeshes[i]) {
auto it = csgrange.begin();
std::advance(it, i);
if (ret == csgrange.end())
ret = it;
}
}
return ret;
}
}
template<class It>
MeshBoolean::cgal::CGALMeshPtr perform_csgmesh_booleans(const Range<It> &csgparts)
{
auto ret = MeshBoolean::cgal::triangle_mesh_to_cgal(indexed_triangle_set{});
if (ret) {
if (ret)
perform_csgmesh_booleans_cgal(ret, csgparts);
return ret;
}
template<class It>
MeshBoolean::mcut::McutMeshPtr perform_csgmesh_booleans_mcut(const Range<It>& csgparts)
{
auto ret = MeshBoolean::mcut::triangle_mesh_to_mcut(indexed_triangle_set{});
if (ret)
perform_csgmesh_booleans_mcut(ret, csgparts);
return ret;
}

315
src/libslic3r/MeshBoolean.cpp
View file

@ -23,6 +23,8 @@
#include <CGAL/property_map.h>
#include <CGAL/boost/graph/copy_face_graph.h>
#include <CGAL/boost/graph/Face_filtered_graph.h>
// BBS: for boolean using mcut
#include "mcut/include/mcut/mcut.h"
namespace Slic3r {
namespace MeshBoolean {
@ -470,5 +472,318 @@ CGALMeshPtr clone(const CGALMesh &m)
}
} // namespace cgal
namespace mcut {
/* BBS: MusangKing
* mcut mesh array format for Boolean Opts calculation
*/
struct McutMesh
{
// variables for mesh data in a format suited for mcut
std::vector<uint32_t> faceSizesArray;
std::vector<uint32_t> faceIndicesArray;
std::vector<double> vertexCoordsArray;
};
void McutMeshDeleter::operator()(McutMesh *ptr) { delete ptr; }
bool empty(const McutMesh &mesh) { return mesh.vertexCoordsArray.empty() || mesh.faceIndicesArray.empty(); }
void triangle_mesh_to_mcut(const TriangleMesh &src_mesh, McutMesh &srcMesh, const Transform3d &src_nm = Transform3d::Identity())
{
// vertices precision convention and copy
srcMesh.vertexCoordsArray.reserve(src_mesh.its.vertices.size() * 3);
for (int i = 0; i < src_mesh.its.vertices.size(); ++i) {
const Vec3d v = src_nm * src_mesh.its.vertices[i].cast<double>();
srcMesh.vertexCoordsArray.push_back(v[0]);
srcMesh.vertexCoordsArray.push_back(v[1]);
srcMesh.vertexCoordsArray.push_back(v[2]);
}
// faces copy
srcMesh.faceIndicesArray.reserve(src_mesh.its.indices.size() * 3);
srcMesh.faceSizesArray.reserve(src_mesh.its.indices.size());
for (int i = 0; i < src_mesh.its.indices.size(); ++i) {
const int &f0 = src_mesh.its.indices[i][0];
const int &f1 = src_mesh.its.indices[i][1];
const int &f2 = src_mesh.its.indices[i][2];
srcMesh.faceIndicesArray.push_back(f0);
srcMesh.faceIndicesArray.push_back(f1);
srcMesh.faceIndicesArray.push_back(f2);
srcMesh.faceSizesArray.push_back((uint32_t) 3);
}
}
McutMeshPtr triangle_mesh_to_mcut(const indexed_triangle_set &M)
{
std::unique_ptr<McutMesh, McutMeshDeleter> out(new McutMesh{});
TriangleMesh trimesh(M);
triangle_mesh_to_mcut(trimesh, *out.get());
return out;
}
TriangleMesh mcut_to_triangle_mesh(const McutMesh &mcutmesh)
{
uint32_t ccVertexCount = mcutmesh.vertexCoordsArray.size() / 3;
auto &ccVertices = mcutmesh.vertexCoordsArray;
auto &ccFaceIndices = mcutmesh.faceIndicesArray;
auto &faceSizes = mcutmesh.faceSizesArray;
uint32_t ccFaceCount = faceSizes.size();
// rearrange vertices/faces and save into result mesh
std::vector<Vec3f> vertices(ccVertexCount);
for (uint32_t i = 0; i < ccVertexCount; i++) {
vertices[i][0] = (float) ccVertices[(uint64_t) i * 3 + 0];
vertices[i][1] = (float) ccVertices[(uint64_t) i * 3 + 1];
vertices[i][2] = (float) ccVertices[(uint64_t) i * 3 + 2];
}
// output faces
int faceVertexOffsetBase = 0;
// for each face in CC
std::vector<Vec3i> faces(ccFaceCount);
for (uint32_t f = 0; f < ccFaceCount; ++f) {
int faceSize = faceSizes.at(f);
// for each vertex in face
for (int v = 0; v < faceSize; v++) { faces[f][v] = ccFaceIndices[(uint64_t) faceVertexOffsetBase + v]; }
faceVertexOffsetBase += faceSize;
}
TriangleMesh out(vertices, faces);
return out;
}
void do_boolean(McutMesh &srcMesh, const McutMesh &cutMesh, const std::string &boolean_opts)
{
// create context
McContext context = MC_NULL_HANDLE;
McResult err = mcCreateContext(&context, static_cast<McFlags>(MC_DEBUG));
// We can either let MCUT compute all possible meshes (including patches etc.), or we can
// constrain the library to compute exactly the boolean op mesh we want. This 'constrained' case
// is done with the following flags.
// NOTE#1: you can extend these flags by bitwise ORing with additional flags (see `McDispatchFlags' in mcut.h)
// NOTE#2: below order of columns MATTERS
const std::map<std::string, McFlags> booleanOpts = {
{"A_NOT_B", MC_DISPATCH_FILTER_FRAGMENT_SEALING_INSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_ABOVE},
{"B_NOT_A", MC_DISPATCH_FILTER_FRAGMENT_SEALING_OUTSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_BELOW},
{"UNION", MC_DISPATCH_FILTER_FRAGMENT_SEALING_OUTSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_ABOVE},
{"INTERSECTION", MC_DISPATCH_FILTER_FRAGMENT_SEALING_INSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_BELOW},
};
std::map<std::string, McFlags>::const_iterator it = booleanOpts.find(boolean_opts);
McFlags boolOpFlags = it->second;
if (srcMesh.vertexCoordsArray.empty() && (boolean_opts == "UNION" || boolean_opts == "B_NOT_A")) {
srcMesh = cutMesh;
return;
}
err = mcDispatch(context,
MC_DISPATCH_VERTEX_ARRAY_DOUBLE | // vertices are in array of doubles
MC_DISPATCH_ENFORCE_GENERAL_POSITION | // perturb if necessary
boolOpFlags, // filter flags which specify the type of output we want
// source mesh
reinterpret_cast<const void *>(srcMesh.vertexCoordsArray.data()), reinterpret_cast<const uint32_t *>(srcMesh.faceIndicesArray.data()),
srcMesh.faceSizesArray.data(), static_cast<uint32_t>(srcMesh.vertexCoordsArray.size() / 3), static_cast<uint32_t>(srcMesh.faceSizesArray.size()),
// cut mesh
reinterpret_cast<const void *>(cutMesh.vertexCoordsArray.data()), cutMesh.faceIndicesArray.data(), cutMesh.faceSizesArray.data(),
static_cast<uint32_t>(cutMesh.vertexCoordsArray.size() / 3), static_cast<uint32_t>(cutMesh.faceSizesArray.size()));
// query the number of available connected component
uint32_t numConnComps;
err = mcGetConnectedComponents(context, MC_CONNECTED_COMPONENT_TYPE_FRAGMENT, 0, NULL, &numConnComps);
std::vector<McConnectedComponent> connectedComponents(numConnComps, MC_NULL_HANDLE);
err = mcGetConnectedComponents(context, MC_CONNECTED_COMPONENT_TYPE_FRAGMENT, (uint32_t) connectedComponents.size(), connectedComponents.data(), NULL);
McutMesh outMesh;
int N_vertices = 0;
// traversal of all connected components
for (int n = 0; n < numConnComps; ++n) {
// query the data of each connected component from MCUT
McConnectedComponent connComp = connectedComponents[n];
// query the vertices
McSize numBytes = 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_VERTEX_DOUBLE, 0, NULL, &numBytes);
uint32_t ccVertexCount = (uint32_t) (numBytes / (sizeof(double) * 3));
std::vector<double> ccVertices((uint64_t) ccVertexCount * 3u, 0);
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_VERTEX_DOUBLE, numBytes, (void *) ccVertices.data(), NULL);
// query the faces
numBytes = 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_FACE_TRIANGULATION, 0, NULL, &numBytes);
std::vector<uint32_t> ccFaceIndices(numBytes / sizeof(uint32_t), 0);
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_FACE_TRIANGULATION, numBytes, ccFaceIndices.data(), NULL);
std::vector<uint32_t> faceSizes(ccFaceIndices.size() / 3, 3);
const uint32_t ccFaceCount = static_cast<uint32_t>(faceSizes.size());
// Here we show, how to know when connected components, pertain particular boolean operations.
McPatchLocation patchLocation = (McPatchLocation) 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_PATCH_LOCATION, sizeof(McPatchLocation), &patchLocation, NULL);
McFragmentLocation fragmentLocation = (McFragmentLocation) 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_FRAGMENT_LOCATION, sizeof(McFragmentLocation), &fragmentLocation, NULL);
outMesh.vertexCoordsArray.insert(outMesh.vertexCoordsArray.end(), ccVertices.begin(), ccVertices.end());
// add offset to face index
for (size_t i = 0; i < ccFaceIndices.size(); i++) {
ccFaceIndices[i] += N_vertices;
}
int faceVertexOffsetBase = 0;
// for each face in CC
std::vector<Vec3i> faces(ccFaceCount);
for (uint32_t f = 0; f < ccFaceCount; ++f) {
bool reverseWindingOrder = (fragmentLocation == MC_FRAGMENT_LOCATION_BELOW) && (patchLocation == MC_PATCH_LOCATION_OUTSIDE);
int faceSize = faceSizes.at(f);
if (reverseWindingOrder) {
std::vector<uint32_t> faceIndex(faceSize);
// for each vertex in face
for (int v = faceSize - 1; v >= 0; v--) { faceIndex[v] = ccFaceIndices[(uint64_t) faceVertexOffsetBase + v]; }
std::copy(faceIndex.begin(), faceIndex.end(), ccFaceIndices.begin() + faceVertexOffsetBase);
}
faceVertexOffsetBase += faceSize;
}
outMesh.faceIndicesArray.insert(outMesh.faceIndicesArray.end(), ccFaceIndices.begin(), ccFaceIndices.end());
outMesh.faceSizesArray.insert(outMesh.faceSizesArray.end(), faceSizes.begin(), faceSizes.end());
N_vertices += ccVertexCount;
}
// free connected component data
err = mcReleaseConnectedComponents(context, (uint32_t) connectedComponents.size(), connectedComponents.data());
// destroy context
err = mcReleaseContext(context);
srcMesh = outMesh;
}
/* BBS: Musang King
* mcut for Mesh Boolean which provides C-style syntax API
*/
std::vector<TriangleMesh> make_boolean(const McutMesh &srcMesh, const McutMesh &cutMesh, const std::string &boolean_opts)
{
// create context
McContext context = MC_NULL_HANDLE;
McResult err = mcCreateContext(&context, static_cast<McFlags>(MC_DEBUG));
// We can either let MCUT compute all possible meshes (including patches etc.), or we can
// constrain the library to compute exactly the boolean op mesh we want. This 'constrained' case
// is done with the following flags.
// NOTE#1: you can extend these flags by bitwise ORing with additional flags (see `McDispatchFlags' in mcut.h)
// NOTE#2: below order of columns MATTERS
const std::map<std::string, McFlags> booleanOpts = {
{"A_NOT_B", MC_DISPATCH_FILTER_FRAGMENT_SEALING_INSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_ABOVE},
{"B_NOT_A", MC_DISPATCH_FILTER_FRAGMENT_SEALING_OUTSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_BELOW},
{"UNION", MC_DISPATCH_FILTER_FRAGMENT_SEALING_OUTSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_ABOVE},
{"INTERSECTION", MC_DISPATCH_FILTER_FRAGMENT_SEALING_INSIDE | MC_DISPATCH_FILTER_FRAGMENT_LOCATION_BELOW},
};
std::map<std::string, McFlags>::const_iterator it = booleanOpts.find(boolean_opts);
McFlags boolOpFlags = it->second;
err = mcDispatch(context,
MC_DISPATCH_VERTEX_ARRAY_DOUBLE | // vertices are in array of doubles
MC_DISPATCH_ENFORCE_GENERAL_POSITION | // perturb if necessary
boolOpFlags, // filter flags which specify the type of output we want
// source mesh
reinterpret_cast<const void *>(srcMesh.vertexCoordsArray.data()), reinterpret_cast<const uint32_t *>(srcMesh.faceIndicesArray.data()),
srcMesh.faceSizesArray.data(), static_cast<uint32_t>(srcMesh.vertexCoordsArray.size() / 3), static_cast<uint32_t>(srcMesh.faceSizesArray.size()),
// cut mesh
reinterpret_cast<const void *>(cutMesh.vertexCoordsArray.data()), cutMesh.faceIndicesArray.data(), cutMesh.faceSizesArray.data(),
static_cast<uint32_t>(cutMesh.vertexCoordsArray.size() / 3), static_cast<uint32_t>(cutMesh.faceSizesArray.size()));
// query the number of available connected component
uint32_t numConnComps;
err = mcGetConnectedComponents(context, MC_CONNECTED_COMPONENT_TYPE_FRAGMENT, 0, NULL, &numConnComps);
std::vector<McConnectedComponent> connectedComponents(numConnComps, MC_NULL_HANDLE);
err = mcGetConnectedComponents(context, MC_CONNECTED_COMPONENT_TYPE_FRAGMENT, (uint32_t) connectedComponents.size(), connectedComponents.data(), NULL);
std::vector<TriangleMesh> outs;
// traversal of all connected components
for (int n = 0; n < numConnComps; ++n) {
// query the data of each connected component from MCUT
McConnectedComponent connComp = connectedComponents[n];
// query the vertices
McSize numBytes = 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_VERTEX_DOUBLE, 0, NULL, &numBytes);
uint32_t ccVertexCount = (uint32_t) (numBytes / (sizeof(double) * 3));
std::vector<double> ccVertices((uint64_t) ccVertexCount * 3u, 0);
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_VERTEX_DOUBLE, numBytes, (void *) ccVertices.data(), NULL);
// query the faces
numBytes = 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_FACE_TRIANGULATION, 0, NULL, &numBytes);
std::vector<uint32_t> ccFaceIndices(numBytes / sizeof(uint32_t), 0);
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_FACE_TRIANGULATION, numBytes, ccFaceIndices.data(), NULL);
std::vector<uint32_t> faceSizes(ccFaceIndices.size() / 3, 3);
const uint32_t ccFaceCount = static_cast<uint32_t>(faceSizes.size());
// Here we show, how to know when connected components, pertain particular boolean operations.
McPatchLocation patchLocation = (McPatchLocation) 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_PATCH_LOCATION, sizeof(McPatchLocation), &patchLocation, NULL);
McFragmentLocation fragmentLocation = (McFragmentLocation) 0;
err = mcGetConnectedComponentData(context, connComp, MC_CONNECTED_COMPONENT_DATA_FRAGMENT_LOCATION, sizeof(McFragmentLocation), &fragmentLocation, NULL);
// rearrange vertices/faces and save into result mesh
std::vector<Vec3f> vertices(ccVertexCount);
for (uint32_t i = 0; i < ccVertexCount; ++i) {
vertices[i][0] = (float) ccVertices[(uint64_t) i * 3 + 0];
vertices[i][1] = (float) ccVertices[(uint64_t) i * 3 + 1];
vertices[i][2] = (float) ccVertices[(uint64_t) i * 3 + 2];
}
// output faces
int faceVertexOffsetBase = 0;
// for each face in CC
std::vector<Vec3i> faces(ccFaceCount);
for (uint32_t f = 0; f < ccFaceCount; ++f) {
bool reverseWindingOrder = (fragmentLocation == MC_FRAGMENT_LOCATION_BELOW) && (patchLocation == MC_PATCH_LOCATION_OUTSIDE);
int faceSize = faceSizes.at(f);
// for each vertex in face
for (int v = (reverseWindingOrder ? (faceSize - 1) : 0); (reverseWindingOrder ? (v >= 0) : (v < faceSize)); v += (reverseWindingOrder ? -1 : 1)) {
faces[f][v] = ccFaceIndices[(uint64_t) faceVertexOffsetBase + v];
}
faceVertexOffsetBase += faceSize;
}
TriangleMesh out(vertices, faces);
outs.emplace_back(out);
}
// free connected component data
err = mcReleaseConnectedComponents(context, (uint32_t) connectedComponents.size(), connectedComponents.data());
// destroy context
err = mcReleaseContext(context);
return outs;
}
void make_boolean(const TriangleMesh &src_mesh, const TriangleMesh &cut_mesh, std::vector<TriangleMesh> &dst_mesh, const std::string &boolean_opts)
{
McutMesh srcMesh, cutMesh;
triangle_mesh_to_mcut(src_mesh, srcMesh);
triangle_mesh_to_mcut(cut_mesh, cutMesh);
dst_mesh = make_boolean(srcMesh, cutMesh, boolean_opts);
}
} // namespace mcut
} // namespace MeshBoolean
} // namespace Slic3r

21
src/libslic3r/MeshBoolean.hpp
View file

@ -72,6 +72,27 @@ bool does_bound_a_volume(const CGALMesh &mesh);
bool empty(const CGALMesh &mesh);
}
namespace mcut {
struct McutMesh;
struct McutMeshDeleter
{
void operator()(McutMesh *ptr);
};
using McutMeshPtr = std::unique_ptr<McutMesh, McutMeshDeleter>;
bool empty(const McutMesh &mesh);
McutMeshPtr triangle_mesh_to_mcut(const indexed_triangle_set &M);
TriangleMesh mcut_to_triangle_mesh(const McutMesh &mcutmesh);
// do boolean and save result to srcMesh
void do_boolean(McutMesh &srcMesh, const McutMesh &cutMesh, const std::string &boolean_opts);
std::vector<TriangleMesh> make_boolean(const McutMesh &srcMesh, const McutMesh &cutMesh, const std::string &boolean_opts);
// do boolean and convert result to TriangleMesh
void make_boolean(const TriangleMesh &src_mesh, const TriangleMesh &cut_mesh, std::vector<TriangleMesh> &dst_mesh, const std::string &boolean_opts);
} // namespace mcut
} // namespace MeshBoolean
} // namespace Slic3r
#endif // libslic3r_MeshBoolean_hpp_

22
src/libslic3r/Model.cpp
View file

@ -1053,6 +1053,28 @@ void ModelObject::assign_new_unique_ids_recursive()
// BBS: production extension
int ModelObject::get_backup_id() const { return m_model ? get_model()->get_object_backup_id(*this) : -1; }
// BBS: Boolean Operations impl. - MusangKing
bool ModelObject::make_boolean(ModelObject *cut_object, const std::string &boolean_opts)
{
// merge meshes into single volume instead of multi-parts object
if (this->volumes.size() != 1) {
// we can't merge meshes if there's not just one volume
return false;
}
std::vector<TriangleMesh> new_meshes;
const TriangleMesh &cut_mesh = cut_object->mesh();
MeshBoolean::mcut::make_boolean(this->mesh(), cut_mesh, new_meshes, boolean_opts);
this->clear_volumes();
int i = 1;
for (TriangleMesh &mesh : new_meshes) {
ModelVolume *vol = this->add_volume(mesh);
vol->name = this->name + "_" + std::to_string(i++);
}
return true;
}
ModelVolume* ModelObject::add_volume(const TriangleMesh &mesh)
{
ModelVolume* v = new ModelVolume(this, mesh);

4
src/libslic3r/Model.hpp
View file

@ -508,6 +508,10 @@ public:
ModelObjectPtrs segment(size_t instance, unsigned int max_extruders, double smoothing_alpha = 0.5, int segment_number = 5);
void split(ModelObjectPtrs* new_objects);
void merge();
// BBS: Boolean opts - Musang King
bool make_boolean(ModelObject *cut_object, const std::string &boolean_opts);
ModelObjectPtrs merge_volumes(std::vector<int>& vol_indeces);//BBS
// Support for non-uniform scaling of instances. If an instance is rotated by angles, which are not multiples of ninety degrees,
// then the scaling in world coordinate system is not representable by the Geometry::Transformation structure.

417
src/mcut/CMakeLists.txt Normal file
View file

@ -0,0 +1,417 @@
#
# Copyright (c) 2021-2022 Floyd M. Chitalu.
# All rights reserved.
#
# NOTE: This file is licensed under GPL-3.0-or-later (default).
# A commercial license can be purchased from Floyd M. Chitalu.
#
# License details:
#
# (A) GNU General Public License ("GPL"); a copy of which you should have
# recieved with this file.
# - see also: <http://www.gnu.org/licenses/>
# (B) Commercial license.
# - email: floyd.m.chitalu@gmail.com
#
# The commercial license options is for users that wish to use MCUT in
# their products for comercial purposes but do not wish to release their
# software products under the GPL license.
#
# Author(s) : Floyd M. Chitalu
#
############################################################################################
#
# You can configure MCUT with the following CMake options:
#
# MCUT_BUILD_AS_SHARED_LIB [default=ON] - Build MCUT as a shared/dynamic library (.so/.dll). -- SET TO OFF
# MCUT_BUILD_DOCUMENTATION [default=OFF] - Build documentation (explicit dependancy on Doxygen)
# MCUT_BUILD_TESTS [default=OFF] - Build the tests (implicit dependancy on GoogleTest)
# MCUT_BUILD_TUTORIALS [default=OFF] - Build tutorials
# MCUT_BUILD_WITH_COMPUTE_HELPER_THREADPOOL [default=ON] - Build as configurable multi-threaded library
#
# This script will define the following CMake cache variables:
#
# MCUT_INCLUDE_DIR - the MCUT include directory
# MCUT_LIB_PATH - path to the MCUT library (i.e. .so/.dll or .a/.lib files)
cmake_minimum_required(VERSION 3.12...3.13 FATAL_ERROR)
project(mcut LANGUAGES CXX C)
get_directory_property(MCUT_PARENT_DIR PARENT_DIRECTORY)
if(NOT MCUT_PARENT_DIR)
set(MCUT_TOPLEVEL_PROJECT ON)
else()
set(MCUT_TOPLEVEL_PROJECT OFF)
endif()
set ( DESCRIPTION "Mesh cutting library." )
if (NOT WIN32 AND NOT CMAKE_BUILD_TYPE)
message(STATUS "No build type selected, default to Release")
set(CMAKE_BUILD_TYPE "Release" CACHE STRING "Choose the type of build, options are: Debug Release RelWithDebInfo MinSizeRel." FORCE)
endif()
set (CMAKE_DEBUG_POSTFIX "d")
#set(CMAKE_LIBRARY_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/bin)
#set(CMAKE_ARCHIVE_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
#set(CMAKE_RUNTIME_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/bin)
list (APPEND CMAKE_MODULE_PATH "${CMAKE_CURRENT_SOURCE_DIR}/cmake")
set(MCUT_MAJOR 1)
set(MCUT_MINOR 1)
set(MCUT_PATCH 0)
set( MCUT_VERSION "${MCUT_MAJOR}.${MCUT_MINOR}.${MCUT_PATCH}" )
message(STATUS "[MCUT] version: ${MCUT_VERSION}")
set (CMAKE_CXX_STANDARD 11)
set (CMAKE_CXX_STANDARD_REQUIRED True)
set (CMAKE_EXPORT_COMPILE_COMMANDS ON)
set (project_API_version_string "${MCUT_VERSION}")
set (project_build_version_string "${MCUT_VERSION}")
set (project_namespace_name MCUT)
# Only do these if this is the main project, and not if it is included through add_subdirectory
if(CMAKE_PROJECT_NAME STREQUAL PROJECT_NAME)
# Let's ensure -std=c++xx instead of -std=g++xx
set(CMAKE_CXX_EXTENSIONS OFF)
# Let's nicely support folders in IDEs
set_property(GLOBAL PROPERTY USE_FOLDERS ON)
# Testing only available if this is the main app
# Note this needs to be done in the main CMakeLists
# since it calls enable_testing, which must be in the
# main CMakeLists.
include(CTest)
if(MCUT_BUILD_DOCUMENTATION)
# Docs only available if this is the main app
find_package(Doxygen)
if(Doxygen_FOUND)
add_subdirectory(docs)
else()
message(STATUS "Doxygen not found, not building docs")
endif()
endif()
endif()
include(CMakePrintHelpers)
#
# User options
#
option(MCUT_BUILD_DOCUMENTATION "Configure to build docs with Doxygen" OFF) # OFF by default
if (MCUT_TOPLEVEL_PROJECT AND NOT MCUT_BUILD_TESTS)
option(MCUT_BUILD_TESTS "Configure to build tests" ON)
endif()
option(MCUT_BUILD_AS_SHARED_LIB "Configure to build MCUT as a shared/dynamic library" OFF)
option(MCUT_BUILD_WITH_COMPUTE_HELPER_THREADPOOL "Configure to build MCUT engine with a shared (amongst contexts) thread-pool" ON)
if (MCUT_TOPLEVEL_PROJECT AND NOT MCUT_BUILD_TUTORIALS)
option(MCUT_BUILD_TUTORIALS "Configure to build MCUT tutorials" ON)
endif()
#
# machine-precision-numbers library targets
#
set (target_name mcut)
#
# MCUT compilation variables/settings
#
set (MCUT_INCLUDE_DIR ${CMAKE_CURRENT_SOURCE_DIR}/include CACHE STRING "The MCUT include directory")
message(STATUS "[MCUT] MCUT_INCLUDE_DIR=${MCUT_INCLUDE_DIR}")
list(APPEND include_dirs ${MCUT_INCLUDE_DIR})
list(APPEND compilation_flags "")
if(MCUT_BUILD_AS_SHARED_LIB)
list(APPEND preprocessor_defs -DMCUT_SHARED_LIB=1)
endif()
find_package(Threads REQUIRED)
list(APPEND extra_libs Threads::Threads)
if (MCUT_BUILD_WITH_COMPUTE_HELPER_THREADPOOL)
list(APPEND preprocessor_defs -DMCUT_WITH_COMPUTE_HELPER_THREADPOOL=1)
endif()
if ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang" OR "${CMAKE_CXX_COMPILER_ID}" STREQUAL "GNU")
list(APPEND compilation_flags -Wall -Wextra)
elseif ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Intel")
list(APPEND compilation_flags w3 -diag-disable:remark)
elseif ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "MSVC")
list(APPEND compilation_flags /W4 /wd26812 /bigobj)
list(APPEND preprocessor_defs -D_CRT_SECURE_NO_WARNINGS)
endif()
message(STATUS "[MCUT] compilation_flags=${compilation_flags}")
message(STATUS "[MCUT] preprocessor_defs=${preprocessor_defs}")
message(STATUS "[MCUT] extra_libs=${extra_libs}")
set ( project_source_files
${CMAKE_CURRENT_SOURCE_DIR}/source/mcut.cpp
${CMAKE_CURRENT_SOURCE_DIR}/source/kernel.cpp
${CMAKE_CURRENT_SOURCE_DIR}/source/hmesh.cpp
${CMAKE_CURRENT_SOURCE_DIR}/source/math.cpp
${CMAKE_CURRENT_SOURCE_DIR}/source/bvh.cpp
${CMAKE_CURRENT_SOURCE_DIR}/source/shewchuk.c
${CMAKE_CURRENT_SOURCE_DIR}/source/frontend.cpp
${CMAKE_CURRENT_SOURCE_DIR}/source/preproc.cpp)
#
# Create MCUT target(s)
#
# This function invokes commands which create a library target (static, shared etc.)
function(create_library_target LIBRARY_TYPE)
message(STATUS "[MCUT] create target: name=${target_name} type=${LIBRARY_TYPE}")
add_library(${target_name} ${LIBRARY_TYPE} ${project_source_files})
target_include_directories(${target_name} PRIVATE ${include_dirs})
target_link_libraries(${target_name} PRIVATE ${extra_libs})
target_compile_options(${target_name} PRIVATE ${compilation_flags})
target_compile_definitions(${target_name} PRIVATE ${preprocessor_defs} )
if (MCUT_BUILD_AS_SHARED_LIB AND WIN32)
# add macro to export .dll symbols
target_compile_definitions(${target_name} PRIVATE -DMCUT_EXPORT_SHARED_LIB_SYMBOLS=1)
endif()
set_property(TARGET ${target_name} PROPERTY VERSION ${project_build_version_string})
set_property(TARGET ${target_name} PROPERTY SOVERSION ${project_API_version_string})
get_target_property(target_type ${target_name} TYPE)
if ("${target_type}" STREQUAL "SHARED")
set_property(TARGET ${target_name} PROPERTY POSITION_INDEPENDENT_CODE ON)
if(MSVC)
set_property(TARGET ${target_name} PROPERTY WINDOWS_EXPORT_ALL_SYMBOLS ON)
endif()
endif()
endfunction()
#
# create target
#
if(MCUT_BUILD_AS_SHARED_LIB)
create_library_target(SHARED)
else()
create_library_target(STATIC)
endif()
set (MCUT_LIB_PATH $<TARGET_FILE:mcut> CACHE STRING "Path to the compiled MCUT library file")
message(STATUS "[MCUT] MCUT_LIB_PATH=${MCUT_LIB_PATH}")
#
# tests, tutorials etc. are dependant on third-party libs (file loaders etc.)
# We enable the code that is used to download those projects here:
#
if(${MCUT_BUILD_TESTS} OR ${MCUT_BUILD_TUTORIALS})
# FetchContent added in CMake 3.11, downloads during the configure step
include(FetchContent)
# FetchContent_MakeAvailable was not added until CMake 3.14; use our shim
if(${CMAKE_VERSION} VERSION_LESS 3.14)
include(cmake/add_FetchContent_MakeAvailable.cmake)
endif()
FetchContent_Populate(
libigl
GIT_REPOSITORY https://github.com/libigl/libigl.git
GIT_TAG v2.3.0
GIT_PROGRESS TRUE
)
#set(libigl_include_dir ${CMAKE_BINARY_DIR}/libigl-src/include)
set(libigl_include_dir ${libigl_SOURCE_DIR}/include)
set(LIBIGL_EIGEN_VERSION 3.3.7 CACHE STRING "Default version of Eigen used by libigl.")
# used by tests & tutorials
#download_project(PROJ eigen
# GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git
# GIT_TAG ${LIBIGL_EIGEN_VERSION}
# ${UPDATE_DISCONNECTED_IF_AVAILABLE}
#)
FetchContent_Declare(
eigen
GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git
GIT_TAG ${LIBIGL_EIGEN_VERSION}
GIT_SHALLOW TRUE
GIT_PROGRESS TRUE
)
set(EIGEN_BUILD_DOC OFF)
# note: To disable eigen tests,
# you should put this code in a add_subdirectory to avoid to change
# BUILD_TESTING for your own project too since variables are directory
# scoped
set(BUILD_TESTING OFF)
set(EIGEN_BUILD_PKGCONFIG OFF)
set( OFF)
FetchContent_MakeAvailable(eigen)
#set(eigen_include_dir ${CMAKE_BINARY_DIR}/eigen-src)
set(eigen_include_dir ${eigen_SOURCE_DIR})
endif()
#
# tests
#
if(MCUT_BUILD_TESTS)
add_subdirectory(tests)
endif()
#
# tutorials
#
if(MCUT_BUILD_TUTORIALS)
add_subdirectory(tutorials)
endif()
#
# documentation
#
if(MCUT_BUILD_DOCUMENTATION)
add_subdirectory(docs)
endif()
########################################################
### PACKAGING ###
### This is a quite INCOMPLETE set of variables that ###
### should be set for the various generators. ###
### Consult the CPack documentations for a full set. ###
########################################################
# https://gitlab.kitware.com/cmake/community/-/wikis/doc/cpack/Component-Install-With-CPack
# https://stackoverflow.com/questions/6003374/what-is-cmake-equivalent-of-configure-prefix-dir-make-all-install
# TODO: package documentation files
if(MCUT_BUILD_AS_SHARED_LIB)
#
# dynamic libs
#
install(TARGETS ${mpn_shared_lib_name}
LIBRARY
DESTINATION lib/shared
COMPONENT dynamic_libraries)
else()
#
# static libs
#
install(TARGETS ${mpn_static_lib_name}
ARCHIVE
DESTINATION lib/static
COMPONENT static_libraries)
endif()
#
# headers
#
install(FILES ${MCUT_INCLUDE_DIR}/mcut/mcut.h ${MCUT_INCLUDE_DIR}/mcut/platform.h
DESTINATION include/mcut
COMPONENT headers)
install(FILES
${CMAKE_CURRENT_SOURCE_DIR}/LICENSE.txt
${CMAKE_CURRENT_SOURCE_DIR}/README.md
DESTINATION ./
COMPONENT text_files)
#
# notify CPack of the names of all of the components in the project
#
set(CPACK_COMPONENTS_ALL static_libraries dynamic_libraries headers text_files) # applications
set(CPACK_COMPONENT_APPLICATIONS_DISPLAY_NAME "MCUT Application")
set(CPACK_COMPONENT_STATIC_LIBRARIES_DISPLAY_NAME "Static Libraries")
set(CPACK_COMPONENT_DYNAMIC_LIBRARIES_DISPLAY_NAME "Dynamics Libraries")
set(CPACK_COMPONENT_HEADERS_DISPLAY_NAME "C++ Headers")
set(CPACK_COMPONENT_APPLICATIONS_DESCRIPTION
"A simple application using MCUT")
set(CPACK_COMPONENT_STATIC_LIBRARIES_DESCRIPTION
"Static libraries used to build programs with MCUT")
set(CPACK_COMPONENT_DYNAMIC_LIBRARIES_DESCRIPTION
"Dynamic libraries used to build programs with MCUT")
set(CPACK_COMPONENT_HEADERS_DESCRIPTION
"C/C++ header files for use with MCUT")
#
# component dependencies
#
set(CPACK_COMPONENT_HEADERS_DEPENDS static_libraries dynamic_libraries)
set(CPACK_COMPONENT_APPLICATIONS_GROUP "Runtime")
set(CPACK_COMPONENT_STATIC_LIBRARIES_GROUP "Development")
set(CPACK_COMPONENT_DYNAMIC_LIBRARIES_GROUP "Development")
set(CPACK_COMPONENT_HEADERS_GROUP "Development")
set(CPACK_COMPONENT_GROUP_DEVELOPMENT_DESCRIPTION
"All of the tools you'll ever need to develop software")
set (CPACK_PACKAGE_NAME "MCUT")
set (CPACK_PACKAGE_VENDOR "Floyd M. Chitalu")
set (CPACK_PACKAGE_DIRECTORY "${CMAKE_CURRENT_BINARY_DIR}")
set (CPACK_PACKAGE_VERSION_MAJOR "${MCUT_MAJOR}")
set (CPACK_PACKAGE_VERSION_MINOR "${MCUT_MINOR}")
set (CPACK_PACKAGE_VERSION_PATCH "${MCUT_PATCH}")
#set (CPACK_PACKAGE_DESCRIPTION "MCUT (pronounced ‘emcut’) is a tool for cutting meshes.")
#set (CPACK_PACKAGE_DESCRIPTION_FILE ${CMAKE_CURRENT_SOURCE_DIR}/DESCRIPTION.txt)
set (CPACK_PACKAGE_DESCRIPTION_SUMMARY "MCUT is a library for cutting meshes to perform tasks like boolean operations and more.")
set (CPACK_PACKAGE_HOMEPAGE_URL "https://cutdigital.github.io/mcut.site/")
set (CPACK_PACKAGE_INSTALL_DIRECTORY ${CPACK_PACKAGE_NAME})
# set (CPACK_PACKAGE_ICON )
set (CPACK_PACKAGE_CHECKSUM SHA256)
#set (CPACK_PROJECT_CONFIG_FILE )
set (CPACK_RESOURCE_FILE_LICENSE ${CMAKE_CURRENT_SOURCE_DIR}/LICENSE.txt) # must also include in install command
set (CPACK_RESOURCE_FILE_README ${CMAKE_CURRENT_SOURCE_DIR}/README.md)
#set (CPACK_RESOURCE_FILE_WELCOME ${CMAKE_CURRENT_SOURCE_DIR}/WELCOME.txt)
if (WIN32)
if (USE_WIX_TOOLSET)
set(CPACK_GENERATOR "WIX") # this need WiX Tooset installed and a path to candle.exe
else ()
set(CPACK_GENERATOR "NSIS") # this needs NSIS installed, and available
endif ()
elseif ( ${CMAKE_SYSTEM_NAME} MATCHES "Darwin" )
set(CPACK_GENERATOR "PackageMake")
else ()
set(CPACK_GENERATOR "TGZ")
endif ()
#set (CPACK_OUTPUT_CONFIG_FILE ) # Defaults to CPackConfig.cmake.
#set (CPACK_PACKAGE_EXECUTABLES )
set (CPACK_STRIP_FILES TRUE)
# set (CPACK_VERBATIM_VARIABLES )
# set (CPACK_SOURCE_PACKAGE_FILE_NAME )
# set (CPACK_SOURCE_STRIP_FILES )
# set (CPACK_SOURCE_GENERATOR )
# set (CPACK_SOURCE_OUTPUT_CONFIG_FILE )
# set (CPACK_SOURCE_IGNORE_FILES )
# set (CPACK_CMAKE_GENERATOR )
# set (CPACK_INSTALL_CMAKE_PROJECTS )
# set (CPACK_INSTALL_CMAKE_PROJECTS )
# set (CPACK_SYSTEM_NAME )
# set (CPACK_PACKAGE_VERSION )
# set (CPACK_TOPLEVEL_TAG )
# set (CPACK_INSTALL_COMMANDS )
# set (CPACK_INSTALLED_DIRECTORIES )
# set ( )
include(CPack)
# eof

674
src/mcut/LICENSE.GPL.txt Normal file
View file

@ -0,0 +1,674 @@
GNU GENERAL PUBLIC LICENSE
Version 3, 29 June 2007
Copyright (C) 2007 Free Software Foundation, Inc. <http://fsf.org/>
Everyone is permitted to copy and distribute verbatim copies
of this license document, but changing it is not allowed.
Preamble
The GNU General Public License is a free, copyleft license for
software and other kinds of works.
The licenses for most software and other practical works are designed
to take away your freedom to share and change the works. By contrast,
the GNU General Public License is intended to guarantee your freedom to
share and change all versions of a program--to make sure it remains free
software for all its users. We, the Free Software Foundation, use the
GNU General Public License for most of our software; it applies also to
any other work released this way by its authors. You can apply it to
your programs, too.
When we speak of free software, we are referring to freedom, not
price. Our General Public Licenses are designed to make sure that you
have the freedom to distribute copies of free software (and charge for
them if you wish), that you receive source code or can get it if you
want it, that you can change the software or use pieces of it in new
free programs, and that you know you can do these things.
To protect your rights, we need to prevent others from denying you
these rights or asking you to surrender the rights. Therefore, you have
certain responsibilities if you distribute copies of the software, or if
you modify it: responsibilities to respect the freedom of others.
For example, if you distribute copies of such a program, whether
gratis or for a fee, you must pass on to the recipients the same
freedoms that you received. You must make sure that they, too, receive
or can get the source code. And you must show them these terms so they
know their rights.
Developers that use the GNU GPL protect your rights with two steps:
(1) assert copyright on the software, and (2) offer you this License
giving you legal permission to copy, distribute and/or modify it.
For the developers' and authors' protection, the GPL clearly explains
that there is no warranty for this free software. For both users' and
authors' sake, the GPL requires that modified versions be marked as
changed, so that their problems will not be attributed erroneously to
authors of previous versions.
Some devices are designed to deny users access to install or run
modified versions of the software inside them, although the manufacturer
can do so. This is fundamentally incompatible with the aim of
protecting users' freedom to change the software. The systematic
pattern of such abuse occurs in the area of products for individuals to
use, which is precisely where it is most unacceptable. Therefore, we
have designed this version of the GPL to prohibit the practice for those
products. If such problems arise substantially in other domains, we
stand ready to extend this provision to those domains in future versions
of the GPL, as needed to protect the freedom of users.
Finally, every program is threatened constantly by software patents.
States should not allow patents to restrict development and use of
software on general-purpose computers, but in those that do, we wish to
avoid the special danger that patents applied to a free program could
make it effectively proprietary. To prevent this, the GPL assures that
patents cannot be used to render the program non-free.
The precise terms and conditions for copying, distribution and
modification follow.
TERMS AND CONDITIONS
0. Definitions.
"This License" refers to version 3 of the GNU General Public License.
"Copyright" also means copyright-like laws that apply to other kinds of
works, such as semiconductor masks.
"The Program" refers to any copyrightable work licensed under this
License. Each licensee is addressed as "you". "Licensees" and
"recipients" may be individuals or organizations.
To "modify" a work means to copy from or adapt all or part of the work
in a fashion requiring copyright permission, other than the making of an
exact copy. The resulting work is called a "modified version" of the
earlier work or a work "based on" the earlier work.
A "covered work" means either the unmodified Program or a work based
on the Program.
To "propagate" a work means to do anything with it that, without
permission, would make you directly or secondarily liable for
infringement under applicable copyright law, except executing it on a
computer or modifying a private copy. Propagation includes copying,
distribution (with or without modification), making available to the
public, and in some countries other activities as well.
To "convey" a work means any kind of propagation that enables other
parties to make or receive copies. Mere interaction with a user through
a computer network, with no transfer of a copy, is not conveying.
An interactive user interface displays "Appropriate Legal Notices"
to the extent that it includes a convenient and prominently visible
feature that (1) displays an appropriate copyright notice, and (2)
tells the user that there is no warranty for the work (except to the
extent that warranties are provided), that licensees may convey the
work under this License, and how to view a copy of this License. If
the interface presents a list of user commands or options, such as a
menu, a prominent item in the list meets this criterion.
1. Source Code.
The "source code" for a work means the preferred form of the work
for making modifications to it. "Object code" means any non-source
form of a work.
A "Standard Interface" means an interface that either is an official
standard defined by a recognized standards body, or, in the case of
interfaces specified for a particular programming language, one that
is widely used among developers working in that language.
The "System Libraries" of an executable work include anything, other
than the work as a whole, that (a) is included in the normal form of
packaging a Major Component, but which is not part of that Major
Component, and (b) serves only to enable use of the work with that
Major Component, or to implement a Standard Interface for which an
implementation is available to the public in source code form. A
"Major Component", in this context, means a major essential component
(kernel, window system, and so on) of the specific operating system
(if any) on which the executable work runs, or a compiler used to
produce the work, or an object code interpreter used to run it.
The "Corresponding Source" for a work in object code form means all
the source code needed to generate, install, and (for an executable
work) run the object code and to modify the work, including scripts to
control those activities. However, it does not include the work's
System Libraries, or general-purpose tools or generally available free
programs which are used unmodified in performing those activities but
which are not part of the work. For example, Corresponding Source
includes interface definition files associated with source files for
the work, and the source code for shared libraries and dynamically
linked subprograms that the work is specifically designed to require,
such as by intimate data communication or control flow between those
subprograms and other parts of the work.
The Corresponding Source need not include anything that users
can regenerate automatically from other parts of the Corresponding
Source.
The Corresponding Source for a work in source code form is that
same work.
2. Basic Permissions.
All rights granted under this License are granted for the term of
copyright on the Program, and are irrevocable provided the stated
conditions are met. This License explicitly affirms your unlimited
permission to run the unmodified Program. The output from running a
covered work is covered by this License only if the output, given its
content, constitutes a covered work. This License acknowledges your
rights of fair use or other equivalent, as provided by copyright law.
You may make, run and propagate covered works that you do not
convey, without conditions so long as your license otherwise remains
in force. You may convey covered works to others for the sole purpose
of having them make modifications exclusively for you, or provide you
with facilities for running those works, provided that you comply with
the terms of this License in conveying all material for which you do
not control copyright. Those thus making or running the covered works
for you must do so exclusively on your behalf, under your direction
and control, on terms that prohibit them from making any copies of
your copyrighted material outside their relationship with you.
Conveying under any other circumstances is permitted solely under
the conditions stated below. Sublicensing is not allowed; section 10
makes it unnecessary.
3. Protecting Users' Legal Rights From Anti-Circumvention Law.
No covered work shall be deemed part of an effective technological
measure under any applicable law fulfilling obligations under article
11 of the WIPO copyright treaty adopted on 20 December 1996, or
similar laws prohibiting or restricting circumvention of such
measures.
When you convey a covered work, you waive any legal power to forbid
circumvention of technological measures to the extent such circumvention
is effected by exercising rights under this License with respect to
the covered work, and you disclaim any intention to limit operation or
modification of the work as a means of enforcing, against the work's
users, your or third parties' legal rights to forbid circumvention of
technological measures.
4. Conveying Verbatim Copies.
You may convey verbatim copies of the Program's source code as you
receive it, in any medium, provided that you conspicuously and
appropriately publish on each copy an appropriate copyright notice;
keep intact all notices stating that this License and any
non-permissive terms added in accord with section 7 apply to the code;
keep intact all notices of the absence of any warranty; and give all
recipients a copy of this License along with the Program.
You may charge any price or no price for each copy that you convey,
and you may offer support or warranty protection for a fee.
5. Conveying Modified Source Versions.
You may convey a work based on the Program, or the modifications to
produce it from the Program, in the form of source code under the
terms of section 4, provided that you also meet all of these conditions:
a) The work must carry prominent notices stating that you modified
it, and giving a relevant date.
b) The work must carry prominent notices stating that it is
released under this License and any conditions added under section
7. This requirement modifies the requirement in section 4 to
"keep intact all notices".
c) You must license the entire work, as a whole, under this
License to anyone who comes into possession of a copy. This
License will therefore apply, along with any applicable section 7
additional terms, to the whole of the work, and all its parts,
regardless of how they are packaged. This License gives no
permission to license the work in any other way, but it does not
invalidate such permission if you have separately received it.
d) If the work has interactive user interfaces, each must display
Appropriate Legal Notices; however, if the Program has interactive
interfaces that do not display Appropriate Legal Notices, your
work need not make them do so.
A compilation of a covered work with other separate and independent
works, which are not by their nature extensions of the covered work,
and which are not combined with it such as to form a larger program,
in or on a volume of a storage or distribution medium, is called an
"aggregate" if the compilation and its resulting copyright are not
used to limit the access or legal rights of the compilation's users
beyond what the individual works permit. Inclusion of a covered work
in an aggregate does not cause this License to apply to the other
parts of the aggregate.
6. Conveying Non-Source Forms.
You may convey a covered work in object code form under the terms
of sections 4 and 5, provided that you also convey the
machine-readable Corresponding Source under the terms of this License,
in one of these ways:
a) Convey the object code in, or embodied in, a physical product
(including a physical distribution medium), accompanied by the
Corresponding Source fixed on a durable physical medium
customarily used for software interchange.
b) Convey the object code in, or embodied in, a physical product
(including a physical distribution medium), accompanied by a
written offer, valid for at least three years and valid for as
long as you offer spare parts or customer support for that product
model, to give anyone who possesses the object code either (1) a
copy of the Corresponding Source for all the software in the
product that is covered by this License, on a durable physical
medium customarily used for software interchange, for a price no
more than your reasonable cost of physically performing this
conveying of source, or (2) access to copy the
Corresponding Source from a network server at no charge.
c) Convey individual copies of the object code with a copy of the
written offer to provide the Corresponding Source. This
alternative is allowed only occasionally and noncommercially, and
only if you received the object code with such an offer, in accord
with subsection 6b.
d) Convey the object code by offering access from a designated
place (gratis or for a charge), and offer equivalent access to the
Corresponding Source in the same way through the same place at no
further charge. You need not require recipients to copy the
Corresponding Source along with the object code. If the place to
copy the object code is a network server, the Corresponding Source
may be on a different server (operated by you or a third party)
that supports equivalent copying facilities, provided you maintain
clear directions next to the object code saying where to find the
Corresponding Source. Regardless of what server hosts the
Corresponding Source, you remain obligated to ensure that it is
available for as long as needed to satisfy these requirements.
e) Convey the object code using peer-to-peer transmission, provided
you inform other peers where the object code and Corresponding
Source of the work are being offered to the general public at no
charge under subsection 6d.
A separable portion of the object code, whose source code is excluded
from the Corresponding Source as a System Library, need not be
included in conveying the object code work.
A "User Product" is either (1) a "consumer product", which means any
tangible personal property which is normally used for personal, family,
or household purposes, or (2) anything designed or sold for incorporation
into a dwelling. In determining whether a product is a consumer product,
doubtful cases shall be resolved in favor of coverage. For a particular
product received by a particular user, "normally used" refers to a
typical or common use of that class of product, regardless of the status
of the particular user or of the way in which the particular user
actually uses, or expects or is expected to use, the product. A product
is a consumer product regardless of whether the product has substantial
commercial, industrial or non-consumer uses, unless such uses represent
the only significant mode of use of the product.
"Installation Information" for a User Product means any methods,
procedures, authorization keys, or other information required to install
and execute modified versions of a covered work in that User Product from
a modified version of its Corresponding Source. The information must
suffice to ensure that the continued functioning of the modified object
code is in no case prevented or interfered with solely because
modification has been made.
If you convey an object code work under this section in, or with, or
specifically for use in, a User Product, and the conveying occurs as
part of a transaction in which the right of possession and use of the
User Product is transferred to the recipient in perpetuity or for a
fixed term (regardless of how the transaction is characterized), the
Corresponding Source conveyed under this section must be accompanied
by the Installation Information. But this requirement does not apply
if neither you nor any third party retains the ability to install
modified object code on the User Product (for example, the work has
been installed in ROM).
The requirement to provide Installation Information does not include a
requirement to continue to provide support service, warranty, or updates
for a work that has been modified or installed by the recipient, or for
the User Product in which it has been modified or installed. Access to a
network may be denied when the modification itself materially and
adversely affects the operation of the network or violates the rules and
protocols for communication across the network.
Corresponding Source conveyed, and Installation Information provided,
in accord with this section must be in a format that is publicly
documented (and with an implementation available to the public in
source code form), and must require no special password or key for
unpacking, reading or copying.
7. Additional Terms.
"Additional permissions" are terms that supplement the terms of this
License by making exceptions from one or more of its conditions.
Additional permissions that are applicable to the entire Program shall
be treated as though they were included in this License, to the extent
that they are valid under applicable law. If additional permissions
apply only to part of the Program, that part may be used separately
under those permissions, but the entire Program remains governed by
this License without regard to the additional permissions.
When you convey a copy of a covered work, you may at your option
remove any additional permissions from that copy, or from any part of
it. (Additional permissions may be written to require their own
removal in certain cases when you modify the work.) You may place
additional permissions on material, added by you to a covered work,
for which you have or can give appropriate copyright permission.
Notwithstanding any other provision of this License, for material you
add to a covered work, you may (if authorized by the copyright holders of
that material) supplement the terms of this License with terms:
a) Disclaiming warranty or limiting liability differently from the
terms of sections 15 and 16 of this License; or
b) Requiring preservation of specified reasonable legal notices or
author attributions in that material or in the Appropriate Legal
Notices displayed by works containing it; or
c) Prohibiting misrepresentation of the origin of that material, or
requiring that modified versions of such material be marked in
reasonable ways as different from the original version; or
d) Limiting the use for publicity purposes of names of licensors or
authors of the material; or
e) Declining to grant rights under trademark law for use of some
trade names, trademarks, or service marks; or
f) Requiring indemnification of licensors and authors of that
material by anyone who conveys the material (or modified versions of
it) with contractual assumptions of liability to the recipient, for
any liability that these contractual assumptions directly impose on
those licensors and authors.
All other non-permissive additional terms are considered "further
restrictions" within the meaning of section 10. If the Program as you
received it, or any part of it, contains a notice stating that it is
governed by this License along with a term that is a further
restriction, you may remove that term. If a license document contains
a further restriction but permits relicensing or conveying under this
License, you may add to a covered work material governed by the terms
of that license document, provided that the further restriction does
not survive such relicensing or conveying.
If you add terms to a covered work in accord with this section, you
must place, in the relevant source files, a statement of the
additional terms that apply to those files, or a notice indicating
where to find the applicable terms.
Additional terms, permissive or non-permissive, may be stated in the
form of a separately written license, or stated as exceptions;
the above requirements apply either way.
8. Termination.
You may not propagate or modify a covered work except as expressly
provided under this License. Any attempt otherwise to propagate or
modify it is void, and will automatically terminate your rights under
this License (including any patent licenses granted under the third
paragraph of section 11).
However, if you cease all violation of this License, then your
license from a particular copyright holder is reinstated (a)
provisionally, unless and until the copyright holder explicitly and
finally terminates your license, and (b) permanently, if the copyright
holder fails to notify you of the violation by some reasonable means
prior to 60 days after the cessation.
Moreover, your license from a particular copyright holder is
reinstated permanently if the copyright holder notifies you of the
violation by some reasonable means, this is the first time you have
received notice of violation of this License (for any work) from that
copyright holder, and you cure the violation prior to 30 days after
your receipt of the notice.
Termination of your rights under this section does not terminate the
licenses of parties who have received copies or rights from you under
this License. If your rights have been terminated and not permanently
reinstated, you do not qualify to receive new licenses for the same
material under section 10.
9. Acceptance Not Required for Having Copies.
You are not required to accept this License in order to receive or
run a copy of the Program. Ancillary propagation of a covered work
occurring solely as a consequence of using peer-to-peer transmission
to receive a copy likewise does not require acceptance. However,
nothing other than this License grants you permission to propagate or
modify any covered work. These actions infringe copyright if you do
not accept this License. Therefore, by modifying or propagating a
covered work, you indicate your acceptance of this License to do so.
10. Automatic Licensing of Downstream Recipients.
Each time you convey a covered work, the recipient automatically
receives a license from the original licensors, to run, modify and
propagate that work, subject to this License. You are not responsible
for enforcing compliance by third parties with this License.
An "entity transaction" is a transaction transferring control of an
organization, or substantially all assets of one, or subdividing an
organization, or merging organizations. If propagation of a covered
work results from an entity transaction, each party to that
transaction who receives a copy of the work also receives whatever
licenses to the work the party's predecessor in interest had or could
give under the previous paragraph, plus a right to possession of the
Corresponding Source of the work from the predecessor in interest, if
the predecessor has it or can get it with reasonable efforts.
You may not impose any further restrictions on the exercise of the
rights granted or affirmed under this License. For example, you may
not impose a license fee, royalty, or other charge for exercise of
rights granted under this License, and you may not initiate litigation
(including a cross-claim or counterclaim in a lawsuit) alleging that
any patent claim is infringed by making, using, selling, offering for
sale, or importing the Program or any portion of it.
11. Patents.
A "contributor" is a copyright holder who authorizes use under this
License of the Program or a work on which the Program is based. The
work thus licensed is called the contributor's "contributor version".
A contributor's "essential patent claims" are all patent claims
owned or controlled by the contributor, whether already acquired or
hereafter acquired, that would be infringed by some manner, permitted
by this License, of making, using, or selling its contributor version,
but do not include claims that would be infringed only as a
consequence of further modification of the contributor version. For
purposes of this definition, "control" includes the right to grant
patent sublicenses in a manner consistent with the requirements of
this License.
Each contributor grants you a non-exclusive, worldwide, royalty-free
patent license under the contributor's essential patent claims, to
make, use, sell, offer for sale, import and otherwise run, modify and
propagate the contents of its contributor version.
In the following three paragraphs, a "patent license" is any express
agreement or commitment, however denominated, not to enforce a patent
(such as an express permission to practice a patent or covenant not to
sue for patent infringement). To "grant" such a patent license to a
party means to make such an agreement or commitment not to enforce a
patent against the party.
If you convey a covered work, knowingly relying on a patent license,
and the Corresponding Source of the work is not available for anyone
to copy, free of charge and under the terms of this License, through a
publicly available network server or other readily accessible means,
then you must either (1) cause the Corresponding Source to be so
available, or (2) arrange to deprive yourself of the benefit of the
patent license for this particular work, or (3) arrange, in a manner
consistent with the requirements of this License, to extend the patent
license to downstream recipients. "Knowingly relying" means you have
actual knowledge that, but for the patent license, your conveying the
covered work in a country, or your recipient's use of the covered work
in a country, would infringe one or more identifiable patents in that
country that you have reason to believe are valid.
If, pursuant to or in connection with a single transaction or
arrangement, you convey, or propagate by procuring conveyance of, a
covered work, and grant a patent license to some of the parties
receiving the covered work authorizing them to use, propagate, modify
or convey a specific copy of the covered work, then the patent license
you grant is automatically extended to all recipients of the covered
work and works based on it.
A patent license is "discriminatory" if it does not include within
the scope of its coverage, prohibits the exercise of, or is
conditioned on the non-exercise of one or more of the rights that are
specifically granted under this License. You may not convey a covered
work if you are a party to an arrangement with a third party that is
in the business of distributing software, under which you make payment
to the third party based on the extent of your activity of conveying
the work, and under which the third party grants, to any of the
parties who would receive the covered work from you, a discriminatory
patent license (a) in connection with copies of the covered work
conveyed by you (or copies made from those copies), or (b) primarily
for and in connection with specific products or compilations that
contain the covered work, unless you entered into that arrangement,
or that patent license was granted, prior to 28 March 2007.
Nothing in this License shall be construed as excluding or limiting
any implied license or other defenses to infringement that may
otherwise be available to you under applicable patent law.
12. No Surrender of Others' Freedom.
If conditions are imposed on you (whether by court order, agreement or
otherwise) that contradict the conditions of this License, they do not
excuse you from the conditions of this License. If you cannot convey a
covered work so as to satisfy simultaneously your obligations under this
License and any other pertinent obligations, then as a consequence you may
not convey it at all. For example, if you agree to terms that obligate you
to collect a royalty for further conveying from those to whom you convey
the Program, the only way you could satisfy both those terms and this
License would be to refrain entirely from conveying the Program.
13. Use with the GNU Affero General Public License.
Notwithstanding any other provision of this License, you have
permission to link or combine any covered work with a work licensed
under version 3 of the GNU Affero General Public License into a single
combined work, and to convey the resulting work. The terms of this
License will continue to apply to the part which is the covered work,
but the special requirements of the GNU Affero General Public License,
section 13, concerning interaction through a network will apply to the
combination as such.
14. Revised Versions of this License.
The Free Software Foundation may publish revised and/or new versions of
the GNU General Public License from time to time. Such new versions will
be similar in spirit to the present version, but may differ in detail to
address new problems or concerns.
Each version is given a distinguishing version number. If the
Program specifies that a certain numbered version of the GNU General
Public License "or any later version" applies to it, you have the
option of following the terms and conditions either of that numbered
version or of any later version published by the Free Software
Foundation. If the Program does not specify a version number of the
GNU General Public License, you may choose any version ever published
by the Free Software Foundation.
If the Program specifies that a proxy can decide which future
versions of the GNU General Public License can be used, that proxy's
public statement of acceptance of a version permanently authorizes you
to choose that version for the Program.
Later license versions may give you additional or different
permissions. However, no additional obligations are imposed on any
author or copyright holder as a result of your choosing to follow a
later version.
15. Disclaimer of Warranty.
THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY
APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT
HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY
OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO,
THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM
IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF
ALL NECESSARY SERVICING, REPAIR OR CORRECTION.
16. Limitation of Liability.
IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING
WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS
THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY
GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE
USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF
DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD
PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS),
EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF
SUCH DAMAGES.
17. Interpretation of Sections 15 and 16.
If the disclaimer of warranty and limitation of liability provided
above cannot be given local legal effect according to their terms,
reviewing courts shall apply local law that most closely approximates
an absolute waiver of all civil liability in connection with the
Program, unless a warranty or assumption of liability accompanies a
copy of the Program in return for a fee.
END OF TERMS AND CONDITIONS
How to Apply These Terms to Your New Programs
If you develop a new program, and you want it to be of the greatest
possible use to the public, the best way to achieve this is to make it
free software which everyone can redistribute and change under these terms.
To do so, attach the following notices to the program. It is safest
to attach them to the start of each source file to most effectively
state the exclusion of warranty; and each file should have at least
the "copyright" line and a pointer to where the full notice is found.
<one line to give the program's name and a brief idea of what it does.>
Copyright (C) <year> <name of author>
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
Also add information on how to contact you by electronic and paper mail.
If the program does terminal interaction, make it output a short
notice like this when it starts in an interactive mode:
<program> Copyright (C) <year> <name of author>
This program comes with ABSOLUTELY NO WARRANTY; for details type `show w'.
This is free software, and you are welcome to redistribute it
under certain conditions; type `show c' for details.
The hypothetical commands `show w' and `show c' should show the appropriate
parts of the General Public License. Of course, your program's commands
might be different; for a GUI interface, you would use an "about box".
You should also get your employer (if you work as a programmer) or school,
if any, to sign a "copyright disclaimer" for the program, if necessary.
For more information on this, and how to apply and follow the GNU GPL, see
<http://www.gnu.org/licenses/>.
The GNU General Public License does not permit incorporating your program
into proprietary programs. If your program is a subroutine library, you
may consider it more useful to permit linking proprietary applications with
the library. If this is what you want to do, use the GNU Lesser General
Public License instead of this License. But first, please read
<http://www.gnu.org/philosophy/why-not-lgpl.html>.

19
src/mcut/LICENSE.txt Normal file
View file

@ -0,0 +1,19 @@
Permission to use this software is governed by two licensing options which are mutually exclusive:
(A) Open source GNU General Public License ("GPL"); a copy of which you should have recieved.
- see also: <http://www.gnu.org/licenses/>
(B) Commercial license, which can be purchased from the owner.
The commercial license option is for users that wish to use MCUT in
their products for comercial purposes but do not wish to release their
software under the GPL.
For more information, email "floyd.m.chitalu@gmail.com".
These options protect the project's commercial value and thus make it possible for the
author to guarantee long term support, maintenance and further development of the
code for the benefit of the project and its users.
This software is also partly based on "CDT" (C++ library for constrained Delaunay triangulation): https://github.com/artem-ogre/CDT
CDT files that were originally under the MPL-2.0 are dual licensed under the MPL-2.0 and the GNU General Public License (GPL) licenses.

112
src/mcut/README.md Normal file
View file

@ -0,0 +1,112 @@
# MCUT Overview
Gist: _A simple and fast C++ library for mesh booleans and more..._
[![Windows](https://github.com/cutdigital/mcut/actions/workflows/windows.yml/badge.svg)](https://github.com/cutdigital/mcut/actions/workflows/windows.yml)
[![MacOS](https://github.com/cutdigital/mcut/actions/workflows/macos.yml/badge.svg)](https://github.com/cutdigital/mcut/actions/workflows/macos.yml) [![Linux](https://github.com/cutdigital/mcut/actions/workflows/linux.yaml/badge.svg)](https://github.com/cutdigital/mcut/actions/workflows/linux.yaml)
This is a software project designed for a broad range of real-world problems relating to 3D modelling and design tasks. Application areas include computer animation, aerospace and automotive engineering, mining, civil and mechanical engineering amongst others.
The project is called "MCUT" (short for 'mesh cutting'), and it provides functionality to perform fast and robust geometry operations, as shown below:
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/github-teaser.png?raw=true">
Figure 1: Generate, slice and perform Booleans without errors.
</p>
The codebase provides a comprehensive tool for ensuring that computer-aided planning tasks for e.g. mine-design, rock strata boring (e.g. underground-tunnel excavations), oil-well drilling and general 3D modelling for animation are achievable with robustness. The tool is developed to take advantage of modern high-performance parallel computing hardware for you, and is demonstrably robust by using precise geometric algorithms that are implemented in C++ and accessed through an intuitive API that resembles the all-familiar C programming language.
Importantly, MCUT is designed with the philosophy that users don't know or don't care about esoteric problems with floating point arithmetic.
# Capabilities
MCUT is a tool for partitioning meshes that represent solids or open surfaces: It is a code library for cutting 3D mesh objects using their geometry to produce crisp fragments at fine scale, which is useful for operations like slicing and boolean operations (union, subtraction and intersection). Supported features include (see images below):
* **Stencilling**: exact cut-outs of the cutting surface
* **Intersection curve access**: geometry representing lines of intersection-contour points
* **Partial cuts**: producing valid results where an open-surface is not necessarily completely cutting through a solid.
* **Concatenation**: merging a solids or open-surfaces with another.
* **Sectioning**: elimination of material/volume on one side of a specified surface (e.g. a plane)
* **Splitting**: partitioning one mesh using another that might be open or solid.
* **Cross-platform**: tested on Windows, Linux (Ubuntu), and macOS
* **Bloat-free**: no external dependencies.
* **Performant**: continuously profiled, and optimized.
* **Numerically robust**: Algorithms rely on robust geometric predicates.
What is being offered is a general solution to the problem of resolving solid- and/or open-mesh intersections. It is a solution that is sought by many companies, researchers, and private individuals for its ability to address extremely difficult problems relating to computational geometry in 3D. A classic application is constructive solid geometry (CSG) i.e. the “boolean operation”, which is shown below, where the resulting meshes/objects are produced with MCUT:
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/teaser2.png?raw=true">
Figure 2: Generate solids / polygons using a robust Boolean engine, where other technologies fail, MCUT solids will be valid.
</p>
# Practical benefits and advantages for users
The expert capabilities of MCUT will allow companies, individuals and researchers-alike to develop robust (and fast) Computer-Aided Design (CAD) and Manufacturing (CAM) tools. For example, these tools could cater to the design of industry-specific structural models like tunnels, drill holes, mechanical instruments and rock-block models. All this alongside the ability to handle general 3D modelling tasks that are typical in industry and academic-fields related to computer graphics (e.g. game-engine level design) and mechanical engineering (e.g. fracture simulation). In essence, users of MCUT are provided with the capability to create robust derivative products and tools for generating (and testing) structural designs in a virtual setting for short- and long-term production operations and feasibility tests/studies.
The following images show more examples of what users can achieve with MCUT:
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/eg-teaser.jpg?raw=true">
Figure 3: Fracture simulation using the Extended Finite Element Method (XFEM) (https://onlinelibrary.wiley.com/doi/abs/10.1111/cgf.13953), where MCUT is used to create fragment geometry by intersecting the simulation domain with propagated cracks.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/image156.png?raw=true">
Figure 4: Intersecting a gear cog with a surface to model the fracturing of steel.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/path1471.png?raw=true">
Figure 5: Merging an engine with the axle shaft to model their connectivity.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/arm-sphere.png?raw=true">
Figure 6: Intersecting a hand model and a sphere, showing how MCUT can also be useful for planning and designing molding processes for e.g. 3D printing.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/arma-bunn.png?raw=true">
Figure 8: Assorted results produced by intersecting the Stanford bunny model and armadillo.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/path1471-2.png?raw=true">
Figure 9: Tunnel excavation of a mountainous terrain for modelling underground construction with a boring machine (represented with cylinder). Note how the input meshes need not be solids.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/image111.png?raw=true">
Figure 10: Creating an Open-Pit mine model on a rough terrain for e.g. pre-planning operations.
</p>
<p align="center">
<img src="https://github.com/cutdigital/mcut.github.io/blob/master/docs/media/repo-teaser/extra-images/path1471-5.png?raw=true">
Figure 11: An example of sectioning with a flat plane, which can be used to eliminate material/volume on either side of this plane or create hollow carve-outs.
</p>
# Source code and test applications
The source code is available for your perusal and evaluation. You can access right here on Github. This is an opportunity for you to trial and experiment with MCUT for your needs. Here is a quick example of how you clone and build the library:
* `git clone https://github.com/cutdigital/mcut.git`
* `mkdir build`
* `cd build`
* `cmake ..` (see `CMakeLists.txt` for available build configuration options)
* run `make -j4` *IF* you are on Linux/MacOS terminal, *ELSE* open the generated `.sln` with e.g. Visual Studio
Next, try out one of the tutorials!
# Licensing
MCUT is available under an Open Source license as well as a commercial license. Users choosing to use MCUT under the free-of-charge Open Source license (e.g. for academic purposes) simply need to comply to its terms, otherwise a commercial license is required. The Open Source license is the "GNU General Public License" (GPL). In cases where the constraints of the Open source license prevent you from using MCUT, a commercial license can be purchased. The library is licensed with an attractively low price which is a one-off sum, requiring no further loyalty fees with guarranteed future updates for free.
These options protect the project's commercial value and thus make it possible for the author to guarantee long term support, maintenance and further development of the code for the benefit of the project and its users.
---
If MCUT helped you please consider adding a star here on GitHub. This means a lot to the author.
_You can also send an [email](floyd.m.chitalu@gmail.com) to the author if you have questions about MCUT_.

7
src/mcut/cmake/add_FetchContent_MakeAvailable.cmake Normal file
View file

@ -0,0 +1,7 @@
macro(FetchContent_MakeAvailable NAME)
FetchContent_GetProperties(${NAME})
if(NOT ${NAME}_POPULATED)
FetchContent_Populate(${NAME})
add_subdirectory(${${NAME}_SOURCE_DIR} ${${NAME}_BINARY_DIR})
endif()
endmacro()

327
src/mcut/include/mcut/internal/bvh.h Normal file
View file

@ -0,0 +1,327 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#ifndef MCUT_BVH_H_
#define MCUT_BVH_H_
#include "mcut/internal/hmesh.h"
#include "mcut/internal/math.h"
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
#include "mcut/internal/tpool.h"
#endif
// OIBVH is over 2-3x faster than the alternative (classic) BVH approaches.
// Our alternative BVH implementations follow: https://www.pbrt.org/chapters/pbrt-2ed-chap4.pdf
#define USE_OIBVH 1
// Expands a 10-bit integer into 30 bits by inserting 2 zeros after each bit.
extern unsigned int expandBits(unsigned int v);
// Calculates a 30-bit Morton code for the given 3D point located within the unit cube [0,1].
extern unsigned int morton3D(float x, float y, float z);
#if defined(USE_OIBVH)
// TODO: just use std::pair
typedef struct
{
int m_left; // node-A ID (implicit index)
int m_right; // node-B ID (implicit index)
} node_pair_t; // collision tree node
// count leading zeros in 32 bit bitfield
extern unsigned int clz(unsigned int x);
// next power of two from x
extern int next_power_of_two(int x);
// check if "x" is a power of two
extern bool is_power_of_two(int x);
// compute log-base-2 of "x"
extern int ilog2(unsigned int x);
// compute index (0...N-1) of the leaf level from the number of leaves
extern int get_leaf_level_from_real_leaf_count(const int t);
// compute tree-level index from implicit index of a node
extern int get_level_from_implicit_idx(const int bvhNodeImplicitIndex);
// compute previous power of two
extern unsigned int flp2(unsigned int x);
// compute size of of Oi-BVH give number of triangles
extern int get_ostensibly_implicit_bvh_size(const int t);
// compute left-most node on a given level
extern int get_level_leftmost_node(const int node_level);
// compute right-most leaf node in tree
extern int get_rightmost_real_leaf(const int bvhLeafLevelIndex, const int num_real_leaf_nodes_in_bvh);
// check if node is a "real node"
extern bool is_real_implicit_tree_node_id(const int bvhNodeImplicitIndex, const int num_real_leaf_nodes_in_bvh);
// get the right most real node on a given tree level
extern int get_level_rightmost_real_node(
const int rightmostRealLeafNodeImplicitIndex,
const int bvhLeafLevelIndex,
const int ancestorLevelIndex);
// compute implicit index of a node's ancestor
extern int get_node_ancestor(
const int nodeImplicitIndex,
const int nodeLevelIndex,
const int ancestorLevelIndex);
// calculate linear memory index of a real node
extern int get_node_mem_index(
const int nodeImplicitIndex,
const int leftmostImplicitIndexOnNodeLevel,
const int bvh_data_base_offset,
const int rightmostRealNodeImplicitIndexOnNodeLevel);
extern void build_oibvh(
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
thread_pool& pool,
#endif
const hmesh_t& mesh,
std::vector<bounding_box_t<vec3>>& bvhAABBs,
std::vector<fd_t>& bvhLeafNodeFaces,
std::vector<bounding_box_t<vec3>>& face_bboxes,
const double& slightEnlargmentEps = double(0.0));
extern void intersectOIBVHs(
std::map<fd_t, std::vector<fd_t>>& ps_face_to_potentially_intersecting_others,
const std::vector<bounding_box_t<vec3>>& srcMeshBvhAABBs,
const std::vector<fd_t>& srcMeshBvhLeafNodeFaces,
const std::vector<bounding_box_t<vec3>>& cutMeshBvhAABBs,
const std::vector<fd_t>& cutMeshBvhLeafNodeFaces);
#else
typedef bounding_box_t<vec3> BBox;
static inline BBox Union(const BBox& a, const BBox& b)
{
BBox out = a;
out.expand(b);
return out;
}
static inline BBox Union(const BBox& a, const vec3& b)
{
BBox out = a;
out.expand(b);
return out;
}
enum class SplitMethod {
SPLIT_SAH = 0,
SPLIT_MIDDLE = 1,
SPLIT_EQUAL_COUNTS = 2,
};
// For each primitive to be stored in the BVH, we store the centroid
// of its bounding box, its complete bounding box, and its index in
// the primitives array
struct BVHPrimitiveInfo {
BVHPrimitiveInfo(int pn, const BBox& b)
: primitiveNumber(pn)
, bounds(b)
{
centroid = b.minimum() * (0.5) + b.maximum() * (0.5);
}
int primitiveNumber;
vec3 centroid;
bounding_box_t<vec3> bounds;
};
// Each BVHBuildNode represents a node of the BVH. All nodes store a BBox, which stores
// the bounds of all of the children beneath the node. Each interior node stores pointers to
// its two children in children. Interior nodes also record the coordinate axis along which
// primitives were sorted for distribution to their two children; this information is used to
// improve the performance of the traversal algorithm. Leaf nodes need to record which
// primitive or primitives are stored in them; the elements of the BVHAccel::primitives
// array from the offset firstPrimOffset up to but not including firstPrimOffset +
// nPrimitives are the primitives in the leaf. (Hence the need for reordering the primi-
// tives array, so that this representation can be used, rather than, for example, storing a
// variable-sized array of primitive indices at each leaf node.)
struct BVHBuildNode {
// The BVHBuildNode constructor only initializes the children pointers; we’ll distinguish
// between leaf and interior nodes by whether their children pointers are NULL or not,
// respectively
BVHBuildNode()
{
children[0] = children[1] = NULL;
}
void InitLeaf(uint32_t first, uint32_t n, const BBox& b)
{
firstPrimOffset = first;
nPrimitives = n;
bounds = b;
}
// The InitInterior() method requires that the two children nodes already have been cre-
// ated, so that their pointers can be passed in. This requirement makes it easy to compute
// the bounds of the interior node, since the children bounds are immediately available.
void InitInterior(uint32_t axis, std::shared_ptr<BVHBuildNode>& c0, std::shared_ptr<BVHBuildNode>& c1)
{
children[0] = c0;
children[1] = c1;
bounds = Union(c0->bounds, c1->bounds);
splitAxis = axis;
nPrimitives = 0;
}
bounding_box_t<vec3> bounds;
std::shared_ptr<BVHBuildNode> children[2];
uint32_t splitAxis, firstPrimOffset, nPrimitives;
};
struct CompareToMid {
CompareToMid(int d, double m)
{
dim = d;
mid = m;
}
int dim;
float mid;
bool operator()(const BVHPrimitiveInfo& a) const
{
return a.centroid[dim] < mid;
}
};
struct ComparePoints {
ComparePoints(int d) { dim = d; }
int dim;
bool operator()(const BVHPrimitiveInfo& a,
const BVHPrimitiveInfo& b) const
{
return a.centroid[dim] < b.centroid[dim];
}
};
struct BucketInfo {
BucketInfo() { count = 0; }
int count;
BBox bounds;
};
struct CompareToBucket {
CompareToBucket(int split, int num, int d, const BBox& b)
: centroidBounds(b)
{
splitBucket = split;
nBuckets = num;
dim = d;
}
// bool operator()(const BVHPrimitiveInfo &p) const;
bool operator()(const BVHPrimitiveInfo& p) const
{
int b = nBuckets * ((p.centroid[dim] - centroidBounds.minimum()[dim]) / (centroidBounds.maximum()[dim] - centroidBounds.minimum()[dim]));
if (b == nBuckets)
b = nBuckets - 1;
return b <= splitBucket;
}
int splitBucket, nBuckets, dim;
const BBox& centroidBounds;
};
// The LinearBVHNode structure stores the information needed to traverse the BVH. In
// addition to the bounding box for each node, for leaf nodes it stores the offset and
// primitive count for the primitives in the node. For interior nodes, it stores the offset to
// the second child as well as which of the coordinate axes the primitives were partitioned
// along when the hierarchy was built; this information is used in the traversal routine below
// to try to visit nodes in front-to-back order along the ray.
struct LinearBVHNode {
BBox bounds;
union {
uint32_t primitivesOffset; // leaf
uint32_t secondChildOffset; // interior
};
uint8_t nPrimitives; // 0 -> interior node
uint8_t axis; // interior node: xyz
uint8_t pad[2]; // ensure 32 byte total size
};
class BoundingVolumeHierarchy {
public:
BoundingVolumeHierarchy();
~BoundingVolumeHierarchy();
// three stages to BVH construction
void buildTree(const hmesh_t& mesh_,
const fixed_precision_number_t& enlargementEps_ = fixed_precision_number_t(0.0),
uint32_t mp_ = 1,
const SplitMethod& sm_ = SplitMethod::SPLIT_MIDDLE);
const BBox& GetPrimitiveBBox(int primitiveIndex) const;
uint32_t flattenBVHTree(std::shared_ptr<BVHBuildNode> node, uint32_t* offset);
// The initial call to recursiveBuild() is given all of the primitives to be stored in
// the tree. It returns a pointer to the root of the tree, which is represented with
// the BVHBuildNode structure.
// responsible for returning a BVH for the subset of primitives represented by the range from
// buildData[start] up to and including buildData[end-1]
std::shared_ptr<BVHBuildNode> recursiveBuild(
std::vector<BVHPrimitiveInfo>& buildData,
uint32_t start,
uint32_t end,
uint32_t* totalNodes,
std::vector<fd_t>& orderedPrims);
int GetNodeCount() const;
const std::shared_ptr<LinearBVHNode>& GetNode(int idx) const;
const fd_t& GetPrimitive(int index) const;
static void intersectBVHTrees(
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
thread_pool& scheduler,
#endif
std::map<fd_t, std::vector<fd_t>>& symmetric_intersecting_pairs,
const BoundingVolumeHierarchy& bvhA,
const BoundingVolumeHierarchy& bvhB,
const uint32_t primitiveOffsetA,
const uint32_t primitiveOffsetB);
private:
const hmesh_t* mesh;
int maxPrimsInNode;
SplitMethod splitMethod;
fixed_precision_number_t enlargementEps; // used to slight enlarge BVH (with bounds of max cut-mesh perturbation magnitude)
std::vector<BVHPrimitiveInfo> buildData;
std::vector<fd_t> primitives; // ordered primitives
std::vector<BBox> primitiveOrderedBBoxes; // unsorted elements correspond to mesh indices
std::vector<std::shared_ptr<LinearBVHNode>> nodes;
};
#endif // #if defined(USE_OIBVH)
#endif // MCUT_BVH_H_

373
src/mcut/include/mcut/internal/cdt/LICENSE.txt Normal file
View file

@ -0,0 +1,373 @@
Mozilla Public License Version 2.0
==================================
1. Definitions
--------------
1.1. "Contributor"
means each individual or legal entity that creates, contributes to
the creation of, or owns Covered Software.
1.2. "Contributor Version"
means the combination of the Contributions of others (if any) used
by a Contributor and that particular Contributor's Contribution.
1.3. "Contribution"
means Covered Software of a particular Contributor.
1.4. "Covered Software"
means Source Code Form to which the initial Contributor has attached
the notice in Exhibit A, the Executable Form of such Source Code
Form, and Modifications of such Source Code Form, in each case
including portions thereof.
1.5. "Incompatible With Secondary Licenses"
means
(a) that the initial Contributor has attached the notice described
in Exhibit B to the Covered Software; or
(b) that the Covered Software was made available under the terms of
version 1.1 or earlier of the License, but not also under the
terms of a Secondary License.
1.6. "Executable Form"
means any form of the work other than Source Code Form.
1.7. "Larger Work"
means a work that combines Covered Software with other material, in
a separate file or files, that is not Covered Software.
1.8. "License"
means this document.
1.9. "Licensable"
means having the right to grant, to the maximum extent possible,
whether at the time of the initial grant or subsequently, any and
all of the rights conveyed by this License.
1.10. "Modifications"
means any of the following:
(a) any file in Source Code Form that results from an addition to,
deletion from, or modification of the contents of Covered
Software; or
(b) any new file in Source Code Form that contains any Covered
Software.
1.11. "Patent Claims" of a Contributor
means any patent claim(s), including without limitation, method,
process, and apparatus claims, in any patent Licensable by such
Contributor that would be infringed, but for the grant of the
License, by the making, using, selling, offering for sale, having
made, import, or transfer of either its Contributions or its
Contributor Version.
1.12. "Secondary License"
means either the GNU General Public License, Version 2.0, the GNU
Lesser General Public License, Version 2.1, the GNU Affero General
Public License, Version 3.0, or any later versions of those
licenses.
1.13. "Source Code Form"
means the form of the work preferred for making modifications.
1.14. "You" (or "Your")
means an individual or a legal entity exercising rights under this
License. For legal entities, "You" includes any entity that
controls, is controlled by, or is under common control with You. For
purposes of this definition, "control" means (a) the power, direct
or indirect, to cause the direction or management of such entity,
whether by contract or otherwise, or (b) ownership of more than
fifty percent (50%) of the outstanding shares or beneficial
ownership of such entity.
2. License Grants and Conditions
--------------------------------
2.1. Grants
Each Contributor hereby grants You a world-wide, royalty-free,
non-exclusive license:
(a) under intellectual property rights (other than patent or trademark)
Licensable by such Contributor to use, reproduce, make available,
modify, display, perform, distribute, and otherwise exploit its
Contributions, either on an unmodified basis, with Modifications, or
as part of a Larger Work; and
(b) under Patent Claims of such Contributor to make, use, sell, offer
for sale, have made, import, and otherwise transfer either its
Contributions or its Contributor Version.
2.2. Effective Date
The licenses granted in Section 2.1 with respect to any Contribution
become effective for each Contribution on the date the Contributor first
distributes such Contribution.
2.3. Limitations on Grant Scope
The licenses granted in this Section 2 are the only rights granted under
this License. No additional rights or licenses will be implied from the
distribution or licensing of Covered Software under this License.
Notwithstanding Section 2.1(b) above, no patent license is granted by a
Contributor:
(a) for any code that a Contributor has removed from Covered Software;
or
(b) for infringements caused by: (i) Your and any other third party's
modifications of Covered Software, or (ii) the combination of its
Contributions with other software (except as part of its Contributor
Version); or
(c) under Patent Claims infringed by Covered Software in the absence of
its Contributions.
This License does not grant any rights in the trademarks, service marks,
or logos of any Contributor (except as may be necessary to comply with
the notice requirements in Section 3.4).
2.4. Subsequent Licenses
No Contributor makes additional grants as a result of Your choice to
distribute the Covered Software under a subsequent version of this
License (see Section 10.2) or under the terms of a Secondary License (if
permitted under the terms of Section 3.3).
2.5. Representation
Each Contributor represents that the Contributor believes its
Contributions are its original creation(s) or it has sufficient rights
to grant the rights to its Contributions conveyed by this License.
2.6. Fair Use
This License is not intended to limit any rights You have under
applicable copyright doctrines of fair use, fair dealing, or other
equivalents.
2.7. Conditions
Sections 3.1, 3.2, 3.3, and 3.4 are conditions of the licenses granted
in Section 2.1.
3. Responsibilities
-------------------
3.1. Distribution of Source Form
All distribution of Covered Software in Source Code Form, including any
Modifications that You create or to which You contribute, must be under
the terms of this License. You must inform recipients that the Source
Code Form of the Covered Software is governed by the terms of this
License, and how they can obtain a copy of this License. You may not
attempt to alter or restrict the recipients' rights in the Source Code
Form.
3.2. Distribution of Executable Form
If You distribute Covered Software in Executable Form then:
(a) such Covered Software must also be made available in Source Code
Form, as described in Section 3.1, and You must inform recipients of
the Executable Form how they can obtain a copy of such Source Code
Form by reasonable means in a timely manner, at a charge no more
than the cost of distribution to the recipient; and
(b) You may distribute such Executable Form under the terms of this
License, or sublicense it under different terms, provided that the
license for the Executable Form does not attempt to limit or alter
the recipients' rights in the Source Code Form under this License.
3.3. Distribution of a Larger Work
You may create and distribute a Larger Work under terms of Your choice,
provided that You also comply with the requirements of this License for
the Covered Software. If the Larger Work is a combination of Covered
Software with a work governed by one or more Secondary Licenses, and the
Covered Software is not Incompatible With Secondary Licenses, this
License permits You to additionally distribute such Covered Software
under the terms of such Secondary License(s), so that the recipient of
the Larger Work may, at their option, further distribute the Covered
Software under the terms of either this License or such Secondary
License(s).
3.4. Notices
You may not remove or alter the substance of any license notices
(including copyright notices, patent notices, disclaimers of warranty,
or limitations of liability) contained within the Source Code Form of
the Covered Software, except that You may alter any license notices to
the extent required to remedy known factual inaccuracies.
3.5. Application of Additional Terms
You may choose to offer, and to charge a fee for, warranty, support,
indemnity or liability obligations to one or more recipients of Covered
Software. However, You may do so only on Your own behalf, and not on
behalf of any Contributor. You must make it absolutely clear that any
such warranty, support, indemnity, or liability obligation is offered by
You alone, and You hereby agree to indemnify every Contributor for any
liability incurred by such Contributor as a result of warranty, support,
indemnity or liability terms You offer. You may include additional
disclaimers of warranty and limitations of liability specific to any
jurisdiction.
4. Inability to Comply Due to Statute or Regulation
---------------------------------------------------
If it is impossible for You to comply with any of the terms of this
License with respect to some or all of the Covered Software due to
statute, judicial order, or regulation then You must: (a) comply with
the terms of this License to the maximum extent possible; and (b)
describe the limitations and the code they affect. Such description must
be placed in a text file included with all distributions of the Covered
Software under this License. Except to the extent prohibited by statute
or regulation, such description must be sufficiently detailed for a
recipient of ordinary skill to be able to understand it.
5. Termination
--------------
5.1. The rights granted under this License will terminate automatically
if You fail to comply with any of its terms. However, if You become
compliant, then the rights granted under this License from a particular
Contributor are reinstated (a) provisionally, unless and until such
Contributor explicitly and finally terminates Your grants, and (b) on an
ongoing basis, if such Contributor fails to notify You of the
non-compliance by some reasonable means prior to 60 days after You have
come back into compliance. Moreover, Your grants from a particular
Contributor are reinstated on an ongoing basis if such Contributor
notifies You of the non-compliance by some reasonable means, this is the
first time You have received notice of non-compliance with this License
from such Contributor, and You become compliant prior to 30 days after
Your receipt of the notice.
5.2. If You initiate litigation against any entity by asserting a patent
infringement claim (excluding declaratory judgment actions,
counter-claims, and cross-claims) alleging that a Contributor Version
directly or indirectly infringes any patent, then the rights granted to
You by any and all Contributors for the Covered Software under Section
2.1 of this License shall terminate.
5.3. In the event of termination under Sections 5.1 or 5.2 above, all
end user license agreements (excluding distributors and resellers) which
have been validly granted by You or Your distributors under this License
prior to termination shall survive termination.
************************************************************************
* *
* 6. Disclaimer of Warranty *
* ------------------------- *
* *
* Covered Software is provided under this License on an "as is" *
* basis, without warranty of any kind, either expressed, implied, or *
* statutory, including, without limitation, warranties that the *
* Covered Software is free of defects, merchantable, fit for a *
* particular purpose or non-infringing. The entire risk as to the *
* quality and performance of the Covered Software is with You. *
* Should any Covered Software prove defective in any respect, You *
* (not any Contributor) assume the cost of any necessary servicing, *
* repair, or correction. This disclaimer of warranty constitutes an *
* essential part of this License. No use of any Covered Software is *
* authorized under this License except under this disclaimer. *
* *
************************************************************************
************************************************************************
* *
* 7. Limitation of Liability *
* -------------------------- *
* *
* Under no circumstances and under no legal theory, whether tort *
* (including negligence), contract, or otherwise, shall any *
* Contributor, or anyone who distributes Covered Software as *
* permitted above, be liable to You for any direct, indirect, *
* special, incidental, or consequential damages of any character *
* including, without limitation, damages for lost profits, loss of *
* goodwill, work stoppage, computer failure or malfunction, or any *
* and all other commercial damages or losses, even if such party *
* shall have been informed of the possibility of such damages. This *
* limitation of liability shall not apply to liability for death or *
* personal injury resulting from such party's negligence to the *
* extent applicable law prohibits such limitation. Some *
* jurisdictions do not allow the exclusion or limitation of *
* incidental or consequential damages, so this exclusion and *
* limitation may not apply to You. *
* *
************************************************************************
8. Litigation
-------------
Any litigation relating to this License may be brought only in the
courts of a jurisdiction where the defendant maintains its principal
place of business and such litigation shall be governed by laws of that
jurisdiction, without reference to its conflict-of-law provisions.
Nothing in this Section shall prevent a party's ability to bring
cross-claims or counter-claims.
9. Miscellaneous
----------------
This License represents the complete agreement concerning the subject
matter hereof. If any provision of this License is held to be
unenforceable, such provision shall be reformed only to the extent
necessary to make it enforceable. Any law or regulation which provides
that the language of a contract shall be construed against the drafter
shall not be used to construe this License against a Contributor.
10. Versions of the License
---------------------------
10.1. New Versions
Mozilla Foundation is the license steward. Except as provided in Section
10.3, no one other than the license steward has the right to modify or
publish new versions of this License. Each version will be given a
distinguishing version number.
10.2. Effect of New Versions
You may distribute the Covered Software under the terms of the version
of the License under which You originally received the Covered Software,
or under the terms of any subsequent version published by the license
steward.
10.3. Modified Versions
If you create software not governed by this License, and you want to
create a new license for such software, you may create and use a
modified version of this License if you rename the license and remove
any references to the name of the license steward (except to note that
such modified license differs from this License).
10.4. Distributing Source Code Form that is Incompatible With Secondary
Licenses
If You choose to distribute Source Code Form that is Incompatible With
Secondary Licenses under the terms of this version of the License, the
notice described in Exhibit B of this License must be attached.
Exhibit A - Source Code Form License Notice
-------------------------------------------
This Source Code Form is subject to the terms of the Mozilla Public
License, v. 2.0. If a copy of the MPL was not distributed with this
file, You can obtain one at http://mozilla.org/MPL/2.0/.
If it is not possible or desirable to put the notice in a particular
file, then You may include the notice in a location (such as a LICENSE
file in a relevant directory) where a recipient would be likely to look
for such a notice.
You may add additional accurate notices of copyright ownership.
Exhibit B - "Incompatible With Secondary Licenses" Notice
---------------------------------------------------------
This Source Code Form is "Incompatible With Secondary Licenses", as
defined by the Mozilla Public License, v. 2.0.

555
src/mcut/include/mcut/internal/cdt/cdt.h Normal file
View file

@ -0,0 +1,555 @@
/* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at https://mozilla.org/MPL/2.0/. */
#ifndef _CONSTRAINED_DELAUNAY_TRIANGULATION_H_
#define _CONSTRAINED_DELAUNAY_TRIANGULATION_H_
#include "mcut/internal/cdt/triangulate.h"
#include "mcut/internal/cdt/utils.h"
#include <algorithm>
#include <cassert>
#include <cstdlib>
#include <iterator>
#include <memory>
#include <stack>
#include <vector>
/// Namespace containing triangulation functionality
namespace cdt {
/** @defgroup API Public API
* Contains API for constrained and conforming Delaunay triangulations
*/
/// @{
/**
* Type used for storing layer depths for triangles
* @note layer_depth_t should support 60K+ layers, which could be to much or
* too little for some use cases. Feel free to re-define this typedef.
*/
typedef unsigned short layer_depth_t;
typedef layer_depth_t boundary_overlap_count_t;
/** @defgroup helpers Helpers
* Helpers for working with cdt::triangulator_t.
*/
/// @{
/**
* Calculate triangles adjacent to vertices (triangles by vertex index)
* @param triangles triangulation
* @param verticesSize total number of vertices to pre-allocate the output
* @return triangles by vertex index (array of size V, where V is the number of vertices; each element is a list of triangles incident to corresponding vertex).
*/
inline std::vector<std::vector<std::uint32_t>> get_vertex_to_triangles_map(
const std::vector<triangle_t>& triangles,
const std::uint32_t verticesSize)
{
std::vector<std::vector<std::uint32_t>> map(verticesSize);
for (std::uint32_t i = 0; i < (std::uint32_t)triangles.size(); ++i) {
const std::uint32_t triangle_index = i;
const triangle_t& triangle = triangles[triangle_index];
const std::array<std::uint32_t, 3>& vertices = triangle.vertices;
for (std::array<std::uint32_t, 3>::const_iterator j = vertices.begin(); j != vertices.end(); ++j) {
const std::uint32_t vertex_index = *j;
map[vertex_index].push_back(i);
}
}
return map;
}
/**
* Information about removed duplicated vertices.
*
* Contains mapping information and removed duplicates indices.
* @note vertices {0,1,2,3,4} where 0 and 3 are the same will produce mapping
* {0,1,2,0,3} (to new vertices {0,1,2,3}) and duplicates {3}
*/
struct duplicates_info_t {
std::vector<std::size_t> mapping; ///< vertex index mapping
std::vector<std::size_t> duplicates; ///< duplicates' indices
};
/*!
* Remove elements in the range [first; last) with indices from the sorted
* unique range [ii_first, ii_last)
*/
template <class ForwardIt, class SortUniqIndsFwdIt>
inline ForwardIt remove_at(
ForwardIt first,
ForwardIt last,
SortUniqIndsFwdIt ii_first,
SortUniqIndsFwdIt ii_last)
{
if (ii_first == ii_last) // no indices-to-remove are given
return last;
typedef typename std::iterator_traits<ForwardIt>::difference_type diff_t;
typedef typename std::iterator_traits<SortUniqIndsFwdIt>::value_type ind_t;
ForwardIt destination = first + static_cast<diff_t>(*ii_first);
while (ii_first != ii_last) {
// advance to an index after a chunk of elements-to-keep
for (ind_t cur = *ii_first++; ii_first != ii_last; ++ii_first) {
const ind_t nxt = *ii_first;
if (nxt - cur > 1)
break;
cur = nxt;
}
// move the chunk of elements-to-keep to new destination
const ForwardIt source_first = first + static_cast<diff_t>(*(ii_first - 1)) + 1;
const ForwardIt source_last = ii_first != ii_last ? first + static_cast<diff_t>(*ii_first) : last;
std::move(source_first, source_last, destination);
destination += source_last - source_first;
}
return destination;
}
/**
* Find duplicates in given custom point-type range
* @note duplicates are points with exactly same X and Y coordinates
* @tparam TVertexIter iterator that dereferences to custom point type
* @tparam TGetVertexCoordX function object getting x coordinate from vertex.
* Getter signature: const TVertexIter::value_type& -> T
* @tparam TGetVertexCoordY function object getting y coordinate from vertex.
* Getter signature: const TVertexIter::value_type& -> T
* @param first beginning of the range of vertices
* @param last end of the range of vertices
* @param get_x_coord getter of X-coordinate
* @param get_y_coord getter of Y-coordinate
* @returns information about vertex duplicates
*/
template <
typename T,
typename TVertexIter,
typename TGetVertexCoordX,
typename TGetVertexCoordY>
duplicates_info_t find_duplicates(
TVertexIter first,
TVertexIter last,
TGetVertexCoordX get_x_coord,
TGetVertexCoordY get_y_coord)
{
std::unordered_map<vec2_<T>, std::size_t> uniqueVerts; // position to index map
const std::size_t verticesSize = std::distance(first, last);
duplicates_info_t di = {
std::vector<std::size_t>(verticesSize),
std::vector<std::size_t>()
};
for (std::size_t iIn = 0, iOut = iIn; iIn < verticesSize; ++iIn, ++first) {
typename std::unordered_map<vec2_<T>, std::size_t>::const_iterator it;
bool isUnique;
// check if the coordinates match [exactly] (in the bitwise sense)
std::tie(it, isUnique) = uniqueVerts.insert(
std::make_pair(
vec2_<T>::make(get_x_coord(*first), get_y_coord(*first)),
iOut));
if (isUnique) {
di.mapping[iIn] = iOut++;
continue;
}
di.mapping[iIn] = it->second; // found a duplicate
di.duplicates.push_back(iIn);
}
return di;
}
/**
* Remove duplicates in-place from vector of custom points
* @tparam TVertex vertex type
* @tparam TAllocator allocator used by input vector of vertices
* @param vertices vertices to remove duplicates from
* @param duplicates information about duplicates
*/
template <typename TVertex, typename TAllocator>
void remove_duplicates(
std::vector<TVertex, TAllocator>& vertices,
const std::vector<std::size_t>& duplicates)
{
vertices.erase(
remove_at(
vertices.begin(),
vertices.end(),
duplicates.begin(),
duplicates.end()),
vertices.end());
}
/**
* Remove duplicated points in-place
*
* @tparam T type of vertex coordinates (e.g., float, double)
* @param[in, out] vertices collection of vertices to remove duplicates from
* @returns information about duplicated vertices that were removed.
*/
template <typename T>
duplicates_info_t remove_duplicates(std::vector<vec2_<T>>& vertices)
{
const duplicates_info_t di = find_duplicates<T>(
vertices.begin(),
vertices.end(),
get_x_coord_vec2d<T>,
get_y_coord_vec2d<T>);
remove_duplicates(vertices, di.duplicates);
return di;
}
/**
* Remap vertex indices in edges (in-place) using given vertex-index mapping.
* @tparam TEdgeIter iterator that dereferences to custom edge type
* @tparam TGetEdgeVertexStart function object getting start vertex index
* from an edge.
* Getter signature: const TEdgeIter::value_type& -> std::uint32_t
* @tparam TGetEdgeVertexEnd function object getting end vertex index from
* an edge. Getter signature: const TEdgeIter::value_type& -> std::uint32_t
* @tparam TMakeEdgeFromStartAndEnd function object that makes new edge from
* start and end vertices
* @param first beginning of the range of edges
* @param last end of the range of edges
* @param mapping vertex-index mapping
* @param getStart getter of edge start vertex index
* @param getEnd getter of edge end vertex index
* @param makeEdge factory for making edge from vetices
*/
template <
typename TEdgeIter,
typename TGetEdgeVertexStart,
typename TGetEdgeVertexEnd,
typename TMakeEdgeFromStartAndEnd>
void remap_edges(
TEdgeIter first,
const TEdgeIter last,
const std::vector<std::size_t>& mapping,
TGetEdgeVertexStart getStart,
TGetEdgeVertexEnd getEnd,
TMakeEdgeFromStartAndEnd makeEdge)
{
for (; first != last; ++first) {
*first = makeEdge(
static_cast<std::uint32_t>(mapping[getStart(*first)]),
static_cast<std::uint32_t>(mapping[getEnd(*first)]));
}
}
/**
* Remap vertex indices in edges (in-place) using given vertex-index mapping.
*
* @note Mapping can be a result of remove_duplicates function
* @param[in,out] edges collection of edges to remap
* @param mapping vertex-index mapping
*/
inline void
remap_edges(std::vector<edge_t>& edges, const std::vector<std::size_t>& mapping)
{
remap_edges(
edges.begin(),
edges.end(),
mapping,
edge_get_v1,
edge_get_v2,
edge_make);
}
/**
* Find point duplicates, remove them from vector (in-place) and remap edges
* (in-place)
* @note Same as a chained call of cdt::find_duplicates, cdt::remove_duplicates,
* and cdt::remap_edges
* @tparam T type of vertex coordinates (e.g., float, double)
* @tparam TVertex type of vertex
* @tparam TGetVertexCoordX function object getting x coordinate from vertex.
* Getter signature: const TVertexIter::value_type& -> T
* @tparam TGetVertexCoordY function object getting y coordinate from vertex.
* Getter signature: const TVertexIter::value_type& -> T
* @tparam TEdgeIter iterator that dereferences to custom edge type
* @tparam TGetEdgeVertexStart function object getting start vertex index
* from an edge.
* Getter signature: const TEdgeIter::value_type& -> std::uint32_t
* @tparam TGetEdgeVertexEnd function object getting end vertex index from
* an edge. Getter signature: const TEdgeIter::value_type& -> std::uint32_t
* @tparam TMakeEdgeFromStartAndEnd function object that makes new edge from
* start and end vertices
* @param[in, out] vertices vertices to remove duplicates from
* @param[in, out] edges collection of edges connecting vertices
* @param get_x_coord getter of X-coordinate
* @param get_y_coord getter of Y-coordinate
* @param edgesFirst beginning of the range of edges
* @param edgesLast end of the range of edges
* @param getStart getter of edge start vertex index
* @param getEnd getter of edge end vertex index
* @param makeEdge factory for making edge from vetices
* @returns information about vertex duplicates
*/
template <
typename T,
typename TVertex,
typename TGetVertexCoordX,
typename TGetVertexCoordY,
typename TVertexAllocator,
typename TEdgeIter,
typename TGetEdgeVertexStart,
typename TGetEdgeVertexEnd,
typename TMakeEdgeFromStartAndEnd>
duplicates_info_t remove_duplicates_and_remap_edges(
std::vector<TVertex, TVertexAllocator>& vertices,
TGetVertexCoordX get_x_coord,
TGetVertexCoordY get_y_coord,
const TEdgeIter edgesFirst,
const TEdgeIter edgesLast,
TGetEdgeVertexStart getStart,
TGetEdgeVertexEnd getEnd,
TMakeEdgeFromStartAndEnd makeEdge)
{
const duplicates_info_t di = find_duplicates<T>(vertices.begin(), vertices.end(), get_x_coord, get_y_coord);
remove_duplicates(vertices, di.duplicates);
remap_edges(edgesFirst, edgesLast, di.mapping, getStart, getEnd, makeEdge);
return di;
}
/**
* Same as a chained call of cdt::remove_duplicates + cdt::remap_edges
*
* @tparam T type of vertex coordinates (e.g., float, double)
* @param[in, out] vertices collection of vertices to remove duplicates from
* @param[in,out] edges collection of edges to remap
*/
template <typename T>
duplicates_info_t remove_duplicates_and_remap_edges(
std::vector<vec2_<T>>& vertices,
std::vector<edge_t>& edges)
{
return remove_duplicates_and_remap_edges<T>(
vertices,
get_x_coord_vec2d<T>,
get_y_coord_vec2d<T>,
edges.begin(),
edges.end(),
edge_get_v1,
edge_get_v2,
edge_make);
}
/**
* Extract all edges of triangles
*
* @param triangles triangles used to extract edges
* @return an unordered set of all edges of triangulation
*/
inline std::unordered_set<edge_t>
extract_edges_from_triangles(const std::vector<triangle_t>& triangles)
{
std::unordered_set<edge_t> edges;
for (std::vector<triangle_t>::const_iterator t = triangles.begin(); t != triangles.end(); ++t) {
edges.insert(edge_t(std::uint32_t(t->vertices[0]), std::uint32_t(t->vertices[1])));
edges.insert(edge_t(std::uint32_t(t->vertices[1]), std::uint32_t(t->vertices[2])));
edges.insert(edge_t(std::uint32_t(t->vertices[2]), std::uint32_t(t->vertices[0])));
}
return edges;
}
/*!
* Converts piece->original_edges mapping to original_edge->pieces
* @param pieceToOriginals maps pieces to original edges
* @return mapping of original edges to pieces
*/
inline std::unordered_map<edge_t, std::vector<edge_t>>
edge_to_pieces_mapping(const std::unordered_map<edge_t, std::vector<edge_t>>& pieceToOriginals)
{
std::unordered_map<edge_t, std::vector<edge_t>> originalToPieces;
typedef std::unordered_map<edge_t, std::vector<edge_t>>::const_iterator Cit;
for (Cit ptoIt = pieceToOriginals.begin(); ptoIt != pieceToOriginals.end();
++ptoIt) {
const edge_t piece = ptoIt->first;
const std::vector<edge_t>& originals = ptoIt->second;
for (std::vector<edge_t>::const_iterator origIt = originals.begin();
origIt != originals.end();
++origIt) {
originalToPieces[*origIt].push_back(piece);
}
}
return originalToPieces;
}
/*!
* Convert edge-to-pieces mapping into edge-to-split-vertices mapping
* @tparam T type of vertex coordinates (e.g., float, double)
* @param edgeToPieces edge-to-pieces mapping
* @param vertices vertex buffer
* @return mapping of edge-to-split-points.
* Split points are sorted from edge's start (v1) to end (v2)
*/
template <typename T>
std::unordered_map<edge_t, std::vector<std::uint32_t>> get_edge_to_split_vertices_map(
const std::unordered_map<edge_t, std::vector<edge_t>>& edgeToPieces,
const std::vector<vec2_<T>>& vertices)
{
typedef std::pair<std::uint32_t, T> VertCoordPair;
struct ComparePred {
bool operator()(const VertCoordPair& a, const VertCoordPair& b) const
{
return a.second < b.second;
}
} comparePred;
std::unordered_map<edge_t, std::vector<std::uint32_t>> edgeToSplitVerts;
for (std::unordered_map<edge_t, std::vector<edge_t>>::const_iterator it = edgeToPieces.begin();
it != edgeToPieces.end();
++it) {
const edge_t& e = it->first;
const T dX = vertices[e.v2()].x() - vertices[e.v1()].x();
const T dY = vertices[e.v2()].y() - vertices[e.v1()].y();
const bool isX = std::abs(dX) >= std::abs(dY); // X-coord longer
const bool isAscending = isX ? dX >= 0 : dY >= 0; // Longer coordinate ascends
const std::vector<edge_t>& pieces = it->second;
std::vector<VertCoordPair> splitVerts;
// size is: 2[ends] + (pieces - 1)[split vertices] = pieces + 1
splitVerts.reserve(pieces.size() + 1);
for (std::vector<edge_t>::const_iterator it = pieces.begin(); it != pieces.end(); ++it) {
const std::array<std::uint32_t, 2> vv = { it->v1(), it->v2() };
for (std::array<std::uint32_t, 2>::const_iterator v = vv.begin(); v != vv.end(); ++v) {
const T c = isX ? vertices[*v].x() : vertices[*v].y();
splitVerts.push_back(std::make_pair(*v, isAscending ? c : -c));
}
}
// sort by longest coordinate
std::sort(splitVerts.begin(), splitVerts.end(), comparePred);
// remove duplicates
splitVerts.erase(
std::unique(splitVerts.begin(), splitVerts.end()),
splitVerts.end());
MCUT_ASSERT(splitVerts.size() > 2); // 2 end points with split vertices
std::pair<edge_t, std::vector<std::uint32_t>> val = std::make_pair(e, std::vector<std::uint32_t>());
val.second.reserve(splitVerts.size());
for (typename std::vector<VertCoordPair>::const_iterator it = splitVerts.begin() + 1;
it != splitVerts.end() - 1;
++it) {
val.second.push_back(it->first);
}
edgeToSplitVerts.insert(val);
}
return edgeToSplitVerts;
}
/// @}
/// @}
} // namespace cdt
//*****************************************************************************
// Implementations of template functionlity
//*****************************************************************************
// hash for vec2_<T>
namespace std {
template <typename T>
struct hash<vec2_<T>> {
size_t operator()(const vec2_<T>& xy) const
{
return std::hash<T>()(xy.x()) ^ std::hash<T>()(xy.y());
}
};
} // namespace std
namespace cdt {
//-----
// API
//-----
/**
* Verify that triangulation topology is consistent.
*
* Checks:
* - for each vertex adjacent triangles contain the vertex
* - each triangle's neighbor in turn has triangle as its neighbor
* - each of triangle's vertices has triangle as adjacent
*
* @tparam T type of vertex coordinates (e.g., float, double)
* @tparam TNearPointLocator class providing locating near point for efficiently
*/
template <typename T, typename TNearPointLocator>
inline bool check_topology(const cdt::triangulator_t<T, TNearPointLocator>& cdt)
{
// Check if vertices' adjacent triangles contain vertex
const std::vector<std::vector<std::uint32_t>> vertTris = cdt.is_finalized()
? get_vertex_to_triangles_map(
cdt.triangles, static_cast<std::uint32_t>(cdt.vertices.size()))
: cdt.vertTris;
for (std::uint32_t iV(0); iV < std::uint32_t(cdt.vertices.size()); ++iV) {
const std::vector<std::uint32_t>& vTris = vertTris[iV];
typedef std::vector<std::uint32_t>::const_iterator TriIndCit;
for (TriIndCit it = vTris.begin(); it != vTris.end(); ++it) {
const std::array<std::uint32_t, 3>& vv = cdt.triangles[*it].vertices;
if (std::find(vv.begin(), vv.end(), iV) == vv.end())
return false;
}
}
// Check if triangle neighbor links are fine
for (std::uint32_t iT(0); iT < std::uint32_t(cdt.triangles.size()); ++iT) {
const triangle_t& t = cdt.triangles[iT];
typedef std::array<std::uint32_t, 3>::const_iterator NCit;
for (NCit it = t.neighbors.begin(); it != t.neighbors.end(); ++it) {
if (*it == null_neighbour)
continue;
const std::array<std::uint32_t, 3>& nn = cdt.triangles[*it].neighbors;
if (std::find(nn.begin(), nn.end(), iT) == nn.end())
return false;
}
}
// Check if triangle's vertices have triangle as adjacent
for (std::uint32_t iT(0); iT < std::uint32_t(cdt.triangles.size()); ++iT) {
const triangle_t& t = cdt.triangles[iT];
typedef std::array<std::uint32_t, 3>::const_iterator VCit;
for (VCit it = t.vertices.begin(); it != t.vertices.end(); ++it) {
const std::vector<std::uint32_t>& tt = vertTris[*it];
if (std::find(tt.begin(), tt.end(), iT) == tt.end())
return false;
}
}
return true;
}
} // namespace cdt
#endif // #ifndef _CONSTRAINED_DELAUNAY_TRIANGULATION_H_

421
src/mcut/include/mcut/internal/cdt/kdtree.h Normal file
View file

@ -0,0 +1,421 @@
/* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at https://mozilla.org/MPL/2.0/. */
#ifndef KDTREE_KDTREE_H
#define KDTREE_KDTREE_H
#include "mcut/internal/cdt/utils.h"
#include "mcut/internal/math.h"
#include <cassert>
#include <limits>
namespace kdt {
struct NodeSplitDirection {
enum Enum {
X,
Y,
};
};
/// Simple tree structure with alternating half splitting nodes
/// @details Simple tree structure
/// - Tree to incrementally add points to the structure.
/// - Get the nearest point to a given input.
/// - Does not check for duplicates, expect unique points.
/// @tparam TCoordType type used for storing point coordinate.
/// @tparam NumVerticesInLeaf The number of points per leaf.
/// @tparam InitialStackDepth initial size of stack depth for nearest query.
/// Should be at least 1.
/// @tparam StackDepthIncrement increment of stack depth for nearest query when
/// stack depth is reached.
template <
typename TCoordType,
size_t NumVerticesInLeaf,
size_t InitialStackDepth,
size_t StackDepthIncrement>
class kdtree_t_ {
public:
typedef TCoordType coord_type;
typedef vec2_<coord_type> point_type;
typedef std::pair<point_type, std::uint32_t> value_type;
typedef std::vector<std::uint32_t> point_data_vec;
typedef point_data_vec::const_iterator pd_cit;
typedef std::array<std::uint32_t, 2> children_type;
/// Stores kd-tree node data
struct Node {
children_type children; ///< two children if not leaf; {0,0} if leaf
point_data_vec data; ///< points' data if leaf
/// Create empty leaf
Node()
{
setChildren(0, 0);
data.reserve(NumVerticesInLeaf);
}
/// Children setter for convenience
void setChildren(const std::uint32_t c1, const std::uint32_t c2)
{
children[0] = c1;
children[1] = c2;
}
/// Check if node is a leaf (has no valid children)
bool isLeaf() const
{
return children[0] == children[1];
}
};
/// Default constructor
kdtree_t_()
: m_rootDir(NodeSplitDirection::X)
, m_min(point_type::make(
-std::numeric_limits<coord_type>::max(),
-std::numeric_limits<coord_type>::max()))
, m_max(point_type::make(
std::numeric_limits<coord_type>::max(),
std::numeric_limits<coord_type>::max()))
, m_isRootBoxInitialized(false)
, m_tasksStack(InitialStackDepth, NearestTask())
{
m_root = addNewNode();
}
/// Constructor with bounding box known in advance
kdtree_t_(const point_type& min, const point_type& max)
: m_rootDir(NodeSplitDirection::X)
, m_min(min)
, m_max(max)
, m_isRootBoxInitialized(true)
, m_tasksStack(InitialStackDepth, NearestTask())
{
m_root = addNewNode();
}
/// Insert a point into kd-tree
/// @note external point-buffer is used to reduce kd-tree's memory footprint
/// @param iPoint index of point in external point-buffer
/// @param points external point-buffer
void
insert(const std::uint32_t& iPoint, const std::vector<point_type>& points)
{
// if point is outside root, extend tree by adding new roots
const point_type& pos = points[iPoint];
while (!isInsideBox(pos, m_min, m_max)) {
extendTree(pos);
}
// now insert the point into the tree
std::uint32_t node = m_root;
point_type min = m_min;
point_type max = m_max;
NodeSplitDirection::Enum dir = m_rootDir;
// below: initialized only to suppress warnings
NodeSplitDirection::Enum newDir(NodeSplitDirection::X);
coord_type mid(0);
point_type newMin, newMax;
while (true) {
if (m_nodes[node].isLeaf()) {
// add point if capacity is not reached
point_data_vec& pd = m_nodes[node].data;
if (pd.size() < NumVerticesInLeaf) {
pd.push_back(iPoint);
return;
}
// initialize bbox first time the root capacity is reached
if (!m_isRootBoxInitialized) {
initializeRootBox(points);
min = m_min;
max = m_max;
}
// split a full leaf node
calcSplitInfo(min, max, dir, mid, newDir, newMin, newMax);
const std::uint32_t c1 = addNewNode(), c2 = addNewNode();
Node& n = m_nodes[node];
n.setChildren(c1, c2);
point_data_vec& c1data = m_nodes[c1].data;
point_data_vec& c2data = m_nodes[c2].data;
// move node's points to children
for (pd_cit it = n.data.begin(); it != n.data.end(); ++it) {
whichChild(points[*it], mid, dir) == 0
? c1data.push_back(*it)
: c2data.push_back(*it);
}
n.data = point_data_vec();
} else {
calcSplitInfo(min, max, dir, mid, newDir, newMin, newMax);
}
// add the point to a child
const std::size_t iChild = whichChild(points[iPoint], mid, dir);
iChild == 0 ? max = newMax : min = newMin;
node = m_nodes[node].children[iChild];
dir = newDir;
}
}
/// Query kd-tree for a nearest neighbor point
/// @note external point-buffer is used to reduce kd-tree's memory footprint
/// @param point query point position
/// @param points external point-buffer
value_type nearest(
const point_type& point,
const std::vector<point_type>& points) const
{
value_type out;
int iTask = -1;
coord_type minDistSq = std::numeric_limits<coord_type>::max();
m_tasksStack[++iTask] = NearestTask(m_root, m_min, m_max, m_rootDir, minDistSq);
while (iTask != -1) {
const NearestTask t = m_tasksStack[iTask--];
if (t.distSq > minDistSq)
continue;
const Node& n = m_nodes[t.node];
if (n.isLeaf()) {
for (pd_cit it = n.data.begin(); it != n.data.end(); ++it) {
const point_type& p = points[*it];
const coord_type distSq = cdt::get_square_distance(point, p);
if (distSq < minDistSq) {
minDistSq = distSq;
out.first = p;
out.second = *it;
}
}
} else {
coord_type mid(0);
NodeSplitDirection::Enum newDir;
point_type newMin, newMax;
calcSplitInfo(t.min, t.max, t.dir, mid, newDir, newMin, newMax);
const coord_type distToMid = t.dir == NodeSplitDirection::X
? (point.x() - mid)
: (point.y() - mid);
const coord_type toMidSq = distToMid * distToMid;
const std::size_t iChild = whichChild(point, mid, t.dir);
if (iTask + 2 >= static_cast<int>(m_tasksStack.size())) {
m_tasksStack.resize(
m_tasksStack.size() + StackDepthIncrement);
}
// node containing point should end up on top of the stack
if (iChild == 0) {
m_tasksStack[++iTask] = NearestTask(
n.children[1], newMin, t.max, newDir, toMidSq);
m_tasksStack[++iTask] = NearestTask(
n.children[0], t.min, newMax, newDir, toMidSq);
} else {
m_tasksStack[++iTask] = NearestTask(
n.children[0], t.min, newMax, newDir, toMidSq);
m_tasksStack[++iTask] = NearestTask(
n.children[1], newMin, t.max, newDir, toMidSq);
}
}
}
return out;
}
private:
/// Add a new node and return it's index in nodes buffer
std::uint32_t addNewNode()
{
const std::uint32_t newNodeIndex = static_cast<std::uint32_t>(m_nodes.size());
m_nodes.push_back(Node());
return newNodeIndex;
}
/// Test which child point belongs to after the split
/// @returns 0 if first child, 1 if second child
std::size_t whichChild(
const point_type& point,
const coord_type& split,
const NodeSplitDirection::Enum dir) const
{
return static_cast<size_t>(
dir == NodeSplitDirection::X ? point.x() > split : point.y() > split);
}
/// Calculate split location, direction, and children boxes
static void calcSplitInfo(
const point_type& min,
const point_type& max,
const NodeSplitDirection::Enum dir,
coord_type& midOut,
NodeSplitDirection::Enum& newDirOut,
point_type& newMinOut,
point_type& newMaxOut)
{
newMaxOut = max;
newMinOut = min;
switch (dir) {
case NodeSplitDirection::X:
midOut = (min.x() + max.x()) / coord_type(2);
newDirOut = NodeSplitDirection::Y;
newMinOut.x() = midOut;
newMaxOut.x() = midOut;
return;
case NodeSplitDirection::Y:
midOut = (min.y() + max.y()) / coord_type(2);
newDirOut = NodeSplitDirection::X;
newMinOut.y() = midOut;
newMaxOut.y() = midOut;
return;
}
}
/// Test if point is inside a box
static bool isInsideBox(
const point_type& p,
const point_type& min,
const point_type& max)
{
return p.x() >= min.x() && p.x() <= max.x() && p.y() >= min.y() && p.y() <= max.y();
}
/// Extend a tree by creating new root with old root and a new node as
/// children
void extendTree(const point_type& point)
{
const std::uint32_t newRoot = addNewNode();
const std::uint32_t newLeaf = addNewNode();
switch (m_rootDir) {
case NodeSplitDirection::X:
m_rootDir = NodeSplitDirection::Y;
point.y() < m_min.y() ? m_nodes[newRoot].setChildren(newLeaf, m_root)
: m_nodes[newRoot].setChildren(m_root, newLeaf);
if (point.y() < m_min.y())
m_min.y() -= m_max.y() - m_min.y();
else if (point.y() > m_max.y())
m_max.y() += m_max.y() - m_min.y();
break;
case NodeSplitDirection::Y:
m_rootDir = NodeSplitDirection::X;
point.x() < m_min.x() ? m_nodes[newRoot].setChildren(newLeaf, m_root)
: m_nodes[newRoot].setChildren(m_root, newLeaf);
if (point.x() < m_min.x())
m_min.x() -= m_max.x() - m_min.x();
else if (point.x() > m_max.x())
m_max.x() += m_max.x() - m_min.x();
break;
}
m_root = newRoot;
}
/// Calculate root's box enclosing all root points
void initializeRootBox(const std::vector<point_type>& points)
{
const point_data_vec& data = m_nodes[m_root].data;
m_min = points[data.front()];
m_max = m_min;
for (pd_cit it = data.begin(); it != data.end(); ++it) {
const point_type& p = points[*it];
m_min = point_type::make(
std::min(m_min.x(), p.x()), std::min(m_min.y(), p.y()));
m_max = point_type::make(
std::max(m_max.x(), p.x()), std::max(m_max.y(), p.y()));
}
// Make sure bounding box does not have a zero size by adding padding:
// zero-size bounding box cannot be extended properly
const TCoordType padding(1);
if (m_min.x() == m_max.x()) {
m_min.x() -= padding;
m_max.x() += padding;
}
if (m_min.y() == m_max.y()) {
m_min.y() -= padding;
m_max.y() += padding;
}
m_isRootBoxInitialized = true;
}
private:
std::uint32_t m_root;
std::vector<Node> m_nodes;
NodeSplitDirection::Enum m_rootDir;
point_type m_min;
point_type m_max;
bool m_isRootBoxInitialized;
// used for nearest query
struct NearestTask {
std::uint32_t node;
point_type min, max;
NodeSplitDirection::Enum dir;
coord_type distSq;
NearestTask()
{
}
NearestTask(
const std::uint32_t node,
const point_type& min,
const point_type& max,
const NodeSplitDirection::Enum dir,
const coord_type distSq)
: node(node)
, min(min)
, max(max)
, dir(dir)
, distSq(distSq)
{
}
};
// allocated in class (not in the 'nearest' method) for better performance
mutable std::vector<NearestTask> m_tasksStack;
};
} // namespace kdt
namespace cdt {
/// KD-tree holding points
template <
typename TCoordType,
size_t NumVerticesInLeaf = 32,
size_t InitialStackDepth = 32,
size_t StackDepthIncrement = 32>
class locator_kdtree_t {
public:
/// Initialize KD-tree with points
void initialize(const std::vector<vec2_<TCoordType>>& points)
{
vec2_<TCoordType> min = points.front();
vec2_<TCoordType> max = min;
for (typename std::vector<vec2_<TCoordType>>::const_iterator it = points.begin();
it != points.end();
++it) {
min = vec2_<TCoordType>::make(std::min(min.x(), it->x()), std::min(min.y(), it->y()));
max = vec2_<TCoordType>::make(std::max(max.x(), it->x()), std::max(max.y(), it->y()));
}
m_tree = kdtree_t(min, max);
for (std::uint32_t i = 0; i < (std::uint32_t)points.size(); ++i) {
m_tree.insert(i, points);
}
}
/// Add point to KD-tree
void add_point(const std::uint32_t i, const std::vector<vec2_<TCoordType>>& points)
{
m_tree.insert(i, points);
}
/// Find nearest point using R-tree
std::uint32_t nearPoint(
const vec2_<TCoordType>& pos,
const std::vector<vec2_<TCoordType>>& points) const
{
return m_tree.nearest(pos, points).second;
}
private:
typedef kdt::kdtree_t_<TCoordType, NumVerticesInLeaf, InitialStackDepth, StackDepthIncrement> kdtree_t;
kdtree_t m_tree;
};
} // namespace cdt
#endif // header guard

2120
src/mcut/include/mcut/internal/cdt/triangulate.h Normal file
View file

File diff suppressed because it is too large Load diff

486
src/mcut/include/mcut/internal/cdt/utils.h Normal file
View file

@ -0,0 +1,486 @@
/* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at https://mozilla.org/MPL/2.0/. */
#ifndef _CDT_UTILITIES_H_
#define _CDT_UTILITIES_H_
#include "mcut/internal/math.h"
#include <cassert>
#include <cmath>
#include <limits>
#include <vector>
#include <array>
#include <functional>
#include <random>
#include <tuple>
#include <unordered_map>
#include <unordered_set>
namespace cdt {
/// X- coordinate getter for vec2d_t
template <typename T>
const T& get_x_coord_vec2d(const vec2_<T>& v)
{
return v.x();
}
/// Y-coordinate getter for vec2d_t
template <typename T>
const T& get_y_coord_vec2d(const vec2_<T>& v)
{
return v.y();
}
/// If two 2D vectors are exactly equal
template <typename T>
bool operator==(const vec2_<T>& lhs, const vec2_<T>& rhs)
{
return lhs.x() == rhs.x() && lhs.y() == rhs.y();
}
/// Constant representing no valid neighbor for a triangle
const static std::uint32_t null_neighbour(std::numeric_limits<std::uint32_t>::max());
/// Constant representing no valid vertex for a triangle
const static std::uint32_t null_vertex(std::numeric_limits<std::uint32_t>::max());
/// 2D bounding box
template <typename T>
struct box2d_t {
vec2_<T> min; ///< min box corner
vec2_<T> max; ///< max box corner
/// Envelop box around a point
void expand_with_point(const vec2_<T>& p)
{
expand_with_point(p.x(), p.y());
}
/// Envelop box around a point with given coordinates
void expand_with_point(const T x, const T y)
{
min.x() = std::min(x, min.x());
max.x() = std::max(x, max.x());
min.y() = std::min(y, min.y());
max.y() = std::max(y, max.y());
}
};
/// Bounding box of a collection of custom 2D points given coordinate getters
template <
typename T,
typename TVertexIter,
typename TGetVertexCoordX,
typename TGetVertexCoordY>
box2d_t<T> expand_with_points(
TVertexIter first,
TVertexIter last,
TGetVertexCoordX get_x_coord,
TGetVertexCoordY get_y_coord)
{
const T max = std::numeric_limits<T>::max();
box2d_t<T> box = { { max, max }, { -max, -max } };
for (; first != last; ++first) {
box.expand_with_point(get_x_coord(*first), get_y_coord(*first));
}
return box;
}
/// Bounding box of a collection of 2D points
template <typename T>
box2d_t<T> expand_with_points(const std::vector<vec2_<T>>& vertices);
/// edge_t connecting two vertices: vertex with smaller index is always first
/// \note: hash edge_t is specialized at the bottom
struct edge_t {
edge_t(std::uint32_t iV1, std::uint32_t iV2)
: m_vertices(
iV1 < iV2 ? std::make_pair(iV1, iV2) : std::make_pair(iV2, iV1))
{
}
inline bool operator==(const edge_t& other) const
{
return m_vertices == other.m_vertices;
}
inline bool operator!=(const edge_t& other) const
{
return !(this->operator==(other));
}
inline std::uint32_t v1() const
{
return m_vertices.first;
}
inline std::uint32_t v2() const
{
return m_vertices.second;
}
inline const std::pair<std::uint32_t, std::uint32_t>& verts() const
{
return m_vertices;
}
private:
std::pair<std::uint32_t, std::uint32_t> m_vertices;
};
/// Get edge first vertex
inline std::uint32_t edge_get_v1(const edge_t& e)
{
return e.v1();
}
/// Get edge second vertex
inline std::uint32_t edge_get_v2(const edge_t& e)
{
return e.v2();
}
/// Get edge second vertex
inline edge_t edge_make(std::uint32_t iV1, std::uint32_t iV2)
{
return edge_t(iV1, iV2);
}
/// triangulator_t triangle (CCW winding)
/* Counter-clockwise winding:
v3
/\
n3/ \n2
/____\
v1 n1 v2 */
struct triangle_t {
std::array<std::uint32_t, 3> vertices;
std::array<std::uint32_t, 3> neighbors;
/**
* Factory method
* @note needed for c++03 compatibility (no uniform initialization
* available)
*/
static triangle_t
make(const std::array<std::uint32_t, 3>& vertices, const std::array<std::uint32_t, 3>& neighbors)
{
triangle_t t = { vertices, neighbors };
return t;
}
};
/// Location of point on a triangle
struct point_to_triangle_location_t {
/// Enum
enum Enum {
INSIDE,
OUTSIDE,
ON_1ST_EDGE,
ON_2ND_EDGE,
ON_3RD_EDGE,
};
};
/// Relative location of point to a line
struct point_to_line_location_t {
/// Enum
enum Enum {
LEFT_SIDE,
RIGHT_SIDE,
COLLINEAR,
};
};
} // namespace cdt
#ifndef CDT_USE_AS_COMPILED_LIBRARY
#include <stdexcept>
namespace cdt {
//*****************************************************************************
// box2d_t
//*****************************************************************************
template <typename T>
box2d_t<T> expand_with_points(const std::vector<vec2_<T>>& vertices)
{
return expand_with_points<T>(
vertices.begin(), vertices.end(), get_x_coord_vec2d<T>, get_y_coord_vec2d<T>);
}
//*****************************************************************************
// Utility functions
//*****************************************************************************
/// Advance vertex or neighbor index counter-clockwise
inline std::uint32_t ccw(std::uint32_t i)
{
return std::uint32_t((i + 1) % 3);
}
/// Advance vertex or neighbor index clockwise
inline std::uint32_t cw(std::uint32_t i)
{
return std::uint32_t((i + 2) % 3);
}
/// Check if location is classified as on any of three edges
inline bool check_on_edge(const point_to_triangle_location_t::Enum location)
{
return location > point_to_triangle_location_t::OUTSIDE;
}
/// Neighbor index from a on-edge location
/// \note Call only if located on the edge!
inline std::uint32_t edge_neighbour(const point_to_triangle_location_t::Enum location)
{
assert(location >= point_to_triangle_location_t::ON_1ST_EDGE);
return static_cast<std::uint32_t>(location - point_to_triangle_location_t::ON_1ST_EDGE);
}
#if 0
/// Orient p against line v1-v2 2D: robust geometric predicate
template <typename T>
T orient2D(const vec2_<T>& p, const vec2_<T>& v1, const vec2_<T>& v2)
{
return orient2d(v1.x(), v1.y(), v2.x(), v2.y(), p.x(), p.y());
}
#endif
/// Classify value of orient2d predicate
template <typename T>
point_to_line_location_t::Enum
classify_orientation(const T orientation, const T orientationTolerance = T(0))
{
if (orientation < -orientationTolerance)
return point_to_line_location_t::RIGHT_SIDE;
if (orientation > orientationTolerance)
return point_to_line_location_t::LEFT_SIDE;
return point_to_line_location_t::COLLINEAR;
}
/// Check if point lies to the left of, to the right of, or on a line
template <typename T>
point_to_line_location_t::Enum locate_point_wrt_line(
const vec2_<T>& p,
const vec2_<T>& v1,
const vec2_<T>& v2,
const T orientationTolerance = T(0))
{
return classify_orientation(orient2d(p, v1, v2), orientationTolerance);
}
/// Check if point a lies inside of, outside of, or on an edge of a triangle
template <typename T>
point_to_triangle_location_t::Enum locate_point_wrt_triangle(
const vec2_<T>& p,
const vec2_<T>& v1,
const vec2_<T>& v2,
const vec2_<T>& v3)
{
point_to_triangle_location_t::Enum result = point_to_triangle_location_t::INSIDE;
point_to_line_location_t::Enum edgeCheck = locate_point_wrt_line(p, v1, v2);
if (edgeCheck == point_to_line_location_t::RIGHT_SIDE)
return point_to_triangle_location_t::OUTSIDE;
if (edgeCheck == point_to_line_location_t::COLLINEAR)
result = point_to_triangle_location_t::ON_1ST_EDGE;
edgeCheck = locate_point_wrt_line(p, v2, v3);
if (edgeCheck == point_to_line_location_t::RIGHT_SIDE)
return point_to_triangle_location_t::OUTSIDE;
if (edgeCheck == point_to_line_location_t::COLLINEAR)
result = point_to_triangle_location_t::ON_2ND_EDGE;
edgeCheck = locate_point_wrt_line(p, v3, v1);
if (edgeCheck == point_to_line_location_t::RIGHT_SIDE)
return point_to_triangle_location_t::OUTSIDE;
if (edgeCheck == point_to_line_location_t::COLLINEAR)
result = point_to_triangle_location_t::ON_3RD_EDGE;
return result;
}
/// Opposed neighbor index from vertex index
inline std::uint32_t get_opposite_neighbour_from_vertex(const std::uint32_t vertIndex)
{
MCUT_ASSERT(vertIndex < 3);
if (vertIndex == std::uint32_t(0))
return std::uint32_t(1);
if (vertIndex == std::uint32_t(1))
return std::uint32_t(2);
if (vertIndex == std::uint32_t(2))
return std::uint32_t(0);
throw std::runtime_error("Invalid vertex index");
}
/// Opposed vertex index from neighbor index
inline std::uint32_t opposite_vertex_from_neighbour(const std::uint32_t neighborIndex)
{
if (neighborIndex == std::uint32_t(0))
return std::uint32_t(2);
if (neighborIndex == std::uint32_t(1))
return std::uint32_t(0);
if (neighborIndex == std::uint32_t(2))
return std::uint32_t(1);
throw std::runtime_error("Invalid neighbor index");
}
/// Index of triangle's neighbor opposed to a vertex
inline std::uint32_t
opposite_triangle_index(const triangle_t& tri, const std::uint32_t iVert)
{
for (std::uint32_t vi = std::uint32_t(0); vi < std::uint32_t(3); ++vi)
if (iVert == tri.vertices[vi])
return get_opposite_neighbour_from_vertex(vi);
throw std::runtime_error("Could not find opposed triangle index");
}
/// Index of triangle's neighbor opposed to an edge
inline std::uint32_t opposite_triangle_index(
const triangle_t& tri,
const std::uint32_t iVedge1,
const std::uint32_t iVedge2)
{
for (std::uint32_t vi = std::uint32_t(0); vi < std::uint32_t(3); ++vi) {
const std::uint32_t iVert = tri.vertices[vi];
if (iVert != iVedge1 && iVert != iVedge2)
return get_opposite_neighbour_from_vertex(vi);
}
throw std::runtime_error("Could not find opposed-to-edge triangle index");
}
/// Index of triangle's vertex opposed to a triangle
inline std::uint32_t
get_opposite_vertex_index(const triangle_t& tri, const std::uint32_t iTopo)
{
for (std::uint32_t ni = std::uint32_t(0); ni < std::uint32_t(3); ++ni)
if (iTopo == tri.neighbors[ni])
return opposite_vertex_from_neighbour(ni);
throw std::runtime_error("Could not find opposed vertex index");
}
/// If triangle has a given neighbor return neighbor-index, throw otherwise
inline std::uint32_t
get_neighbour_index(const triangle_t& tri, std::uint32_t iTnbr)
{
for (std::uint32_t ni = std::uint32_t(0); ni < std::uint32_t(3); ++ni)
if (iTnbr == tri.neighbors[ni])
return ni;
throw std::runtime_error("Could not find neighbor triangle index");
}
/// If triangle has a given vertex return vertex-index, throw otherwise
inline std::uint32_t get_vertex_index(const triangle_t& tri, const std::uint32_t iV)
{
for (std::uint32_t i = std::uint32_t(0); i < std::uint32_t(3); ++i)
if (iV == tri.vertices[i])
return i;
throw std::runtime_error("Could not find vertex index in triangle");
}
/// Given triangle and a vertex find opposed triangle
inline std::uint32_t
get_opposite_triangle_index(const triangle_t& tri, const std::uint32_t iVert)
{
return tri.neighbors[opposite_triangle_index(tri, iVert)];
}
/// Given two triangles, return vertex of first triangle opposed to the second
inline std::uint32_t
get_opposed_vertex_index(const triangle_t& tri, std::uint32_t iTopo)
{
return tri.vertices[get_opposite_vertex_index(tri, iTopo)];
}
/// Test if point lies in a circumscribed circle of a triangle
template <typename T>
bool check_is_in_circumcircle(
const vec2_<T>& p,
const vec2_<T>& v1,
const vec2_<T>& v2,
const vec2_<T>& v3)
{
const double p_[2] = { static_cast<double>(p.x()), static_cast<double>(p.y()) };
const double v1_[2] = { static_cast<double>(v1.x()), static_cast<double>(v1.y()) };
const double v2_[2] = { static_cast<double>(v2.x()), static_cast<double>(v2.y()) };
const double v3_[2] = { static_cast<double>(v3.x()), static_cast<double>(v3.y()) };
return ::incircle(v1_, v2_, v3_, p_) > T(0);
}
/// Test if two vertices share at least one common triangle
inline bool check_vertices_share_edge(const std::vector<std::uint32_t>& aTris, const std::vector<std::uint32_t>& bTris)
{
for (std::vector<std::uint32_t>::const_iterator it = aTris.begin(); it != aTris.end(); ++it)
if (std::find(bTris.begin(), bTris.end(), *it) != bTris.end())
return true;
return false;
}
template <typename T>
T get_square_distance(const T ax, const T ay, const T bx, const T by)
{
const T dx = bx - ax;
const T dy = by - ay;
return dx * dx + dy * dy;
}
template <typename T>
T distance(const T ax, const T ay, const T bx, const T by)
{
return std::sqrt(get_square_distance(ax, ay, bx, by));
}
template <typename T>
T distance(const vec2_<T>& a, const vec2_<T>& b)
{
return distance(a.x(), a.y(), b.x(), b.y());
}
template <typename T>
T get_square_distance(const vec2_<T>& a, const vec2_<T>& b)
{
return get_square_distance(a.x(), a.y(), b.x(), b.y());
}
} // namespace cdt
#endif
//*****************************************************************************
// Specialize hash functions
//*****************************************************************************
namespace std {
/// edge_t hasher
template <>
struct hash<cdt::edge_t> {
/// Hash operator
std::size_t operator()(const cdt::edge_t& e) const
{
return get_hashed_edge_index(e);
}
private:
static void combine_hash_values(std::size_t& seed, const std::uint32_t& key)
{
seed ^= std::hash<std::uint32_t>()(key) + 0x9e3779b9 + (seed << 6) + (seed >> 2);
}
static std::size_t get_hashed_edge_index(const cdt::edge_t& e)
{
const std::pair<std::uint32_t, std::uint32_t>& vv = e.verts();
std::size_t seed1(0);
combine_hash_values(seed1, vv.first);
combine_hash_values(seed1, vv.second);
std::size_t seed2(0);
combine_hash_values(seed2, vv.second);
combine_hash_values(seed2, vv.first);
return std::min(seed1, seed2);
}
};
} // namespace std/boost
#endif // header guard

751
src/mcut/include/mcut/internal/frontend.h Normal file
View file

@ -0,0 +1,751 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
/**
* @file mcut.h
* @author Floyd M. Chitalu
* @date 11 July 2022
*
* @brief API-function implementations.
*
* NOTE: This header file defines the pre- and post-cutting processing of mesh
* data, which includes any intermediate correctons/modifications to the user's
* input meshes like 'polygon partitioning'.
*
*/
#ifndef _FRONTEND_H_
#define _FRONTEND_H_
#include "mcut/mcut.h"
#include <future>
#include <map>
#include <memory>
#include <string>
#include "mcut/internal/tpool.h"
#include "mcut/internal/kernel.h"
/*
std::invalid_argument: related to the input parameters
std::runtime_error: system runtime error e.g. out of memory
std::logic_error: a bug caught through an assertion failure
std::exception: unknown error source e.g. probably another bug
*/
#define CATCH_POSSIBLE_EXCEPTIONS(logstr) \
catch (std::invalid_argument & e0) \
{ \
logstr = e0.what(); \
return_value = McResult::MC_INVALID_VALUE; \
} \
catch (std::runtime_error & e1) \
{ \
logstr = e1.what(); \
return_value = McResult::MC_INVALID_OPERATION; \
} \
catch (std::logic_error & e2) \
{ \
logstr = e2.what(); \
return_value = McResult::MC_RESULT_MAX_ENUM; \
} \
catch (std::exception & e3) \
{ \
logstr = e3.what(); \
return_value = McResult::MC_RESULT_MAX_ENUM; \
}
extern thread_local std::string per_thread_api_log_str; // frontend.cpp
extern "C" void create_context_impl(
McContext* pContext, McFlags flags, uint32_t num_helper_threads) noexcept(false);
extern "C" void debug_message_callback_impl(
McContext context,
pfn_mcDebugOutput_CALLBACK cb,
const McVoid* userParam) noexcept(false);
extern "C" void get_debug_message_log_impl(McContext context,
McUint32 count, McSize bufSize,
McDebugSource* sources, McDebugType* types, McDebugSeverity* severities,
McSize* lengths, McChar* messageLog, McUint32& numFetched);
extern "C" void debug_message_control_impl(
McContext context,
McDebugSource source,
McDebugType type,
McDebugSeverity severity,
bool enabled) noexcept(false);
extern "C" void get_info_impl(
const McContext context,
McFlags info,
McSize bytes,
McVoid* pMem,
McSize* pNumBytes) noexcept(false);
extern "C" void create_user_event_impl(McEvent* event, McContext context);
extern "C" void set_user_event_status_impl(McEvent event, McInt32 execution_status);
extern "C" void get_event_info_impl(
const McEvent event,
McFlags info,
McSize bytes,
McVoid* pMem,
McSize* pNumBytes) noexcept(false);
extern "C" void set_event_callback_impl(
McEvent eventHandle,
pfn_McEvent_CALLBACK eventCallback,
McVoid* data);
extern "C" void wait_for_events_impl(
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McResult& runtimeStatusFromAllPrecedingEvents) noexcept(false);
extern "C" void dispatch_impl(
McContext context,
McFlags flags,
const McVoid* pSrcMeshVertices,
const uint32_t* pSrcMeshFaceIndices,
const uint32_t* pSrcMeshFaceSizes,
uint32_t numSrcMeshVertices,
uint32_t numSrcMeshFaces,
const McVoid* pCutMeshVertices,
const uint32_t* pCutMeshFaceIndices,
const uint32_t* pCutMeshFaceSizes,
uint32_t numCutMeshVertices,
uint32_t numCutMeshFaces,
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McEvent* pEvent) noexcept(false);
extern "C" void get_connected_components_impl(
const McContext context,
const McConnectedComponentType connectedComponentType,
const uint32_t numEntries,
McConnectedComponent* pConnComps,
uint32_t* numConnComps,
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McEvent* pEvent) noexcept(false);
extern "C" void get_connected_component_data_impl(
const McContext context,
const McConnectedComponent connCompId,
McFlags flags,
McSize bytes,
McVoid* pMem,
McSize* pNumBytes,
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McEvent* pEvent) noexcept(false);
extern "C" void release_connected_components_impl(
const McContext context,
uint32_t numConnComps,
const McConnectedComponent* pConnComps) noexcept(false);
extern "C" void release_context_impl(
McContext context) noexcept(false);
extern "C" void release_events_impl(uint32_t numEvents, const McEvent* pEvents);
// base struct from which other structs represent connected components inherit
struct connected_component_t {
virtual ~connected_component_t() {};
McConnectedComponentType type = (McConnectedComponentType)0;
McConnectedComponent m_user_handle = MC_NULL_HANDLE;
// array_mesh_t indexArrayMesh;
// hmesh_t mesh;
std::shared_ptr<output_mesh_info_t> kernel_hmesh_data;
//
std::shared_ptr< //
std::unordered_map< //
fd_t /*child face*/,
fd_t /*parent face in the [user-provided] source mesh*/
> //
>
source_hmesh_child_to_usermesh_birth_face; // fpPartitionChildFaceToCorrespondingInputSrcMeshFace
std::shared_ptr< //
std::unordered_map< //
fd_t /*child face*/,
fd_t /*parent face in the [user-provided] cut mesh*/
>>
cut_hmesh_child_to_usermesh_birth_face; // fpPartitionChildFaceToCorrespondingInputCutMeshFace
// descriptors and coordinates of new vertices that are added into an input mesh (source mesh or cut mesh)
// in order to carry out partitioning
std::shared_ptr<std::unordered_map<vd_t, vec3>> source_hmesh_new_poly_partition_vertices; // addedFpPartitioningVerticesOnCorrespondingInputSrcMesh
std::shared_ptr<std::unordered_map<vd_t, vec3>> cut_hmesh_new_poly_partition_vertices; // addedFpPartitioningVerticesOnCorrespondingInputCutMesh
uint32_t internal_sourcemesh_vertex_count; // init from source_hmesh.number_of_vertices()
uint32_t client_sourcemesh_vertex_count; // init from numSrcMeshVertices
uint32_t internal_sourcemesh_face_count; // init from source_hmesh.number_of_faces()
uint32_t client_sourcemesh_face_count; // init from source_hmesh_face_count OR numSrcMeshFaces
// Stores the contiguous array of unsigned integers that define
// a triangulation of all [non-triangle faces] of the connected component.
// This vector is only populated if client invokes mcGetConnectedComponnentData
// with flag MC_CONNECTED_COMPONENT_DATA_FACE_TRIANGULATION and has the effect of
// triangulating every non-triangle face in the connected component.
std::vector<uint32_t> cdt_index_cache;
bool cdt_index_cache_initialized = false;
// stores the mapping between a CDT triangle in the connected component and
// the original "birth-face" in an input mesh (source mesh or cut mesh)
std::vector<uint32_t> cdt_face_map_cache;
bool cdt_face_map_cache_initialized = false;
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
// Stores the number of vertices per face of CC. This is an optimization
// because there is a possibility that face-sizes may (at-minimum) be queried
// twice by user. The first case is during the populating (i.e. second) call to the API
// mcGetConnectedComponentData(..., MC_CONNECTED_COMPONENT_DATA_FACE_SIZE, ...);
// The second case is during the populating (i.e. second) call to the API
// mcGetConnectedComponentData(..., MC_CONNECTED_COMPONENT_DATA_FACE, ...);.
// The key detail here is that the second case requires knowledge of the
// number of vertices in each face in order to know how to schedule parallel
// work with prefix-sums etc.. Thus, the optimization is useful only if
// building MCUT with multi-threading
std::vector<uint32_t> face_sizes_cache;
bool face_sizes_cache_initialized = false;
// see documentation of face_sizes_cache above
// Similar concepts but applied to MC_CONNECTED_COMPONENT_DATA_FACE_ADJACENT_FACE
std::vector<uint32_t> face_adjacent_faces_size_cache;
bool face_adjacent_faces_size_cache_initialized = false;
#endif // #if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
};
// struct representing a fragment
struct fragment_cc_t : public connected_component_t {
McFragmentLocation fragmentLocation = (McFragmentLocation)0;
McFragmentSealType srcMeshSealType = (McFragmentSealType)0;
McPatchLocation patchLocation = (McPatchLocation)0;
};
// struct representing a patch
struct patch_cc_t : public connected_component_t {
McPatchLocation patchLocation = (McPatchLocation)0;
};
// struct representing a seam
struct seam_cc_t : public connected_component_t {
McSeamOrigin origin = (McSeamOrigin)0;
};
// struct representing an input (user provided mesh)
struct input_cc_t : public connected_component_t {
McInputOrigin origin = (McInputOrigin)0;
};
struct event_t {
std::future<void> m_future; // used to wait on event
// used to synchronise access to variables associated with the callback
// function.
// This also also allows us to overcome the edgecase that mcSetEventCallback
// is called after the task associated with an event has been completed,
// in which case the new callback will be invoked immediately.
// see "set_callback_data()" below
std::mutex m_callback_mutex;
struct {
// optional user callback, which is invoked when associated task is finished
pfn_McEvent_CALLBACK m_fn_ptr;
// pointer passed to user provided callback function
McVoid* m_data_ptr;
// atomic boolean flag indicating whether the callback associated with event
// object has been called
std::atomic<bool> m_invoked;
} m_callback_info;
// atomic boolean flag indicating whether the task associated with event
// object has completed running
std::atomic<bool> m_finished;
McEvent m_user_handle; // handle used by client app to reference this event object
// the Manager thread which was assigned the task of managing the task associated with this event object.
std::atomic<uint32_t> m_responsible_thread_id;
std::atomic<int32_t> m_runtime_exec_status; // API return code associated with respective task (for user to query)
std::atomic<size_t> m_timestamp_submit;
std::atomic<size_t> m_timestamp_start;
std::atomic<size_t> m_timestamp_end;
std::atomic<uint32_t> m_command_exec_status;
bool m_profiling_enabled;
McCommandType m_command_type;
// A callable object that also holds an std::future. Its purpose is to emulate
// the internal representation of an API task, where this time the task is actually
// a user command since this pointer is define ONLY for user events.
// The std::future object of the packaged task is used to initialize m_future
// when this event is a user event. This our internal mechanism allowing for
// API command to be able to effectively wait on user events.
std::unique_ptr<std::packaged_task<void()>> m_user_API_command_task_emulator;
McContext m_context;
event_t()
: m_user_handle(MC_NULL_HANDLE)
, m_responsible_thread_id(UINT32_MAX)
, m_runtime_exec_status(MC_NO_ERROR)
, m_timestamp_submit(0)
, m_timestamp_start(0)
, m_timestamp_end(0)
, m_command_exec_status(MC_RESULT_MAX_ENUM)
, m_profiling_enabled(true)
, m_command_type(MC_COMMAND_UKNOWN)
, m_user_API_command_task_emulator(nullptr)
, m_context(nullptr)
{
log_msg("[MCUT] Create event " << this);
m_callback_info.m_fn_ptr = nullptr;
m_callback_info.m_data_ptr = nullptr;
m_finished.store(false);
m_callback_info.m_invoked.store(true); // so that we do not call a null pointer/needless invoke the callback in the destructor
}
~event_t()
{
if (m_callback_info.m_invoked.load() == false && m_callback_info.m_fn_ptr != nullptr && m_runtime_exec_status.load() == MC_NO_ERROR) {
MCUT_ASSERT(m_user_handle != MC_NULL_HANDLE);
(*(m_callback_info.m_fn_ptr))(m_user_handle, m_callback_info.m_data_ptr);
}
log_msg("[MCUT] Destroy event " << this << "(" << m_user_handle << ")");
}
inline std::size_t get_time_since_epoch()
{
return std::chrono::duration_cast<std::chrono::nanoseconds>(std::chrono::system_clock::now().time_since_epoch()).count();
}
inline void log_submit_time()
{
if (m_profiling_enabled) {
this->m_timestamp_submit.store(get_time_since_epoch());
}
// TODO: use specific acquire-release semantics
// see e.g.: https://stackoverflow.com/questions/13632344/understanding-c11-memory-fences
m_command_exec_status.store(McEventCommandExecStatus::MC_SUBMITTED);
}
inline void log_start_time()
{
if (m_profiling_enabled) {
this->m_timestamp_start.store(get_time_since_epoch());
}
m_command_exec_status = McEventCommandExecStatus::MC_RUNNING;
}
// logs time and sets m_command_exec_status to MC_COMPLETE
inline void log_end_time()
{
if (m_profiling_enabled) {
this->m_timestamp_end.store(get_time_since_epoch());
}
m_command_exec_status.store(McEventCommandExecStatus::MC_COMPLETE);
}
// thread-safe function to set the callback function for an event object
void set_callback_data(McEvent handle, pfn_McEvent_CALLBACK fn_ptr, McVoid* data_ptr)
{
std::lock_guard<std::mutex> lock(m_callback_mutex); // exclusive access
m_user_handle = handle;
m_callback_info.m_fn_ptr = fn_ptr;
m_callback_info.m_data_ptr = data_ptr;
m_callback_info.m_invoked.store(false);
if (m_finished.load() == true) { // see mutex documentation
// immediately invoke the callback
(*(m_callback_info.m_fn_ptr))(m_user_handle, m_callback_info.m_data_ptr);
m_callback_info.m_invoked.store(true);
}
}
// update the status of the event object to "finished"
void notify_task_complete(McResult exec_status)
{
m_finished = true;
m_runtime_exec_status = exec_status;
std::lock_guard<std::mutex> lock(m_callback_mutex);
if (m_callback_info.m_invoked.load() == false && m_callback_info.m_fn_ptr != nullptr) {
MCUT_ASSERT(m_user_handle != MC_NULL_HANDLE);
(*(m_callback_info.m_fn_ptr))(m_user_handle, m_callback_info.m_data_ptr);
m_callback_info.m_invoked.store(true);
}
}
};
// init in frontened.cpp
extern threadsafe_list<std::shared_ptr<event_t>> g_events;
extern std::atomic<std::uintptr_t> g_objects_counter; // a counter that is used to assign a unique value to a McObject handle that will be returned to the user
extern std::once_flag g_objects_counter_init_flag;
// our custome deleter function for std::unique_ptr variable of an array type
template <typename Derived>
void fn_delete_cc(connected_component_t* p)
{
log_msg("[MCUT] Destroy connected component " << p->m_user_handle);
delete dynamic_cast<Derived*>(p);
}
// struct defining the state of a context object
struct context_t {
private:
std::atomic<bool> m_done;
std::vector<thread_safe_queue<function_wrapper>> m_queues;
// Master/Manager thread(s) which are responsible for running the API calls
// When a user of MCUT calls one of the APIs (e.g. mcEnqueueDispatch) the task
// of actually executing everything related to that task will be handled by
// one Manager thread. This manager thread itself will be involved in
// computing some/all part of the respective task (think of it as the "main"
// thread insofar as the API task is concerned). Some API tasks contain code
// sections that run in parallel, which is where the Manager thread will also
// submit tasks to the shared compute threadpool ("m_compute_threadpool").
// NOTE: must be declared after "thread_pool_terminate" and "work_queues"
std::vector<std::thread> m_api_threads;
// join_threads m_joiner;
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
// A pool of threads that is shared by manager threads to execute e.g. parallel
// section of MCUT API tasks (see e.g. frontend.cpp and kernel.cpp)
std::unique_ptr<thread_pool> m_compute_threadpool;
#endif
// The state and flag variable current used to configure the next dispatch call
McFlags m_flags = (McFlags)0;
void api_thread_main(uint32_t thread_id)
{
log_msg("[MCUT] Launch API thread " << std::this_thread::get_id() << " (" << thread_id << ")");
do {
function_wrapper task;
// We try_pop() first in case the task "producer" (API) thread
// already invoked cond_var.notify_one() of "m_queues[thread_id]""
// BEFORE current thread first-entered this function.
if (!m_queues[thread_id].try_pop(task)) {
m_queues[thread_id].wait_and_pop(task);
}
if (m_done) {
break;
}
task();
} while (true);
log_msg("[MCUT] Shutdown API thread " << std::this_thread::get_id() << " (" << thread_id << ")");
}
public:
context_t(McContext handle, McFlags flags
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
,
uint32_t num_compute_threads
#endif
)
: m_done(false)
// , m_joiner(m_api_threads)
, m_flags(flags)
, m_user_handle(handle)
, dbgCallbackBitfieldSource(0)
, dbgCallbackBitfieldType(0)
, dbgCallbackBitfieldSeverity(0)
{
log_msg("\n[MCUT] Create context " << m_user_handle);
try {
const uint32_t manager_thread_count = (flags & MC_OUT_OF_ORDER_EXEC_MODE_ENABLE) ? 2 : 1;
m_queues = std::vector<thread_safe_queue<function_wrapper>>(manager_thread_count);
for (uint32_t i = 0; i < manager_thread_count; ++i) {
m_queues[i].set_done_ptr(&m_done);
m_api_threads.push_back(std::thread(&context_t::api_thread_main, this, i));
}
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
// create the pool of compute threads. These are the worker threads that
// can be tasked with work from any manager-thread. Thus, manager threads
// share the available/user-specified compute threads.
m_compute_threadpool = std::unique_ptr<thread_pool>(new thread_pool(num_compute_threads, manager_thread_count));
#endif
} catch (...) {
shutdown();
log_msg("[MCUT] Destroy context due to exception" << m_user_handle);
throw;
}
}
~context_t()
{
shutdown();
log_msg("[MCUT] Destroy context " << m_user_handle);
}
void shutdown()
{
m_done = true;
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
m_compute_threadpool.reset();
#endif
for (int i = (int)m_api_threads.size() - 1; i >= 0; --i) {
m_queues[i].disrupt_wait_for_data();
if (m_api_threads[i].joinable()) {
m_api_threads[i].join();
}
}
}
McContext m_user_handle;
// returns the flags which determine the runtime configuration of this context
const McFlags& get_flags() const
{
return this->m_flags;
}
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
thread_pool& get_shared_compute_threadpool()
{
return m_compute_threadpool.get()[0];
}
#endif
template <typename FunctionType>
McEvent prepare_and_submit_API_task(McCommandType cmdType, uint32_t numEventsInWaitlist, const McEvent* pEventWaitList, FunctionType api_fn)
{
//
// create the event object associated with the enqueued task
//
std::shared_ptr<event_t> event_ptr = std::shared_ptr<event_t>(new event_t);
MCUT_ASSERT(event_ptr != nullptr);
g_events.push_front(event_ptr);
event_ptr->m_user_handle = reinterpret_cast<McEvent>(g_objects_counter++);
event_ptr->m_profiling_enabled = (this->m_flags & MC_PROFILING_ENABLE) != 0;
event_ptr->m_command_type = cmdType;
event_ptr->log_submit_time();
// List of events the enqueued task depends on
//
// local copy that will be captured by-value (user permitted to re-use pEventWaitList)
const std::vector<McEvent> event_waitlist(pEventWaitList, pEventWaitList + numEventsInWaitlist);
//
// Determine which manager thread to assign the task
//
// the id of manager thread that will be assigned the current task
uint32_t responsible_thread_id = UINT32_MAX;
for (std::vector<McEvent>::const_iterator waitlist_iter = event_waitlist.cbegin(); waitlist_iter != event_waitlist.cend(); ++waitlist_iter) {
const McEvent& parent_task_event_handle = *waitlist_iter;
const std::shared_ptr<event_t> parent_task_event_ptr = g_events.find_first_if([=](std::shared_ptr<event_t> e) { return e->m_user_handle == parent_task_event_handle; });
MCUT_ASSERT(parent_task_event_ptr != nullptr);
const bool parent_task_is_not_finished = parent_task_event_ptr->m_finished.load() == false;
if (parent_task_is_not_finished && parent_task_event_ptr->m_command_type != McCommandType::MC_COMMAND_USER) {
// id of manager thread, which was assigned the parent task
responsible_thread_id = parent_task_event_ptr->m_responsible_thread_id.load();
MCUT_ASSERT(responsible_thread_id != UINT32_MAX);
MCUT_ASSERT(responsible_thread_id < (uint32_t)m_api_threads.size());
break;
}
}
const bool have_responsible_thread = responsible_thread_id != UINT32_MAX;
if (!have_responsible_thread) {
uint32_t thread_with_empty_queue = UINT32_MAX;
for (uint32_t i = 0; i < (uint32_t)m_api_threads.size(); ++i) {
if (m_queues[i].empty() == true) {
thread_with_empty_queue = i;
break;
}
}
if (thread_with_empty_queue != UINT32_MAX) {
responsible_thread_id = thread_with_empty_queue;
} else { // all threads have work to do
responsible_thread_id = 0; // just pick thread 0
}
}
//
// Package-up the task as a synchronised operation that will wait for
// other tasks in the event_waitlist, compute the operation, and finally update
// the respective event state with the completion status.
//
std::weak_ptr<event_t> event_weak_ptr(event_ptr);
std::packaged_task<void()> task(
[=]() {
McResult runtime_status_from_all_preceding_events = McResult::MC_NO_ERROR;
if (!event_waitlist.empty()) {
wait_for_events_impl((uint32_t)event_waitlist.size(), &event_waitlist[0], runtime_status_from_all_preceding_events); // block until events are done
}
// if any previous event failed then we cannot proceed with this task.
// i.e. no-Op
if (runtime_status_from_all_preceding_events != McResult::MC_NO_ERROR) {
return;
}
MCUT_ASSERT(!event_weak_ptr.expired());
{
std::shared_ptr<event_t> event = event_weak_ptr.lock();
MCUT_ASSERT(event != nullptr);
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
event->log_start_time();
try {
api_fn(); // execute the API function.
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str); // exceptions may be thrown due to runtime errors, which must be reported back to user
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s (EventID=%p)\n", __FUNCTION__, per_thread_api_log_str.c_str(), event == nullptr ? (McEvent)0 : event->m_user_handle);
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
event->notify_task_complete(return_value); // updated event state to indicate task completion (lock-based)
event->log_end_time();
}
}
});
event_ptr->m_future = task.get_future(); // the future we can later wait on via mcWaitForEVents
event_ptr->m_responsible_thread_id = responsible_thread_id;
m_queues[responsible_thread_id].push(std::move(task)); // enqueue task to be executed when responsible API thread is free
return event_ptr->m_user_handle;
}
// the current set of connected components associated with context
threadsafe_list<std::shared_ptr<connected_component_t>> connected_components;
// McFlags dispatchFlags = (McFlags)0;
// client/user debugging variable
// ------------------------------
// function pointer to user-define callback function for status/erro reporting
pfn_mcDebugOutput_CALLBACK debugCallback = nullptr;
// user provided data for callback
const McVoid* debugCallbackUserParam = nullptr;
std::mutex debugCallbackMutex;
// controller for permmited messages based on the source of message
std::atomic<McFlags> dbgCallbackBitfieldSource;
// controller for permmited messages based on the type of message
std::atomic<McFlags> dbgCallbackBitfieldType;
// controller for permmited messages based on the severity of message
std::atomic<McFlags> dbgCallbackBitfieldSeverity;
bool dbgCallbackAllowAsyncCalls = true;
void set_debug_callback_data(pfn_mcDebugOutput_CALLBACK cb, const McVoid* data_ptr)
{
std::lock_guard<std::mutex> lguard(debugCallbackMutex);
debugCallback = cb;
debugCallbackUserParam = data_ptr;
}
struct debug_log_msg_t {
McDebugSource source;
McDebugType type;
McDebugSeverity severity;
std::string str;
};
std::vector<debug_log_msg_t> m_debug_logs;
// function to invoke the user-provided debug call back
void dbg_cb(McDebugSource source,
McDebugType type,
unsigned int id, // unused
McDebugSeverity severity,
const std::string& message)
{
if (this->m_flags & McContextCreationFlags::MC_DEBUG) // information logged only during debug mode
{
std::unique_lock<std::mutex> ulock(debugCallbackMutex, std::defer_lock);
if (!dbgCallbackAllowAsyncCalls) {
ulock.lock();
} // otherwise the callback will be invoked asynchronously
// can we log this type of message? (based on user preferences via mcDebugMessageControl)
const bool canLog = ((uint32_t)source & dbgCallbackBitfieldSource.load(std::memory_order_acquire)) && //
((uint32_t)type & dbgCallbackBitfieldType.load(std::memory_order_acquire)) && //
((uint32_t)severity & dbgCallbackBitfieldSeverity.load(std::memory_order_acquire));
if (canLog) {
if (debugCallback != nullptr) { // user gave us a callback function pointer
(*debugCallback)(source, type, id, severity, message.length(), message.c_str(), debugCallbackUserParam);
} else // write to the internal log
{
m_debug_logs.emplace_back(debug_log_msg_t());
debug_log_msg_t& dbg_log = m_debug_logs.back();
dbg_log.source = source;
dbg_log.type = type;
dbg_log.severity = severity;
dbg_log.str = message;
}
}
}
}
};
// list of contexts created by client/user
extern "C" threadsafe_list<std::shared_ptr<context_t>> g_contexts;
#endif // #ifndef _FRONTEND_H_

639
src/mcut/include/mcut/internal/hmesh.h Normal file
View file

@ -0,0 +1,639 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#ifndef MCUT_HALFEDGE_MESH_H_
#define MCUT_HALFEDGE_MESH_H_
#include "mcut/internal/math.h"
#include "mcut/internal/utils.h"
#include <algorithm>
#include <limits>
#include <map>
#include <memory>
#include <vector>
template <typename T>
class descriptor_t_ {
public:
typedef unsigned int index_type;
descriptor_t_() { }
virtual ~descriptor_t_() { }
explicit descriptor_t_(index_type i = (std::numeric_limits<index_type>::max)())
: m_value(i)
{
}
operator index_type() const { return m_value; }
void reset() { m_value = (std::numeric_limits<index_type>::max)(); }
bool is_valid() const
{
index_type inf = (std::numeric_limits<index_type>::max)();
return m_value != inf;
}
descriptor_t_& operator=(const index_type& _rhs)
{
m_value = _rhs;
return *this;
}
descriptor_t_& operator=(const T& _rhs) const
{
m_value = _rhs.m_value;
return *this;
}
bool operator==(const T& _rhs) const
{
return m_value == _rhs.m_value;
}
bool operator!=(const T& _rhs) const
{
return m_value != _rhs.m_value;
}
bool operator<(const T& _rhs) const
{
return m_value < _rhs.m_value;
}
descriptor_t_& operator++()
{
++m_value;
return *this;
}
descriptor_t_& operator--()
{
--m_value;
return *this;
}
descriptor_t_ operator++(int)
{
descriptor_t_ tmp(*this);
++m_value;
return tmp;
}
descriptor_t_ operator--(int)
{
descriptor_t_ tmp(*this);
--m_value;
return tmp;
}
descriptor_t_& operator+=(std::ptrdiff_t n)
{
m_value = (unsigned int)(m_value + n);
return *this;
}
protected:
unsigned int m_value;
};
class halfedge_descriptor_t : public descriptor_t_<halfedge_descriptor_t> {
public:
halfedge_descriptor_t()
: descriptor_t_<halfedge_descriptor_t>(std::numeric_limits<index_type>::max())
{
}
explicit halfedge_descriptor_t(descriptor_t_<halfedge_descriptor_t>::index_type idx)
: descriptor_t_<halfedge_descriptor_t>(idx)
{
}
virtual ~halfedge_descriptor_t()
{
}
};
class edge_descriptor_t : public descriptor_t_<edge_descriptor_t> {
public:
edge_descriptor_t()
: descriptor_t_<edge_descriptor_t>((std::numeric_limits<index_type>::max)())
{
}
explicit edge_descriptor_t(descriptor_t_<edge_descriptor_t>::index_type idx)
: descriptor_t_<edge_descriptor_t>(idx)
{
}
virtual ~edge_descriptor_t()
{
}
};
class face_descriptor_t : public descriptor_t_<face_descriptor_t> {
public:
face_descriptor_t()
: descriptor_t_<face_descriptor_t>((std::numeric_limits<index_type>::max)())
{
}
explicit face_descriptor_t(descriptor_t_<face_descriptor_t>::index_type idx)
: descriptor_t_<face_descriptor_t>(idx)
{
}
virtual ~face_descriptor_t()
{
}
};
class vertex_descriptor_t : public descriptor_t_<vertex_descriptor_t> {
public:
vertex_descriptor_t()
: descriptor_t_<vertex_descriptor_t>((std::numeric_limits<index_type>::max)())
{
}
explicit vertex_descriptor_t(descriptor_t_<vertex_descriptor_t>::index_type idx)
: descriptor_t_<vertex_descriptor_t>(idx)
{
}
virtual ~vertex_descriptor_t()
{
}
};
template <typename T>
struct id_ {
typedef T type;
};
template <typename V>
class array_iterator_t;
struct halfedge_data_t : id_<halfedge_descriptor_t> {
vertex_descriptor_t t; // target vertex
face_descriptor_t f; // face
halfedge_descriptor_t o; // opposite halfedge
halfedge_descriptor_t n; // next halfedge
halfedge_descriptor_t p; // previous halfedge
edge_descriptor_t e; // edge
halfedge_data_t()
//: o(null_halfedge()), n(null_halfedge()), p(null_halfedge()), t(null_vertex()), e(null_edge()), f(null_face())
{
}
};
struct edge_data_t : id_<edge_descriptor_t> {
halfedge_descriptor_t h; // primary halfedge (even idx)
};
struct face_data_t : id_<face_descriptor_t> {
std::vector<halfedge_descriptor_t> m_halfedges;
};
struct vertex_data_t : id_<vertex_descriptor_t> {
vec3 p; // geometry coordinates
//std::vector<face_descriptor_t> m_faces; // ... incident to vertex // TODO: this is not needed (can be inferred from "m_halfedges")
std::vector<halfedge_descriptor_t> m_halfedges; // ... which point to vertex (note: can be used to infer edges too)
};
typedef std::vector<vertex_data_t> vertex_array_t;
typedef std::vector<edge_data_t> edge_array_t;
typedef std::vector<halfedge_data_t> halfedge_array_t;
typedef std::vector<face_data_t> face_array_t;
typedef array_iterator_t<face_array_t> face_array_iterator_t;
typedef array_iterator_t<vertex_array_t> vertex_array_iterator_t;
typedef array_iterator_t<edge_array_t> edge_array_iterator_t;
typedef array_iterator_t<halfedge_array_t> halfedge_array_iterator_t;
/*
Internal mesh data structure used for cutting meshes
Memory Management
Memory management is semi-automatic. Memory grows as more elements are added to the structure but does not shrink when elements are removed.
When you add elements and the capacity of the underlying vector is exhausted, the vector reallocates memory.
As descriptors are basically indices, they refer to the same element after a reallocation.
When you remove an element it is only marked as removed.
Internally it is put in a free list, and when you add elements to the surface mesh, they are taken from the free list in case it is not empty.
For all elements there is a function to obtain the number of used elements, as well as the number of used [and] removed elements.
For vertices the functions are hmesh_t::number_of_vertices() and hmesh_t::number_of_internal_vertices(), respectively.
The first function is slightly different from the free function num_vertices(const G&) of the BGL package.
Iterators such as hmesh_t::vertex_iterator_t only enumerate elements that are not marked as deleted.
*/
class hmesh_t {
public:
hmesh_t();
~hmesh_t();
// static member functions
// -----------------------
static vertex_descriptor_t null_vertex();
static halfedge_descriptor_t null_halfedge();
static edge_descriptor_t null_edge();
static face_descriptor_t null_face();
// regular member functions
// ------------------------
// excluding removed elements
int number_of_vertices() const;
int number_of_edges() const;
int number_of_halfedges() const;
int number_of_faces() const;
vertex_descriptor_t source(const halfedge_descriptor_t& h) const;
vertex_descriptor_t target(const halfedge_descriptor_t& h) const;
halfedge_descriptor_t opposite(const halfedge_descriptor_t& h) const;
halfedge_descriptor_t prev(const halfedge_descriptor_t& h) const;
halfedge_descriptor_t next(const halfedge_descriptor_t& h) const;
void set_next(const halfedge_descriptor_t& h, const halfedge_descriptor_t& nxt);
void set_previous(const halfedge_descriptor_t& h, const halfedge_descriptor_t& prev);
edge_descriptor_t edge(const halfedge_descriptor_t& h) const;
face_descriptor_t face(const halfedge_descriptor_t& h) const;
vertex_descriptor_t vertex(const edge_descriptor_t e, const int v) const;
bool is_border(const halfedge_descriptor_t h);
bool is_border(const edge_descriptor_t e);
halfedge_descriptor_t halfedge(const edge_descriptor_t e, const int i) const;
// finds a halfedge between two vertices. Returns a default constructed halfedge descriptor, if source and target are not connected.
halfedge_descriptor_t halfedge(const vertex_descriptor_t s, const vertex_descriptor_t t, bool strict_check = false) const;
// finds an edge between two vertices. Returns a default constructed halfedge descriptor, if source and target are not connected.
edge_descriptor_t edge(const vertex_descriptor_t s, const vertex_descriptor_t t, bool strict_check = false) const;
vertex_descriptor_t add_vertex(const vec3& point);
vertex_descriptor_t add_vertex(const double& x, const double& y, const double& z);
// adds an edges into the mesh data structure, creating incident halfedges, and returns the
// halfedge whole target is "v1"
halfedge_descriptor_t add_edge(const vertex_descriptor_t v0, const vertex_descriptor_t v1);
face_descriptor_t add_face(const std::vector<vertex_descriptor_t>& vi);
// checks whether adding this face will violate 2-manifoldness (i.e. halfedge
// construction rules) which would lead to creating a non-manifold edge
// (one that is referenced by more than 2 faces which is illegal).
bool is_insertable(const std::vector<vertex_descriptor_t> &vi) const;
// also disassociates (not remove) any halfedges(s) and vertices incident to face
void remove_face(const face_descriptor_t f);
// also disassociates (not remove) the halfedges(s) and vertex incident to this halfedge
void remove_halfedge(halfedge_descriptor_t h);
// also disassociates (not remove) any face(s) incident to edge via its halfedges, and also disassociates the halfedges
void remove_edge(const edge_descriptor_t e, bool remove_halfedges = true);
void remove_vertex(const vertex_descriptor_t v);
void remove_elements();
void reset();
int number_of_internal_faces() const;
int number_of_internal_edges() const;
int number_of_internal_halfedges() const;
int number_of_internal_vertices() const;
int number_of_vertices_removed() const;
int number_of_edges_removed() const;
int number_of_halfedges_removed() const;
int number_of_faces_removed() const;
bool is_removed(face_descriptor_t f) const;
bool is_removed(edge_descriptor_t e) const;
bool is_removed(halfedge_descriptor_t h) const;
bool is_removed(vertex_descriptor_t v) const;
void reserve_for_additional_vertices(std::uint32_t n);
void reserve_for_additional_edges(std::uint32_t n);
void reserve_for_additional_halfedges(std::uint32_t n);
void reserve_for_additional_faces(std::uint32_t n);
void reserve_for_additional_elements(std::uint32_t additional_vertices);
///
template <typename I = int>
I get_removed_elements(id_<I>)
{
return I(); // unused
}
const std::vector<vertex_descriptor_t>& get_removed_elements(id_<array_iterator_t<vertex_array_t>>) const;
const std::vector<edge_descriptor_t>& get_removed_elements(id_<array_iterator_t<edge_array_t>>) const;
const std::vector<halfedge_descriptor_t>& get_removed_elements(id_<array_iterator_t<halfedge_array_t>>) const;
const std::vector<face_descriptor_t>& get_removed_elements(id_<array_iterator_t<face_array_t>>) const;
//
template <typename I = int>
I elements_begin_(id_<I>)
{
return I(); // unused
}
const vertex_array_iterator_t elements_begin_(id_<vertex_array_iterator_t>, bool account_for_removed_elems = true) const;
const edge_array_iterator_t elements_begin_(id_<edge_array_iterator_t>, bool account_for_removed_elems = true) const;
const halfedge_array_iterator_t elements_begin_(id_<halfedge_array_iterator_t>, bool account_for_removed_elems = true) const;
const face_array_iterator_t elements_begin_(id_<face_array_iterator_t>, bool account_for_removed_elems = true) const;
// returns the number of removed mesh elements (vertices, edges, faces or halfedges) between [start, end)
template <typename I>
uint32_t count_removed_elements_in_range(const array_iterator_t<I>& start, const array_iterator_t<I>& end) const
{
const long long N = (uint32_t)(end - start); // length including removed elements
MCUT_ASSERT(N >= 0);
if (N == 0) {
return 0;
}
// raw starting ptr offset
const uint32_t start_ = (std::uint32_t)(start - elements_begin_(id_<array_iterator_t<I>> {}, false));
uint32_t n = 0;
for (auto elem_descr : get_removed_elements(id_<array_iterator_t<I>> {})) {
const uint32_t descr = (uint32_t)elem_descr;
if (descr >= start_ && (descr <= (start_ + (uint32_t)(N - 1)))) {
++n;
}
}
return n;
}
const vec3& vertex(const vertex_descriptor_t& vd) const;
// returns vector of halfedges which point to vertex (i.e. "v" is their target)
const std::vector<halfedge_descriptor_t>& get_halfedges_around_vertex(const vertex_descriptor_t v) const;
std::vector<vertex_descriptor_t> get_vertices_around_face(const face_descriptor_t f, uint32_t prepend_offset = 0) const;
void get_vertices_around_face(std::vector<vertex_descriptor_t>& vertex_descriptors, const face_descriptor_t f, uint32_t prepend_offset=0) const;
std::vector<vertex_descriptor_t> get_vertices_around_vertex(const vertex_descriptor_t v) const;
void get_vertices_around_vertex(std::vector<vertex_descriptor_t>& vertices_around_vertex, const vertex_descriptor_t v) const;
uint32_t get_num_vertices_around_face(const face_descriptor_t f) const;
const std::vector<halfedge_descriptor_t>& get_halfedges_around_face(const face_descriptor_t f) const;
const std::vector<face_descriptor_t> get_faces_around_face(const face_descriptor_t f, const std::vector<halfedge_descriptor_t>* halfedges_around_face_ = nullptr) const;
void get_faces_around_face( std::vector<face_descriptor_t>& faces_around_face, const face_descriptor_t f, const std::vector<halfedge_descriptor_t>* halfedges_around_face_ = nullptr) const;
uint32_t get_num_faces_around_face(const face_descriptor_t f, const std::vector<halfedge_descriptor_t>* halfedges_around_face_ = nullptr) const;
// iterators
// ---------
vertex_array_iterator_t vertices_begin(bool account_for_removed_elems = true) const;
vertex_array_iterator_t vertices_end() const;
edge_array_iterator_t edges_begin(bool account_for_removed_elems = true) const;
edge_array_iterator_t edges_end() const;
halfedge_array_iterator_t halfedges_begin(bool account_for_removed_elems = true) const;
halfedge_array_iterator_t halfedges_end() const;
face_array_iterator_t faces_begin(bool account_for_removed_elems = true) const;
face_array_iterator_t faces_end() const;
const std::vector<vertex_descriptor_t>& get_removed_vertices() const;
const std::vector<edge_descriptor_t>& get_removed_edges() const;
const std::vector<halfedge_descriptor_t>& get_removed_halfedges() const;
const std::vector<face_descriptor_t>& get_removed_faces() const;
private:
// member variables
// ----------------
std::vector<vertex_data_t> m_vertices;
std::vector<edge_data_t> m_edges;
std::vector<halfedge_data_t> m_halfedges;
std::vector<face_data_t> m_faces;
// NOTE: I use std::vector because we'll have very few (typically zero)
// elements removed at a given time. In fact removal only happens during
// input-mesh face-partitioning to resolve floating polygons, which is
// rare. Maybe in the future things change...
std::vector<face_descriptor_t> m_faces_removed;
std::vector<edge_descriptor_t> m_edges_removed;
std::vector<halfedge_descriptor_t> m_halfedges_removed;
std::vector<vertex_descriptor_t> m_vertices_removed;
}; // class hmesh_t {
typedef vertex_descriptor_t vd_t;
typedef halfedge_descriptor_t hd_t;
typedef edge_descriptor_t ed_t;
typedef face_descriptor_t fd_t;
void write_off(const char* fpath, const hmesh_t& mesh);
void read_off(hmesh_t& mesh, const char* fpath);
template <typename V = face_array_t>
class array_iterator_t : public V::const_iterator {
const hmesh_t* mesh_ptr;
typedef typename V::const_iterator std_iterator_base_class;
typedef typename V::value_type::type element_descriptor_type;
typename V::value_type* operator->() = delete;
public:
array_iterator_t()
: V::const_iterator()
, mesh_ptr(nullptr) {};
array_iterator_t(typename V::const_iterator it_, const hmesh_t* const mesh)
: V::const_iterator(it_)
, mesh_ptr(mesh)
{
}
const hmesh_t* get_mesh_ptr() const
{
return mesh_ptr;
}
typename V::value_type::type operator*()
{
size_t raw_index = (*this) - cbegin<>(false);
element_descriptor_type d((std::uint32_t)raw_index);
return d;
}
// prefix increment (++i)
// increment pointer to the next valid element (i.e. we skip removed elements).
array_iterator_t<V>& operator++()
{
bool cur_elem_is_removed = false;
bool reached_end = false;
do {
V::const_iterator::operator++();
reached_end = (*this) == cend<>();
cur_elem_is_removed = false;
if (!reached_end) {
const std::size_t diff = ((*this) - cbegin<array_iterator_t<V>>(false));
element_descriptor_type raw_descriptor((std::uint32_t)diff); // std::distance(cbegin<array_iterator_t<V>>(false), (*this)); // O(1) ??
cur_elem_is_removed = mesh_ptr->is_removed(raw_descriptor);
if (!cur_elem_is_removed) {
break;
}
}
// keep iterating until the value pointed to after the (++i) operator is a valid element
// i.e. one that is not marked removed!
} while (cur_elem_is_removed && !reached_end);
return (*this);
}
// we provide this overide to ensure that stl functions like std::advance, work properly
// by accounting for removed elements
array_iterator_t<V>& operator+=(std::ptrdiff_t n)
{
V::const_iterator::operator+=(n); // raw ptr shift (i.e. ignoring that there may be removed elements)
bool cur_elem_is_removed = false;
bool reached_end = (*this) == cend<>();
cur_elem_is_removed = mesh_ptr->is_removed(*(*this));
while (!reached_end && cur_elem_is_removed) {
V::const_iterator::operator++(); //++(*this);
size_t raw_descriptor = *(*this); // (*this) - cbegin<array_iterator_t<V>>(false); //std::distance(cbegin<array_iterator_t<V>>(false), (*this)); // O(1) ??
cur_elem_is_removed = mesh_ptr->is_removed(element_descriptor_type((std::uint32_t)raw_descriptor));
if (!cur_elem_is_removed) {
break;
}
reached_end = (*this) == cend<>();
}
return *this;
}
// The following are helper functions which are specialised (via type-deduction)
// for the type of mesh elements that *this* iterator walks over in "mesh_ptr"
// e.g. faces. These functions are used to determine when *this* iterator has
// reached the end of the respective std::map data structure over which we are
// iterating.
template <typename I = array_iterator_t<V>>
I cend()
{
return cend(id_<I>()); // https://stackoverflow.com/questions/3052579/explicit-specialization-in-non-namespace-scope
}
template <typename I = array_iterator_t<V>>
I cbegin(bool account_for_removed_elems)
{
return cbegin(account_for_removed_elems, id_<I>());
}
private:
// postfix increment (i++)
// NOTE: This overide is private to simplify implementation, and we don't need it
array_iterator_t<V> operator++(int)
{
MCUT_ASSERT(false);
return cend<>();
}
template <typename I = array_iterator_t<V>>
I cend(id_<I>)
{
return I(); // unused stub
}
template <typename I = array_iterator_t<V>>
I cbegin(bool account_for_removed_elems, id_<I>)
{
return I(account_for_removed_elems); // unused stub
}
vertex_array_iterator_t cbegin(bool account_for_removed_elems, id_<vertex_array_iterator_t> = {});
vertex_array_iterator_t cend(id_<vertex_array_iterator_t>);
edge_array_iterator_t cbegin(bool account_for_removed_elems, id_<edge_array_iterator_t> = {});
edge_array_iterator_t cend(id_<edge_array_iterator_t>);
halfedge_array_iterator_t cbegin(bool account_for_removed_elems, id_<halfedge_array_iterator_t> = {});
halfedge_array_iterator_t cend(id_<halfedge_array_iterator_t>);
face_array_iterator_t cbegin(bool account_for_removed_elems, id_<face_array_iterator_t> = {});
face_array_iterator_t cend(id_<face_array_iterator_t>);
}; // class array_iterator_t : public V::const_iterator
namespace std {
#if 1
template <>
inline typename edge_array_iterator_t::difference_type distance(
edge_array_iterator_t first,
edge_array_iterator_t last)
{
MCUT_ASSERT(first.get_mesh_ptr() == last.get_mesh_ptr());
edge_array_iterator_t it = first;
edge_array_iterator_t::difference_type dist = last - first;
uint32_t r = it.get_mesh_ptr()->count_removed_elements_in_range(first, last);
if (r > 0) {
dist = dist - r;
}
MCUT_ASSERT(dist >= 0);
return dist;
}
#endif
#if 0
template <>
void advance(
hmesh_t::array_iterator_t<hmesh_t::edge_array_t> &iter,
typename std::iterator_traits<hmesh_t::array_iterator_t<hmesh_t::edge_array_t>>::difference_type n);
#endif
template <>
struct hash<vertex_descriptor_t> {
std::size_t operator()(const vertex_descriptor_t& k) const
{
return std::hash<typename vertex_descriptor_t::index_type>()(static_cast<typename vertex_descriptor_t::index_type>(k));
}
};
template <>
struct hash<edge_descriptor_t> {
std::size_t operator()(const edge_descriptor_t& k) const
{
return std::hash<typename edge_descriptor_t::index_type>()(static_cast<typename edge_descriptor_t::index_type>(k));
}
};
template <>
struct hash<halfedge_descriptor_t> {
std::size_t operator()(const halfedge_descriptor_t& k) const
{
return std::hash<typename halfedge_descriptor_t::index_type>()(static_cast<typename halfedge_descriptor_t::index_type>(k));
}
};
template <>
struct hash<face_descriptor_t> {
std::size_t operator()(const face_descriptor_t& k) const
{
return std::hash<typename face_descriptor_t::index_type>()(static_cast<typename face_descriptor_t::index_type>(k));
}
};
}
#endif // #ifndef MCUT_HALFEDGE_MESH_H_

233
src/mcut/include/mcut/internal/kernel.h Normal file
View file

@ -0,0 +1,233 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#ifndef MCUT_KERNEL_H
#define MCUT_KERNEL_H
#include <mcut/internal/bvh.h>
#include <mcut/internal/hmesh.h>
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
#include <mcut/internal/tpool.h>
#endif
#include <map>
#include <unordered_map>
#include <vector>
//
// final execution states (i.e. did anything go wrong..?)
//
enum class status_t {
SUCCESS = 0,
// mesh is malformed
// * vertices less than 3
// * no faces
// * non-manifold
// * contains more than one connected component
INVALID_SRC_MESH = -1, // TODO: these error flags should be generated in mcut.cpp not the kernel
INVALID_CUT_MESH = -2,
// EDGE_EDGE_INTERSECTION = -3, // Found an edge-edge intersection.
// FACE_VERTEX_INTERSECTION = -4, // Found an face-vertex intersection.
// there exists no edge in the input mesh which intersects cut-surface polygon
INVALID_MESH_INTERSECTION = -3,
// The bounding volume heirarchies of the input mesh and cut surface do not overlap
// INVALID_BVH_INTERSECTION = -6,
/*
MCUT is formulated for inputs in general position. Here the notion of general position is defined with
respect to the orientation predicate (as evaluated on the intersecting polygons). Thus, a set of points
is in general position if no three points are collinear and also no four points are coplanar.
MCUT uses the "GENERAL_POSITION_VIOLATION" flag to inform of when to use perturbation (of the
cut-mesh) so as to bring the input into general position. In such cases, the idea is to solve the cutting
problem not on the given input, but on a nearby input. The nearby input is obtained by perturbing the given
input. The perturbed input will then be in general position and, since it is near the original input,
the result for the perturbed input will hopefully still be useful. This is justified by the fact that
the task of MCUT is not to decide whether the input is in general position but rather to make perturbation
on the input (if) necessary within the available precision of the computing device.
*/
GENERAL_POSITION_VIOLATION = -4,
/*
TODO: add documentation
*/
DETECTED_FLOATING_POLYGON = -5
};
//
// Position of a cut surface patch with respect to the input mesh
//
enum class cm_patch_location_t : unsigned char {
INSIDE, // + : The patch is located inside the input mesh volume (i.e. it is used to seal holes)
OUTSIDE, // - : The patch is located outside the input mesh volume (boolean union).
UNDEFINED // ~ : The notion of INSIDE/OUTSIDE is not applicable because input mesh is non-watertight
};
//
// Position of a connected component (CC) relative to cut-surface
//
enum class sm_frag_location_t : unsigned char {
ABOVE, // + : The CC is on positive side of the cut-surface (normal direction)
BELOW, // - : The CC is on negative side of the cut-surface (normal direction)
UNDEFINED // ~ : The notion of ABOVE/BELOW is not applicable because the CC has [partially] cut
};
//
// The winding order of the polygons of a cut surface patch
//
enum class cm_patch_winding_order_t : unsigned char {
DEFAULT, // + : The polygons of the patch have the [same] winding order as the cut-surface (e.g. CCW)
REVERSE, // - : The polygons of the patch have the [opposite] winding order as the cut-surface (e.g. CW)
};
struct floating_polygon_info_t {
// normal of polygon
vec3 polygon_normal;
int polygon_normal_largest_component = -1;
// the positions of the vertices of the floating polygon (order implies connectivity i.e. two points next to each other share a vertex)
std::vector<vec3> polygon_vertices;
};
//
// settings for how to execute the function "dispatch(...)"
//
struct input_t {
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
thread_pool* scheduler = nullptr;
#endif
/*const*/ std::shared_ptr<hmesh_t> src_mesh = nullptr;
/*const*/ std::shared_ptr<hmesh_t> cut_mesh = nullptr;
// NOTE: we use std::map because it is beneficial that keys are sorted when
// extracting edge-face intersection pairs
const std::map<fd_t, std::vector<fd_t>>* ps_face_to_potentially_intersecting_others = nullptr;
#if defined(USE_OIBVH)
const std::vector<bounding_box_t<vec3>>* source_hmesh_face_aabb_array_ptr = nullptr;
const std::vector<bounding_box_t<vec3>>* cut_hmesh_face_aabb_array_ptr = nullptr;
#else
BoundingVolumeHierarchy* source_hmesh_BVH;
BoundingVolumeHierarchy* cut_hmesh_BVH;
#endif
bool verbose = true;
// bool keep_partially_sealed_connected_components = false;
bool require_looped_cutpaths = false; // ... i.e. bail on partial cuts (any!)
bool populate_vertex_maps = false; // compute data relating vertices in cc to original input mesh
bool populate_face_maps = false; // compute data relating face in cc to original input mesh
bool enforce_general_position = false;
// counts how many times we have perturbed the cut-mesh to enforce general-position
int general_position_enforcement_count = 0;
// NOTE TO SELF: if the user simply wants seams, then kernel should not have to proceed to stitching!!!
bool keep_srcmesh_seam = false;
bool keep_cutmesh_seam = false;
//
bool keep_unsealed_fragments = false;
//
bool keep_inside_patches = false;
bool keep_outside_patches = false;
//
bool keep_fragments_below_cutmesh = false;
bool keep_fragments_above_cutmesh = false;
//
bool keep_fragments_partially_cut = false; // TODO: remove
bool keep_fragments_sealed_inside = false;
bool keep_fragments_sealed_outside = false;
// bool include_fragment_sealed_partial = false; // See: variable above "keep_partially_sealed_connected_components"
bool keep_fragments_sealed_inside_exhaustive = false; // TODO remove
bool keep_fragments_sealed_outside_exhaustive = false; // TODO remove
// NOTE TO SELF: if the user simply wants patches, then kernel should not have to proceed to stitching!!!
};
struct output_mesh_data_maps_t {
// std::map<
// vd_t, // vertex descriptor in connected component
// vd_t // vertex descriptor in input-mesh e.g. source mesh or cut mesh. ("null_vertex()" if vertex is an intersection point)
// >
std::vector<vd_t> vertex_map;
// std::map<
// fd_t, // face descriptor in connected component
// fd_t // face descriptor in input-mesh e.g. source mesh or cut mesh. (new polygons resulting from clipping are mapped to the same input mesh face)
// >
std::vector<fd_t> face_map;
};
struct output_mesh_info_t {
std::shared_ptr<hmesh_t> mesh;
std::vector<vd_t> seam_vertices;
output_mesh_data_maps_t data_maps;
};
//
// the output returned from the function "dispatch"
//
struct output_t {
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
std::atomic<status_t> status;
#else
status_t status = status_t::SUCCESS;
#endif
logger_t logger;
// fragments
std::map<sm_frag_location_t, std::map<cm_patch_location_t, std::vector<std::shared_ptr<output_mesh_info_t>>>> connected_components;
std::map<sm_frag_location_t, std::vector<std::shared_ptr<output_mesh_info_t>>> unsealed_cc; // connected components before hole-filling
// patches
std::map<cm_patch_winding_order_t, std::vector<std::shared_ptr<output_mesh_info_t>>> inside_patches; // .. between neigbouring connected ccsponents (cs-sealing patches)
std::map<cm_patch_winding_order_t, std::vector<std::shared_ptr<output_mesh_info_t>>> outside_patches;
// the input meshes which also include the edges that define the cut path
// NOTE: not always defined (depending on the arising cutpath configurations)
std::shared_ptr<output_mesh_info_t> seamed_src_mesh;
std::shared_ptr<output_mesh_info_t> seamed_cut_mesh;
// floating polygon handling
std::map<
// the face of the origin-mesh on which floating polygon(s) are discovered
// NOTE: this is a descriptor into the polygon soup. Thus, we'll need to
// subtract the number of source-mesh faces if this face belongs to the cut
// mesh.
fd_t,
// info about floating polygons contained on ps-face
std::vector<floating_polygon_info_t>>
detected_floating_polygons;
};
// internal main
void dispatch(output_t& out, const input_t& in);
int find_connected_components(
#if defined(MCUT_WITH_COMPUTE_HELPER_THREADPOOL)
thread_pool& scheduler,
#endif
std::vector<int>& fccmap, const hmesh_t& mesh, std::vector<int>& cc_to_vertex_count,
std::vector<int>& cc_to_face_count);
// return true if point p lies on the plane of every three vertices of f
bool point_on_face_plane(const hmesh_t& m, const fd_t& f, const vec3& p, int& fv_count);
//
// returns string equivalent value (e.g. for printing)
//
std::string to_string(const sm_frag_location_t&);
std::string to_string(const cm_patch_location_t&);
std::string to_string(const status_t&);
std::string to_string(const cm_patch_winding_order_t&);
#endif // #ifndef MCUT_KERNEL_H

695
src/mcut/include/mcut/internal/math.h Normal file
View file

@ -0,0 +1,695 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#ifndef MCUT_MATH_H_
#define MCUT_MATH_H_
#include <cmath>
#include <iostream>
#include <limits>
#include <memory>
#include <vector>
#include "mcut/internal/utils.h"
enum sign_t {
ON_NEGATIVE_SIDE = -1, // left
ON_ORIENTED_BOUNDARY = 0, // on boundary
ON_POSITIVE_SIDE = 1, // right
//
NEGATIVE = ON_NEGATIVE_SIDE,
ZERO = ON_ORIENTED_BOUNDARY,
POSITIVE = ON_POSITIVE_SIDE,
};
template <typename T = double>
class vec2_ {
public:
typedef T element_type;
vec2_()
: m_x(0.0)
, m_y(0.0)
{
}
vec2_(const T& value)
: m_x(value)
, m_y(value)
{
}
vec2_(const T& x, const T& y)
: m_x(x)
, m_y(y)
{
}
virtual ~vec2_()
{
}
static vec2_ make(const T x, const T y)
{
return vec2_<T>(x, y);
}
static int cardinality()
{
return 2;
}
bool operator==(const vec2_<T>& other) const
{
return this->m_x == other.x() && this->m_y == other.y();
}
const T& operator[](int index) const
{
MCUT_ASSERT(index >= 0 && index <= 1);
if (index == 0) {
return m_x;
} else {
return m_y;
}
}
T& operator[](int index)
{
MCUT_ASSERT(index >= 0 && index <= 1);
if (index == 0) {
return m_x;
}
return m_y;
}
const vec2_ operator-(const vec2_& other) const
{
return vec2_(m_x - other.m_x, m_y - other.m_y);
}
const vec2_ operator+(const vec2_& other) const
{
return vec2_(m_x + other.m_x, m_y + other.m_y);
}
const vec2_ operator/(const T& number) const
{
return vec2_(m_x / number, m_y / number);
}
const vec2_ operator*(const T& number) const
{
return vec2_(m_x * number, m_y * number);
}
const T& x() const
{
return m_x;
}
const T& y() const
{
return m_y;
}
T& x()
{
return m_x;
}
T& y()
{
return m_y;
}
protected:
T m_x, m_y;
}; // vec2_
typedef vec2_<> vec2;
template <typename T = double>
class vec3_ : public vec2_<T> {
public:
vec3_()
: vec2_<T>(0.0, 0.0)
, m_z(0.0)
{
}
vec3_(const T& value)
: vec2_<T>(value, value)
, m_z(value)
{
}
vec3_(const T& x, const T& y, const T& z)
: vec2_<T>(x, y)
, m_z(z)
{
}
~vec3_()
{
}
static int cardinality()
{
return 3;
}
const T& operator[](int index) const
{
MCUT_ASSERT(index >= 0 && index <= 2);
if (index == 0) {
return this->m_x;
} else if (index == 1) {
return this->m_y;
} else {
return this->m_z;
}
}
T& operator[](int index)
{
MCUT_ASSERT(index >= 0 && index <= 2);
if (index == 0) {
return this->m_x;
} else if (index == 1) {
return this->m_y;
} else {
return this->m_z;
}
}
vec3_ operator-(const vec3_& other) const
{
return vec3_(this->m_x - other.m_x, this->m_y - other.m_y, this->m_z - other.m_z);
}
vec3_ operator+(const vec3_& other) const
{
return vec3_(this->m_x + other.m_x, this->m_y + other.m_y, this->m_z + other.m_z);
}
const vec3_ operator/(const T& number) const
{
return vec3_(this->m_x / number, this->m_y / number, this->m_z / number);
}
const vec3_ operator*(const T& number) const
{
return vec3_(this->m_x * number, this->m_y * number, this->m_z * number);
}
const T& z() const
{
return m_z;
}
protected:
T m_z;
}; // vec3_
typedef vec3_<> vec3;
template <typename T = int>
class matrix_t {
public:
matrix_t()
: m_rows(-1)
, m_cols(-1)
{
}
matrix_t(unsigned int rows, unsigned int cols)
: m_rows(rows)
, m_cols(cols)
, m_entries(std::vector<T>((size_t)rows * cols, T(0.0))) // zeroes
{
}
T& operator()(unsigned int row, unsigned int col)
{
return m_entries[(size_t)row * m_cols + col];
}
T operator()(unsigned int row, unsigned int col) const
{
return m_entries[(size_t)row * m_cols + col];
}
// multiply
matrix_t<T> operator*(const matrix_t<T>& other) const
{
matrix_t<T> result(this->rows(), other.cols());
for (int i = 0; i < this->rows(); ++i) {
for (int j = 0; j < other.cols(); ++j) {
for (int k = 0; k < this->cols(); ++k) {
result(i, j) += (*this)(i, k) * other(k, j);
}
}
}
return result;
}
matrix_t<T> operator*(const double& s) const
{
matrix_t<T> result(m_rows, m_cols);
for (int i = 0; i < m_rows; ++i) {
for (int j = 0; j < m_cols; ++j) {
result(i, j) = (*this)(i, j) * s;
}
}
return result;
}
matrix_t<T> operator/(const double& s) const
{
matrix_t<T> result(m_rows, m_cols);
for (int i = 0; i < m_rows; ++i) {
for (int j = 0; j < m_cols; ++j) {
result(i, j) = (*this)(i, j) / s;
}
}
return result;
}
matrix_t<T> operator-(const matrix_t<T>& m) const
{
MCUT_ASSERT(m.rows() == this->rows());
MCUT_ASSERT(m.cols() == this->cols());
matrix_t<T> result(m_rows, m_cols);
for (int i = 0; i < m_rows; ++i) {
for (int j = 0; j < m_cols; ++j) {
result(i, j) = (*this)(i, j) - m(i, j);
}
}
return result;
}
// 2x3 matrix times 3x1 vector
vec2 operator*(const vec3& v) const
{
MCUT_ASSERT(this->cols() == vec3::cardinality());
vec2 result(double(0.0));
MCUT_ASSERT(this->rows() == vec2::cardinality());
for (int col = 0; col < this->cols(); ++col) {
for (int row = 0; row < vec2::cardinality(); ++row) {
result[row] = result[row] + ((*this)(row, col) * v[col]);
}
}
#if 0
// columns
const Vec c0((*this)(0, 0), (*this)(1, 0), (*this)(2, 0));
const Vec c1((*this)(0, 1), (*this)(1, 1), (*this)(2, 1));
result = c0 * v[0] + c1 * v[1];
if (this->rows() == 3 && (this->cols() == 3)
{
const Vec c2((*this)(0, 2), (*this)(1, 2), (*this)(2, 2));
result = result + c2 * v[2];
}
#endif
return result;
}
inline int rows() const
{
return m_rows;
}
inline int cols() const
{
return m_cols;
}
private:
int m_rows;
int m_cols;
std::vector<T> m_entries;
};
extern double square_root(const double& number);
extern double absolute_value(const double& number);
extern sign_t sign(const double& number);
extern std::ostream& operator<<(std::ostream& os, const vec3& v);
template <typename U>
std::ostream& operator<<(std::ostream& os, const matrix_t<U>& m)
{
for (int i = 0; i < m.rows(); i++) {
for (int j = 0; j < m.cols(); j++) {
os << m(i, j) << ", ";
}
os << "\n";
}
return os;
}
template <typename T>
const T& min(const T& a, const T& b)
{
return ((b < a) ? b : a);
}
template <typename T>
const T& max(const T& a, const T& b)
{
return ((a < b) ? b : a);
}
extern bool operator==(const vec3& a, const vec3& b);
template <typename T>
vec2_<T> compwise_min(const vec2_<T>& a, const vec2_<T>& b)
{
return vec2_<T>(min(a.x(), b.x()), min(a.y(), b.y()));
}
template <typename T>
vec2_<T> compwise_max(const vec2_<T>& a, const vec2_<T>& b)
{
return vec2_<T>(max(a.x(), b.x()), max(a.y(), b.y()));
}
template <typename T>
vec3_<T> compwise_min(const vec3_<T>& a, const vec3_<T>& b)
{
return vec3_<T>(min(a.x(), b.x()), min(a.y(), b.y()), min(a.z(), b.z()));
}
template <typename T>
vec3_<T> compwise_max(const vec3_<T>& a, const vec3_<T>& b)
{
return vec3_<T>(max(a.x(), b.x()), max(a.y(), b.y()), max(a.z(), b.z()));
}
extern vec3 cross_product(const vec3& a, const vec3& b);
template <typename vector_type>
double dot_product(const vector_type& a, const vector_type& b)
{
double out(0.0);
for (int i = 0; i < vector_type::cardinality(); ++i) {
out += (a[i] * b[i]);
}
return out;
}
// compute a*b^T, which is a matrix
template <typename vector_type>
matrix_t<typename vector_type::element_type> outer_product(const vector_type& a, const vector_type& b)
{
matrix_t<typename vector_type::element_type> out(vector_type::cardinality(), vector_type::cardinality());
const vector_type c0 = a * b[0]; // colmuns
const vector_type c1 = a * b[1];
if (vector_type::cardinality() == 3) {
const vector_type c2 = a * b[2];
out(0, 0) = c0[0];
out(1, 0) = c0[1];
out(2, 0) = c0[2];
out(0, 1) = c1[0];
out(1, 1) = c1[1];
out(2, 1) = c1[2];
out(0, 2) = c2[0];
out(1, 2) = c2[1];
out(2, 2) = c2[2];
} else {
out(0, 0) = c0[0];
out(1, 0) = c0[1];
out(0, 1) = c1[0];
out(1, 1) = c1[1];
}
return out;
}
template <typename vector_type>
double squared_length(const vector_type& v)
{
return dot_product(v, v);
}
template <typename vector_type>
double length(const vector_type& v)
{
return square_root(squared_length(v));
}
template <typename vector_type>
vector_type normalize(const vector_type& v)
{
return v / length(v);
}
double orient2d(const vec2& pa, const vec2& pb, const vec2& pc);
double orient3d(const vec3& pa, const vec3& pb, const vec3& pc,
const vec3& pd);
// Compute a polygon's plane coefficients (i.e. normal and d parameters).
// The computed normal is not normalized. This function returns the largest component of the normal.
int compute_polygon_plane_coefficients(vec3& normal, double& d_coeff,
const vec3* polygon_vertices, const int polygon_vertex_count);
// Compute the intersection point between a line (not a segment) and a plane defined by a polygon.
//
// Parameters:
// 'p' : output intersection point (computed if line does indeed intersect the plane)
// 'q' : first point defining your line
// 'r' : second point defining your line
// 'polygon_vertices' : the vertices of the polygon defineing the plane (assumed to not be degenerate)
// 'polygon_vertex_count' : number of olygon vertices
// 'polygon_normal_max_comp' : largest component of polygon normal.
// 'polygon_plane_normal' : normal of the given polygon
// 'polygon_plane_d_coeff' : the distance coefficient of the plane equation corresponding to the polygon's plane
//
// Return values:
// '0': line is parallel to plane or polygon is degenerate (within available precision)
// '1': an intersection exists.
// 'p': q and r lie in the plane (technically they are parallel to the plane too).
char compute_line_plane_intersection(vec3& p, // intersection point
const vec3& q, const vec3& r, const vec3* polygon_vertices,
const int polygon_vertex_count, const int polygon_normal_max_comp,
const vec3& polygon_plane_normal);
// Test if a line segment intersects with a plane, and yeild the intersection point if so.
//
// Return values:
// 'p': The segment lies wholly within the plane.
// 'q': The(first) q endpoint is on the plane (but not 'p').
// 'r' : The(second) r endpoint is on the plane (but not 'p').
// '0' : The segment lies strictly to one side or the other of the plane.
// '1': The segment intersects the plane, and none of {p, q, r} hold.
char compute_segment_plane_intersection(vec3& p, const vec3& normal, const double& d_coeff,
const vec3& q, const vec3& r);
// Similar to "compute_segment_plane_intersection" but simply checks the [type] of intersection using
// exact arithmetic
//
// Return values:
// 'p': The segment lies wholly within the plane.
// 'q': The(first) q endpoint is on the plane (but not 'p').
// 'r' : The(second) r endpoint is on the plane (but not 'p').
// '0' : The segment lies strictly to one side or the other of the plane.
// '1': The segment intersects the plane, and none of {p, q, r} hold.
char compute_segment_plane_intersection_type(const vec3& q, const vec3& r,
const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal, const int polygon_normal_largest_component);
// Test if a point 'q' (in 2D) lies inside or outside a given polygon (count the number ray crossings).
//
// Return values:
// 'i': q is strictly interior
// 'o': q is strictly exterior (outside).
// 'e': q is on an edge, but not an endpoint.
// 'v': q is a vertex.
char compute_point_in_polygon_test(const vec2& q, const std::vector<vec2>& polygon_vertices);
// Test if a point 'q' (in 3D) lies inside or outside a given polygon (count the number ray crossings).
//
// Return values:
// 'i': q is strictly interior
// 'o': q is strictly exterior (outside).
// 'e': q is on an edge, but not an endpoint.
// 'v': q is a vertex.
char compute_point_in_polygon_test(const vec3& p, const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal, const int polygon_normal_largest_component);
// project a 3d polygon to 3d by eliminating the largest component of its normal
void project_to_2d(std::vector<vec2>& out, const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal, const int polygon_normal_largest_component);
void project_to_2d(std::vector<vec2>& out, const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal);
bool coplaner(const vec3& pa, const vec3& pb, const vec3& pc,
const vec3& pd);
bool collinear(const vec2& a, const vec2& b, const vec2& c, double& predResult);
bool collinear(const vec2& a, const vec2& b, const vec2& c);
/*
Compute the intersection of two line segments. Can also be used to calculate where the respective lines intersect.
Parameters:
'a' and 'b': end points of first segment
'c' and 'd': end points of second segment
'p': the intersection point
's': the parameter for parametric equation of segment a,b (0..1)
't': the parameter for parametric equation of segment c,d (0..1)
Return values:
'e': The segments collinearly overlap, sharing a point; 'e' stands for 'edge.'
'v': An endpoint of one segment is on the other segment, but 'e' doesn't hold; 'v' stands for 'vertex.'
'1': The segments intersect properly (i.e., they share a point and neither 'v' nor 'e' holds); '1' stands for TRUE.
'0': The segments do not intersect (i.e., they share no points); '0' stands for FALSE
*/
char compute_segment_intersection(const vec2& a, const vec2& b, const vec2& c, const vec2& d,
vec2& p, double& s, double& t);
template <typename vector_type>
struct bounding_box_t {
vector_type m_minimum;
vector_type m_maximum;
bounding_box_t(const vector_type& minimum, const vector_type& maximum)
{
m_minimum = minimum;
m_maximum = maximum;
}
bounding_box_t()
{
m_minimum = vector_type(std::numeric_limits<double>::max());
m_maximum = vector_type(-std::numeric_limits<double>::max());
}
inline const vector_type& minimum() const
{
return m_minimum;
}
inline const vector_type& maximum() const
{
return m_maximum;
}
inline void expand(const vector_type& point)
{
m_maximum = compwise_max(m_maximum, point);
m_minimum = compwise_min(m_minimum, point);
}
inline void expand(const bounding_box_t<vector_type>& bbox)
{
m_maximum = compwise_max(m_maximum, bbox.maximum());
m_minimum = compwise_min(m_minimum, bbox.minimum());
}
inline void enlarge(const typename vector_type::element_type& eps_)
{
m_maximum = m_maximum + eps_;
m_minimum = m_minimum - eps_;
}
float SurfaceArea() const
{
vector_type d = m_maximum - m_minimum;
return typename vector_type::element_type(2.0) * (d.x() * d.y() + d.x() * d.z() + d.y() * d.z());
}
int MaximumExtent() const
{
vector_type diag = m_maximum - m_minimum;
if (diag.x() > diag.y() && diag.x() > diag.z())
return 0;
else if (diag.y() > diag.z())
return 1;
else
return 2;
}
};
template <typename T>
inline bool intersect_bounding_boxes(const bounding_box_t<vec3_<T>>& a, const bounding_box_t<vec3_<T>>& b)
{
const vec3_<T>& amin = a.minimum();
const vec3_<T>& amax = a.maximum();
const vec3_<T>& bmin = b.minimum();
const vec3_<T>& bmax = b.maximum();
return (amin.x() <= bmax.x() && amax.x() >= bmin.x()) && //
(amin.y() <= bmax.y() && amax.y() >= bmin.y()) && //
(amin.z() <= bmax.z() && amax.z() >= bmin.z());
}
bool point_in_bounding_box(const vec2& point, const bounding_box_t<vec2>& bbox);
bool point_in_bounding_box(const vec3& point, const bounding_box_t<vec3>& bbox);
template <typename vector_type>
void make_bbox(bounding_box_t<vector_type>& bbox, const vector_type* vertices, const int num_vertices)
{
MCUT_ASSERT(vertices != nullptr);
MCUT_ASSERT(num_vertices >= 3);
for (int i = 0; i < num_vertices; ++i) {
const vector_type& vertex = vertices[i];
bbox.expand(vertex);
}
}
// Shewchuk predicates : shewchuk.c
extern "C" {
// void exactinit();
double orient2d(const double* pa, const double* pb, const double* pc);
double orient3d(const double* pa, const double* pb, const double* pc, const double* pd);
double orient3dfast(const double* pa, const double* pb, const double* pc, const double* pd);
double incircle(const double* pa, const double* pb, const double* pc, const double* pd);
double insphere(const double* pa, const double* pb, const double* pc, const double* pd, const double* pe);
}
#endif // MCUT_MATH_H_

55
src/mcut/include/mcut/internal/preproc.h Normal file
View file

@ -0,0 +1,55 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
/**
* @file mcut.h
* @author Floyd M. Chitalu
* @date 22 July 2022
*
* @brief API-function implementations.
*
* NOTE: This header file defines the frontend implementation of mcDispatch,
* which handles mesh preparation (BVH building, traversal, polygon partitioning
* etc).
*
*/
#ifndef _FRONTEND_INTERSECT_H_
#define _FRONTEND_INTERSECT_H_
#include "mcut/internal/frontend.h"
extern "C" void preproc(
std::shared_ptr<context_t> context_uptr,
McFlags dispatchFlags,
const void* pSrcMeshVertices,
const uint32_t* pSrcMeshFaceIndices,
const uint32_t* pSrcMeshFaceSizes,
uint32_t numSrcMeshVertices,
uint32_t numSrcMeshFaces,
const void* pCutMeshVertices,
const uint32_t* pCutMeshFaceIndices,
const uint32_t* pCutMeshFaceSizes,
uint32_t numCutMeshVertices,
uint32_t numCutMeshFaces) noexcept(false);
#endif // #ifndef _FRONTEND_INTERSECT_H_

87
src/mcut/include/mcut/internal/timer.h Normal file
View file

@ -0,0 +1,87 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#ifndef MCUT_TIMER_H_
#define MCUT_TIMER_H_
#include <chrono>
#include <map>
#include <memory>
#include <stack>
#include <sstream>
#include "mcut/internal/utils.h"
#include "mcut/internal/tpool.h"
//#define PROFILING_BUILD
class mini_timer {
std::chrono::time_point<std::chrono::steady_clock> m_start;
const std::string m_name;
bool m_valid = true;
public:
mini_timer(const std::string& name)
: m_start(std::chrono::steady_clock::now())
, m_name(name)
{
}
~mini_timer()
{
if (m_valid) {
const std::chrono::time_point<std::chrono::steady_clock> now = std::chrono::steady_clock::now();
const std::chrono::milliseconds elapsed = std::chrono::duration_cast<std::chrono::milliseconds>(now - m_start);
unsigned long long elapsed_ = elapsed.count();
log_msg("[MCUT][PROF:" << std::this_thread::get_id() << "]: \"" << m_name << "\" ("<< elapsed_ << "ms)");
}
}
void set_invalid()
{
m_valid = false;
}
};
extern thread_local std::stack<std::unique_ptr<mini_timer>> g_thrd_loc_timerstack;
#if defined(PROFILING_BUILD)
#define TIMESTACK_PUSH(name) \
g_thrd_loc_timerstack.push(std::unique_ptr<mini_timer>(new mini_timer(name)))
#define TIMESTACK_POP() \
g_thrd_loc_timerstack.pop()
#define TIMESTACK_RESET() \
while (!g_thrd_loc_timerstack.empty()) { \
g_thrd_loc_timerstack.top()->set_invalid(); \
g_thrd_loc_timerstack.pop(); \
}
#define SCOPED_TIMER(name) \
mini_timer _1mt(name)
#else
#define SCOPED_TIMER(name)
#define TIMESTACK_PUSH(name)
#define TIMESTACK_POP()
#define TIMESTACK_RESET()
#endif
#endif

1010
src/mcut/include/mcut/internal/tpool.h Normal file
View file

File diff suppressed because it is too large Load diff

260
src/mcut/include/mcut/internal/utils.h Normal file
View file

@ -0,0 +1,260 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#ifndef MCUT_UTILS_H_
#define MCUT_UTILS_H_
// check if c++11 is supported
#if __cplusplus >= 201103L || (defined(_MSC_VER) && _MSC_VER >= 1900)
#define MCUT_CXX11_IS_SUPPORTED
#elif !defined(__cplusplus) && !defined(_MSC_VER)
typedef char couldnt_parse_cxx_standard[-1]; ///< Error: couldn't parse standard
#endif
#if defined(_WIN64) || defined(_WIN32)
#define MCUT_BUILD_WINDOWS 1
#elif defined(__APPLE__)
#define MCUT_BUILD_APPLE 1
#elif defined(__linux__) || defined(__unix__)
#define MCUT_BUILD_LINUX 1
#endif // #if defined(_WIN64) || defined(_WIN32)
#define MCUT_MAKE_STRING__(s) #s
#define MCUT_MAKE_STRING_(s) MCUT_MAKE_STRING__(s)
#ifndef NDEBUG
#define MCUT_DEBUG_BUILD 1
#endif
// debug macros
#if defined(MCUT_DEBUG_BUILD)
//
// MCUT_DEBUG_BREAKPOINT
//
#if defined(MCUT_BUILD_WINDOWS)
#define MCUT_DEBUG_BREAKPOINT_() __debugbreak()
#else // #if defined(MCUT_BUILD_WINDOWS)
#define MCUT_DEBUG_BREAKPOINT_() std::abort()
#endif // #if defined(MCUT_BUILD_WINDOWS)
//
// MCUT_ASSERT
//
#define MCUT_ASSERT(a) \
do { \
if (false == (a)) { \
std::fprintf(stderr, \
"Assertion failed: %s, " \
"%d at \'%s\'\n", \
__FILE__, \
__LINE__, \
MCUT_MAKE_STRING_(a)); \
MCUT_DEBUG_BREAKPOINT_(); \
} \
} while (0)
#define DEBUG_CODE_MASK(code) code
#else // #if defined(MCUT_DEBUG_BUILD)
//
// MCUT_ASSERT
//
#define MCUT_ASSERT(a) // do nothing
#define DEBUG_CODE_MASK(code) // do nothing
#endif // #if defined(MCUT_DEBUG_BUILD)
#include <fstream>
#include <iostream>
#include <sstream>
#if MCUT_BUILD_WINDOWS
#define EXCEPTION_THROWN throw()
#else
#define EXCEPTION_THROWN
#endif
#define PEDANTIC_SUBSCRIPT_ACCESS 1
#if defined(PEDANTIC_SUBSCRIPT_ACCESS)
#define SAFE_ACCESS(var, i) var.at(i)
#else
#define SAFE_ACCESS(var, i) var[i]
#endif
static inline int wrap_integer(int x, const int lo, const int hi)
{
const int range_size = hi - lo + 1;
if (x < lo) {
x += range_size * ((lo - x) / range_size + 1);
}
return lo + (x - lo) % range_size;
}
class logger_t {
std::stringstream m_buffer;
bool m_verbose;
std::string m_prepend;
std::string m_reason_for_failure;
public:
typedef std::ostream& (*ManipFn)(std::ostream&);
typedef std::ios_base& (*FlagsFn)(std::ios_base&);
logger_t()
: m_buffer()
, m_verbose(false)
, m_prepend()
, m_reason_for_failure()
{
}
logger_t(const logger_t& other) = delete;
logger_t& operator=(const logger_t& other) = delete;
~logger_t()
{
}
std::string get_log_string()
{
return m_buffer.str();
}
void set_reason_for_failure(const std::string& msg)
{
if (m_reason_for_failure.empty()) // NOTE
m_reason_for_failure = msg;
}
std::string get_reason_for_failure()
{
std::string s(m_reason_for_failure); // copy
return s;
}
inline bool verbose()
{
return m_verbose;
}
inline void set_verbose(bool b)
{
m_verbose = b;
}
inline void reset()
{
m_prepend.clear();
}
inline void indent()
{
if (!verbose()) {
return;
}
m_prepend.append(" ");
}
inline void unindent()
{
if (!verbose()) {
return;
}
m_prepend.pop_back();
m_prepend.pop_back();
}
template <class T> // int, double, strings, etc
inline logger_t& operator<<(const T& output)
{
if (verbose()) {
m_buffer << output;
}
return *this;
}
inline logger_t& operator<<(ManipFn manip) /// endl, flush, setw, setfill, etc.
{
if (verbose()) {
manip(m_buffer);
if (manip == static_cast<ManipFn>(std::flush) || manip == static_cast<ManipFn>(std::endl)) {
this->flush();
}
}
return *this;
}
inline logger_t& operator<<(FlagsFn manip) /// setiosflags, resetiosflags
{
if (verbose()) {
manip(m_buffer);
}
return *this;
}
inline void flush()
{
if (!(verbose())) {
return;
}
#if 0 // dump log to terminal [immediately]
std::cout << m_prepend << "::" << m_buffer.str();
m_buffer.str(std::string());
m_buffer.clear();
#endif
}
};
template <typename T>
struct pair : std::pair<T, T> {
pair(const T a, const T b)
: std::pair<T, T>(a < b ? a : b, a < b ? b : a)
{
}
};
template <typename T>
pair<T> make_pair(const T a, const T b)
{
return pair<T>(a, b);
}
// Threadsafe logging to console which prevents std::cerr from mixing strings when
// concatenating with the operator<< multiple time per string, across multiple
// threads.
#define log_msg(msg_str) \
{ \
std::stringstream ss; \
ss << msg_str << std::endl; \
std::cerr << ss.str(); \
}
// used to marked/label unused function parameters to prevent warnings
#define UNUSED(x) [&x] {}()
#endif // MCUT_UTILS_H_

1265
src/mcut/include/mcut/mcut.h Normal file
View file

File diff suppressed because it is too large Load diff

109
src/mcut/include/mcut/platform.h Normal file
View file

@ -0,0 +1,109 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
/**
* @file platform.h
* @author Floyd M. Chitalu
* @date 1 Jan 2021
* @brief File containing platform-specific types and definitions.
*
* This header file defines platform specific directives and integral type
* declarations.
* Platforms should define these so that MCUT users call MCUT commands
* with the same calling conventions that the MCUT implementation expects.
*
* MCAPI_ATTR - Placed before the return type in function declarations.
* Useful for C++11 and GCC/Clang-style function attribute syntax.
* MCAPI_CALL - Placed after the return type in function declarations.
* Useful for MSVC-style calling convention syntax.
* MCAPI_PTR - Placed between the '(' and '*' in function pointer types.
*
* Function declaration: MCAPI_ATTR void MCAPI_CALL mcFunction(void);
* Function pointer type: typedef void (MCAPI_PTR *PFN_mcFunction)(void);
*
*/
#ifndef MC_PLATFORM_H_
#define MC_PLATFORM_H_
#ifdef __cplusplus
extern "C"
{
#endif // __cplusplus
/** @file */
#if defined(_WIN32)
#ifdef MCUT_SHARED_LIB
#if defined(MCUT_EXPORT_SHARED_LIB_SYMBOLS)
//** Symbol visibilty */
#define MCAPI_ATTR __declspec(dllexport)
#else
//** Symbol visibilty */
#define MCAPI_ATTR __declspec(dllimport)
#endif
#else // MCUT_SHARED_LIB
//** Symbol visibilty */
#define MCAPI_ATTR
#endif // MCUT_SHARED_LIB
//** Function calling convention */
#define MCAPI_CALL __stdcall
//** Function pointer-type declaration helper */
#define MCAPI_PTR MCAPI_CALL
#else // #if defined(_WIN32)
//** Symbol visibilty */
#define MCAPI_ATTR __attribute__((visibility("default")))
//** Function calling convention */
#define MCAPI_CALL
//** Function pointer-type declaration helper */
#define MCAPI_PTR
#endif // #if defined(_WIN32)
#include <stddef.h> // standard type definitions
#if !defined(MC_NO_STDINT_H)
#if defined(_MSC_VER) && (_MSC_VER < 1600)
typedef signed __int8 int8_t;
typedef unsigned __int8 uint8_t;
typedef signed __int16 int16_t;
typedef unsigned __int16 uint16_t;
typedef signed __int32 int32_t;
typedef unsigned __int32 uint32_t;
typedef signed __int64 int64_t;
typedef unsigned __int64 uint64_t;
#else
#include <stdint.h> // integer types
#endif
#endif // !defined(MC_NO_STDINT_H)
#ifdef __cplusplus
} // extern "C"
#endif // __cplusplus
#endif

1148
src/mcut/source/bvh.cpp Normal file
View file

File diff suppressed because it is too large Load diff

3058
src/mcut/source/frontend.cpp Normal file
View file

File diff suppressed because it is too large Load diff

1411
src/mcut/source/hmesh.cpp Normal file
View file

File diff suppressed because it is too large Load diff

10634
src/mcut/source/kernel.cpp Normal file
View file

File diff suppressed because it is too large Load diff

903
src/mcut/source/math.cpp Normal file
View file

@ -0,0 +1,903 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#include "mcut/internal/math.h"
#include <algorithm> // std::sort
#include <cstdlib>
#include <tuple> // std::make_tuple std::get<>
double square_root(const double& number)
{
#if !defined(MCUT_WITH_ARBITRARY_PRECISION_NUMBERS)
return std::sqrt(number);
#else
arbitrary_precision_number_t out(number);
mpfr_sqrt(out.get_mpfr_handle(), number.get_mpfr_handle(), arbitrary_precision_number_t::get_default_rounding_mode());
return out;
#endif // #if !defined(MCUT_WITH_ARBITRARY_PRECISION_NUMBERS)
}
double absolute_value(const double& number)
{
#if !defined(MCUT_WITH_ARBITRARY_PRECISION_NUMBERS)
return std::fabs(number);
#else
double out(number);
mpfr_abs(out.get_mpfr_handle(), number.get_mpfr_handle(), arbitrary_precision_number_t::get_default_rounding_mode());
return out;
#endif // #if defined(MCUT_WITH_ARBITRARY_PRECISION_NUMBERS)
}
sign_t sign(const double& number)
{
#if !defined(MCUT_WITH_ARBITRARY_PRECISION_NUMBERS)
int s = (double(0) < number) - (number < double(0));
sign_t result = sign_t::ZERO;
if (s > 0) {
result = sign_t::POSITIVE;
} else if (s < 0) {
result = sign_t::NEGATIVE;
}
return result;
#else
double out(number);
int s = mpfr_sgn(number.get_mpfr_handle());
sign_t result = sign_t::ZERO;
if (s > 0) {
result = sign_t::POSITIVE;
} else if (s < 0) {
result = sign_t::NEGATIVE;
}
return result;
#endif // #if defined(MCUT_WITH_ARBITRARY_PRECISION_NUMBERS)
}
std::ostream& operator<<(std::ostream& os, const vec3& v)
{
return os << static_cast<double>(v.x()) << ", " << static_cast<double>(v.y()) << ", " << static_cast<double>(v.z());
}
bool operator==(const vec3& a, const vec3& b)
{
return (a.x() == b.x()) && (a.y() == b.y()) && (a.z() == b.z());
}
vec3 cross_product(const vec3& a, const vec3& b)
{
return vec3(
a.y() * b.z() - a.z() * b.y(),
a.z() * b.x() - a.x() * b.z(),
a.x() * b.y() - a.y() * b.x());
}
double orient2d(const vec2& pa, const vec2& pb, const vec2& pc)
{
const double pa_[2] = { static_cast<double>(pa.x()), static_cast<double>(pa.y()) };
const double pb_[2] = { static_cast<double>(pb.x()), static_cast<double>(pb.y()) };
const double pc_[2] = { static_cast<double>(pc.x()), static_cast<double>(pc.y()) };
return ::orient2d(pa_, pb_, pc_); // shewchuk predicate
}
double orient3d(const vec3& pa, const vec3& pb, const vec3& pc,
const vec3& pd)
{
const double pa_[3] = { static_cast<double>(pa.x()), static_cast<double>(pa.y()), static_cast<double>(pa.z()) };
const double pb_[3] = { static_cast<double>(pb.x()), static_cast<double>(pb.y()), static_cast<double>(pb.z()) };
const double pc_[3] = { static_cast<double>(pc.x()), static_cast<double>(pc.y()), static_cast<double>(pc.z()) };
const double pd_[3] = { static_cast<double>(pd.x()), static_cast<double>(pd.y()), static_cast<double>(pd.z()) };
return ::orient3d(pa_, pb_, pc_, pd_); // shewchuk predicate
}
#if 0
void polygon_normal(vec3& normal, const vec3* vertices, const int num_vertices)
{
normal = vec3(0.0);
for (int i = 0; i < num_vertices; ++i) {
normal = normal + cross_product(vertices[i] - vertices[0], vertices[(i + 1) % num_vertices] - vertices[0]);
}
normal = normalize(normal);
}
#endif
int compute_polygon_plane_coefficients(vec3& normal, double& d_coeff,
const vec3* polygon_vertices, const int polygon_vertex_count)
{
// compute polygon normal (http://cs.haifa.ac.il/~gordon/plane.pdf)
normal = vec3(0.0);
for (int i = 0; i < polygon_vertex_count - 1; ++i) {
normal = normal + cross_product(polygon_vertices[i] - polygon_vertices[0], polygon_vertices[(i + 1) % polygon_vertex_count] - polygon_vertices[0]);
}
// In our calculations we need the normal be be of unit length so that the d-coeff
// represents the distance of the plane from the origin
// normal = normalize(normal);
d_coeff = dot_product(polygon_vertices[0], normal);
double largest_component(0.0);
double tmp(0.0);
int largest_component_idx = 0;
for (int i = 0; i < 3; ++i) {
tmp = absolute_value(normal[i]);
if (tmp > largest_component) {
largest_component = tmp;
largest_component_idx = i;
}
}
return largest_component_idx;
}
// Intersect a line segment with a plane
//
// Return values:
// 'p': The segment lies wholly within the plane.
// 'q': The(first) q endpoint is on the plane (but not 'p').
// 'r' : The(second) r endpoint is on the plane (but not 'p').
// '0' : The segment lies strictly to one side or the other of the plane.
// '1': The segment intersects the plane, and none of {p, q, r} hold.
char compute_segment_plane_intersection(vec3& p, const vec3& normal, const double& d_coeff,
const vec3& q, const vec3& r)
{
// vec3 planeNormal;
// double planeDCoeff;
// planeNormalLargestComponent = compute_polygon_plane_coefficients(planeNormal, planeDCoeff, polygonVertices,
// polygonVertexCount);
double num = d_coeff - dot_product(q, normal);
const vec3 rq = (r - q);
double denom = dot_product(rq, normal);
if (denom == double(0.0) /* Segment is parallel to plane.*/) {
if (num == double(0.0)) { // 'q' is on plane.
return 'p'; // The segment lies wholly within the plane
} else {
return '0';
}
}
double t = num / denom;
for (int i = 0; i < 3; ++i) {
p[i] = q[i] + t * (r[i] - q[i]);
}
if ((double(0.0) < t) && (t < double(1.0))) {
return '1'; // The segment intersects the plane, and none of {p, q, r} hold
} else if (num == double(0.0)) // t==0
{
return 'q'; // The (first) q endpoint is on the plane (but not 'p').
} else if (num == denom) // t==1
{
return 'r'; // The (second) r endpoint is on the plane (but not 'p').
} else {
return '0'; // The segment lies strictly to one side or the other of the plane
}
}
bool determine_three_noncollinear_vertices(int& i, int& j, int& k, const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal, const int polygon_normal_largest_component)
{
const int polygon_vertex_count = (int)polygon_vertices.size();
MCUT_ASSERT(polygon_vertex_count >= 3);
std::vector<vec2> x;
project_to_2d(x, polygon_vertices, polygon_normal, polygon_normal_largest_component);
MCUT_ASSERT(x.size() == (size_t)polygon_vertex_count);
/*
NOTE: We cannot just use _any_/the first result of "colinear(x[i], x[j], x[k])" which returns true since
any three points that are _nearly_ colinear (in the limit of floating point precision)
will be found to be not colinear when using the exact predicate "orient2d".
This would have implications "further down the line", where e.g. the three nearly colinear points
are determined to be non-colinear (in the exact sense) but using them to evaluate
operations like segment-plane intersection would then give a false positive.
To overcome this, must find the three vertices of the polygon, which maximise non-colinearity.
NOTE: if the polygon has 3 vertices and they are indeed nearly colinear then answer is determined by the
exact predicate (i.e. i j k = 0 1 2)
*/
// get any three vertices that are not collinear
i = 0;
j = i + 1;
k = j + 1;
std::vector<std::tuple<int, int, int, double>> non_colinear_triplets;
for (; i < polygon_vertex_count; ++i) {
for (; j < polygon_vertex_count; ++j) {
for (; k < polygon_vertex_count; ++k) {
double predRes;
if (!collinear(x[i], x[j], x[k], predRes)) {
non_colinear_triplets.emplace_back(std::make_tuple(i, j, k, predRes));
}
}
}
}
std::sort(non_colinear_triplets.begin(), non_colinear_triplets.end(),
[](const std::tuple<int, int, int, double>& a, const std::tuple<int, int, int, double>& b) {
return std::fabs(std::get<3>(a)) > std::fabs(std::get<3>(b));
});
std::tuple<int, int, int, double> best_triplet = non_colinear_triplets.front(); // maximising non-colinearity
i = std::get<0>(best_triplet);
j = std::get<1>(best_triplet);
k = std::get<2>(best_triplet);
return !non_colinear_triplets.empty(); // need at least one non-colinear triplet
}
char compute_segment_plane_intersection_type(const vec3& q, const vec3& r,
const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal,
const int polygon_normal_largest_component)
{
// TODO: we could also return i,j and k so that "determine_three_noncollinear_vertices" is not called multiple times,
// which we do to determine the type of intersection and the actual intersection point
const int polygon_vertex_count = (int)polygon_vertices.size();
// ... any three vertices that are not collinear
int i = 0;
int j = 1;
int k = 2;
if (polygon_vertex_count > 3) { // case where we'd have the possibility of noncollinearity
bool b = determine_three_noncollinear_vertices(i, j, k, polygon_vertices, polygon_normal,
polygon_normal_largest_component);
if (!b) {
return '0'; // all polygon points are collinear
}
}
double qRes = orient3d(polygon_vertices[i], polygon_vertices[j], polygon_vertices[k], q);
double rRes = orient3d(polygon_vertices[i], polygon_vertices[j], polygon_vertices[k], r);
if (qRes == double(0.0) && rRes == double(0.0)) {
return 'p';
} else if (qRes == double(0.0)) {
return 'q';
} else if (rRes == double(0.0)) {
return 'r';
} else if ((rRes < double(0.0) && qRes < double(0.0)) || (rRes > double(0.0) && qRes > double(0.0))) {
return '0';
} else {
return '1';
}
}
char compute_segment_line_plane_intersection_type(const vec3& q, const vec3& r,
const std::vector<vec3>& polygon_vertices,
const int polygon_normal_max_comp,
const vec3& polygon_plane_normal)
{
const int polygon_vertex_count = (int)polygon_vertices.size();
// ... any three vertices that are not collinear
int i = 0;
int j = 1;
int k = 2;
if (polygon_vertex_count > 3) { // case where we'd have the possibility of noncollinearity
bool b = determine_three_noncollinear_vertices(i, j, k, polygon_vertices, polygon_plane_normal,
polygon_normal_max_comp);
if (!b) {
return '0'; // all polygon points are collinear
}
}
double qRes = orient3d(polygon_vertices[i], polygon_vertices[j], polygon_vertices[k], q);
double rRes = orient3d(polygon_vertices[i], polygon_vertices[j], polygon_vertices[k], r);
if (qRes == double(0.0) && rRes == double(0.0)) {
// both points used to define line lie on plane therefore we have an in-plane intersection
// or the polygon is a degenerate triangle
return 'p';
} else {
if ((rRes < double(0.0) && qRes < double(0.0)) || (rRes > double(0.0) && qRes > double(0.0))) { // both points used to define line lie on same side of plane
// check if line is parallel to plane
// const double num = polygon_plane_d_coeff - dot_product(q, polygon_plane_normal);
const vec3 rq = (r - q);
const double denom = dot_product(rq, polygon_plane_normal);
if (denom == double(0.0) /* Segment is parallel to plane.*/) {
// MCUT_ASSERT(num != 0.0); // implies 'q' is on plane (see: "compute_segment_plane_intersection(...)")
// but we have already established that q and r are on same side.
return '0';
}
}
}
// q and r are on difference sides of the plane, therefore we have an intersection
return '1';
}
// Compute the intersection point between a line (not a segment) and a plane defined by a polygon.
//
// Parameters:
// 'p' : output intersection point (computed if line does indeed intersect the plane)
// 'q' : first point defining your line
// 'r' : second point defining your line
// 'polygon_vertices' : the vertices of the polygon defineing the plane (assumed to not be degenerate)
// 'polygon_vertex_count' : number of olygon vertices
// 'polygon_normal_max_comp' : largest component of polygon normal.
// 'polygon_plane_normal' : normal of the given polygon
// 'polygon_plane_d_coeff' : the distance coefficient of the plane equation corresponding to the polygon's plane
//
// Return values:
// '0': line is parallel to plane (or polygon is degenerate ... within available precision)
// '1': an intersection exists.
// 'p': q and r lie in the plane (technically they are parallel to the plane too but we need to report this because it
// violates GP).
char compute_line_plane_intersection(vec3& p, // intersection point
const vec3& q, const vec3& r,
const std::vector<vec3>& polygon_vertices, const int polygon_normal_max_comp,
const vec3& polygon_plane_normal,
const double& polygon_plane_d_coeff)
{
const int polygon_vertex_count = (int)polygon_vertices.size();
// ... any three vertices that are not collinear
int i = 0;
int j = 1;
int k = 2;
if (polygon_vertex_count > 3) { // case where we'd have the possibility of noncollinearity
bool b = determine_three_noncollinear_vertices(i, j, k, polygon_vertices, polygon_plane_normal,
polygon_normal_max_comp);
if (!b) {
return '0'; // all polygon points are collinear
}
}
double qRes = orient3d(polygon_vertices[i], polygon_vertices[j], polygon_vertices[k], q);
double rRes = orient3d(polygon_vertices[i], polygon_vertices[j], polygon_vertices[k], r);
if (qRes == double(0.) && rRes == double(0.)) {
return 'p'; // both points used to define line lie on plane therefore we have an in-plane intersection
} else {
const double num = polygon_plane_d_coeff - dot_product(q, polygon_plane_normal);
const vec3 rq = (r - q);
const double denom = dot_product(rq, polygon_plane_normal);
if ((rRes < double(0.) && qRes < double(0.)) || (rRes > double(0.) && qRes > double(0.))) { // both q an r are on same side of plane
if (denom == double(0.) /* line is parallel to plane.*/) {
return '0';
}
}
// q and r are on difference sides of the plane, therefore we have an intersection
// compute the intersection point
const double t = num / denom;
for (int it = 0; it < 3; ++it) {
p[it] = q[it] + t * (r[it] - q[it]);
}
return '1';
}
}
// Count the number ray crossings to determine if a point 'q' lies inside or outside a given polygon.
//
// Return values:
// 'i': q is strictly interior
// 'o': q is strictly exterior (outside).
// 'e': q is on an edge, but not an endpoint.
// 'v': q is a vertex.
char compute_point_in_polygon_test(const vec2& q, const std::vector<vec2>& polygon_vertices)
{
const int polygon_vertex_count = (int)polygon_vertices.size();
#if 0 // Algorithm 1 in :
// https://pdf.sciencedirectassets.com/271512/1-s2.0-S0925772100X00545/1-s2.0-S0925772101000128/main.pdf?X-Amz-Security-Token=IQoJb3JpZ2luX2VjEAEaCXVzLWVhc3QtMSJHMEUCIBr6Fu%2F%2FtscErn%2Fl4pn2UGNA45sAw9vggQscK7Tnl0ssAiEAzfnyy4B4%2BXkZp8xcEk7utDBrxgmGyH7pyqk0efKOUEoq%2BgMIKhAEGgwwNTkwMDM1NDY4NjUiDHQT4kOfk4%2B6kBB0xCrXAzd8SWlnELJbM5l93HSvv4nUgtO85%2FGyx5%2BOoHmYoUuVwJjCXChmLJmlz%2BxYUQE%2Fr8vQa2hPlUEPfiTVGgoHaK8NkSMP6LRs%2F3WjyZ9OxzHbqSwZixNceW34OAkq0E1QgYLdGVjpPKxNK1haXDfBTMN6bF2IOU9dKi9wTv3uTlep0kHDa1BNNl6M6yZk5QlF2bPF9XmNjjZCpFQLhr%2BPoHo%2Bx4xy39aH8hCkkTqGdy2KrKGN6lv0%2FduIaItyZfqalYS%2BwW6cll2F5G11g0tSu7yKu6mug94LUTzsRmzD0UhzeGl2WV6Ev2qhw26mwFEKgnTMqGst8XAHjFjjAQyMzXNHCQqNBRIevHIzVWuUY4QZMylSRsodo0dfwwCFhzl0%2BJA1ZXb0%2BoyB8c11meQBO8FpMSshngNcbiAYUtIOFcHNCTCUMk0JPOZ%2FxvANsopnivbrPARL71zU4PaTujg5jfC2zoO6ZDUL8E6Vn%2BNtfb%2BuQV7DwtIH51Bv%2F%2F1m6h5mjbpRiGYcH%2F5SDD%2Fk%2BpHfKRQzKu7jJc%2B0XO0bQvoLSrAte0Qk10PwOvDl5jMIMdmxTBDDiDGErRykYxMQhq5EwjyiWPXzM3ll9uK59Vy0bAEj5Qemj5by1jCDht6IBjqlAV4okAPQO5wWdWojCYeKvluKfXCvudrUxrLhsxb7%2BTZNMODTG%2Fn%2Fbw875Yr6fOQ42u40pOsnPB%2FTP3cWwjiB%2BEHzDqN8AhCVQyoedw7QrU3OBWlSl6lB%2BGLAVqrhcizgFUiM2nj3HaVP2m7S%2FXqpv%2FoWlEXt4gR8iI9XsIlh6L6SBE22FqbsU5ewCxXaqip19VXhAGnlvjTihXUg6yZGWhExHj%2BKcA%3D%3D&X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Date=20210814T103958Z&X-Amz-SignedHeaders=host&X-Amz-Expires=300&X-Amz-Credential=ASIAQ3PHCVTY3PPQSW2L%2F20210814%2Fus-east-1%2Fs3%2Faws4_request&X-Amz-Signature=94403dce5c12ede9507e9612c000db5772607e59a2134d30d4e5b44212056dc7&hash=cada083b165f9786e6042e21d3d22147ec08b1e2c8fa2572ceeb2445cb5e730a&host=68042c943591013ac2b2430a89b270f6af2c76d8dfd086a07176afe7c76c2c61&pii=S0925772101000128&tid=spdf-f9e7c4b7-808d-4a8e-a474-6430ff687798&sid=8dddff5b32fff2499908bf6501965aec3d0fgxrqb&type=client
const double zero(0.0);
int i = 0;
int k = 0;
double f = zero;
double u1 = 0.0;
double v1 = 0.0;
double u2 = 0.0;
double v2 = 0.0;
const int n = polygon_vertex_count;
for (i = 0; i < n; ++i)
{
v1 = polygon_vertices[i].y() - q.y();
v2 = polygon_vertices[(i + 1) % n].y() - q.y();
if ((v1 < zero && v2 < zero) || (v1 > zero && v2 > zero))
{
continue;
}
u1 = polygon_vertices[i].x() - q.x();
u2 = polygon_vertices[(i + 1) % n].x() - q.x();
if (v2 > zero && v1 <= zero)
{
f = u1 * v2 - u2 * v1;
if (f > zero)
{
k = k + 1;
}
else if (f == zero)
{
return 'e'; //-1; // boundary
}
}
else if (v1 > zero && v2 <= zero)
{
f = u1 * v2 - u2 * v1;
if (f < zero)
{
k = k + 1;
}
else if (f == zero)
{
return 'e'; //-1; // boundary
}
}
else if (v2 == zero && v1 < zero)
{
f = u1 * v2 - u2 * v1;
if (f == zero)
{
return 'e'; //-1; // // boundary
}
}
else if (v1 == zero && v2 < zero)
{
f = u1 * v2 - u2 * v1;
if (f == zero)
{
return 'e'; //-1; // boundary
}
}
else if (v1 == zero && v2 == zero)
{
if (u2 <= zero && u1 >= zero)
{
return 'e'; //-1; // boundary
}
else if (u1 <= zero && u2 >= zero)
{
return 'e'; //-1; // boundary
}
}
}
if (k % 2 == 0)
{
return 'o'; //0; // outside
}
else
{
return 'i'; //1; // inside
}
#else // http://www.science.smith.edu/~jorourke/books/compgeom.html
std::vector<vec2> vertices(polygon_vertex_count, vec2());
// Shift so that q is the origin. Note this destroys the polygon.
// This is done for pedagogical clarity.
for (int i = 0; i < polygon_vertex_count; i++) {
for (int d = 0; d < 2; d++) {
const double& a = polygon_vertices[i][d];
const double& b = q[d];
double& c = vertices[i][d];
c = a - b;
}
}
int Rcross = 0; /* number of right edge/ray crossings */
int Lcross = 0; /* number ofleft edge/ray crossings */
/* For each edge e = (i<>lj), see if crosses ray. */
for (int i = 0; i < polygon_vertex_count; i++) {
/* First check if q = (0, 0) is a vertex. */
if (vertices[i].x() == double(0.) && vertices[i].y() == double(0.)) {
return 'v';
}
int il = (i + polygon_vertex_count - 1) % polygon_vertex_count;
// Check if e straddles x axis, with bias above/below.
// Rstrad is TRUE iff one endpoint of e is strictly above the x axis and the other is not (i.e., the other is on
// or below)
bool Rstrad = (vertices[i].y() > double(0.)) != (vertices[il].y() > double(0.));
bool Lstrad = (vertices[i].y() < double(0.)) != (vertices[il].y() < double(0.));
if (Rstrad || Lstrad) {
/* Compute intersection of e with x axis. */
// The computation of x is needed whenever either of these straddle variables is TRUE, which
// only excludes edges passing through q = (0, 0) (and incidentally protects against division by 0).
double x = (vertices[i].x() * vertices[il].y() - vertices[il].x() * vertices[i].y()) / (vertices[il].y() - vertices[i].y());
if (Rstrad && x > double(0.)) {
Rcross++;
}
if (Lstrad && x < double(0.)) {
Lcross++;
}
} /* end straddle computation*/
} // end for
/* q on an edge if L/Rcross counts are not the same parity.*/
if ((Rcross % 2) != (Lcross % 2)) {
return 'e';
}
/* q inside iff an odd number of crossings. */
if ((Rcross % 2) == 1) {
return 'i';
} else {
return 'o';
}
#endif
}
// given a normal vector (not necessarily normalized) and its largest component, calculate
// a matrix P that will project any 3D vertex to 2D by removing the 3D component that
// corresponds to the largest component of the normal..
// https://math.stackexchange.com/questions/180418/calculate-rotation-matrix-to-align-vector-a-to-vector-b-in-3d/2672702#2672702
matrix_t<double> calculate_projection_matrix(const vec3& polygon_normal,
const int polygon_normal_largest_component)
{
MCUT_ASSERT(squared_length(polygon_normal) > double(0.0));
MCUT_ASSERT(polygon_normal_largest_component >= 0);
MCUT_ASSERT(polygon_normal_largest_component <= 2);
// unit length normal vector of polygon
const vec3 a = normalize(polygon_normal);
// unit length basis vector corresponding to the largest component of the normal vector
const vec3 b = [&]() {
vec3 x(double(0.0));
const sign_t s = sign(polygon_normal[polygon_normal_largest_component]);
MCUT_ASSERT(s != sign_t::ZERO); // implies that the normal vector has a magnitude of zero
// The largest component of the normal is the one we will "remove"
// NOTE: we multiple by the sign here to ensure that "a_plus_b" below is not zero when a == b
x[polygon_normal_largest_component] = double((1.0) * static_cast<long double>(s));
return x;
}();
matrix_t<double> I(3, 3); // 3x3 identity
I(0, 0) = 1.0;
I(1, 1) = 1.0;
I(2, 2) = 1.0;
matrix_t<double> R = I;
// NOTE: While this will map vectors a to b, it can add a lot of unnecessary "twist".
// For example, if a=b=e_{z} this formula will produce a 180-degree rotation about the z-axis rather
// than the identity one might expect.
if ((a[0] != b[0]) || (a[1] != b[1]) || (a[2] != b[2])) // a != b
{
const vec3 a_plus_b = a + b;
const double a_dot_b = dot_product(a, b);
// this will never be zero because we set 'b' as the canonical basis vector
// that has the largest projection onto the normal
MCUT_ASSERT(a_dot_b != double(0.0));
const matrix_t<double> outer = outer_product(a_plus_b, a_plus_b);
// compute the 3x3 rotation matrix R to orient the polygon into the canonical axes "i" and "j",
// where "i" and "j" != "polygon_normal_largest_component"
R = ((outer / a_dot_b) * 2.0) - I; // rotation
}
// compute the 2x3 selector matrix K to select the polygon-vertex components that correspond
// to canonical axes "i" and "j".
matrix_t<double> K = matrix_t<double>(2, 3);
if (polygon_normal_largest_component == 0) // project by removing x-component
{
// 1st row
K(0, 0) = 0.0; // col 0
K(0, 1) = 1.0; // col 1
K(0, 2) = 0.0; // col 2
// 2nd row
K(1, 0) = 0.0; // col 0
K(1, 1) = 0.0; // col 1
K(1, 2) = 1.0; // col 2
} else if (polygon_normal_largest_component == 1) // project by removing y-component
{
// 1st row
K(0, 0) = 1.0; // col 0
K(0, 1) = 0.0; // col 1
K(0, 2) = 0.0; // col 2
// 2nd row
K(1, 0) = 0.0; // col 0
K(1, 1) = 0.0; // col 1
K(1, 2) = 1.0; // col 2
} else if (polygon_normal_largest_component == 2) // project by removing z-component
{
// 1st row
K(0, 0) = 1.0; // col 0
K(0, 1) = 0.0; // col 1
K(0, 2) = 0.0; // col 2
// 2nd row
K(1, 0) = 0.0; // col 0
K(1, 1) = 1.0; // col 1
K(1, 2) = 0.0; // col 2
}
return K * R;
}
// https://answers.unity.com/questions/1522620/converting-a-3d-polygon-into-a-2d-polygon.html
void project_to_2d(std::vector<vec2>& out, const std::vector<vec3>& polygon_vertices, const vec3& polygon_normal)
{
const uint32_t N = polygon_vertices.size();
out.resize(N);
const vec3 normal = normalize(polygon_normal);
// first unit vector on plane (rotation by 90 degrees)
vec3 u(-normal.y(), normal.x(), normal.z());
vec3 v = normalize(cross_product(u, normal));
for(uint32_t i =0; i < N; ++i)
{
const vec3 point = polygon_vertices[i];
const vec2 projected(
dot_product(point, u),
dot_product(point, v)
);
out[i] = projected;
}
}
void project_to_2d(std::vector<vec2>& out, const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal, const int polygon_normal_largest_component)
{
const int polygon_vertex_count = (int)polygon_vertices.size();
out.clear();
out.resize(polygon_vertex_count);
#if 1
// 3x3 matrix for projecting a point to 2D
matrix_t<double> P = calculate_projection_matrix(polygon_normal, polygon_normal_largest_component);
for (int i = 0; i < polygon_vertex_count; ++i) { // for each vertex
const vec3& x = polygon_vertices[i];
out[i] = P * x; // vertex in xz plane
}
#else // This code is not reliable because it shadow-projects a polygons which can skew computations
for (int i = 0; i < polygon_vertex_count; ++i) { // for each vertex
vec2& Tp = out[i];
int k = 0;
for (int j = 0; j < 3; j++) { // for each component
if (j != polygon_plane_normal_largest_component) { /* skip largest coordinate */
Tp[k] = polygon_vertices[i][j];
k++;
}
}
}
#endif
}
// TODO: update this function to use "project_to_2d" for projection step
char compute_point_in_polygon_test(const vec3& p, const std::vector<vec3>& polygon_vertices,
const vec3& polygon_normal, const int polygon_normal_largest_component)
{
const int polygon_vertex_count = (int)polygon_vertices.size();
/* Project out coordinate m in both p and the triangular face */
vec2 pp; /*projected p */
#if 0
int k = 0;
for (int j = 0; j < 3; j++)
{ // for each component
if (j != polygon_plane_normal_largest_component)
{ /* skip largest coordinate */
pp[k] = p[j];
k++;
}
}
#endif
const matrix_t<double> P = calculate_projection_matrix(polygon_normal, polygon_normal_largest_component);
pp = P * p;
std::vector<vec2> polygon_vertices2d(polygon_vertex_count, vec2());
for (int i = 0; i < polygon_vertex_count; ++i) { // for each vertex
const vec3& x = polygon_vertices[i];
polygon_vertices2d[i] = P * x; // vertex in xz plane
}
#if 0
for (int i = 0; i < polygon_vertex_count; ++i)
{ // for each vertex
vec2 &Tp = polygon_vertices2d[i];
k = 0;
for (int j = 0; j < 3; j++)
{ // for each component
if (j != polygon_plane_normal_largest_component)
{ /* skip largest coordinate */
Tp[k] = polygon_vertices[i][j];
k++;
}
}
}
#endif
return compute_point_in_polygon_test(pp, polygon_vertices2d);
}
inline bool Between(const vec2& a, const vec2& b, const vec2& c)
{
vec2 ba, ca;
/* If ab not vertical check betweenness on x; else on y. */
if (a[0] != b[0])
return ((a[0] <= c[0]) && (c[0] <= b[0])) || //
((a[0] >= c[0]) && (c[0] >= b[0]));
else
return ((a[1] <= c[1]) && (c[1] <= b[1])) || //
((a[1] >= c[1]) && (c[1] >= b[1]));
}
bool coplaner(const vec3& pa, const vec3& pb, const vec3& pc,
const vec3& pd)
{
const double val = orient3d(pa, pb, pc, pd);
// typedef std::numeric_limits<double> dbl;
// double d = 3.14159265358979;
// std::cout.precision(dbl::max_digits10);
// std::cout << "value=" << (double)val << std::endl;
// NOTE: thresholds are chosen based on benchmark meshes that are used for testing.
// It is extremely difficult to get this right because of intermediate conversions
// between exact and fixed precision representations during cutting.
// Thus, general advise is for user to ensure that the input polygons are really
// co-planer. It might be possible for MCUT to help here (see eps used during poly
// partitioning).
return absolute_value(val) <= double(4e-7);
}
bool collinear(const vec2& a, const vec2& b, const vec2& c, double& predResult)
{
predResult = orient2d(a, b, c);
return predResult == double(0.);
}
bool collinear(const vec2& a, const vec2& b, const vec2& c)
{
return orient2d(a, b, c) == double(0.);
}
char Parallellnt(const vec2& a, const vec2& b, const vec2& c, const vec2& d, vec2& p)
{
if (!collinear(a, b, c)) {
return '0';
}
if (Between(a, b, c)) {
p = c;
return 'e';
}
if (Between(a, b, d)) {
p = d;
return 'e';
}
if (Between(c, d, a)) {
p = a;
return 'e';
}
if (Between(c, d, b)) {
p = b;
return 'e';
}
return '0';
}
char compute_segment_intersection(
const vec2& a, const vec2& b, const vec2& c,
const vec2& d, vec2& p, double& s, double& t)
{
// double s, t; /* The two parameters of the parametric eqns. */
double num, denom; /* Numerator and denominator of equations. */
char code = '?'; /* Return char characterizing intersection.*/
denom = a[0] * (d[1] - c[1]) + //
b[0] * (c[1] - d[1]) + //
d[0] * (b[1] - a[1]) + //
c[0] * (a[1] - b[1]);
/* If denom is zero, then segments are parallel: handle separately. */
if (denom == double(0.0)) {
return Parallellnt(a, b, c, d, p);
}
num = a[0] * (d[1] - c[1]) + //
c[0] * (a[1] - d[1]) + //
d[0] * (c[1] - a[1]);
if ((num == double(0.0)) || (num == denom)) {
code = 'v';
}
s = num / denom;
num = -(a[0] * (c[1] - b[1]) + //
b[0] * (a[1] - c[1]) + //
c[0] * (b[1] - a[1]));
if ((num == double(0.0)) || (num == denom)) {
code = 'v';
}
t = num / denom;
if ((double(0.0) < s) && (s < double(1.0)) && (double(0.0) < t) && (t < double(1.0))) {
code = '1';
} else if ((double(0.0) > s) || (s > double(1.0)) || (double(0.0) > t) || (t > double(1.0))) {
code = '0';
}
p[0] = a[0] + s * (b[0] - a[0]);
p[1] = a[1] + s * (b[1] - a[1]);
return code;
}
inline bool point_in_bounding_box(const vec2& point, const bounding_box_t<vec2>& bbox)
{
if ((point.x() < bbox.m_minimum.x() || point.x() > bbox.m_maximum.x()) || //
(point.y() < bbox.m_minimum.y() || point.y() > bbox.m_maximum.y())) {
return false;
} else {
return true;
}
}
inline bool point_in_bounding_box(const vec3& point, const bounding_box_t<vec3>& bbox)
{
if ((point.x() < bbox.m_minimum.x() || point.x() > bbox.m_maximum.x()) || //
(point.y() < bbox.m_minimum.y() || point.y() > bbox.m_maximum.y()) || //
(point.z() < bbox.m_minimum.z() || point.z() > bbox.m_maximum.z())) { //
return false;
} else {
return true;
}
}

764
src/mcut/source/mcut.cpp Normal file
View file

@ -0,0 +1,764 @@
/**
* Copyright (c) 2021-2022 Floyd M. Chitalu.
* All rights reserved.
*
* NOTE: This file is licensed under GPL-3.0-or-later (default).
* A commercial license can be purchased from Floyd M. Chitalu.
*
* License details:
*
* (A) GNU General Public License ("GPL"); a copy of which you should have
* recieved with this file.
* - see also: <http://www.gnu.org/licenses/>
* (B) Commercial license.
* - email: floyd.m.chitalu@gmail.com
*
* The commercial license options is for users that wish to use MCUT in
* their products for comercial purposes but do not wish to release their
* software products under the GPL license.
*
* Author(s) : Floyd M. Chitalu
*/
#include "mcut/mcut.h"
#include "mcut/internal/frontend.h"
#include "mcut/internal/timer.h"
#include "mcut/internal/utils.h"
#include <exception>
#include <stdexcept>
#if defined(MCUT_BUILD_WINDOWS)
#pragma warning(disable : 26812)
#endif
MCAPI_ATTR McResult MCAPI_CALL mcCreateContext(McContext* pOutContext, McFlags contextFlags)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (pOutContext == nullptr) {
per_thread_api_log_str = "context ptr undef (NULL)";
return_value = McResult::MC_INVALID_VALUE;
} else {
try {
// no helper threads
// only manager threads (2 managers if context is created with MC_OUT_OF_ORDER_EXEC_MODE_ENABLE, otherwise 1 manager)
create_context_impl(pOutContext, contextFlags, 0);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (return_value != McResult::MC_NO_ERROR) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcCreateContextWithHelpers(McContext* pOutContext, McFlags contextFlags, uint32_t helperThreadCount)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (pOutContext == nullptr) {
per_thread_api_log_str = "context ptr undef (NULL)";
return_value = McResult::MC_INVALID_VALUE;
} else {
try {
create_context_impl(pOutContext, contextFlags, helperThreadCount);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (return_value != McResult::MC_NO_ERROR) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcDebugMessageCallback(McContext pContext, pfn_mcDebugOutput_CALLBACK cb, const McVoid* userParam)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (pContext == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (cb == nullptr) {
per_thread_api_log_str = "callback function ptr (param1) undef (NULL)";
} else {
try {
debug_message_callback_impl(pContext, cb, userParam);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcGetDebugMessageLog(
McContext context,
McUint32 count, McSize bufSize,
McDebugSource* sources, McDebugType* types, McDebugSeverity* severities,
McSize* lengths, McChar* messageLog, McUint32* numFetched)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (
count == 0) {
per_thread_api_log_str = "count must be > 0";
} else if (bufSize == 0) {
per_thread_api_log_str = "bufSize must be > 0";
} else if (numFetched == nullptr) {
per_thread_api_log_str = "numFetched undef (NULL)";
} else {
try {
get_debug_message_log_impl(context,
count, bufSize,
sources, types, severities,
lengths, messageLog, *numFetched);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcDebugMessageControl(McContext pContext, McDebugSource source, McDebugType type, McDebugSeverity severity, bool enabled)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (pContext == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (
false == (source == McDebugSource::MC_DEBUG_SOURCE_ALL || //
source == McDebugSource::MC_DEBUG_SOURCE_KERNEL || //
source == McDebugSource::MC_DEBUG_SOURCE_ALL)) {
per_thread_api_log_str = "debug source val (param1) invalid";
} else if (
false == (type == McDebugType::MC_DEBUG_TYPE_ALL || //
type == McDebugType::MC_DEBUG_TYPE_DEPRECATED_BEHAVIOR || //
type == McDebugType::MC_DEBUG_TYPE_ERROR || //
type == McDebugType::MC_DEBUG_TYPE_OTHER)) {
per_thread_api_log_str = "debug type val (param2) invalid";
} else if (
false == (severity == McDebugSeverity::MC_DEBUG_SEVERITY_HIGH || //
severity == McDebugSeverity::MC_DEBUG_SEVERITY_MEDIUM || //
severity == McDebugSeverity::MC_DEBUG_SEVERITY_LOW || //
severity == McDebugSeverity::MC_DEBUG_SEVERITY_NOTIFICATION || //
severity == McDebugSeverity::MC_DEBUG_SEVERITY_ALL)) {
per_thread_api_log_str = "debug severity val (param3) invalid";
} else {
try {
debug_message_control_impl(pContext, source, type, severity, enabled);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcGetInfo(const McContext context, McFlags info, McSize bytes, McVoid* pMem, McSize* pNumBytes)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (bytes != 0 && pMem == nullptr) {
per_thread_api_log_str = "invalid specification (param2 & param3)";
} else if (false == (info == MC_CONTEXT_FLAGS || info == MC_MAX_DEBUG_MESSAGE_LENGTH)) // check all possible values
{
per_thread_api_log_str = "invalid info flag val (param1)";
} else if ((info == MC_CONTEXT_FLAGS) && (pMem != nullptr && bytes != sizeof(McFlags))) {
per_thread_api_log_str = "invalid byte size (param2)"; // leads to e.g. "out of bounds" memory access during memcpy
} else {
try {
get_info_impl(context, info, bytes, pMem, pNumBytes);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcCreateUserEvent(
McEvent* event,
McContext context)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (event == nullptr) {
per_thread_api_log_str = "event ptr (param0) undef (NULL)";
} else if (context == nullptr) {
per_thread_api_log_str = "context handle undefined (NULL)";
} else {
try {
create_user_event_impl(event, context);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcSetUserEventStatus(
McEvent event,
McInt32 execution_status)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (event == nullptr) {
per_thread_api_log_str = "event ptr (param0) undef (NULL)";
} else {
try {
set_user_event_status_impl(event, execution_status);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcGetEventInfo(const McEvent event, McFlags info, McSize bytes, McVoid* pMem, McSize* pNumBytes)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (event == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (bytes != 0 && pMem == nullptr) {
per_thread_api_log_str = "invalid specification (param2 & param3)";
} else if ((info == MC_EVENT_RUNTIME_EXECUTION_STATUS) && (pMem != nullptr && bytes != sizeof(McResult))) {
per_thread_api_log_str = "invalid byte size (param2)"; // leads to e.g. "out of bounds" memory access during memcpy
} else if ((info == MC_EVENT_COMMAND_EXECUTION_STATUS) && (pMem != nullptr && bytes != sizeof(McFlags))) {
per_thread_api_log_str = "invalid byte size (param2)"; // leads to e.g. "out of bounds" memory access during memcpy
} else if ((info == MC_EVENT_TIMESTAMP_SUBMIT || info == MC_EVENT_TIMESTAMP_START || info == MC_EVENT_TIMESTAMP_END) && (pMem != nullptr && bytes != sizeof(McSize))) {
per_thread_api_log_str = "invalid byte size (param2)"; // leads to e.g. "out of bounds" memory access during memcpy
} else {
try {
get_event_info_impl(event, info, bytes, pMem, pNumBytes);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcWaitForEvents(
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (pEventWaitList == nullptr && numEventsInWaitlist > 0) {
per_thread_api_log_str = "invalid event waitlist ptr (NULL)";
} else if (pEventWaitList != nullptr && numEventsInWaitlist == 0) {
per_thread_api_log_str = "invalid event waitlist size (zero)";
} else {
try {
McResult waitliststatus = MC_NO_ERROR;
wait_for_events_impl(numEventsInWaitlist, pEventWaitList, waitliststatus);
if (waitliststatus != McResult::MC_NO_ERROR) {
per_thread_api_log_str = "event in waitlist has an error";
}
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcSetEventCallback(
McEvent eventHandle,
pfn_McEvent_CALLBACK eventCallback,
McVoid* data)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (eventHandle == nullptr) {
per_thread_api_log_str = "invalid event ptr (NULL)";
}
if (eventCallback == nullptr) {
per_thread_api_log_str = "invalid event callback function ptr (NULL)";
} else {
try {
set_event_callback_impl(eventHandle, eventCallback, data);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcEnqueueDispatch(
const McContext context,
McFlags dispatchFlags,
const McVoid* pSrcMeshVertices,
const uint32_t* pSrcMeshFaceIndices,
const uint32_t* pSrcMeshFaceSizes,
uint32_t numSrcMeshVertices,
uint32_t numSrcMeshFaces,
const McVoid* pCutMeshVertices,
const uint32_t* pCutMeshFaceIndices,
const uint32_t* pCutMeshFaceSizes,
uint32_t numCutMeshVertices,
uint32_t numCutMeshFaces,
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McEvent* pEvent)
{
TIMESTACK_RESET(); // reset tracking vars
SCOPED_TIMER(__FUNCTION__);
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (dispatchFlags == 0) {
per_thread_api_log_str = "dispatch flags unspecified";
} else if ((dispatchFlags & MC_DISPATCH_REQUIRE_THROUGH_CUTS) && //
(dispatchFlags & MC_DISPATCH_FILTER_FRAGMENT_LOCATION_UNDEFINED)) {
// The user states that she does not want a partial cut but yet also states that she
// wants to keep fragments with partial cuts. These two options are mutually exclusive!
per_thread_api_log_str = "use of mutually-exclusive flags: MC_DISPATCH_REQUIRE_THROUGH_CUTS & MC_DISPATCH_FILTER_FRAGMENT_LOCATION_UNDEFINED";
} else if ((dispatchFlags & MC_DISPATCH_VERTEX_ARRAY_FLOAT) == 0 && (dispatchFlags & MC_DISPATCH_VERTEX_ARRAY_DOUBLE) == 0) {
per_thread_api_log_str = "dispatch vertex aray type unspecified";
} else if (pSrcMeshVertices == nullptr) {
per_thread_api_log_str = "source-mesh vertex-position array ptr undef (NULL)";
} else if (numSrcMeshVertices < 3) {
per_thread_api_log_str = "invalid source-mesh vertex count";
} else if (pSrcMeshFaceIndices == nullptr) {
per_thread_api_log_str = "source-mesh face-index array ptr undef (NULL)";
} /*else if (pSrcMeshFaceSizes == nullptr) {
per_thread_api_log_str = "source-mesh face-size array ptr undef (NULL)";
}*/
else if (numSrcMeshFaces < 1) {
per_thread_api_log_str = "invalid source-mesh vertex count";
} else if (pCutMeshVertices == nullptr) {
per_thread_api_log_str = "cut-mesh vertex-position array ptr undef (NULL)";
} else if (numCutMeshVertices < 3) {
per_thread_api_log_str = "invalid cut-mesh vertex count";
} else if (pCutMeshFaceIndices == nullptr) {
per_thread_api_log_str = "cut-mesh face-index array ptr undef (NULL)";
} /*else if (pCutMeshFaceSizes == nullptr) {
per_thread_api_log_str = "cut-mesh face-size array ptr undef (NULL)";
} */
else if (numCutMeshFaces < 1) {
per_thread_api_log_str = "invalid cut-mesh vertex count";
} else if (pEventWaitList == nullptr && numEventsInWaitlist > 0) {
per_thread_api_log_str = "invalid event waitlist ptr (NULL)";
} else if (pEventWaitList != nullptr && numEventsInWaitlist == 0) {
per_thread_api_log_str = "invalid event waitlist size (zero)";
} else if (pEventWaitList == nullptr && numEventsInWaitlist == 0 && pEvent == nullptr) {
per_thread_api_log_str = "invalid event ptr (zero)";
} else {
try {
dispatch_impl(
context,
dispatchFlags,
pSrcMeshVertices,
pSrcMeshFaceIndices,
pSrcMeshFaceSizes,
numSrcMeshVertices,
numSrcMeshFaces,
pCutMeshVertices,
pCutMeshFaceIndices,
pCutMeshFaceSizes,
numCutMeshVertices,
numCutMeshFaces,
numEventsInWaitlist,
pEventWaitList,
pEvent);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcDispatch(
const McContext context,
McFlags dispatchFlags,
const McVoid* pSrcMeshVertices,
const uint32_t* pSrcMeshFaceIndices,
const uint32_t* pSrcMeshFaceSizes,
uint32_t numSrcMeshVertices,
uint32_t numSrcMeshFaces,
const McVoid* pCutMeshVertices,
const uint32_t* pCutMeshFaceIndices,
const uint32_t* pCutMeshFaceSizes,
uint32_t numCutMeshVertices,
uint32_t numCutMeshFaces)
{
McEvent event = MC_NULL_HANDLE;
McResult return_value = mcEnqueueDispatch(
context,
dispatchFlags,
pSrcMeshVertices,
pSrcMeshFaceIndices,
pSrcMeshFaceSizes,
numSrcMeshVertices,
numSrcMeshFaces,
pCutMeshVertices,
pCutMeshFaceIndices,
pCutMeshFaceSizes,
numCutMeshVertices,
numCutMeshFaces,
0,
nullptr,
&event);
if (return_value == MC_NO_ERROR) { // API parameter checks are fine
if (event != MC_NULL_HANDLE) // event must exist to wait on and query
{
McResult waitliststatus = MC_NO_ERROR;
wait_for_events_impl(1, &event, waitliststatus); // block until event of mcEnqueueDispatch is completed!
if (waitliststatus != McResult::MC_NO_ERROR) {
return_value = waitliststatus;
}
release_events_impl(1, &event); // destroy
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcEnqueueGetConnectedComponents(
const McContext context,
const McConnectedComponentType connectedComponentType,
const uint32_t numEntries,
McConnectedComponent* pConnComps,
uint32_t* numConnComps,
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McEvent* pEvent)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (connectedComponentType == 0) {
per_thread_api_log_str = "invalid type-parameter (param1) (0)";
} else if (numConnComps == nullptr && pConnComps == nullptr) {
per_thread_api_log_str = "output parameters undef (param3 & param4)";
} else if (pEventWaitList == nullptr && numEventsInWaitlist > 0) {
per_thread_api_log_str = "invalid event waitlist ptr (NULL)";
} else if (pEventWaitList != nullptr && numEventsInWaitlist == 0) {
per_thread_api_log_str = "invalid event waitlist size (zero)";
} else if (pEventWaitList == nullptr && numEventsInWaitlist == 0 && pEvent == nullptr) {
per_thread_api_log_str = "invalid event ptr (zero)";
} else {
try {
get_connected_components_impl(context, connectedComponentType, numEntries, pConnComps, numConnComps, numEventsInWaitlist, pEventWaitList, pEvent);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcGetConnectedComponents(
const McContext context,
const McConnectedComponentType connectedComponentType,
const uint32_t numEntries,
McConnectedComponent* pConnComps,
uint32_t* numConnComps)
{
McEvent event = MC_NULL_HANDLE;
McResult return_value = mcEnqueueGetConnectedComponents(context, connectedComponentType, numEntries, pConnComps, numConnComps, 0, nullptr, &event);
if (event != MC_NULL_HANDLE) // event must exist to wait on and query
{
McResult waitliststatus = MC_NO_ERROR;
wait_for_events_impl(1, &event, waitliststatus); // block until event of mcEnqueueDispatch is completed!
if (waitliststatus != McResult::MC_NO_ERROR) {
return_value = waitliststatus;
}
release_events_impl(1, &event); // destroy
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcEnqueueGetConnectedComponentData(
const McContext context,
const McConnectedComponent connCompId,
McFlags queryFlags,
McSize bytes,
McVoid* pMem,
McSize* pNumBytes,
uint32_t numEventsInWaitlist,
const McEvent* pEventWaitList,
McEvent* pEvent)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
}
if (connCompId == nullptr) {
per_thread_api_log_str = "connected component ptr (param1) undef (NULL)";
} else if (queryFlags == 0) {
per_thread_api_log_str = "flags (param1) undef (0)";
} else if (bytes != 0 && pMem == nullptr) {
per_thread_api_log_str = "null parameter (param3 & param4)";
} else if (pEventWaitList == nullptr && numEventsInWaitlist > 0) {
per_thread_api_log_str = "invalid event waitlist ptr (NULL)";
} else if (pEventWaitList != nullptr && numEventsInWaitlist == 0) {
per_thread_api_log_str = "invalid event waitlist size (zero)";
} else if (pEventWaitList == nullptr && numEventsInWaitlist == 0 && pEvent == nullptr) {
per_thread_api_log_str = "invalid event ptr (zero)";
} else {
try {
get_connected_component_data_impl(context, connCompId, queryFlags, bytes, pMem, pNumBytes, numEventsInWaitlist, pEventWaitList, pEvent);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcGetConnectedComponentData(
const McContext context,
const McConnectedComponent connCompId,
McFlags queryFlags,
McSize bytes,
McVoid* pMem,
McSize* pNumBytes)
{
McEvent event = MC_NULL_HANDLE;
McResult return_value = mcEnqueueGetConnectedComponentData(context, connCompId, queryFlags, bytes, pMem, pNumBytes, 0, nullptr, &event);
if (event != MC_NULL_HANDLE) // event must exist to wait on and query
{
McResult waitliststatus = MC_NO_ERROR;
wait_for_events_impl(1, &event, waitliststatus); // block until event of mcEnqueueDispatch is completed!
if (waitliststatus != McResult::MC_NO_ERROR) {
return_value = waitliststatus;
}
release_events_impl(1, &event); // destroy
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcReleaseEvents(
uint32_t numEvents,
const McEvent* pEvents)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (numEvents > 0 && pEvents == NULL) {
per_thread_api_log_str = "invalid pointer to events";
} else if (numEvents == 0 && pEvents != NULL) {
per_thread_api_log_str = "number of events not set";
} else {
try {
release_events_impl(numEvents, pEvents);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcReleaseConnectedComponents(
const McContext context,
uint32_t numConnComps,
const McConnectedComponent* pConnComps)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else if (numConnComps > 0 && pConnComps == NULL) {
per_thread_api_log_str = "invalid pointer to connected components";
} else if (numConnComps == 0 && pConnComps != NULL) {
per_thread_api_log_str = "number of connected components not set";
} else {
try {
release_connected_components_impl(context, numConnComps, pConnComps);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}
MCAPI_ATTR McResult MCAPI_CALL mcReleaseContext(const McContext context)
{
McResult return_value = McResult::MC_NO_ERROR;
per_thread_api_log_str.clear();
if (context == nullptr) {
per_thread_api_log_str = "context ptr (param0) undef (NULL)";
} else {
try {
release_context_impl(context);
}
CATCH_POSSIBLE_EXCEPTIONS(per_thread_api_log_str);
}
if (!per_thread_api_log_str.empty()) {
std::fprintf(stderr, "%s(...) -> %s\n", __FUNCTION__, per_thread_api_log_str.c_str());
if (return_value == McResult::MC_NO_ERROR) // i.e. problem with basic local parameter checks
{
return_value = McResult::MC_INVALID_VALUE;
}
}
return return_value;
}

2198
src/mcut/source/preproc.cpp Normal file
View file

File diff suppressed because it is too large Load diff

8682
src/mcut/source/shewchuk.c Normal file
View file

File diff suppressed because it is too large Load diff

6
src/slic3r/GUI/Plater.cpp
View file

@ -9535,6 +9535,7 @@ void Plater::export_core_3mf()
export_3mf(path_u8, SaveStrategy::Silence);
}
#define USE_CGAL_BOOLEAN 0
void Plater::export_stl(bool extended, bool selection_only)
{
if (p->model.objects.empty()) { return; }
@ -9561,8 +9562,13 @@ void Plater::export_stl(bool extended, bool selection_only)
if (csg::check_csgmesh_booleans(Range{ std::begin(csgmesh), std::end(csgmesh) }) == csgmesh.end()) {
try {
#if USE_CGAL_BOOLEAN
auto meshPtr = csg::perform_csgmesh_booleans(Range{ std::begin(csgmesh), std::end(csgmesh) });
mesh = MeshBoolean::cgal::cgal_to_triangle_mesh(*meshPtr);
#else
MeshBoolean::mcut::McutMeshPtr meshPtr = csg::perform_csgmesh_booleans_mcut(Range{std::begin(csgmesh), std::end(csgmesh)});
mesh = MeshBoolean::mcut::mcut_to_triangle_mesh(*meshPtr);
#endif
} catch (...) {}
}