Developer Documentation
Sqrt3T.hh
Go to the documentation of this file.
1 /* ========================================================================= *
2  * *
3  * OpenMesh *
4  * Copyright (c) 2001-2015, RWTH-Aachen University *
5  * Department of Computer Graphics and Multimedia *
6  * All rights reserved. *
7  * www.openmesh.org *
8  * *
9  *---------------------------------------------------------------------------*
10  * This file is part of OpenMesh. *
11  *---------------------------------------------------------------------------*
12  * *
13  * Redistribution and use in source and binary forms, with or without *
14  * modification, are permitted provided that the following conditions *
15  * are met: *
16  * *
17  * 1. Redistributions of source code must retain the above copyright notice, *
18  * this list of conditions and the following disclaimer. *
19  * *
20  * 2. Redistributions in binary form must reproduce the above copyright *
21  * notice, this list of conditions and the following disclaimer in the *
22  * documentation and/or other materials provided with the distribution. *
23  * *
24  * 3. Neither the name of the copyright holder nor the names of its *
25  * contributors may be used to endorse or promote products derived from *
26  * this software without specific prior written permission. *
27  * *
28  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS *
29  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED *
30  * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A *
31  * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER *
32  * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, *
33  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, *
34  * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR *
35  * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF *
36  * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING *
37  * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS *
38  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
39  * *
40  * ========================================================================= */
41 
42 /*===========================================================================*\
43  * *
44  * $Revision$ *
45  * $Date$ *
46  * *
47 \*===========================================================================*/
48 
53 //=============================================================================
54 //
55 // CLASS Sqrt3T
56 //
57 //=============================================================================
58 
59 #ifndef OPENMESH_SUBDIVIDER_UNIFORM_SQRT3T_HH
60 #define OPENMESH_SUBDIVIDER_UNIFORM_SQRT3T_HH
61 
62 
63 //== INCLUDES =================================================================
64 
65 #include <OpenMesh/Core/Mesh/Handles.hh>
66 #include <OpenMesh/Core/System/config.hh>
68 #if defined(_DEBUG) || defined(DEBUG)
69 // Makes life lot easier, when playing/messing around with low-level topology
70 // changing methods of OpenMesh
71 # include <OpenMesh/Tools/Utils/MeshCheckerT.hh>
72 # define ASSERT_CONSISTENCY( T, m ) \
73  assert(OpenMesh::Utils::MeshCheckerT<T>(m).check())
74 #else
75 # define ASSERT_CONSISTENCY( T, m )
76 #endif
77 // -------------------- STL
78 #include <vector>
79 #if defined(OM_CC_MIPS)
80 # include <math.h>
81 #else
82 # include <cmath>
83 #endif
84 
85 
86 //== NAMESPACE ================================================================
87 
88 namespace OpenMesh { // BEGIN_NS_OPENMESH
89 namespace Subdivider { // BEGIN_NS_DECIMATER
90 namespace Uniform { // BEGIN_NS_DECIMATER
91 
92 
93 //== CLASS DEFINITION =========================================================
94 
95 
102 template <typename MeshType, typename RealType = float>
103 class Sqrt3T : public SubdividerT< MeshType, RealType >
104 {
105 public:
106 
107  typedef RealType real_t;
108  typedef MeshType mesh_t;
110 
111  typedef std::pair< real_t, real_t > weight_t;
112  typedef std::vector< std::pair<real_t,real_t> > weights_t;
113 
114 public:
115 
116 
117  Sqrt3T(void) : parent_t(), _1over3( real_t(1.0/3.0) ), _1over27( real_t(1.0/27.0) )
118  { init_weights(); }
119 
120  Sqrt3T(MeshType &_m) : parent_t(_m), _1over3( real_t(1.0/3.0) ), _1over27( real_t(1.0/27.0) )
121  { init_weights(); }
122 
123  virtual ~Sqrt3T() {}
124 
125 
126 public:
127 
128 
129  const char *name() const { return "Uniform Sqrt3"; }
130 
131 
133  void init_weights(size_t _max_valence=50)
134  {
135  weights_.resize(_max_valence);
136  std::generate(weights_.begin(), weights_.end(), compute_weight());
137  }
138 
139 
140 protected:
141 
142 
143  bool prepare( MeshType& _m )
144  {
145  _m.request_edge_status();
146  _m.add_property( vp_pos_ );
147  _m.add_property( ep_nv_ );
148  _m.add_property( mp_gen_ );
149  _m.property( mp_gen_ ) = 0;
150 
151  return _m.has_edge_status() && vp_pos_.is_valid()
152  && ep_nv_.is_valid() && mp_gen_.is_valid();
153  }
154 
155 
156  bool cleanup( MeshType& _m )
157  {
158  _m.release_edge_status();
159  _m.remove_property( vp_pos_ );
160  _m.remove_property( ep_nv_ );
161  _m.remove_property( mp_gen_ );
162  return true;
163  }
164 
165  bool subdivide( MeshType& _m, size_t _n , const bool _update_points = true)
166  {
167 
169 
170  typename MeshType::VertexIter vit;
171  typename MeshType::VertexVertexIter vvit;
172  typename MeshType::EdgeIter eit;
173  typename MeshType::FaceIter fit;
174  typename MeshType::FaceVertexIter fvit;
175  typename MeshType::VertexHandle vh;
176  typename MeshType::HalfedgeHandle heh;
177  typename MeshType::Point pos(0,0,0), zero(0,0,0);
178  size_t &gen = _m.property( mp_gen_ );
179 
180  for (size_t l=0; l<_n; ++l)
181  {
182  // tag existing edges
183  for (eit=_m.edges_begin(); eit != _m.edges_end();++eit)
184  {
185  _m.status( *eit ).set_tagged( true );
186  if ( (gen%2) && _m.is_boundary(*eit) )
187  compute_new_boundary_points( _m, *eit ); // *) creates new vertices
188  }
189 
190  // do relaxation of old vertices, but store new pos in property vp_pos_
191 
192  for (vit=_m.vertices_begin(); vit!=_m.vertices_end(); ++vit)
193  {
194  if ( _m.is_boundary(*vit) )
195  {
196  if ( gen%2 )
197  {
198  heh = _m.halfedge_handle(*vit);
199  if (heh.is_valid()) // skip isolated newly inserted vertices *)
200  {
201  typename OpenMesh::HalfedgeHandle
202  prev_heh = _m.prev_halfedge_handle(heh);
203 
204  assert( _m.is_boundary(heh ) );
205  assert( _m.is_boundary(prev_heh) );
206 
207  pos = _m.point(_m.to_vertex_handle(heh));
208  pos += _m.point(_m.from_vertex_handle(prev_heh));
209  pos *= real_t(4.0);
210 
211  pos += real_t(19.0) * _m.point( *vit );
212  pos *= _1over27;
213 
214  _m.property( vp_pos_, *vit ) = pos;
215  }
216  }
217  else
218  _m.property( vp_pos_, *vit ) = _m.point( *vit );
219  }
220  else
221  {
222  size_t valence=0;
223 
224  pos = zero;
225  for ( vvit = _m.vv_iter(*vit); vvit.is_valid(); ++vvit)
226  {
227  pos += _m.point( *vvit );
228  ++valence;
229  }
230  pos *= weights_[ valence ].second;
231  pos += weights_[ valence ].first * _m.point(*vit);
232  _m.property( vp_pos_, *vit ) = pos;
233  }
234  }
235 
236  // insert new vertices, but store pos in vp_pos_
237  typename MeshType::FaceIter fend = _m.faces_end();
238  for (fit = _m.faces_begin();fit != fend; ++fit)
239  {
240  if ( (gen%2) && _m.is_boundary(*fit))
241  {
242  boundary_split( _m, *fit );
243  }
244  else
245  {
246  fvit = _m.fv_iter( *fit );
247  pos = _m.point( *fvit);
248  pos += _m.point(*(++fvit));
249  pos += _m.point(*(++fvit));
250  pos *= _1over3;
251  vh = _m.add_vertex( zero );
252  _m.property( vp_pos_, vh ) = pos;
253  _m.split( *fit, vh );
254  }
255  }
256 
257  // commit new positions (now iterating over all vertices)
258  for (vit=_m.vertices_begin();vit != _m.vertices_end(); ++vit)
259  _m.set_point(*vit, _m.property( vp_pos_, *vit ) );
260 
261  // flip old edges
262  for (eit=_m.edges_begin(); eit != _m.edges_end(); ++eit)
263  if ( _m.status( *eit ).tagged() && !_m.is_boundary( *eit ) )
264  _m.flip(*eit);
265 
266  // Now we have an consistent mesh!
267  ASSERT_CONSISTENCY( MeshType, _m );
268 
269  // increase generation by one
270  ++gen;
271  }
272  return true;
273  }
274 
275 private:
276 
280  {
281  compute_weight() : valence(-1) { }
282  weight_t operator() (void)
283  {
284 #if !defined(OM_CC_MIPS)
285  using std::cos;
286 #endif
287  if (++valence)
288  {
289  real_t alpha = real_t( (4.0-2.0*cos(2.0*M_PI / real_t(valence)) )/9.0 );
290  return weight_t( real_t(1)-alpha, alpha/real_t(valence) );
291  }
292  return weight_t(real_t(0.0), real_t(0.0) );
293  }
294  int valence;
295  };
296 
297 private:
298 
299  // Pre-compute location of new boundary points for odd generations
300  // and store them in the edge property ep_nv_;
301  void compute_new_boundary_points( MeshType& _m,
302  const typename MeshType::EdgeHandle& _eh)
303  {
304  assert( _m.is_boundary(_eh) );
305 
306  typename MeshType::HalfedgeHandle heh;
307  typename MeshType::VertexHandle vh1, vh2, vh3, vh4, vhl, vhr;
308  typename MeshType::Point zero(0,0,0), P1, P2, P3, P4;
309 
310  /*
311  // *---------*---------*
312  // / \ / \ / \
313  // / \ / \ / \
314  // / \ / \ / \
315  // / \ / \ / \
316  // *---------*--#---#--*---------*
317  //
318  // ^ ^ ^ ^ ^ ^
319  // P1 P2 pl pr P3 P4
320  */
321  // get halfedge pointing from P3 to P2 (outer boundary halfedge)
322 
323  heh = _m.halfedge_handle(_eh,
324  _m.is_boundary(_m.halfedge_handle(_eh,1)));
325 
326  assert( _m.is_boundary( _m.next_halfedge_handle( heh ) ) );
327  assert( _m.is_boundary( _m.prev_halfedge_handle( heh ) ) );
328 
329  vh1 = _m.to_vertex_handle( _m.next_halfedge_handle( heh ) );
330  vh2 = _m.to_vertex_handle( heh );
331  vh3 = _m.from_vertex_handle( heh );
332  vh4 = _m.from_vertex_handle( _m.prev_halfedge_handle( heh ));
333 
334  P1 = _m.point(vh1);
335  P2 = _m.point(vh2);
336  P3 = _m.point(vh3);
337  P4 = _m.point(vh4);
338 
339  vhl = _m.add_vertex(zero);
340  vhr = _m.add_vertex(zero);
341 
342  _m.property(vp_pos_, vhl ) = (P1 + real_t(16.0f) * P2 + real_t(10.0f) * P3) * _1over27;
343  _m.property(vp_pos_, vhr ) = ( real_t(10.0f) * P2 + real_t(16.0f) * P3 + P4) * _1over27;
344  _m.property(ep_nv_, _eh).first = vhl;
345  _m.property(ep_nv_, _eh).second = vhr;
346  }
347 
348 
349  void boundary_split( MeshType& _m, const typename MeshType::FaceHandle& _fh )
350  {
351  assert( _m.is_boundary(_fh) );
352 
353  typename MeshType::VertexHandle vhl, vhr;
354  typename MeshType::FaceEdgeIter fe_it;
355  typename MeshType::HalfedgeHandle heh;
356 
357  // find boundary edge
358  for( fe_it=_m.fe_iter( _fh ); fe_it.is_valid() && !_m.is_boundary( *fe_it ); ++fe_it ) {};
359 
360  // use precomputed, already inserted but not linked vertices
361  vhl = _m.property(ep_nv_, *fe_it).first;
362  vhr = _m.property(ep_nv_, *fe_it).second;
363 
364  /*
365  // *---------*---------*
366  // / \ / \ / \
367  // / \ / \ / \
368  // / \ / \ / \
369  // / \ / \ / \
370  // *---------*--#---#--*---------*
371  //
372  // ^ ^ ^ ^ ^ ^
373  // P1 P2 pl pr P3 P4
374  */
375  // get halfedge pointing from P2 to P3 (inner boundary halfedge)
376 
377  heh = _m.halfedge_handle(*fe_it,
378  _m.is_boundary(_m.halfedge_handle(*fe_it,0)));
379 
380  typename MeshType::HalfedgeHandle pl_P3;
381 
382  // split P2->P3 (heh) in P2->pl (heh) and pl->P3
383  boundary_split( _m, heh, vhl ); // split edge
384  pl_P3 = _m.next_halfedge_handle( heh ); // store next halfedge handle
385  boundary_split( _m, heh ); // split face
386 
387  // split pl->P3 in pl->pr and pr->P3
388  boundary_split( _m, pl_P3, vhr );
389  boundary_split( _m, pl_P3 );
390 
391  assert( _m.is_boundary( vhl ) && _m.halfedge_handle(vhl).is_valid() );
392  assert( _m.is_boundary( vhr ) && _m.halfedge_handle(vhr).is_valid() );
393  }
394 
395  void boundary_split(MeshType& _m,
396  const typename MeshType::HalfedgeHandle& _heh,
397  const typename MeshType::VertexHandle& _vh)
398  {
399  assert( _m.is_boundary( _m.edge_handle(_heh) ) );
400 
401  typename MeshType::HalfedgeHandle
402  heh(_heh),
403  opp_heh( _m.opposite_halfedge_handle(_heh) ),
404  new_heh, opp_new_heh;
405  typename MeshType::VertexHandle to_vh(_m.to_vertex_handle(heh));
406  typename MeshType::HalfedgeHandle t_heh;
407 
408  /*
409  * P5
410  * *
411  * /|\
412  * / \
413  * / \
414  * / \
415  * / \
416  * /_ heh new \
417  * *-----\*-----\*\-----*
418  * ^ ^ t_heh
419  * _vh to_vh
420  *
421  * P1 P2 P3 P4
422  */
423  // Re-Setting Handles
424 
425  // find halfedge point from P4 to P3
426  for(t_heh = heh;
427  _m.next_halfedge_handle(t_heh) != opp_heh;
428  t_heh = _m.opposite_halfedge_handle(_m.next_halfedge_handle(t_heh)))
429  {}
430 
431  assert( _m.is_boundary( t_heh ) );
432 
433  new_heh = _m.new_edge( _vh, to_vh );
434  opp_new_heh = _m.opposite_halfedge_handle(new_heh);
435 
436  // update halfedge connectivity
437 
438  _m.set_next_halfedge_handle(t_heh, opp_new_heh); // P4-P3 -> P3-P2
439  // P2-P3 -> P3-P5
440  _m.set_next_halfedge_handle(new_heh, _m.next_halfedge_handle(heh));
441  _m.set_next_halfedge_handle(heh, new_heh); // P1-P2 -> P2-P3
442  _m.set_next_halfedge_handle(opp_new_heh, opp_heh); // P3-P2 -> P2-P1
443 
444  // both opposite halfedges point to same face
445  _m.set_face_handle(opp_new_heh, _m.face_handle(opp_heh));
446 
447  // let heh finally point to new inserted vertex
448  _m.set_vertex_handle(heh, _vh);
449 
450  // let heh and new_heh point to same face
451  _m.set_face_handle(new_heh, _m.face_handle(heh));
452 
453  // let opp_new_heh be the new outgoing halfedge for to_vh
454  // (replaces for opp_heh)
455  _m.set_halfedge_handle( to_vh, opp_new_heh );
456 
457  // let opp_heh be the outgoing halfedge for _vh
458  _m.set_halfedge_handle( _vh, opp_heh );
459  }
460 
461  void boundary_split( MeshType& _m,
462  const typename MeshType::HalfedgeHandle& _heh)
463  {
464  assert( _m.is_boundary( _m.opposite_halfedge_handle( _heh ) ) );
465 
466  typename MeshType::HalfedgeHandle
467  heh(_heh),
468  n_heh(_m.next_halfedge_handle(heh));
469 
470  typename MeshType::VertexHandle
471  to_vh(_m.to_vertex_handle(heh));
472 
473  typename MeshType::HalfedgeHandle
474  heh2(_m.new_edge(to_vh,
475  _m.to_vertex_handle(_m.next_halfedge_handle(n_heh)))),
476  heh3(_m.opposite_halfedge_handle(heh2));
477 
478  typename MeshType::FaceHandle
479  new_fh(_m.new_face()),
480  fh(_m.face_handle(heh));
481 
482  // Relink (half)edges
483 
484 #define set_next_heh set_next_halfedge_handle
485 #define next_heh next_halfedge_handle
486 
487  _m.set_face_handle(heh, new_fh);
488  _m.set_face_handle(heh2, new_fh);
489  _m.set_next_heh(heh2, _m.next_heh(_m.next_heh(n_heh)));
490  _m.set_next_heh(heh, heh2);
491  _m.set_face_handle( _m.next_heh(heh2), new_fh);
492 
493  // _m.set_face_handle( _m.next_heh(_m.next_heh(heh2)), new_fh);
494 
495  _m.set_next_heh(heh3, n_heh);
496  _m.set_next_heh(_m.next_halfedge_handle(n_heh), heh3);
497  _m.set_face_handle(heh3, fh);
498  // _m.set_face_handle(n_heh, fh);
499 
500  _m.set_halfedge_handle( fh, n_heh);
501  _m.set_halfedge_handle(new_fh, heh);
502 
503 #undef set_next_halfedge_handle
504 #undef next_halfedge_handle
505 
506  }
507 
508 private:
509 
510  weights_t weights_;
512  OpenMesh::EPropHandleT< std::pair< typename MeshType::VertexHandle,
513  typename MeshType::VertexHandle> > ep_nv_;
515 
516  const real_t _1over3;
517  const real_t _1over27;
518 };
519 
520 
521 //=============================================================================
522 } // END_NS_UNIFORM
523 } // END_NS_SUBDIVIDER
524 } // END_NS_OPENMESH
525 //=============================================================================
526 #endif // OPENMESH_SUBDIVIDER_UNIFORM_SQRT3T_HH
527 //=============================================================================
bool cleanup(MeshType &_m)
Cleanup mesh after usage, e.g. remove added properties.
Definition: Sqrt3T.hh:156
bool subdivide(MeshType &_m, size_t _n, const bool _update_points=true)
Subdivide mesh _m _n times.
Definition: Sqrt3T.hh:165
bool is_valid() const
The handle is valid iff the index is not equal to -1.
Definition: Handles.hh:77
void init_weights(size_t _max_valence=50)
Pre-compute weights.
Definition: Sqrt3T.hh:133
Handle for a halfedge entity.
Definition: Handles.hh:132
bool prepare(MeshType &_m)
Prepare mesh, e.g. add properties.
Definition: Sqrt3T.hh:143
const char * name() const
Return name of subdivision algorithm.
Definition: Sqrt3T.hh:129