Loading [MathJax]/extensions/TeX/AMSsymbols.js
DGtal 2.0.0
dgtalCalculus.cpp
1
25#include <iostream>
26#include <DGtal/base/Common.h>
27#include <DGtal/helpers/StdDefs.h>
28#include <DGtal/helpers/Shortcuts.h>
29#include <DGtal/helpers/ShortcutsGeometry.h>
30#include <DGtal/shapes/SurfaceMesh.h>
31#include <DGtal/geometry/surfaces/DigitalSurfaceRegularization.h>
32#include <DGtal/dec/PolygonalCalculus.h>
33
34#include <polyscope/polyscope.h>
35#include <polyscope/surface_mesh.h>
36#include <polyscope/point_cloud.h>
37#include <Eigen/Dense>
38#include <Eigen/Sparse>
39
40using namespace DGtal;
41using namespace Z3i;
42
43// Using standard 3D digital space.
46// The following typedefs are useful
50typedef SurfMesh::Face Face;
52
53//Polyscope global
54polyscope::SurfaceMesh *psMesh;
56std::vector<double> phiV;
57float scale = 0.1;
58
59//Restriction of an ambient scalar function to vertices
60double phiVertex(const Vertex v)
61{
62 return cos(scale*(surfmesh.position(v)[0]))*sin(scale*surfmesh.position(v)[1]);
63}
64
65//Restriction of an ambient scalar function to vertices
67{
68 auto vertices = surfmesh.incidentVertices(f);
69 auto nf = vertices.size();
70 Eigen::VectorXd ph(nf);
71 size_t cpt=0;
72 for(auto v: vertices)
73 {
74 ph(cpt) = phiVertex(v);
75 ++cpt;
76 }
77 return ph;
78}
79
80//Vertex valued function for polyscope
81void initPhi()
82{
83 phiV.clear();
84 for(auto i = 0; i < surfmesh.nbVertices(); ++i)
85 phiV.push_back(phiVertex(i));
86 psMesh->addVertexScalarQuantity("Phi", phiV);
87}
88
89void initQuantities()
90{
92
93 std::vector<PolygonalCalculus<SH3::RealPoint,SH3::RealVector>::Vector> gradients;
94 std::vector<PolygonalCalculus<SH3::RealPoint,SH3::RealVector>::Vector> cogradients;
95 std::vector<PolygonalCalculus<SH3::RealPoint,SH3::RealVector>::Real3dVector> normals;
96 std::vector<PolygonalCalculus<SH3::RealPoint,SH3::RealVector>::Real3dVector> vectorArea;
97 std::vector<PolygonalCalculus<SH3::RealPoint,SH3::RealVector>::Real3dPoint> centroids;
98 std::vector<double> faceArea;
99
100 for(auto f=0; f < surfmesh.nbFaces(); ++f)
101 {
104 gradients.push_back( grad );
106 cogradients.push_back( cograd );
107 normals.push_back(calculus.faceNormalAsDGtalVector(f));
108
109 auto vA = calculus.vectorArea(f);
110 vectorArea.push_back({vA(0) , vA(1), vA(2)});
111
112 faceArea.push_back( calculus.faceArea(f));
113
114 centroids.push_back( calculus.centroidAsDGtalPoint(f) );
115 }
116
117 psMesh->addFaceVectorQuantity("Gradients", gradients);
118 psMesh->addFaceVectorQuantity("co-Gradients", cogradients);
119 psMesh->addFaceVectorQuantity("Normals", normals);
120 psMesh->addFaceScalarQuantity("Face area", faceArea);
121 psMesh->addFaceVectorQuantity("Vector area", vectorArea);
122
123 polyscope::registerPointCloud("Centroids", centroids);
124}
125
126
127void myCallback()
128{
129 ImGui::SliderFloat("Phi scale", &scale, 0., 1.);
130 if (ImGui::Button("Init phi"))
131 initPhi();
132
133 if (ImGui::Button("Compute quantities"))
135
136}
137
138int main()
139{
141
142 auto h=.3 ; //gridstep
143 params( "polynomial", "goursat" )( "gridstep", h );
144 auto implicit_shape = SH3::makeImplicitShape3D ( params );
145 auto digitized_shape = SH3::makeDigitizedImplicitShape3D( implicit_shape, params );
146 auto K = SH3::getKSpace( params );
147 auto binary_image = SH3::makeBinaryImage( digitized_shape, params );
149 SH3::Cell2Index c2i;
150 auto primalSurface = SH3::makePrimalSurfaceMesh(c2i, surface);
151
152 // Convert faces to appropriate indexed format
153 std::vector<std::vector<SH3::SurfaceMesh::Vertex>> faces;
154 for(auto face= 0 ; face < primalSurface->nbFaces(); ++face)
155 faces.push_back(primalSurface->incidentVertices( face ));
156
157 //Recasting to vector of vertices
158 auto positions = primalSurface->positions();
159
160 surfmesh = SurfMesh(positions.begin(),
161 positions.end(),
162 faces.begin(),
163 faces.end());
164
165 // Initialize polyscope
166 polyscope::init();
167
168 psMesh = polyscope::registerSurfaceMesh("digital surface", positions, faces);
169
170 // Set the callback function
171 polyscope::state::userCallback = myCallback;
172 polyscope::show();
173 return EXIT_SUCCESS;
174}
Implements differential operators on polygonal surfaces from degoes2020discrete.
LinAlg::DenseVector Vector
Type of Vector.
Aim: This class is used to simplify shape and surface creation. With it, you can create new shapes an...
Aim: This class is used to simplify shape and surface creation. With it, you can create new shapes an...
Definition Shortcuts.h:106
static KSpace getKSpace(const Point &low, const Point &up, Parameters params=parametersKSpace())
Definition Shortcuts.h:333
static CountedPtr< DigitizedImplicitShape3D > makeDigitizedImplicitShape3D(CountedPtr< ImplicitShape3D > shape, Parameters params=parametersDigitizedImplicitShape3D())
Definition Shortcuts.h:524
std::map< Cell, IdxVertex > Cell2Index
Definition Shortcuts.h:190
static CountedPtr< DigitalSurface > makeDigitalSurface(CountedPtr< TPointPredicate > bimage, const KSpace &K, const Parameters &params=parametersDigitalSurface())
Definition Shortcuts.h:1210
static Parameters defaultParameters()
Definition Shortcuts.h:204
static CountedPtr< SurfaceMesh > makePrimalSurfaceMesh(Cell2Index &c2i, CountedPtr< ::DGtal::DigitalSurface< TContainer > > aSurface)
Definition Shortcuts.h:2373
static CountedPtr< BinaryImage > makeBinaryImage(Domain shapeDomain)
Definition Shortcuts.h:562
static CountedPtr< ImplicitShape3D > makeImplicitShape3D(const Parameters &params=parametersImplicitShape3D())
Definition Shortcuts.h:283
SurfaceMesh< RealPoint, RealVector > SurfMesh
float scale
double phiVertex(const Vertex v)
void initPhi()
void myCallback()
polyscope::SurfaceMesh * psMesh
SurfMesh surfmesh
void initQuantities()
PolyCalculus * calculus
CountedPtr< SH3::DigitalSurface > surface
CountedPtr< SH3::BinaryImage > binary_image
SMesh::Vertices Vertices
Z3i this namespace gathers the standard of types for 3D imagery.
Space::RealPoint RealPoint
Definition StdDefs.h:170
DGtal is the top-level namespace which contains all DGtal functions and types.
std::pair< typename graph_traits< DGtal::DigitalSurface< TDigitalSurfaceContainer > >::vertex_iterator, typename graph_traits< DGtal::DigitalSurface< TDigitalSurfaceContainer > >::vertex_iterator > vertices(const DGtal::DigitalSurface< TDigitalSurfaceContainer > &digSurf)
Shortcuts< KSpace > SH3
int main()
Definition testBits.cpp:56
KSpace K
ShortcutsGeometry< Z3i::KSpace > SHG3
TriMesh::Face Face
TriMesh::Vertex Vertex