diff --git a/Data/Testing/shrink_plane.vtp b/Data/Testing/shrink_plane.vtp new file mode 100644 index 0000000000000000000000000000000000000000..ae39dc501b655c1f525c02a8243607fe2f99a1d8 --- /dev/null +++ b/Data/Testing/shrink_plane.vtp @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:dab096e6f0f3fcc8fd02fcec8195d9194d3e5c162bb5d5a401e33b2d0e2ece6c +size 16766 diff --git a/vespa/Delaunay/CMakeLists.txt b/vespa/Delaunay/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..689690e821c711f6926496a9eebe25277d499383 --- /dev/null +++ b/vespa/Delaunay/CMakeLists.txt @@ -0,0 +1,6 @@ +set(vtkcgaldelaunay_files + vtkCGALDelaunay2 +) +vtk_module_add_module(vtkCGALDelaunay + CLASSES ${vtkcgaldelaunay_files} +) diff --git a/vespa/Delaunay/Testing/CMakeLists.txt b/vespa/Delaunay/Testing/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..06696843a5ec976069b357a2e5daf577ab17f162 --- /dev/null +++ b/vespa/Delaunay/Testing/CMakeLists.txt @@ -0,0 +1,41 @@ +if (TARGET VTK::vtkpython) + # import + add_test(NAME "import_Delaunay" + COMMAND + "$" + "${CMAKE_CURRENT_LIST_DIR}/import_Delaunay.py") + set_property(TEST "import_Delaunay" APPEND + PROPERTY + ENVIRONMENT "PYTHONPATH=${CMAKE_BINARY_DIR}/${python_destination}") + if (WIN32) + set(test_path "$ENV{PATH};${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_BINDIR}") + string(REPLACE ";" "\;" test_path "${test_path}") + set_property(TEST "import_Delaunay" APPEND + PROPERTY + ENVIRONMENT "PATH=${test_path}") + endif () + + # execute + add_test(NAME "execute_Delaunay" + COMMAND + "$" + "${CMAKE_CURRENT_LIST_DIR}/execute_Delaunay.py") + set_property(TEST "execute_Delaunay" APPEND + PROPERTY + ENVIRONMENT "PYTHONPATH=${CMAKE_BINARY_DIR}/${python_destination}") + if (WIN32) + set(test_path "$ENV{PATH};${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_BINDIR}") + string(REPLACE ";" "\;" test_path "${test_path}") + set_property(TEST "execute_Delaunay" APPEND + PROPERTY + ENVIRONMENT "PATH=${test_path}") + endif () +endif () + +vtk_add_test_cxx(vtkCGALDelaunayCxxTests no_data_tests + NO_DATA NO_VALID NO_OUTPUT + TestDelaunayExecution.cxx + + ${PROJECT_SOURCE_DIR}/Data/Testing/ +) +vtk_test_cxx_executable(vtkCGALDelaunayCxxTests no_data_tests) diff --git a/vespa/Delaunay/Testing/TestDelaunayExecution.cxx b/vespa/Delaunay/Testing/TestDelaunayExecution.cxx new file mode 100644 index 0000000000000000000000000000000000000000..2750f7ea247b6aaad3a6ef4b4b2b8a0fe2f92ff9 --- /dev/null +++ b/vespa/Delaunay/Testing/TestDelaunayExecution.cxx @@ -0,0 +1,33 @@ +#include + +#include +#include +#include +#include + +#include "vtkCGALDelaunay2.h" + +int TestDelaunayExecution(int, char* argv[]) +{ + // Open data + + vtkNew reader2; + std::string cfname2(argv[1]); + cfname2 += "shrink_plane.vtp"; + reader2->SetFileName(cfname2.c_str()); + + // Remesh + + vtkNew rm2; + rm2->SetInputConnection(reader2->GetOutputPort()); + + // Save result + + vtkNew writer; + + writer->SetInputConnection(rm2->GetOutputPort()); + writer->SetFileName("delaunay2_remesh.vtp"); + writer->Write(); + + return 0; +} diff --git a/vespa/Delaunay/Testing/execute_Delaunay.py b/vespa/Delaunay/Testing/execute_Delaunay.py new file mode 100644 index 0000000000000000000000000000000000000000..950f1aeb2984edd9a7c74be0f9b1f0529c37ef58 --- /dev/null +++ b/vespa/Delaunay/Testing/execute_Delaunay.py @@ -0,0 +1,11 @@ +from vtk import vtkSphereSource, vtkShrinkPolyData +from vespa import vtkCGALDelaunay + +sp = vtkSphereSource() + +pd = vtkShrinkPolyData() +pd.SetInputConnection(sp.GetOutputPort()) + +d2 = vtkCGALDelaunay.vtkCGALDelaunay2() +d2.SetInputConnection(pd.GetOutputPort()) +d2.Update() diff --git a/vespa/Delaunay/Testing/import_Delaunay.py b/vespa/Delaunay/Testing/import_Delaunay.py new file mode 100644 index 0000000000000000000000000000000000000000..30a07955fa2c2e3ed5cb042c26dc8b3df254f2ce --- /dev/null +++ b/vespa/Delaunay/Testing/import_Delaunay.py @@ -0,0 +1,4 @@ +from vespa import vtkCGALDelaunay + +d2 = vtkCGALDelaunay.vtkCGALDelaunay2() +help(d2) diff --git a/vespa/Delaunay/vtk.module b/vespa/Delaunay/vtk.module new file mode 100644 index 0000000000000000000000000000000000000000..3ad2b5e03e8cb084052a5e630bacadf308ff69a3 --- /dev/null +++ b/vespa/Delaunay/vtk.module @@ -0,0 +1,20 @@ +NAME + vtkCGALDelaunay +DESCRIPTION + "This module contains filters related to the Delaunay module of CGAL." +GROUPS + Meshing +DEPENDS + VTK::CommonCore + VTK::CommonDataModel + VTK::CommonExecutionModel + VTK::FiltersCore + VTK::FiltersExtraction + VTK::FiltersGeometry + CGAL::CGAL +TEST_DEPENDS + VTK::CommonSystem + VTK::FiltersSources + VTK::IOInfovis + VTK::IOXML + VTK::TestingCore diff --git a/vespa/Delaunay/vtkCGALDelaunay2.cxx b/vespa/Delaunay/vtkCGALDelaunay2.cxx new file mode 100644 index 0000000000000000000000000000000000000000..b2eca9fc3ad7974ccac4aef4f49ea92341c8f0c1 --- /dev/null +++ b/vespa/Delaunay/vtkCGALDelaunay2.cxx @@ -0,0 +1,179 @@ +#include "vtkCGALDelaunay2.h" + +// VTK related includes +#include "vtkCellArrayIterator.h" +#include "vtkDataSet.h" +#include "vtkIdList.h" +#include "vtkInformationVector.h" +#include "vtkObjectFactory.h" + +// CGAL related includes +#include +#include + +vtkStandardNewMacro(vtkCGALDelaunay2); + +// TODO May try to use ProjectionTraits_3 to handle open 3D surfaces +// Look at perf then +// caution, a sphere won't work: intersection +// caution, infinit loop on some tests +using CDT2 = CGAL::Constrained_Delaunay_triangulation_2; + +//------------------------------------------------------------------------------ +void vtkCGALDelaunay2::PrintSelf(ostream& os, vtkIndent indent) +{ + this->Superclass::PrintSelf(os, indent); +} + +//------------------------------------------------------------------------------ +int vtkCGALDelaunay2::RequestData( + vtkInformation*, vtkInformationVector** inputVector, vtkInformationVector* outputVector) +{ + // Get the input and output data objects. + vtkPolyData* input = vtkPolyData::GetData(inputVector[0]); + vtkPolyData* output = vtkPolyData::GetData(outputVector); + + // Create the surface mesh for CGAL + // -------------------------------- + + vtkPoints* vtkPts = input->GetPoints(); + vtkIdType nbPts = input->GetNumberOfPoints(); + vtkDataArray* ptsArr = vtkPts->GetData(); + const auto pointRange = vtk::DataArrayTupleRange<3>(ptsArr); + + // manually handle the plannar coordinate + // should be along the x, y or z axis + double rangeVal[3]; + for (int i = 0; i < 3; i++) + { + double range[2]; + ptsArr->GetRange(range, i); + rangeVal[i] = range[1] - range[0]; + } + + if (rangeVal[0] && rangeVal[1] && rangeVal[2]) + { + vtkErrorMacro("This dataset is 3D"); + } + int d1 = 0, d2 = 1, d3 = 2; // z null + if (!rangeVal[0]) + { + d1 = 2; + d2 = 0; + d3 = 1; + } + if (!rangeVal[1]) + { + d1 = 0; + d2 = 2; + d3 = 1; + } + + std::vector> pts; + pts.reserve(nbPts); + for (const auto pt : pointRange) + { + pts.emplace_back(CDT2::Point(pt[d1], pt[d2]), true); + } + + // CGAL Processing + // --------------- + + CDT2 delaunay; + try + { + // Add constraints (lines and polys) + vtkCellArray* polys = input->GetPolys(); + auto polysIt = vtk::TakeSmartPointer(polys->NewIterator()); + // each poly + for (polysIt->GoToFirstCell(); !polysIt->IsDoneWithTraversal(); polysIt->GoToNextCell()) + { + vtkIdList* p = polysIt->GetCurrentCell(); + std::list poly; + // each segment of poly + for (vtkIdType i = 0; i < p->GetNumberOfIds(); i++) + { + auto point = std::get<0>(pts[p->GetId(i)]); + poly.emplace_back(point); + std::get<1>(pts[p->GetId(i)]) = false; + } + delaunay.insert_constraint(poly.begin(), poly.end(), true); + } + + vtkCellArray* lines = input->GetLines(); + auto linesIt = vtk::TakeSmartPointer(lines->NewIterator()); + // each line + for (linesIt->GoToFirstCell(); !linesIt->IsDoneWithTraversal(); linesIt->GoToNextCell()) + { + // each segment of line + vtkIdList* l = linesIt->GetCurrentCell(); + std::list line; + for (vtkIdType i = 1; i < l->GetNumberOfIds(); i++) + { + if (!std::get<1>(pts[l->GetId(i)])) + { + std::cerr << "invalid line: " << l->GetId(i) << std::endl; + continue; + } + auto point = std::get<0>(pts[l->GetId(i)]); + line.emplace_back(point); + std::get<1>(pts[l->GetId(i)]) = false; + } + delaunay.insert_constraint(line.begin(), line.end()); + } + + // Add points + for (auto point : pts) + { + if (std::get<1>(point)) + { + delaunay.push_back(CDT2::Point(std::get<0>(point))); + } + } + } + catch (std::exception& e) + { + vtkErrorMacro("CGAL Exception: " << e.what()); + return 0; + } + + // VTK Output + // ---------- + + vtkNew outPts; + const vtkIdType outNPts = delaunay.number_of_vertices(); + outPts->Allocate(outNPts); + std::map vmap; + + for (auto vertex : delaunay.finite_vertex_handles()) + { + double coords[3]; + coords[d1] = vertex->point()[0]; + coords[d2] = vertex->point()[1]; + coords[d3] = rangeVal[d3]; + vtkIdType id = outPts->InsertNextPoint(coords); + vmap[vertex->point()] = id; + } + outPts->Squeeze(); + + // cells + vtkNew cells; + cells->AllocateEstimate(delaunay.number_of_faces(), 3); + + for (auto face : delaunay.finite_face_handles()) + { + vtkNew ids; + ids->InsertNextId(vmap[face->vertex(0)->point()]); + ids->InsertNextId(vmap[face->vertex(1)->point()]); + ids->InsertNextId(vmap[face->vertex(2)->point()]); + + cells->InsertNextCell(ids); + } + cells->Squeeze(); + + // VTK dataset + output->SetPoints(outPts); + output->SetPolys(cells); + + return 1; +} diff --git a/vespa/Delaunay/vtkCGALDelaunay2.h b/vespa/Delaunay/vtkCGALDelaunay2.h new file mode 100644 index 0000000000000000000000000000000000000000..55f74d65d2686bd0853a62e3bb43c45ea37a6c05 --- /dev/null +++ b/vespa/Delaunay/vtkCGALDelaunay2.h @@ -0,0 +1,42 @@ +/** + * @class vtkCGALDelaunay2 + * @brief remesh using the CGAL delaunay + * + * vtkCGALDelaunay2 allows to create plannar delaunay meshes + * from a set of planar points, edges and polygons. + * From now on, the input mesh needs to be planar along x, y or z. + * Constraints should not overlap each others. + */ + +#ifndef vtkCGALDelaunay2_h +#define vtkCGALDelaunay2_h + +#include "vtkPolyDataAlgorithm.h" + +// CGAL includes +#include +#include + +using CGAL_Kernel = CGAL::Simple_cartesian; + +#include "vtkCGALDelaunayModule.h" // For export macro + +class VTKCGALDELAUNAY_EXPORT vtkCGALDelaunay2 : public vtkPolyDataAlgorithm +{ +public: + static vtkCGALDelaunay2* New(); + vtkTypeMacro(vtkCGALDelaunay2, vtkPolyDataAlgorithm); + void PrintSelf(ostream& os, vtkIndent indent) override; + +protected: + vtkCGALDelaunay2() = default; + ~vtkCGALDelaunay2() override = default; + + int RequestData(vtkInformation*, vtkInformationVector**, vtkInformationVector*) override; + +private: + vtkCGALDelaunay2(const vtkCGALDelaunay2&) = delete; + void operator=(const vtkCGALDelaunay2&) = delete; +}; + +#endif diff --git a/vespa/PolygonMeshProcessing/Testing/CMakeLists.txt b/vespa/PolygonMeshProcessing/Testing/CMakeLists.txt index 1666c1574b93323d6a60ce145382abf2f5bc0d17..8bf63c81f34675e99993501dcfa8f0c312ec3e2e 100644 --- a/vespa/PolygonMeshProcessing/Testing/CMakeLists.txt +++ b/vespa/PolygonMeshProcessing/Testing/CMakeLists.txt @@ -32,7 +32,7 @@ if (TARGET VTK::vtkpython) endif () endif () -vtk_add_test_cxx(vtkCGALCxxTests no_data_tests +vtk_add_test_cxx(vtkCGALPMPCxxTests no_data_tests NO_DATA NO_VALID NO_OUTPUT TestPMPInstance.cxx TestPMPBooleanExecution.cxx @@ -44,4 +44,4 @@ vtk_add_test_cxx(vtkCGALCxxTests no_data_tests ${PROJECT_SOURCE_DIR}/Data/Testing/ ) -vtk_test_cxx_executable(vtkCGALCxxTests no_data_tests) +vtk_test_cxx_executable(vtkCGALPMPCxxTests no_data_tests) diff --git a/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.cxx b/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.cxx index c115936cba54c3bff9c5e8cb343ef3fe00892121..caf1f8e4b3c43dd286c8050adeaaf071db84d1fd 100644 --- a/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.cxx +++ b/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.cxx @@ -38,26 +38,18 @@ std::unique_ptr vtkCGALPolyDataAlgorithm::toCGAL(vtkPolyData* vtkMesh } // Cells - std::array tri; - auto cit = vtk::TakeSmartPointer(vtkMesh->NewCellIterator()); for (cit->InitTraversal(); !cit->IsDoneWithTraversal(); cit->GoToNextCell()) { - // Sanity check - if (cit->GetCellType() != VTK_TRIANGLE) - { - vtkIdType id = cit->GetCellId(); - vtkErrorMacro("Cell " << id << " is not a triangle. Abort."); - return 0; - } - - // Add the triangle + // Add the cell vtkIdList* ids = cit->GetPointIds(); - for (vtkIdType i = 0; i < 3; i++) + vtkIdType nbIds = cit->GetNumberOfPoints(); + std::vector cell(nbIds); + for (vtkIdType i = 0; i < nbIds; i++) { - tri[i] = surfaceVertices[ids->GetId(i)]; + cell[i] = surfaceVertices[ids->GetId(i)]; } - CGAL::Euler::add_face(tri, cgalMesh->surface); + CGAL::Euler::add_face(cell, cgalMesh->surface); } return cgalMesh; diff --git a/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.h b/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.h index e30b1801d150508bc4e754f2cbe6be3d97d5c88d..abb0a4bfe8a93ca4be9d6d9fdcefe7ef66b7107c 100644 --- a/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.h +++ b/vespa/PolygonMeshProcessing/vtkCGALPolyDataAlgorithm.h @@ -17,14 +17,14 @@ #include #include -#include "vtkCGALPMPModule.h" // For export macro - using CGAL_Kernel = CGAL::Simple_cartesian; using CGAL_Surface = CGAL::Surface_mesh; using Graph_Verts = boost::graph_traits::vertex_descriptor; using Graph_Faces = boost::graph_traits::face_descriptor; using Graph_Coord = boost::property_map::type; +#include "vtkCGALPMPModule.h" // For export macro + // Container for CGAL related info struct CGAL_Mesh {