Developer Documentation
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
RulesT.cc
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 // Rules - IMPLEMENTATION
56 //
57 //=============================================================================
58 
59 
60 #define OPENMESH_SUBDIVIDER_ADAPTIVE_RULEST_CC
61 
62 
63 //== INCLUDES =================================================================
64 
65 #include <OpenMesh/Core/System/config.h>
66 #include <OpenMesh/Core/IO/MeshIO.hh>
67 #include "RulesT.hh"
68 // --------------------
69 #if defined(OM_CC_MIPS)
70 # include <math.h>
71 #else
72 # include <cmath>
73 #endif
74 
75 #if defined(OM_CC_MSVC)
76 # pragma warning(disable:4244)
77 #endif
78 
79 //== NAMESPACE ================================================================
80 
81 namespace OpenMesh { // BEGIN_NS_OPENMESH
82 namespace Subdivider { // BEGIN_NS_DECIMATER
83 namespace Adaptive { // BEGIN_NS_ADAPTIVE
84 
85 
86 //== IMPLEMENTATION ==========================================================
87 
88 #define MOBJ Base::mesh_.data
89 #define FH face_handle
90 #define VH vertex_handle
91 #define EH edge_handle
92 #define HEH halfedge_handle
93 #define NHEH next_halfedge_handle
94 #define PHEH prev_halfedge_handle
95 #define OHEH opposite_halfedge_handle
96 #define TVH to_vertex_handle
97 #define FVH from_vertex_handle
98 
99 // ------------------------------------------------------------------ Tvv3 ----
100 
101 
102 template<class M>
103 void
104 Tvv3<M>::raise(typename M::FaceHandle& _fh, state_t _target_state)
105 {
106  if (MOBJ(_fh).state() < _target_state)
107  {
108  this->update(_fh, _target_state);
109 
110  typename M::VertexVertexIter vv_it;
111  typename M::FaceVertexIter fv_it;
112  typename M::VertexHandle vh;
113  typename M::Point position(0.0, 0.0, 0.0);
114  typename M::Point face_position;
115  const typename M::Point zero_point(0.0, 0.0, 0.0);
116  std::vector<typename M::VertexHandle> vertex_vector;
117 
118 
119  // raise all adjacent vertices to level x-1
120  for (fv_it = Base::mesh_.fv_iter(_fh); fv_it.is_valid(); ++fv_it) {
121 
122  vertex_vector.push_back(*fv_it);
123  }
124 
125  while(!vertex_vector.empty()) {
126 
127  vh = vertex_vector.back();
128  vertex_vector.pop_back();
129 
130  if (_target_state > 1)
131  Base::prev_rule()->raise(vh, _target_state - 1);
132  }
133 
134  face_position = MOBJ(_fh).position(_target_state - 1);
135 
136  typename M::EdgeHandle eh;
137  std::vector<typename M::EdgeHandle> edge_vector;
138 
139  // interior face
140  if (!Base::mesh_.is_boundary(_fh) || MOBJ(_fh).final()) {
141 
142  // insert new vertex
143  vh = Base::mesh_.new_vertex();
144 
145  Base::mesh_.split(_fh, vh);
146 
147  typename M::Scalar valence(0.0);
148 
149  // calculate display position for new vertex
150  for (vv_it = Base::mesh_.vv_iter(vh); vv_it.is_valid(); ++vv_it)
151  {
152  position += Base::mesh_.point(*vv_it);
153  valence += 1.0;
154  }
155 
156  position /= valence;
157 
158  // set attributes for new vertex
159  Base::mesh_.set_point(vh, position);
160  MOBJ(vh).set_position(_target_state, zero_point);
161  MOBJ(vh).set_state(_target_state);
162  MOBJ(vh).set_not_final();
163 
164  typename M::VertexOHalfedgeIter voh_it;
165  // check for edge flipping
166  for (voh_it = Base::mesh_.voh_iter(vh); voh_it.is_valid(); ++voh_it) {
167 
168  if (Base::mesh_.FH(*voh_it).is_valid()) {
169 
170  MOBJ(Base::mesh_.FH(*voh_it)).set_state(_target_state);
171  MOBJ(Base::mesh_.FH(*voh_it)).set_not_final();
172  MOBJ(Base::mesh_.FH(*voh_it)).set_position(_target_state - 1, face_position);
173 
174 
175  for (state_t j = 0; j < _target_state; ++j) {
176  MOBJ(Base::mesh_.FH(*voh_it)).set_position(j, MOBJ(_fh).position(j));
177  }
178 
179  if (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))).is_valid()) {
180 
181  if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it)))).state() == _target_state) {
182 
183  if (Base::mesh_.is_flip_ok(Base::mesh_.EH(Base::mesh_.NHEH(*voh_it)))) {
184 
185  edge_vector.push_back(Base::mesh_.EH(Base::mesh_.NHEH(*voh_it)));
186  }
187  }
188  }
189  }
190  }
191  }
192 
193  // boundary face
194  else {
195 
196  typename M::VertexHandle vh1 = Base::mesh_.new_vertex(),
197  vh2 = Base::mesh_.new_vertex();
198 
199  typename M::HalfedgeHandle hh2 = Base::mesh_.HEH(_fh),
200  hh1, hh3;
201 
202  while (!Base::mesh_.is_boundary(Base::mesh_.OHEH(hh2)))
203  hh2 = Base::mesh_.NHEH(hh2);
204 
205  eh = Base::mesh_.EH(hh2);
206 
207  hh2 = Base::mesh_.NHEH(hh2);
208  hh1 = Base::mesh_.NHEH(hh2);
209 
210  assert(Base::mesh_.is_boundary(eh));
211 
212  Base::mesh_.split(eh, vh1);
213 
214  eh = Base::mesh_.EH(Base::mesh_.PHEH(hh2));
215 
216  assert(Base::mesh_.is_boundary(eh));
217 
218  Base::mesh_.split(eh, vh2);
219 
220  hh3 = Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(hh1)));
221 
222  typename M::VertexHandle vh0(Base::mesh_.TVH(hh1)),
223  vh3(Base::mesh_.FVH(hh2));
224 
225  // set display position and attributes for new vertices
226  Base::mesh_.set_point(vh1, (Base::mesh_.point(vh0) * 2.0 + Base::mesh_.point(vh3)) / 3.0);
227 
228  MOBJ(vh1).set_position(_target_state, zero_point);
229  MOBJ(vh1).set_state(_target_state);
230  MOBJ(vh1).set_not_final();
231 
232  MOBJ(vh0).set_position(_target_state, MOBJ(vh0).position(_target_state - 1) * 3.0);
233  MOBJ(vh0).set_state(_target_state);
234  MOBJ(vh0).set_not_final();
235 
236  // set display position and attributes for old vertices
237  Base::mesh_.set_point(vh2, (Base::mesh_.point(vh3) * 2.0 + Base::mesh_.point(vh0)) / 3.0);
238  MOBJ(vh2).set_position(_target_state, zero_point);
239  MOBJ(vh2).set_state(_target_state);
240  MOBJ(vh2).set_not_final();
241 
242  MOBJ(vh3).set_position(_target_state, MOBJ(vh3).position(_target_state - 1) * 3.0);
243  MOBJ(vh3).set_state(_target_state);
244  MOBJ(vh3).set_not_final();
245 
246  // init 3 faces
247  MOBJ(Base::mesh_.FH(hh1)).set_state(_target_state);
248  MOBJ(Base::mesh_.FH(hh1)).set_not_final();
249  MOBJ(Base::mesh_.FH(hh1)).set_position(_target_state - 1, face_position);
250 
251  MOBJ(Base::mesh_.FH(hh2)).set_state(_target_state);
252  MOBJ(Base::mesh_.FH(hh2)).set_not_final();
253  MOBJ(Base::mesh_.FH(hh2)).set_position(_target_state - 1, face_position);
254 
255  MOBJ(Base::mesh_.FH(hh3)).set_state(_target_state);
256  MOBJ(Base::mesh_.FH(hh3)).set_final();
257  MOBJ(Base::mesh_.FH(hh3)).set_position(_target_state - 1, face_position);
258 
259 
260  for (state_t j = 0; j < _target_state; ++j) {
261  MOBJ(Base::mesh_.FH(hh1)).set_position(j, MOBJ(_fh).position(j));
262  }
263 
264  for (state_t j = 0; j < _target_state; ++j) {
265 
266  MOBJ(Base::mesh_.FH(hh2)).set_position(j, MOBJ(_fh).position(j));
267  }
268 
269  for (state_t j = 0; j < _target_state; ++j) {
270 
271  MOBJ(Base::mesh_.FH(hh3)).set_position(j, MOBJ(_fh).position(j));
272  }
273 
274  // check for edge flipping
275  if (Base::mesh_.FH(Base::mesh_.OHEH(hh1)).is_valid()) {
276 
277  if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(hh1))).state() == _target_state) {
278 
279  if (Base::mesh_.is_flip_ok(Base::mesh_.EH(hh1))) {
280 
281  edge_vector.push_back(Base::mesh_.EH(hh1));
282  }
283  }
284  }
285 
286  if (Base::mesh_.FH(Base::mesh_.OHEH(hh2)).is_valid()) {
287 
288  if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(hh2))).state() == _target_state) {
289 
290  if (Base::mesh_.is_flip_ok(Base::mesh_.EH(hh2))) {
291 
292  edge_vector.push_back(Base::mesh_.EH(hh2));
293  }
294  }
295  }
296  }
297 
298  // flip edges
299  while (!edge_vector.empty()) {
300 
301  eh = edge_vector.back();
302  edge_vector.pop_back();
303 
304  assert(Base::mesh_.is_flip_ok(eh));
305 
306  Base::mesh_.flip(eh);
307 
308  MOBJ(Base::mesh_.FH(Base::mesh_.HEH(eh, 0))).set_final();
309  MOBJ(Base::mesh_.FH(Base::mesh_.HEH(eh, 1))).set_final();
310 
311  MOBJ(Base::mesh_.FH(Base::mesh_.HEH(eh, 0))).set_state(_target_state);
312  MOBJ(Base::mesh_.FH(Base::mesh_.HEH(eh, 1))).set_state(_target_state);
313 
314  MOBJ(Base::mesh_.FH(Base::mesh_.HEH(eh, 0))).set_position(_target_state, face_position);
315  MOBJ(Base::mesh_.FH(Base::mesh_.HEH(eh, 1))).set_position(_target_state, face_position);
316  }
317  }
318 }
319 
320 
321 template<class M>
322 void Tvv3<M>::raise(typename M::VertexHandle& _vh, state_t _target_state) {
323 
324  if (MOBJ(_vh).state() < _target_state) {
325 
326  this->update(_vh, _target_state);
327 
328  // multiply old position by 3
329  MOBJ(_vh).set_position(_target_state, MOBJ(_vh).position(_target_state - 1) * 3.0);
330 
331  MOBJ(_vh).inc_state();
332 
333  assert(MOBJ(_vh).state() == _target_state);
334  }
335 }
336 
337 
338 // ------------------------------------------------------------------ Tvv4 ----
339 
340 
341 template<class M>
342 void
343 Tvv4<M>::raise(typename M::FaceHandle& _fh, state_t _target_state)
344 {
345 
346  if (MOBJ(_fh).state() < _target_state) {
347 
348  this->update(_fh, _target_state);
349 
350  typename M::FaceVertexIter fv_it;
351  typename M::VertexHandle temp_vh;
352  typename M::Point face_position;
353  const typename M::Point zero_point(0.0, 0.0, 0.0);
354  std::vector<typename M::VertexHandle> vertex_vector;
355  std::vector<typename M::HalfedgeHandle> halfedge_vector;
356 
357  // raise all adjacent vertices to level x-1
358  for (fv_it = Base::mesh_.fv_iter(_fh); fv_it.is_valid(); ++fv_it) {
359 
360  vertex_vector.push_back(*fv_it);
361  }
362 
363  while(!vertex_vector.empty()) {
364 
365  temp_vh = vertex_vector.back();
366  vertex_vector.pop_back();
367 
368  if (_target_state > 1) {
369  Base::prev_rule()->raise(temp_vh, _target_state - 1);
370  }
371  }
372 
373  face_position = MOBJ(_fh).position(_target_state - 1);
374 
375  typename M::HalfedgeHandle hh[3];
376  typename M::VertexHandle vh[3];
377  typename M::VertexHandle new_vh[3];
378  typename M::FaceHandle fh[4];
379  typename M::EdgeHandle eh;
380  typename M::HalfedgeHandle temp_hh;
381 
382  // normal (final) face
383  if (MOBJ(_fh).final()) {
384 
385  // define three halfedge handles around the face
386  hh[0] = Base::mesh_.HEH(_fh);
387  hh[1] = Base::mesh_.NHEH(hh[0]);
388  hh[2] = Base::mesh_.NHEH(hh[1]);
389 
390  assert(hh[0] == Base::mesh_.NHEH(hh[2]));
391 
392  vh[0] = Base::mesh_.TVH(hh[0]);
393  vh[1] = Base::mesh_.TVH(hh[1]);
394  vh[2] = Base::mesh_.TVH(hh[2]);
395 
396  new_vh[0] = Base::mesh_.add_vertex(zero_point);
397  new_vh[1] = Base::mesh_.add_vertex(zero_point);
398  new_vh[2] = Base::mesh_.add_vertex(zero_point);
399 
400  // split three edges
401  split_edge(hh[0], new_vh[0], _target_state);
402  eh = Base::mesh_.EH(Base::mesh_.PHEH(hh[2]));
403  split_edge(hh[1], new_vh[1], _target_state);
404  split_edge(hh[2], new_vh[2], _target_state);
405 
406  assert(Base::mesh_.FVH(hh[2]) == vh[1]);
407  assert(Base::mesh_.FVH(hh[1]) == vh[0]);
408  assert(Base::mesh_.FVH(hh[0]) == vh[2]);
409 
410  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[0])).is_valid())
411  {
412  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[0])));
413  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
414  halfedge_vector.push_back(temp_hh);
415  }
416 
417  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[1])).is_valid()) {
418 
419  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[1])));
420  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
421  halfedge_vector.push_back(temp_hh);
422  }
423 
424  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[2])).is_valid()) {
425 
426  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[2])));
427  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
428  halfedge_vector.push_back(temp_hh);
429  }
430  }
431 
432  // splitted face, check for type
433  else {
434 
435  // define red halfedge handle
436  typename M::HalfedgeHandle red_hh(MOBJ(_fh).red_halfedge());
437 
438  if (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh))).is_valid()
439  && Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)))).is_valid()
440  && MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh)))).red_halfedge() == red_hh
441  && MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh))))).red_halfedge() == red_hh)
442  {
443 
444  // three times divided face
445  vh[0] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)))));
446  vh[1] = Base::mesh_.TVH(red_hh);
447  vh[2] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh))));
448 
449  new_vh[0] = Base::mesh_.FVH(red_hh);
450  new_vh[1] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)));
451  new_vh[2] = Base::mesh_.TVH(Base::mesh_.NHEH(red_hh));
452 
453  hh[0] = Base::mesh_.PHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh)));
454  hh[1] = Base::mesh_.PHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh))));
455  hh[2] = Base::mesh_.NHEH(red_hh);
456 
457  eh = Base::mesh_.EH(red_hh);
458  }
459 
460  else
461  {
462 
463  if ((Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh))).is_valid() &&
464  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh)))).red_halfedge()
465  == red_hh )
466  || (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)))).is_valid()
467  && MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh))))).red_halfedge() == red_hh))
468  {
469 
470  // double divided face
471  if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh)))).red_halfedge() == red_hh)
472  {
473  // first case
474  vh[0] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)));
475  vh[1] = Base::mesh_.TVH(red_hh);
476  vh[2] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh))));
477 
478  new_vh[0] = Base::mesh_.FVH(red_hh);
479  new_vh[1] = Base::mesh_.add_vertex(zero_point);
480  new_vh[2] = Base::mesh_.TVH(Base::mesh_.NHEH(red_hh));
481 
482  hh[0] = Base::mesh_.PHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(red_hh)));
483  hh[1] = Base::mesh_.PHEH(Base::mesh_.OHEH(red_hh));
484  hh[2] = Base::mesh_.NHEH(red_hh);
485 
486  // split one edge
487  eh = Base::mesh_.EH(red_hh);
488 
489  split_edge(hh[1], new_vh[1], _target_state);
490 
491  assert(Base::mesh_.FVH(hh[2]) == vh[1]);
492  assert(Base::mesh_.FVH(hh[1]) == vh[0]);
493  assert(Base::mesh_.FVH(hh[0]) == vh[2]);
494 
495  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[1])).is_valid())
496  {
497  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[1])));
498  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
499  halfedge_vector.push_back(temp_hh);
500  }
501  }
502  else
503  {
504 
505  // second case
506  vh[0] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)))));
507  vh[1] = Base::mesh_.TVH(red_hh);
508  vh[2] = Base::mesh_.TVH(Base::mesh_.NHEH(red_hh));
509 
510  new_vh[0] = Base::mesh_.FVH(red_hh);
511  new_vh[1] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)));
512  new_vh[2] = Base::mesh_.add_vertex(zero_point);
513 
514  hh[0] = Base::mesh_.PHEH(red_hh);
515  hh[1] = Base::mesh_.PHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh))));
516  hh[2] = Base::mesh_.NHEH(red_hh);
517 
518  // split one edge
519  eh = Base::mesh_.EH(red_hh);
520 
521  split_edge(hh[2], new_vh[2], _target_state);
522 
523  assert(Base::mesh_.FVH(hh[2]) == vh[1]);
524  assert(Base::mesh_.FVH(hh[1]) == vh[0]);
525  assert(Base::mesh_.FVH(hh[0]) == vh[2]);
526 
527  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[2])).is_valid()) {
528 
529  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[2])));
530  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
531  halfedge_vector.push_back(temp_hh);
532  }
533  }
534  }
535 
536  else {
537 
538  // one time divided face
539  vh[0] = Base::mesh_.TVH(Base::mesh_.NHEH(Base::mesh_.OHEH(red_hh)));
540  vh[1] = Base::mesh_.TVH(red_hh);
541  vh[2] = Base::mesh_.TVH(Base::mesh_.NHEH(red_hh));
542 
543  new_vh[0] = Base::mesh_.FVH(red_hh);
544  new_vh[1] = Base::mesh_.add_vertex(zero_point);
545  new_vh[2] = Base::mesh_.add_vertex(zero_point);
546 
547  hh[0] = Base::mesh_.PHEH(red_hh);
548  hh[1] = Base::mesh_.PHEH(Base::mesh_.OHEH(red_hh));
549  hh[2] = Base::mesh_.NHEH(red_hh);
550 
551  // split two edges
552  eh = Base::mesh_.EH(red_hh);
553 
554  split_edge(hh[1], new_vh[1], _target_state);
555  split_edge(hh[2], new_vh[2], _target_state);
556 
557  assert(Base::mesh_.FVH(hh[2]) == vh[1]);
558  assert(Base::mesh_.FVH(hh[1]) == vh[0]);
559  assert(Base::mesh_.FVH(hh[0]) == vh[2]);
560 
561  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[1])).is_valid()) {
562 
563  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[1])));
564  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
565  halfedge_vector.push_back(temp_hh);
566  }
567 
568  if (Base::mesh_.FH(Base::mesh_.OHEH(hh[2])).is_valid()) {
569 
570  temp_hh = Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(hh[2])));
571  if (MOBJ(Base::mesh_.FH(temp_hh)).red_halfedge() != temp_hh)
572  halfedge_vector.push_back(temp_hh);
573  }
574  }
575  }
576  }
577 
578  // continue here for all cases
579 
580  // flip edge
581  if (Base::mesh_.is_flip_ok(eh)) {
582 
583  Base::mesh_.flip(eh);
584  }
585 
586  // search new faces
587  fh[0] = Base::mesh_.FH(hh[0]);
588  fh[1] = Base::mesh_.FH(hh[1]);
589  fh[2] = Base::mesh_.FH(hh[2]);
590  fh[3] = Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(hh[0])));
591 
592  assert(_fh == fh[0] || _fh == fh[1] || _fh == fh[2] || _fh == fh[3]);
593 
594  // init new faces
595  for (int i = 0; i <= 3; ++i) {
596 
597  MOBJ(fh[i]).set_state(_target_state);
598  MOBJ(fh[i]).set_final();
599  MOBJ(fh[i]).set_position(_target_state, face_position);
600  MOBJ(fh[i]).set_red_halfedge(Base::mesh_.InvalidHalfedgeHandle);
601  }
602 
603  // init new vertices and edges
604  for (int i = 0; i <= 2; ++i) {
605 
606  MOBJ(new_vh[i]).set_position(_target_state, zero_point);
607  MOBJ(new_vh[i]).set_state(_target_state);
608  MOBJ(new_vh[i]).set_not_final();
609 
610  Base::mesh_.set_point(new_vh[i], (Base::mesh_.point(vh[i]) + Base::mesh_.point(vh[(i + 2) % 3])) * 0.5);
611 
612  MOBJ(Base::mesh_.EH(hh[i])).set_state(_target_state);
613  MOBJ(Base::mesh_.EH(hh[i])).set_position(_target_state, zero_point);
614  MOBJ(Base::mesh_.EH(hh[i])).set_final();
615 
616  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh[i]))).set_state(_target_state);
617  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh[i]))).set_position(_target_state, zero_point);
618  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh[i]))).set_final();
619 
620  MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh[i]))).set_state(_target_state);
621  MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh[i]))).set_position(_target_state, zero_point);
622  MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh[i]))).set_final();
623  }
624 
625  // check, if opposite triangle needs splitting
626  while (!halfedge_vector.empty()) {
627 
628  temp_hh = halfedge_vector.back();
629  halfedge_vector.pop_back();
630 
631  check_edge(temp_hh, _target_state);
632  }
633 
634  assert(MOBJ(fh[0]).state() == _target_state);
635  assert(MOBJ(fh[1]).state() == _target_state);
636  assert(MOBJ(fh[2]).state() == _target_state);
637  assert(MOBJ(fh[3]).state() == _target_state);
638  }
639 }
640 
641 
642 template<class M>
643 void
644 Tvv4<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
645 {
646 
647  if (MOBJ(_vh).state() < _target_state)
648  {
649 
650  this->update(_vh, _target_state);
651 
652  // multiply old position by 4
653  MOBJ(_vh).set_position(_target_state, MOBJ(_vh).position(_target_state - 1) * 4.0);
654 
655  MOBJ(_vh).inc_state();
656  }
657 }
658 
659 
660 template<class M>
661 void
662 Tvv4<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state)
663 {
664  if (MOBJ(_eh).state() < _target_state)
665  {
666  this->update(_eh, _target_state);
667 
668  typename M::FaceHandle fh(Base::mesh_.FH(Base::mesh_.HEH(_eh, 0)));
669 
670  if (!fh.is_valid())
671  fh=Base::mesh_.FH(Base::mesh_.HEH(_eh, 1));
672 
673  raise(fh, _target_state);
674 
675  assert(MOBJ(_eh).state() == _target_state);
676  }
677 }
678 
679 #ifndef DOXY_IGNORE_THIS
680 
681 template<class M>
682 void
683 Tvv4<M>::split_edge(typename M::HalfedgeHandle &_hh,
684  typename M::VertexHandle &_vh,
685  state_t _target_state)
686  {
687  typename M::HalfedgeHandle temp_hh;
688 
689  if (Base::mesh_.FH(Base::mesh_.OHEH(_hh)).is_valid())
690  {
691  if (!MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).final())
692  {
693  if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).red_halfedge().is_valid())
694  {
695  temp_hh = MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).red_halfedge();
696  }
697  else
698  {
699  // two cases for divided, but not visited face
700  if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(_hh))))).state()
701  == MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).state())
702  {
703  temp_hh = Base::mesh_.PHEH(Base::mesh_.OHEH(_hh));
704  }
705 
706  else if (MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(_hh))))).state()
707  == MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).state())
708  {
709  temp_hh = Base::mesh_.NHEH(Base::mesh_.OHEH(_hh));
710  }
711  }
712  }
713  else
714  temp_hh = Base::mesh_.InvalidHalfedgeHandle;
715  }
716  else
717  temp_hh = Base::mesh_.InvalidHalfedgeHandle;
718 
719  // split edge
720  Base::mesh_.split(Base::mesh_.EH(_hh), _vh);
721 
722  if (Base::mesh_.FVH(_hh) == _vh)
723  {
724  MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(_hh))))).set_state(MOBJ(Base::mesh_.EH(_hh)).state());
725  _hh = Base::mesh_.PHEH(Base::mesh_.OHEH(Base::mesh_.PHEH(_hh)));
726  }
727 
728  if (Base::mesh_.FH(Base::mesh_.OHEH(_hh)).is_valid()) {
729 
730  MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(Base::mesh_.OHEH(_hh)))).set_not_final();
731  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).set_state(_target_state-1);
732  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(_hh)))))).set_state(_target_state-1);
733 
734  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).set_not_final();
735  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(_hh)))))).set_not_final();
736 
737  MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(Base::mesh_.OHEH(_hh)))).set_state(_target_state);
738 
739  if (temp_hh.is_valid()) {
740 
741  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(_hh))).set_red_halfedge(temp_hh);
742  MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(_hh)))))).set_red_halfedge(temp_hh);
743  }
744  else {
745 
746  typename M::FaceHandle
747  fh1(Base::mesh_.FH(Base::mesh_.OHEH(_hh))),
748  fh2(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(_hh))))));
749 
750  MOBJ(fh1).set_red_halfedge(Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(_hh))));
751  MOBJ(fh2).set_red_halfedge(Base::mesh_.OHEH(Base::mesh_.PHEH(Base::mesh_.OHEH(_hh))));
752 
753  const typename M::Point zero_point(0.0, 0.0, 0.0);
754 
755  MOBJ(fh1).set_position(_target_state - 1, zero_point);
756  MOBJ(fh2).set_position(_target_state - 1, zero_point);
757  }
758  }
759 
760  // init edges
761  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(_hh))))).set_state(_target_state - 1);
762  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(Base::mesh_.OHEH(Base::mesh_.NHEH(_hh))))).set_final();
763 
764  MOBJ(Base::mesh_.EH(_hh)).set_state(_target_state - 1);
765  MOBJ(Base::mesh_.EH(_hh)).set_final();
766 }
767 
768 
769 template<class M>
770 void Tvv4<M>::check_edge(const typename M::HalfedgeHandle& _hh,
771  state_t _target_state)
772 {
773  typename M::FaceHandle fh1(Base::mesh_.FH(_hh)),
774  fh2(Base::mesh_.FH(Base::mesh_.OHEH(_hh)));
775 
776  assert(fh1.is_valid());
777  assert(fh2.is_valid());
778 
779  typename M::HalfedgeHandle red_hh(MOBJ(fh1).red_halfedge());
780 
781  if (!MOBJ(fh1).final()) {
782 
783  assert (MOBJ(fh1).final() == MOBJ(fh2).final());
784  assert (!MOBJ(fh1).final());
785  assert (MOBJ(fh1).red_halfedge() == MOBJ(fh2).red_halfedge());
786 
787  const typename M::Point zero_point(0.0, 0.0, 0.0);
788 
789  MOBJ(fh1).set_position(_target_state - 1, zero_point);
790  MOBJ(fh2).set_position(_target_state - 1, zero_point);
791 
792  assert(red_hh.is_valid());
793 
794  if (!red_hh.is_valid()) {
795 
796  MOBJ(fh1).set_state(_target_state - 1);
797  MOBJ(fh2).set_state(_target_state - 1);
798 
799  MOBJ(fh1).set_red_halfedge(_hh);
800  MOBJ(fh2).set_red_halfedge(_hh);
801 
802  MOBJ(Base::mesh_.EH(_hh)).set_not_final();
803  MOBJ(Base::mesh_.EH(_hh)).set_state(_target_state - 1);
804  }
805 
806  else {
807 
808  MOBJ(Base::mesh_.EH(_hh)).set_not_final();
809  MOBJ(Base::mesh_.EH(_hh)).set_state(_target_state - 1);
810 
811  raise(fh1, _target_state);
812 
813  assert(MOBJ(fh1).state() == _target_state);
814  }
815  }
816 }
817 
818 
819 // -------------------------------------------------------------------- VF ----
820 
821 
822 template<class M>
823 void VF<M>::raise(typename M::FaceHandle& _fh, state_t _target_state)
824 {
825  if (MOBJ(_fh).state() < _target_state) {
826 
827  this->update(_fh, _target_state);
828 
829  // raise all neighbor vertices to level x-1
830  typename M::FaceVertexIter fv_it;
831  typename M::VertexHandle vh;
832  std::vector<typename M::VertexHandle> vertex_vector;
833 
834  if (_target_state > 1) {
835 
836  for (fv_it = Base::mesh_.fv_iter(_fh); fv_it.is_valid(); ++fv_it) {
837 
838  vertex_vector.push_back(*fv_it);
839  }
840 
841  while (!vertex_vector.empty()) {
842 
843  vh = vertex_vector.back();
844  vertex_vector.pop_back();
845 
846  Base::prev_rule()->raise(vh, _target_state - 1);
847  }
848  }
849 
850  // calculate new position
851  typename M::Point position(0.0, 0.0, 0.0);
852  typename M::Scalar valence(0.0);
853 
854  for (fv_it = Base::mesh_.fv_iter(_fh); fv_it.is_valid(); ++fv_it) {
855 
856  valence += 1.0;
857  position += Base::mesh_.data(*fv_it).position(_target_state - 1);
858  }
859 
860  position /= valence;
861 
862  // boundary rule
863  if (Base::number() == Base::subdiv_rule()->number() + 1 &&
864  Base::mesh_.is_boundary(_fh) &&
865  !MOBJ(_fh).final())
866  position *= 0.5;
867 
868  MOBJ(_fh).set_position(_target_state, position);
869  MOBJ(_fh).inc_state();
870 
871  assert(_target_state == MOBJ(_fh).state());
872  }
873 }
874 
875 
876 // -------------------------------------------------------------------- FF ----
877 
878 
879 template<class M>
880 void FF<M>::raise(typename M::FaceHandle& _fh, state_t _target_state) {
881 
882  if (MOBJ(_fh).state() < _target_state) {
883 
884  this->update(_fh, _target_state);
885 
886  // raise all neighbor faces to level x-1
887  typename M::FaceFaceIter ff_it;
888  typename M::FaceHandle fh;
889  std::vector<typename M::FaceHandle> face_vector;
890 
891  if (_target_state > 1) {
892 
893  for (ff_it = Base::mesh_.ff_iter(_fh); ff_it.is_valid(); ++ff_it) {
894 
895  face_vector.push_back(*ff_it);
896  }
897 
898  while (!face_vector.empty()) {
899 
900  fh = face_vector.back();
901  face_vector.pop_back();
902 
903  Base::prev_rule()->raise(fh, _target_state - 1);
904  }
905 
906  for (ff_it = Base::mesh_.ff_iter(_fh); ff_it.is_valid(); ++ff_it) {
907 
908  face_vector.push_back(*ff_it);
909  }
910 
911  while (!face_vector.empty()) {
912 
913  fh = face_vector.back();
914  face_vector.pop_back();
915 
916  while (MOBJ(fh).state() < _target_state - 1)
917  Base::prev_rule()->raise(fh, _target_state - 1);
918  }
919  }
920 
921  // calculate new position
922  typename M::Point position(0.0, 0.0, 0.0);
923  typename M::Scalar valence(0.0);
924 
925  for (ff_it = Base::mesh_.ff_iter(_fh); ff_it.is_valid(); ++ff_it) {
926 
927  valence += 1.0;
928 
929  position += Base::mesh_.data(*ff_it).position(_target_state - 1);
930  }
931 
932  position /= valence;
933 
934  MOBJ(_fh).set_position(_target_state, position);
935  MOBJ(_fh).inc_state();
936  }
937 }
938 
939 
940 // ------------------------------------------------------------------- FFc ----
941 
942 
943 template<class M>
944 void FFc<M>::raise(typename M::FaceHandle& _fh, state_t _target_state)
945 {
946  if (MOBJ(_fh).state() < _target_state) {
947 
948  this->update(_fh, _target_state);
949 
950  // raise all neighbor faces to level x-1
951  typename M::FaceFaceIter ff_it(Base::mesh_.ff_iter(_fh));
952  typename M::FaceHandle fh;
953  std::vector<typename M::FaceHandle> face_vector;
954 
955  if (_target_state > 1)
956  {
957  for (; ff_it.is_valid(); ++ff_it)
958  face_vector.push_back(*ff_it);
959 
960  while (!face_vector.empty())
961  {
962  fh = face_vector.back();
963  face_vector.pop_back();
964  Base::prev_rule()->raise(fh, _target_state - 1);
965  }
966 
967  for (ff_it = Base::mesh_.ff_iter(_fh); ff_it.is_valid(); ++ff_it)
968  face_vector.push_back(*ff_it);
969 
970  while (!face_vector.empty()) {
971 
972  fh = face_vector.back();
973  face_vector.pop_back();
974 
975  while (MOBJ(fh).state() < _target_state - 1)
976  Base::prev_rule()->raise(fh, _target_state - 1);
977  }
978  }
979 
980  // calculate new position
981  typename M::Point position(0.0, 0.0, 0.0);
982  typename M::Scalar valence(0.0);
983 
984  for (ff_it = Base::mesh_.ff_iter(_fh); ff_it.is_valid(); ++ff_it)
985  {
986  valence += 1.0;
987  position += Base::mesh_.data(*ff_it).position(_target_state - 1);
988  }
989 
990  position /= valence;
991 
992  // choose coefficient c
993  typename M::Scalar c = Base::coeff();
994 
995  position *= (1.0 - c);
996  position += MOBJ(_fh).position(_target_state - 1) * c;
997 
998  MOBJ(_fh).set_position(_target_state, position);
999  MOBJ(_fh).inc_state();
1000  }
1001 }
1002 
1003 
1004 // -------------------------------------------------------------------- FV ----
1005 
1006 
1007 template<class M>
1008 void FV<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
1009 {
1010 
1011  if (MOBJ(_vh).state() < _target_state) {
1012 
1013  this->update(_vh, _target_state);
1014 
1015  // raise all neighbor vertices to level x-1
1016  typename M::VertexFaceIter vf_it(Base::mesh_.vf_iter(_vh));
1017  typename M::FaceHandle fh;
1018  std::vector<typename M::FaceHandle> face_vector;
1019 
1020  if (_target_state > 1) {
1021 
1022  for (; vf_it.is_valid(); ++vf_it) {
1023 
1024  face_vector.push_back(*vf_it);
1025  }
1026 
1027  while (!face_vector.empty()) {
1028 
1029  fh = face_vector.back();
1030  face_vector.pop_back();
1031 
1032  Base::prev_rule()->raise(fh, _target_state - 1);
1033  }
1034 
1035  for (vf_it = Base::mesh_.vf_iter(_vh); vf_it.is_valid(); ++vf_it) {
1036 
1037  face_vector.push_back(*vf_it);
1038  }
1039 
1040  while (!face_vector.empty()) {
1041 
1042  fh = face_vector.back();
1043  face_vector.pop_back();
1044 
1045  while (MOBJ(fh).state() < _target_state - 1)
1046  Base::prev_rule()->raise(fh, _target_state - 1);
1047  }
1048  }
1049 
1050  // calculate new position
1051  typename M::Point position(0.0, 0.0, 0.0);
1052  typename M::Scalar valence(0.0);
1053 
1054  for (vf_it = Base::mesh_.vf_iter(_vh); vf_it.is_valid(); ++vf_it) {
1055 
1056  valence += 1.0;
1057  position += Base::mesh_.data(*vf_it).position(_target_state - 1);
1058  }
1059 
1060  position /= valence;
1061 
1062  MOBJ(_vh).set_position(_target_state, position);
1063  MOBJ(_vh).inc_state();
1064 
1065  if (Base::number() == Base::n_rules() - 1) {
1066 
1067  Base::mesh_.set_point(_vh, position);
1068  MOBJ(_vh).set_final();
1069  }
1070  }
1071 }
1072 
1073 
1074 // ------------------------------------------------------------------- FVc ----
1075 
1076 
1077 template<class M>
1078 void FVc<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
1079 {
1080  if (MOBJ(_vh).state() < _target_state) {
1081 
1082  this->update(_vh, _target_state);
1083 
1084  typename M::VertexOHalfedgeIter voh_it;
1085  typename M::FaceHandle fh;
1086  std::vector<typename M::FaceHandle> face_vector;
1087  int valence(0);
1088 
1089  face_vector.clear();
1090 
1091  // raise all neighbour faces to level x-1
1092  if (_target_state > 1) {
1093 
1094  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it) {
1095 
1096  if (Base::mesh_.FH(*voh_it).is_valid()) {
1097 
1098  face_vector.push_back(Base::mesh_.FH(*voh_it));
1099 
1100  if (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))).is_valid()) {
1101 
1102  face_vector.push_back(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))));
1103  }
1104  }
1105  }
1106 
1107  while (!face_vector.empty()) {
1108 
1109  fh = face_vector.back();
1110  face_vector.pop_back();
1111 
1112  Base::prev_rule()->raise(fh, _target_state - 1);
1113  }
1114 
1115  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it) {
1116 
1117  if (Base::mesh_.FH(*voh_it).is_valid()) {
1118 
1119  face_vector.push_back(Base::mesh_.FH(*voh_it));
1120 
1121  if (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))).is_valid()) {
1122 
1123  face_vector.push_back(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))));
1124  }
1125  }
1126  }
1127 
1128  while (!face_vector.empty()) {
1129 
1130  fh = face_vector.back();
1131  face_vector.pop_back();
1132 
1133  while (MOBJ(fh).state() < _target_state - 1)
1134  Base::prev_rule()->raise(fh, _target_state - 1);
1135  }
1136 
1137  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it) {
1138 
1139  if (Base::mesh_.FH(*voh_it).is_valid()) {
1140 
1141  face_vector.push_back(Base::mesh_.FH(*voh_it));
1142 
1143  if (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))).is_valid()) {
1144 
1145  face_vector.push_back(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))));
1146  }
1147  }
1148  }
1149 
1150  while (!face_vector.empty()) {
1151 
1152  fh = face_vector.back();
1153  face_vector.pop_back();
1154 
1155  while (MOBJ(fh).state() < _target_state - 1)
1156  Base::prev_rule()->raise(fh, _target_state - 1);
1157  }
1158  }
1159 
1160  // calculate new position
1161  typename M::Point position(0.0, 0.0, 0.0);
1162  typename M::Scalar c;
1163 #if 0
1164  const typename M::Scalar _2pi(2.0*M_PI);
1165  const typename M::Scalar _2over3(2.0/3.0);
1166 
1167  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it; ++voh_it)
1168  {
1169  ++valence;
1170  }
1171 
1172  // choose coefficient c
1173  c = _2over3 * ( cos( _2pi / valence) + 1.0);
1174 #else
1175  valence = Base::mesh_.valence(_vh);
1176  c = typename M::Scalar(coeff(valence));
1177 #endif
1178 
1179 
1180  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it) {
1181 
1182  fh = Base::mesh_.FH(*voh_it);
1183  if (fh.is_valid())
1184  Base::prev_rule()->raise(fh, _target_state - 1);
1185 
1186  fh = Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it)));
1187  if (fh.is_valid())
1188  Base::prev_rule()->raise(fh, _target_state - 1);
1189 
1190  if (Base::mesh_.FH(*voh_it).is_valid()) {
1191 
1192  if (Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it))).is_valid()) {
1193 
1194  position += MOBJ(Base::mesh_.FH(*voh_it)).position(_target_state - 1) * c;
1195 
1196  position += MOBJ(Base::mesh_.FH(Base::mesh_.OHEH(Base::mesh_.NHEH(*voh_it)))).position(_target_state - 1) * ( typename M::Scalar(1.0) - c);
1197  }
1198  else {
1199 
1200  position += MOBJ(Base::mesh_.FH(*voh_it)).position(_target_state - 1);
1201  }
1202  }
1203 
1204  else {
1205 
1206  --valence;
1207  }
1208  }
1209 
1210  position /= typename M::Scalar(valence);
1211 
1212  MOBJ(_vh).set_position(_target_state, position);
1213  MOBJ(_vh).inc_state();
1214 
1215  assert(MOBJ(_vh).state() == _target_state);
1216 
1217  // check if last rule
1218  if (Base::number() == Base::n_rules() - 1) {
1219 
1220  Base::mesh_.set_point(_vh, position);
1221  MOBJ(_vh).set_final();
1222  }
1223  }
1224 }
1225 
1226 template<class M>
1227 std::vector<double> FVc<M>::coeffs_;
1228 
1229 template <class M>
1230 void FVc<M>::init_coeffs(size_t _max_valence)
1231 {
1232  if ( coeffs_.size() == _max_valence+1 )
1233  return;
1234 
1235  if ( coeffs_.size() < _max_valence+1 )
1236  {
1237  const double _2pi(2.0*M_PI);
1238  const double _2over3(2.0/3.0);
1239 
1240  if (coeffs_.empty())
1241  coeffs_.push_back(0.0); // dummy for valence 0
1242 
1243  for(size_t v=coeffs_.size(); v <= _max_valence; ++v)
1244  coeffs_.push_back(_2over3 * ( cos( _2pi / v) + 1.0));
1245  }
1246 }
1247 
1248 
1249 // -------------------------------------------------------------------- VV ----
1250 
1251 
1252 template<class M>
1253 void VV<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
1254 {
1255  if (MOBJ(_vh).state() < _target_state)
1256  {
1257  this->update(_vh, _target_state);
1258 
1259  // raise all neighbor vertices to level x-1
1260  typename M::VertexVertexIter vv_it(Base::mesh_.vv_iter(_vh));
1261  typename M::VertexHandle vh;
1262  std::vector<typename M::VertexHandle> vertex_vector;
1263 
1264  if (_target_state > 1) {
1265 
1266  for (; vv_it.is_valid(); ++vv_it) {
1267 
1268  vertex_vector.push_back(*vv_it);
1269  }
1270 
1271  while (!vertex_vector.empty()) {
1272 
1273  vh = vertex_vector.back();
1274  vertex_vector.pop_back();
1275 
1276  Base::prev_rule()->raise(vh, _target_state - 1);
1277  }
1278 
1279  for (; vv_it.is_valid(); ++vv_it) {
1280 
1281  vertex_vector.push_back(*vv_it);
1282  }
1283 
1284  while (!vertex_vector.empty()) {
1285 
1286  vh = vertex_vector.back();
1287  vertex_vector.pop_back();
1288 
1289  Base::prev_rule()->raise(vh, _target_state - 1);
1290  }
1291  }
1292 
1293  // calculate new position
1294  typename M::Point position(0.0, 0.0, 0.0);
1295  typename M::Scalar valence(0.0);
1296 
1297  for (vv_it = Base::mesh_.vv_iter(_vh); vv_it.is_valid(); ++vv_it) {
1298 
1299  valence += 1.0;
1300  position += Base::mesh_.data(*vv_it).position(_target_state - 1);
1301  }
1302 
1303  position /= valence;
1304 
1305  MOBJ(_vh).set_position(_target_state, position);
1306  MOBJ(_vh).inc_state();
1307 
1308  // check if last rule
1309  if (Base::number() == Base::n_rules() - 1) {
1310 
1311  Base::mesh_.set_point(_vh, position);
1312  MOBJ(_vh).set_final();
1313  }
1314  }
1315 }
1316 
1317 
1318 // ------------------------------------------------------------------- VVc ----
1319 
1320 
1321 template<class M>
1322 void VVc<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
1323 {
1324  if (MOBJ(_vh).state() < _target_state) {
1325 
1326  this->update(_vh, _target_state);
1327 
1328  // raise all neighbor vertices to level x-1
1329  typename M::VertexVertexIter vv_it(Base::mesh_.vv_iter(_vh));
1330  typename M::VertexHandle vh;
1331  std::vector<typename M::VertexHandle> vertex_vector;
1332 
1333  if (_target_state > 1) {
1334 
1335  for (; vv_it.is_valid(); ++vv_it) {
1336 
1337  vertex_vector.push_back(*vv_it);
1338  }
1339 
1340  while (!vertex_vector.empty()) {
1341 
1342  vh = vertex_vector.back();
1343  vertex_vector.pop_back();
1344 
1345  Base::prev_rule()->raise(vh, _target_state - 1);
1346  }
1347 
1348  for (; vv_it.is_valid(); ++vv_it) {
1349 
1350  vertex_vector.push_back(*vv_it);
1351  }
1352 
1353  while (!vertex_vector.empty()) {
1354 
1355  vh = vertex_vector.back();
1356  vertex_vector.pop_back();
1357 
1358  Base::prev_rule()->raise(vh, _target_state - 1);
1359  }
1360  }
1361 
1362  // calculate new position
1363  typename M::Point position(0.0, 0.0, 0.0);
1364  typename M::Scalar valence(0.0);
1365  typename M::Scalar c;
1366 
1367  for (vv_it = Base::mesh_.vv_iter(_vh); vv_it.is_valid(); ++vv_it)
1368  {
1369  valence += 1.0;
1370  position += Base::mesh_.data(*vv_it).position(_target_state - 1);
1371  }
1372 
1373  position /= valence;
1374 
1375  // choose coefficient c
1376  c = Base::coeff();
1377 
1378  position *= (1.0 - c);
1379  position += MOBJ(_vh).position(_target_state - 1) * c;
1380 
1381  MOBJ(_vh).set_position(_target_state, position);
1382  MOBJ(_vh).inc_state();
1383 
1384  if (Base::number() == Base::n_rules() - 1)
1385  {
1386  Base::mesh_.set_point(_vh, position);
1387  MOBJ(_vh).set_final();
1388  }
1389  }
1390 }
1391 
1392 
1393 // -------------------------------------------------------------------- VE ----
1394 
1395 
1396 template<class M>
1397 void VE<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state)
1398 {
1399  if (MOBJ(_eh).state() < _target_state) {
1400 
1401  this->update(_eh, _target_state);
1402 
1403  // raise all neighbour vertices to level x-1
1404  typename M::VertexHandle vh;
1405  typename M::HalfedgeHandle hh1(Base::mesh_.HEH(_eh, 0)),
1406  hh2(Base::mesh_.HEH(_eh, 1));
1407 
1408  if (_target_state > 1) {
1409 
1410  vh = Base::mesh_.TVH(hh1);
1411 
1412  Base::prev_rule()->raise(vh, _target_state - 1);
1413 
1414  vh = Base::mesh_.TVH(hh2);
1415 
1416  Base::prev_rule()->raise(vh, _target_state - 1);
1417  }
1418 
1419  // calculate new position
1420  typename M::Point position(0.0, 0.0, 0.0);
1421  const typename M::Scalar valence(2.0);
1422 
1423  position += MOBJ(Base::mesh_.TVH(hh1)).position(_target_state - 1);
1424  position += MOBJ(Base::mesh_.TVH(hh2)).position(_target_state - 1);
1425 
1426  position /= valence;
1427 
1428  MOBJ(_eh).set_position(_target_state, position);
1429  MOBJ(_eh).inc_state();
1430  }
1431 }
1432 
1433 
1434 // ------------------------------------------------------------------- VdE ----
1435 
1436 
1437 template<class M>
1438 void VdE<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state)
1439 {
1440  if (MOBJ(_eh).state() < _target_state)
1441  {
1442  this->update(_eh, _target_state);
1443 
1444  // raise all neighbor vertices to level x-1
1445  typename M::VertexHandle vh;
1446  typename M::HalfedgeHandle hh1(Base::mesh_.HEH(_eh, 0)),
1447  hh2(Base::mesh_.HEH(_eh, 1));
1448  typename M::FaceHandle fh1, fh2;
1449 
1450  if (_target_state > 1) {
1451 
1452  fh1 = Base::mesh_.FH(hh1);
1453  fh2 = Base::mesh_.FH(hh2);
1454 
1455  if (fh1.is_valid()) {
1456 
1457  Base::prev_rule()->raise(fh1, _target_state - 1);
1458 
1459  vh = Base::mesh_.TVH(Base::mesh_.NHEH(hh1));
1460  Base::prev_rule()->raise(vh, _target_state - 1);
1461  }
1462 
1463  if (fh2.is_valid()) {
1464 
1465  Base::prev_rule()->raise(fh2, _target_state - 1);
1466 
1467  vh = Base::mesh_.TVH(Base::mesh_.NHEH(hh2));
1468  Base::prev_rule()->raise(vh, _target_state - 1);
1469  }
1470 
1471  vh = Base::mesh_.TVH(hh1);
1472  Base::prev_rule()->raise(vh, _target_state - 1);
1473 
1474  vh = Base::mesh_.TVH(hh2);
1475  Base::prev_rule()->raise(vh, _target_state - 1);
1476  }
1477 
1478  // calculate new position
1479  typename M::Point position(0.0, 0.0, 0.0);
1480  typename M::Scalar valence(2.0);
1481 
1482  position += MOBJ(Base::mesh_.TVH(hh1)).position(_target_state - 1);
1483  position += MOBJ(Base::mesh_.TVH(hh2)).position(_target_state - 1);
1484 
1485  if (fh1.is_valid()) {
1486 
1487  position += MOBJ(Base::mesh_.TVH(Base::mesh_.NHEH(hh1))).position(_target_state - 1);
1488  valence += 1.0;
1489  }
1490 
1491  if (fh2.is_valid()) {
1492 
1493  position += MOBJ(Base::mesh_.TVH(Base::mesh_.NHEH(hh2))).position(_target_state - 1);
1494  valence += 1.0;
1495  }
1496 
1497  if (Base::number() == Base::subdiv_rule()->Base::number() + 1)
1498  valence = 4.0;
1499 
1500  position /= valence;
1501 
1502  MOBJ(_eh).set_position(_target_state, position);
1503  MOBJ(_eh).inc_state();
1504  }
1505 }
1506 
1507 
1508 // ------------------------------------------------------------------ VdEc ----
1509 
1510 
1511 template<class M>
1512 void
1513 VdEc<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state)
1514 {
1515  if (MOBJ(_eh).state() < _target_state)
1516  {
1517  this->update(_eh, _target_state);
1518 
1519  // raise all neighbor vertices to level x-1
1520  typename M::VertexHandle vh;
1521  typename M::HalfedgeHandle hh1(Base::mesh_.HEH(_eh, 0)),
1522  hh2(Base::mesh_.HEH(_eh, 1));
1523  std::vector<typename M::VertexHandle> vertex_vector;
1524  typename M::FaceHandle fh1, fh2;
1525 
1526  if (_target_state > 1) {
1527 
1528  fh1 = Base::mesh_.FH(Base::mesh_.HEH(_eh, 0));
1529  fh2 = Base::mesh_.FH(Base::mesh_.HEH(_eh, 1));
1530 
1531  Base::prev_rule()->raise(fh1, _target_state - 1);
1532  Base::prev_rule()->raise(fh2, _target_state - 1);
1533 
1534  vertex_vector.push_back(Base::mesh_.TVH(hh1));
1535  vertex_vector.push_back(Base::mesh_.TVH(hh2));
1536 
1537  vertex_vector.push_back(Base::mesh_.TVH(Base::mesh_.NHEH(hh1)));
1538  vertex_vector.push_back(Base::mesh_.TVH(Base::mesh_.NHEH(hh2)));
1539 
1540  while (!vertex_vector.empty()) {
1541 
1542  vh = vertex_vector.back();
1543  vertex_vector.pop_back();
1544 
1545  Base::prev_rule()->raise(vh, _target_state - 1);
1546  }
1547 
1548  vertex_vector.push_back(Base::mesh_.TVH(hh1));
1549  vertex_vector.push_back(Base::mesh_.TVH(hh2));
1550 
1551  vertex_vector.push_back(Base::mesh_.TVH(Base::mesh_.NHEH(hh1)));
1552  vertex_vector.push_back(Base::mesh_.TVH(Base::mesh_.NHEH(hh2)));
1553 
1554  while (!vertex_vector.empty()) {
1555 
1556  vh = vertex_vector.back();
1557  vertex_vector.pop_back();
1558 
1559  Base::prev_rule()->raise(vh, _target_state - 1);
1560  }
1561  }
1562 
1563  // calculate new position
1564  typename M::Point position(0.0, 0.0, 0.0);
1565  const typename M::Scalar valence(4.0);
1566  typename M::Scalar c;
1567 
1568  // choose coefficient c
1569  c = Base::coeff();
1570 
1571  position += MOBJ(Base::mesh_.TVH(hh1)).position(_target_state - 1) * c;
1572  position += MOBJ(Base::mesh_.TVH(hh2)).position(_target_state - 1) * c;
1573  position += MOBJ(Base::mesh_.TVH(Base::mesh_.NHEH(hh1))).position(_target_state - 1) * (0.5 - c);
1574  position += MOBJ(Base::mesh_.TVH(Base::mesh_.NHEH(hh2))).position(_target_state - 1) * (0.5 - c);
1575 
1576  position /= valence;
1577 
1578  MOBJ(_eh).set_position(_target_state, position);
1579  MOBJ(_eh).inc_state();
1580  }
1581 }
1582 
1583 
1584 // -------------------------------------------------------------------- EV ----
1585 
1586 
1587 template<class M>
1588 void EV<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
1589 {
1590  if (MOBJ(_vh).state() < _target_state) {
1591 
1592  this->update(_vh, _target_state);
1593 
1594  // raise all neighbor vertices to level x-1
1595  typename M::VertexEdgeIter ve_it(Base::mesh_.ve_iter(_vh));
1596  typename M::EdgeHandle eh;
1597  std::vector<typename M::EdgeHandle> edge_vector;
1598 
1599  if (_target_state > 1) {
1600 
1601  for (; ve_it.is_valid(); ++ve_it) {
1602 
1603  edge_vector.push_back(*ve_it);
1604  }
1605 
1606  while (!edge_vector.empty()) {
1607 
1608  eh = edge_vector.back();
1609  edge_vector.pop_back();
1610 
1611  Base::prev_rule()->raise(eh, _target_state - 1);
1612  }
1613 
1614  for (ve_it = Base::mesh_.ve_iter(_vh); ve_it.is_valid(); ++ve_it) {
1615 
1616  edge_vector.push_back(*ve_it);
1617  }
1618 
1619  while (!edge_vector.empty()) {
1620 
1621  eh = edge_vector.back();
1622  edge_vector.pop_back();
1623 
1624  while (MOBJ(eh).state() < _target_state - 1)
1625  Base::prev_rule()->raise(eh, _target_state - 1);
1626  }
1627  }
1628 
1629  // calculate new position
1630  typename M::Point position(0.0, 0.0, 0.0);
1631  typename M::Scalar valence(0.0);
1632 
1633  for (ve_it = Base::mesh_.ve_iter(_vh); ve_it.is_valid(); ++ve_it) {
1634 
1635  if (Base::mesh_.data(*ve_it).final()) {
1636 
1637  valence += 1.0;
1638 
1639  position += Base::mesh_.data(*ve_it).position(_target_state - 1);
1640  }
1641  }
1642 
1643  position /= valence;
1644 
1645  MOBJ(_vh).set_position(_target_state, position);
1646  MOBJ(_vh).inc_state();
1647 
1648  // check if last rule
1649  if (Base::number() == Base::n_rules() - 1) {
1650 
1651  Base::mesh_.set_point(_vh, position);
1652  MOBJ(_vh).set_final();
1653  }
1654  }
1655 }
1656 
1657 
1658 // ------------------------------------------------------------------- EVc ----
1659 
1660 template<class M>
1661 std::vector<double> EVc<M>::coeffs_;
1662 
1663 template<class M>
1664 void EVc<M>::raise(typename M::VertexHandle& _vh, state_t _target_state)
1665 {
1666  if (MOBJ(_vh).state() < _target_state)
1667  {
1668  this->update(_vh, _target_state);
1669 
1670  // raise all neighbour vertices to level x-1
1671  typename M::VertexOHalfedgeIter voh_it;
1672  typename M::EdgeHandle eh;
1673  typename M::FaceHandle fh;
1674  std::vector<typename M::EdgeHandle> edge_vector;
1675  std::vector<typename M::FaceHandle> face_vector;
1676 
1677  if (_target_state > 1) {
1678 
1679  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it) {
1680 
1681  face_vector.push_back(Base::mesh_.FH(*voh_it));
1682  }
1683 
1684  while (!face_vector.empty()) {
1685 
1686  fh = face_vector.back();
1687  face_vector.pop_back();
1688 
1689  if (fh.is_valid())
1690  Base::prev_rule()->raise(fh, _target_state - 1);
1691  }
1692 
1693  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it) {
1694 
1695  edge_vector.push_back(Base::mesh_.EH(*voh_it));
1696 
1697  edge_vector.push_back(Base::mesh_.EH(Base::mesh_.NHEH(*voh_it)));
1698  }
1699 
1700  while (!edge_vector.empty()) {
1701 
1702  eh = edge_vector.back();
1703  edge_vector.pop_back();
1704 
1705  while (MOBJ(eh).state() < _target_state - 1)
1706  Base::prev_rule()->raise(eh, _target_state - 1);
1707  }
1708  }
1709 
1710 
1711  // calculate new position
1712  typename M::Point position(0.0, 0.0, 0.0);
1713  typename M::Scalar c;
1714  typename M::Point zero_point(0.0, 0.0, 0.0);
1715  int valence(0);
1716 
1717  valence = Base::mesh_.valence(_vh);
1718  c = coeff( valence );
1719 
1720  for (voh_it = Base::mesh_.voh_iter(_vh); voh_it.is_valid(); ++voh_it)
1721  {
1722  if (MOBJ(Base::mesh_.EH(*voh_it)).final())
1723  {
1724  position += MOBJ(Base::mesh_.EH(*voh_it)).position(_target_state-1)*c;
1725 
1726  if ( Base::mesh_.FH(*voh_it).is_valid() &&
1727  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(*voh_it))).final() &&
1728  MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(*voh_it))).position(_target_state - 1) != zero_point)
1729  {
1730  position += MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(*voh_it))).position(_target_state-1) * (1.0-c);
1731  }
1732  else {
1733  position += MOBJ(Base::mesh_.EH(*voh_it)).position(_target_state - 1) * (1.0 - c);
1734  }
1735  }
1736  else {
1737  --valence;
1738  }
1739  }
1740 
1741  position /= valence;
1742 
1743  MOBJ(_vh).set_position(_target_state, position);
1744  MOBJ(_vh).inc_state();
1745 
1746  // check if last rule
1747  if (Base::number() == Base::n_rules() - 1) {
1748 
1749  Base::mesh_.set_point(_vh, position);
1750  MOBJ(_vh).set_final();
1751  }
1752  }
1753 }
1754 
1755 
1756 template <class M>
1757 void
1758 EVc<M>::init_coeffs(size_t _max_valence)
1759 {
1760  if ( coeffs_.size() == _max_valence+1 ) // equal? do nothing
1761  return;
1762 
1763  if (coeffs_.size() < _max_valence+1) // less than? add additional valences
1764  {
1765  const double _2pi = 2.0*M_PI;
1766 
1767  if (coeffs_.empty())
1768  coeffs_.push_back(0.0); // dummy for invalid valences 0,1,2
1769 
1770  for(size_t v=coeffs_.size(); v <= _max_valence; ++v)
1771  {
1772  // ( 3/2 + cos ( 2 PI / valence ) )� / 2 - 1
1773  double c = 1.5 + cos( _2pi / v );
1774  c = c * c * 0.5 - 1.0;
1775  coeffs_.push_back(c);
1776  }
1777  }
1778 }
1779 
1780 
1781 // -------------------------------------------------------------------- EF ----
1782 
1783 template<class M>
1784 void
1785 EF<M>::raise(typename M::FaceHandle& _fh, state_t _target_state) {
1786 
1787  if (MOBJ(_fh).state() < _target_state) {
1788 
1789  this->update(_fh, _target_state);
1790 
1791  // raise all neighbour edges to level x-1
1792  typename M::FaceEdgeIter fe_it(Base::mesh_.fe_iter(_fh));
1793  typename M::EdgeHandle eh;
1794  std::vector<typename M::EdgeHandle> edge_vector;
1795 
1796  if (_target_state > 1) {
1797 
1798  for (; fe_it.is_valid(); ++fe_it) {
1799 
1800  edge_vector.push_back(*fe_it);
1801  }
1802 
1803  while (!edge_vector.empty()) {
1804 
1805  eh = edge_vector.back();
1806  edge_vector.pop_back();
1807 
1808  Base::prev_rule()->raise(eh, _target_state - 1);
1809  }
1810 
1811  for (fe_it = Base::mesh_.fe_iter(_fh); fe_it.is_valid(); ++fe_it) {
1812 
1813  edge_vector.push_back(*fe_it);
1814  }
1815 
1816  while (!edge_vector.empty()) {
1817 
1818  eh = edge_vector.back();
1819  edge_vector.pop_back();
1820 
1821  while (MOBJ(eh).state() < _target_state - 1)
1822  Base::prev_rule()->raise(eh, _target_state - 1);
1823  }
1824  }
1825 
1826  // calculate new position
1827  typename M::Point position(0.0, 0.0, 0.0);
1828  typename M::Scalar valence(0.0);
1829 
1830  for (fe_it = Base::mesh_.fe_iter(_fh); fe_it.is_valid(); ++fe_it) {
1831 
1832  if (Base::mesh_.data(*fe_it).final()) {
1833 
1834  valence += 1.0;
1835 
1836  position += Base::mesh_.data(*fe_it).position(_target_state - 1);
1837  }
1838  }
1839 
1840  assert (valence == 3.0);
1841 
1842  position /= valence;
1843 
1844  MOBJ(_fh).set_position(_target_state, position);
1845  MOBJ(_fh).inc_state();
1846  }
1847 }
1848 
1849 
1850 // -------------------------------------------------------------------- FE ----
1851 
1852 
1853 template<class M>
1854 void
1855 FE<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state) {
1856 
1857  if (MOBJ(_eh).state() < _target_state) {
1858 
1859  this->update(_eh, _target_state);
1860 
1861  // raise all neighbor faces to level x-1
1862  typename M::FaceHandle fh;
1863 
1864  if (_target_state > 1) {
1865 
1866  fh = Base::mesh_.FH(Base::mesh_.HEH(_eh, 0));
1867  Base::prev_rule()->raise(fh, _target_state - 1);
1868 
1869  fh = Base::mesh_.FH(Base::mesh_.HEH(_eh, 1));
1870  Base::prev_rule()->raise(fh, _target_state - 1);
1871 
1872  fh = Base::mesh_.FH(Base::mesh_.HEH(_eh, 0));
1873  Base::prev_rule()->raise(fh, _target_state - 1);
1874 
1875  fh = Base::mesh_.FH(Base::mesh_.HEH(_eh, 1));
1876  Base::prev_rule()->raise(fh, _target_state - 1);
1877  }
1878 
1879  // calculate new position
1880  typename M::Point position(0.0, 0.0, 0.0);
1881  const typename M::Scalar valence(2.0);
1882 
1883  position += MOBJ(Base::mesh_.FH(Base::mesh_.HEH(_eh, 0))).position(_target_state - 1);
1884 
1885  position += MOBJ(Base::mesh_.FH(Base::mesh_.HEH(_eh, 1))).position(_target_state - 1);
1886 
1887  position /= valence;
1888 
1889  MOBJ(_eh).set_position(_target_state, position);
1890  MOBJ(_eh).inc_state();
1891  }
1892 }
1893 
1894 
1895 // ------------------------------------------------------------------- EdE ----
1896 
1897 
1898 template<class M>
1899 void
1900 EdE<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state) {
1901 
1902  if (MOBJ(_eh).state() < _target_state) {
1903 
1904  this->update(_eh, _target_state);
1905 
1906  // raise all neighbor faces and edges to level x-1
1907  typename M::HalfedgeHandle hh1, hh2;
1908  typename M::FaceHandle fh;
1909  typename M::EdgeHandle eh;
1910 
1911  hh1 = Base::mesh_.HEH(_eh, 0);
1912  hh2 = Base::mesh_.HEH(_eh, 1);
1913 
1914  if (_target_state > 1) {
1915 
1916  fh = Base::mesh_.FH(hh1);
1917  Base::prev_rule()->raise(fh, _target_state - 1);
1918 
1919  fh = Base::mesh_.FH(hh2);
1920  Base::prev_rule()->raise(fh, _target_state - 1);
1921 
1922  eh = Base::mesh_.EH(Base::mesh_.NHEH(hh1));
1923  Base::prev_rule()->raise(eh, _target_state - 1);
1924 
1925  eh = Base::mesh_.EH(Base::mesh_.PHEH(hh1));
1926  Base::prev_rule()->raise(eh, _target_state - 1);
1927 
1928  eh = Base::mesh_.EH(Base::mesh_.NHEH(hh2));
1929  Base::prev_rule()->raise(eh, _target_state - 1);
1930 
1931  eh = Base::mesh_.EH(Base::mesh_.PHEH(hh2));
1932  Base::prev_rule()->raise(eh, _target_state - 1);
1933  }
1934 
1935  // calculate new position
1936  typename M::Point position(0.0, 0.0, 0.0);
1937  const typename M::Scalar valence(4.0);
1938 
1939  position += MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh1))).position(_target_state - 1);
1940  position += MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh1))).position(_target_state - 1);
1941  position += MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh2))).position(_target_state - 1);
1942  position += MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh2))).position(_target_state - 1);
1943 
1944  position /= valence;
1945 
1946  MOBJ(_eh).set_position(_target_state, position);
1947  MOBJ(_eh).inc_state();
1948  }
1949 }
1950 
1951 
1952 // ------------------------------------------------------------------ EdEc ----
1953 
1954 
1955 template<class M>
1956 void
1957 EdEc<M>::raise(typename M::EdgeHandle& _eh, state_t _target_state)
1958 {
1959  if (MOBJ(_eh).state() < _target_state) {
1960 
1961  this->update(_eh, _target_state);
1962 
1963  // raise all neighbor faces and edges to level x-1
1964  typename M::HalfedgeHandle hh1, hh2;
1965  typename M::FaceHandle fh;
1966  typename M::EdgeHandle eh;
1967 
1968  hh1 = Base::mesh_.HEH(_eh, 0);
1969  hh2 = Base::mesh_.HEH(_eh, 1);
1970 
1971  if (_target_state > 1) {
1972 
1973  fh = Base::mesh_.FH(hh1);
1974  Base::prev_rule()->raise(fh, _target_state - 1);
1975 
1976  fh = Base::mesh_.FH(hh2);
1977  Base::prev_rule()->raise(fh, _target_state - 1);
1978 
1979  eh = Base::mesh_.EH(Base::mesh_.NHEH(hh1));
1980  Base::prev_rule()->raise(eh, _target_state - 1);
1981 
1982  eh = Base::mesh_.EH(Base::mesh_.PHEH(hh1));
1983  Base::prev_rule()->raise(eh, _target_state - 1);
1984 
1985  eh = Base::mesh_.EH(Base::mesh_.NHEH(hh2));
1986  Base::prev_rule()->raise(eh, _target_state - 1);
1987 
1988  eh = Base::mesh_.EH(Base::mesh_.PHEH(hh2));
1989  Base::prev_rule()->raise(eh, _target_state - 1);
1990  }
1991 
1992  // calculate new position
1993  typename M::Point position(0.0, 0.0, 0.0);
1994  const typename M::Scalar valence(4.0);
1995  typename M::Scalar c;
1996 
1997  position += MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh1))).position(_target_state - 1);
1998  position += MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh1))).position(_target_state - 1);
1999  position += MOBJ(Base::mesh_.EH(Base::mesh_.NHEH(hh2))).position(_target_state - 1);
2000  position += MOBJ(Base::mesh_.EH(Base::mesh_.PHEH(hh2))).position(_target_state - 1);
2001 
2002  position /= valence;
2003 
2004  // choose coefficient c
2005  c = Base::coeff();
2006 
2007  position *= (1.0 - c);
2008 
2009  position += MOBJ(_eh).position(_target_state - 1) * c;
2010 
2011  MOBJ(_eh).set_position(_target_state, position);
2012  MOBJ(_eh).inc_state();
2013  }
2014 }
2015 
2016 #endif // DOXY_IGNORE_THIS
2017 
2018 #undef FH
2019 #undef VH
2020 #undef EH
2021 #undef HEH
2022 #undef M
2023 #undef MOBJ
2024 
2025 //=============================================================================
2026 } // END_NS_ADAPTIVE
2027 } // END_NS_SUBDIVIDER
2028 } // END_NS_OPENMESH
2029 //=============================================================================
OpenMesh::Subdivider::Adaptive::FF::raise
void raise(typename M::FaceHandle &_fh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::FV::raise
void raise(typename M::VertexHandle &_vh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::VdEc::raise
void raise(typename M::EdgeHandle &_eh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::VV::raise
void raise(typename M::VertexHandle &_vh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::FE::raise
void raise(typename M::EdgeHandle &_eh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::VdE::raise
void raise(typename M::EdgeHandle &_eh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::EF::raise
void raise(typename M::FaceHandle &_fh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::Tvv4::raise
void raise(typename M::FaceHandle &_fh, state_t _target_state)
Raise item to target state _target_state.
Definition: RulesT.cc:343
OpenMesh::Subdivider::Adaptive::FVc::raise
void raise(typename M::VertexHandle &_vh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::VVc::raise
void raise(typename M::VertexHandle &_vh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::VF::raise
void raise(typename M::FaceHandle &_fh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::EVc::raise
void raise(typename M::VertexHandle &_vh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::VE::raise
void raise(typename M::EdgeHandle &_eh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::EdE::raise
void raise(typename M::EdgeHandle &_eh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::Tvv3::raise
void raise(typename M::FaceHandle &_fh, state_t _target_state)
Raise item to target state _target_state.
Definition: RulesT.cc:104
OpenMesh::Subdivider::Adaptive::EdEc::raise
void raise(typename M::EdgeHandle &_eh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::FFc::raise
void raise(typename M::FaceHandle &_fh, state_t _target_state)
Raise item to target state _target_state.
OpenMesh::Subdivider::Adaptive::state_t
CompositeTraits::state_t state_t
Definition: CompositeTraits.hh:255
OpenMesh::Subdivider::Adaptive::EV::raise
void raise(typename M::VertexHandle &_vh, state_t _target_state)
Raise item to target state _target_state.