cut_mesh.cpp 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327
  1. #include <igl/cut_mesh.h>
  2. #include <igl/vertex_triangle_adjacency.h>
  3. #include <igl/triangle_triangle_adjacency.h>
  4. #include <igl/is_border_vertex.h>
  5. #include <igl/HalfEdgeIterator.h>
  6. #include <set>
  7. namespace igl {
  8. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  9. class MeshCutterMini
  10. {
  11. public:
  12. // Input
  13. //mesh
  14. const Eigen::PlainObjectBase<DerivedV> &V;
  15. const Eigen::PlainObjectBase<DerivedF> &F;
  16. const Eigen::PlainObjectBase<DerivedTT> &TT;
  17. const Eigen::PlainObjectBase<DerivedTT> &TTi;
  18. const std::vector<std::vector<VFType> >& VF;
  19. const std::vector<std::vector<VFType> >& VFi;
  20. const std::vector<bool> &V_border; // bool
  21. //edges to cut
  22. const Eigen::PlainObjectBase<DerivedC> &Handle_Seams; // 3 bool
  23. // total number of scalar variables
  24. int num_scalar_variables;
  25. // per face indexes of vertex in the solver
  26. #warning "Constructing Eigen::PlainObjectBase directly is deprecated"
  27. Eigen::PlainObjectBase<DerivedF> HandleS_Index;
  28. // per vertex variable indexes
  29. std::vector<std::vector<int> > HandleV_Integer;
  30. IGL_INLINE MeshCutterMini(const Eigen::PlainObjectBase<DerivedV> &_V,
  31. const Eigen::PlainObjectBase<DerivedF> &_F,
  32. const Eigen::PlainObjectBase<DerivedTT> &_TT,
  33. const Eigen::PlainObjectBase<DerivedTT> &_TTi,
  34. const std::vector<std::vector<VFType> > &_VF,
  35. const std::vector<std::vector<VFType> > &_VFi,
  36. const std::vector<bool> &_V_border,
  37. const Eigen::PlainObjectBase<DerivedC> &_Handle_Seams);
  38. // vertex to variable mapping
  39. // initialize the mapping for a given sampled mesh
  40. IGL_INLINE void InitMappingSeam();
  41. private:
  42. IGL_INLINE void FirstPos(const int v, int &f, int &edge);
  43. IGL_INLINE int AddNewIndex(const int v0);
  44. IGL_INLINE bool IsSeam(const int f0, const int f1);
  45. // find initial position of the pos to
  46. // assing face to vert inxex correctly
  47. IGL_INLINE void FindInitialPos(const int vert, int &edge, int &face);
  48. // intialize the mapping given an initial pos
  49. // whih must be initialized with FindInitialPos
  50. IGL_INLINE void MapIndexes(const int vert, const int edge_init, const int f_init);
  51. // intialize the mapping for a given vertex
  52. IGL_INLINE void InitMappingSeam(const int vert);
  53. };
  54. }
  55. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  56. IGL_INLINE igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  57. MeshCutterMini(const Eigen::PlainObjectBase<DerivedV> &_V,
  58. const Eigen::PlainObjectBase<DerivedF> &_F,
  59. const Eigen::PlainObjectBase<DerivedTT> &_TT,
  60. const Eigen::PlainObjectBase<DerivedTT> &_TTi,
  61. const std::vector<std::vector<VFType> > &_VF,
  62. const std::vector<std::vector<VFType> > &_VFi,
  63. const std::vector<bool> &_V_border,
  64. const Eigen::PlainObjectBase<DerivedC> &_Handle_Seams):
  65. V(_V),
  66. F(_F),
  67. TT(_TT),
  68. TTi(_TTi),
  69. VF(_VF),
  70. VFi(_VFi),
  71. V_border(_V_border),
  72. Handle_Seams(_Handle_Seams)
  73. {
  74. num_scalar_variables=0;
  75. HandleS_Index.setConstant(F.rows(),3,-1);
  76. HandleV_Integer.resize(V.rows());
  77. }
  78. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  79. IGL_INLINE void igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  80. FirstPos(const int v, int &f, int &edge)
  81. {
  82. f = VF[v][0]; // f=v->cVFp();
  83. edge = VFi[v][0]; // edge=v->cVFi();
  84. }
  85. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  86. IGL_INLINE int igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  87. AddNewIndex(const int v0)
  88. {
  89. num_scalar_variables++;
  90. HandleV_Integer[v0].push_back(num_scalar_variables);
  91. return num_scalar_variables;
  92. }
  93. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  94. IGL_INLINE bool igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  95. IsSeam(const int f0, const int f1)
  96. {
  97. for (int i=0;i<3;i++)
  98. {
  99. int f_clos = TT(f0,i);
  100. if (f_clos == -1)
  101. continue; ///border
  102. if (f_clos == f1)
  103. return(Handle_Seams(f0,i));
  104. }
  105. assert(0);
  106. return false;
  107. }
  108. ///find initial position of the pos to
  109. // assing face to vert inxex correctly
  110. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  111. IGL_INLINE void igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  112. FindInitialPos(const int vert,
  113. int &edge,
  114. int &face)
  115. {
  116. int f_init;
  117. int edge_init;
  118. FirstPos(vert,f_init,edge_init); // todo manually the function
  119. igl::HalfEdgeIterator<DerivedF> VFI(F,TT,TTi,f_init,edge_init);
  120. bool vertexB = V_border[vert];
  121. bool possible_split=false;
  122. bool complete_turn=false;
  123. do
  124. {
  125. int curr_f = VFI.Fi();
  126. int curr_edge=VFI.Ei();
  127. VFI.NextFE();
  128. int next_f=VFI.Fi();
  129. ///test if I've just crossed a border
  130. bool on_border=(TT(curr_f,curr_edge)==-1);
  131. //bool mismatch=false;
  132. bool seam=false;
  133. ///or if I've just crossed a seam
  134. ///if I'm on a border I MUST start from the one next t othe border
  135. if (!vertexB)
  136. //seam=curr_f->IsSeam(next_f);
  137. seam=IsSeam(curr_f,next_f);
  138. // if (vertexB)
  139. // assert(!Handle_Singular(vert));
  140. // ;
  141. //assert(!vert->IsSingular());
  142. possible_split=((on_border)||(seam));
  143. complete_turn = next_f == f_init;
  144. } while ((!possible_split)&&(!complete_turn));
  145. face=VFI.Fi();
  146. edge=VFI.Ei();
  147. }
  148. ///intialize the mapping given an initial pos
  149. ///whih must be initialized with FindInitialPos
  150. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  151. IGL_INLINE void igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  152. MapIndexes(const int vert,
  153. const int edge_init,
  154. const int f_init)
  155. {
  156. ///check that is not on border..
  157. ///in such case maybe it's non manyfold
  158. ///insert an initial index
  159. int curr_index=AddNewIndex(vert);
  160. ///and initialize the jumping pos
  161. igl::HalfEdgeIterator<DerivedF> VFI(F,TT,TTi,f_init,edge_init);
  162. bool complete_turn=false;
  163. do
  164. {
  165. int curr_f = VFI.Fi();
  166. int curr_edge = VFI.Ei();
  167. ///assing the current index
  168. HandleS_Index(curr_f,curr_edge) = curr_index;
  169. VFI.NextFE();
  170. int next_f = VFI.Fi();
  171. ///test if I've finiseh with the face exploration
  172. complete_turn = (next_f==f_init);
  173. ///or if I've just crossed a mismatch
  174. if (!complete_turn)
  175. {
  176. bool seam=false;
  177. //seam=curr_f->IsSeam(next_f);
  178. seam=IsSeam(curr_f,next_f);
  179. if (seam)
  180. {
  181. ///then add a new index
  182. curr_index=AddNewIndex(vert);
  183. }
  184. }
  185. } while (!complete_turn);
  186. }
  187. ///initialize the mapping for a given vertex
  188. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  189. IGL_INLINE void igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  190. InitMappingSeam(const int vert)
  191. {
  192. ///first rotate until find the first pos after a mismatch
  193. ///or a border or return to the first position...
  194. int f_init = VF[vert][0];
  195. int indexE = VFi[vert][0];
  196. igl::HalfEdgeIterator<DerivedF> VFI(F,TT,TTi,f_init,indexE);
  197. int edge_init;
  198. int face_init;
  199. FindInitialPos(vert,edge_init,face_init);
  200. MapIndexes(vert,edge_init,face_init);
  201. }
  202. ///vertex to variable mapping
  203. ///initialize the mapping for a given sampled mesh
  204. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  205. IGL_INLINE void igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC>::
  206. InitMappingSeam()
  207. {
  208. num_scalar_variables=-1;
  209. for (unsigned int i=0;i<V.rows();i++)
  210. InitMappingSeam(i);
  211. for (unsigned int j=0;j<V.rows();j++)
  212. assert(HandleV_Integer[j].size()>0);
  213. }
  214. template <typename DerivedV, typename DerivedF, typename VFType, typename DerivedTT, typename DerivedC>
  215. IGL_INLINE void igl::cut_mesh(
  216. const Eigen::PlainObjectBase<DerivedV> &V,
  217. const Eigen::PlainObjectBase<DerivedF> &F,
  218. const std::vector<std::vector<VFType> >& VF,
  219. const std::vector<std::vector<VFType> >& VFi,
  220. const Eigen::PlainObjectBase<DerivedTT>& TT,
  221. const Eigen::PlainObjectBase<DerivedTT>& TTi,
  222. const std::vector<bool> &V_border,
  223. const Eigen::PlainObjectBase<DerivedC> &cuts,
  224. Eigen::PlainObjectBase<DerivedV> &Vcut,
  225. Eigen::PlainObjectBase<DerivedF> &Fcut)
  226. {
  227. //finding the cuts is done, now we need to actually generate a cut mesh
  228. igl::MeshCutterMini<DerivedV, DerivedF, VFType, DerivedTT, DerivedC> mc(V, F, TT, TTi, VF, VFi, V_border, cuts);
  229. mc.InitMappingSeam();
  230. Fcut = mc.HandleS_Index;
  231. //we have the faces, we need the vertices;
  232. int newNumV = Fcut.maxCoeff()+1;
  233. Vcut.setZero(newNumV,3);
  234. for (int vi=0; vi<V.rows(); ++vi)
  235. for (int i=0; i<mc.HandleV_Integer[vi].size();++i)
  236. Vcut.row(mc.HandleV_Integer[vi][i]) = V.row(vi);
  237. //ugly hack to fix some problematic cases (border vertex that is also on the boundary of the hole
  238. for (int fi =0; fi<Fcut.rows(); ++fi)
  239. for (int k=0; k<3; ++k)
  240. if (Fcut(fi,k)==-1)
  241. {
  242. //we need to add a vertex
  243. Fcut(fi,k) = newNumV;
  244. newNumV ++;
  245. Vcut.conservativeResize(newNumV, Eigen::NoChange);
  246. Vcut.row(newNumV-1) = V.row(F(fi,k));
  247. }
  248. }
  249. //Wrapper of the above with only vertices and faces as mesh input
  250. template <typename DerivedV, typename DerivedF, typename DerivedC>
  251. IGL_INLINE void igl::cut_mesh(
  252. const Eigen::PlainObjectBase<DerivedV> &V,
  253. const Eigen::PlainObjectBase<DerivedF> &F,
  254. const Eigen::PlainObjectBase<DerivedC> &cuts,
  255. Eigen::PlainObjectBase<DerivedV> &Vcut,
  256. Eigen::PlainObjectBase<DerivedF> &Fcut)
  257. {
  258. std::vector<std::vector<int> > VF, VFi;
  259. igl::vertex_triangle_adjacency(V,F,VF,VFi);
  260. Eigen::MatrixXd Vt = V;
  261. Eigen::MatrixXi Ft = F;
  262. Eigen::MatrixXi TT, TTi;
  263. igl::triangle_triangle_adjacency(Ft,TT,TTi);
  264. std::vector<bool> V_border = igl::is_border_vertex(V,F);
  265. igl::cut_mesh(V, F, VF, VFi, TT, TTi, V_border, cuts, Vcut, Fcut);
  266. }
  267. #ifdef IGL_STATIC_LIBRARY
  268. // Explicit template specialization
  269. template void igl::cut_mesh<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, int, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, std::vector<bool, std::allocator<bool> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
  270. template void igl::cut_mesh<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
  271. template void igl::cut_mesh<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);
  272. template void igl::cut_mesh<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
  273. #endif