Developer Documentation
Loading...
Searching...
No Matches
BaseRemesherT_impl.hh
1/*===========================================================================*\
2* *
3* OpenFlipper *
4 * Copyright (c) 2001-2015, RWTH-Aachen University *
5 * Department of Computer Graphics and Multimedia *
6 * All rights reserved. *
7 * www.openflipper.org *
8 * *
9 *---------------------------------------------------------------------------*
10 * This file is part of OpenFlipper. *
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// CLASS BaseRemesherT - IMPLEMENTATION
45//
46//=============================================================================
47
48#define BASE_REMESHERT_C
49
50
51//== INCLUDES =================================================================
52
53
54#include "BaseRemesherT.hh"
55
56// Remshing
57#include "DiffGeoT.hh"
58
59// OpenMesh
60#include <OpenMesh/Core/Utils/Property.hh>
61#include <OpenMesh/Core/System/omstream.hh>
62
63// Selection Stuff
64#include <MeshTools/MeshSelectionT.hh>
65
66
67//== NAMESPACES ===============================================================
68
69namespace Remeshing {
70
71
72//== IMPLEMENTATION ==========================================================
73
74
75template <class Mesh>
76BaseRemesherT<Mesh>::
77BaseRemesherT(Mesh& _mesh, ProgressEmitter* _progress)
78 : mesh_(_mesh), refmesh_(0), bsp_(0), nothing_selected_(true),progress_(_progress) {
79
80}
81
82
83//-----------------------------------------------------------------------------
84
85
86template <class Mesh>
87BaseRemesherT<Mesh>::
88~BaseRemesherT() {
89
90 delete bsp_;
91 delete refmesh_;
92}
93
94
95//-----------------------------------------------------------------------------
96
97
98template <class Mesh>
99void
100BaseRemesherT<Mesh>::
101init_reference()
102{
103 typename Mesh::VIter v_it, v_end;
104 typename Mesh::FIter f_it, f_end;
105 typename Mesh::CFVIter fv_it;
106 typename Mesh::VHandle v0, v1, v2;
107
108
109 // clear old stuff first
110 if (bsp_) delete_reference();
111
112 // tag non-locked vertices + one ring
113 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
114 v_it!=v_end; ++v_it)
115 mesh_.status(*v_it).set_tagged(!mesh_.status(*v_it).locked());
116
117 MeshSelection::growVertexSelection(&mesh_);
118
119 // create reference mesh
120 refmesh_ = new Mesh();
121
122 OpenMesh::VPropHandleT<typename Mesh::VHandle> vhandles;
123 mesh_.add_property(vhandles);
124
125 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
126 v_it!=v_end; ++v_it)
127 if (mesh_.status(*v_it).tagged())
128 mesh_.property(vhandles, *v_it) = refmesh_->add_vertex(mesh_.point(*v_it));
129
130 for (f_it=mesh_.faces_begin(), f_end=mesh_.faces_end();
131 f_it!=f_end; ++f_it)
132 {
133 fv_it = mesh_.cfv_iter(*f_it);
134 v0 = *fv_it;
135 v1 = *(++fv_it);
136 v2 = *(++fv_it);
137
138 if (mesh_.status(v0).tagged() &&
139 mesh_.status(v1).tagged() &&
140 mesh_.status(v2).tagged())
141 {
142 refmesh_->add_face( mesh_.property(vhandles, v0),
143 mesh_.property(vhandles, v1),
144 mesh_.property(vhandles, v2) );
145 }
146 }
147
148 MeshSelection::clearVertexSelection(&mesh_);
149 mesh_.remove_property(vhandles);
150
151 // need vertex normals
152 refmesh_->request_face_normals();
153 refmesh_->request_vertex_normals();
154 refmesh_->update_normals(/*Mesh::ANGLE_WEIGHTED*/);
155
156
157 // build BSP
158 bsp_ = new BSP(*refmesh_);
159 bsp_->reserve(refmesh_->n_faces());
160 for (f_it=refmesh_->faces_begin(), f_end=refmesh_->faces_end();
161 f_it!=f_end; ++f_it)
162 bsp_->push_back(*f_it);
163 bsp_->build(10, 100);
164}
165
166
167//-----------------------------------------------------------------------------
168
169
170template <class Mesh>
171void
172BaseRemesherT<Mesh>::
173delete_reference()
174{
175 delete bsp_; bsp_ = 0;
176 delete refmesh_; refmesh_ = 0;
177}
178
179
180//-----------------------------------------------------------------------------
181
182
183template <class Mesh>
184void
185BaseRemesherT<Mesh>::
186project_to_reference(typename Mesh::VertexHandle _vh) const
187{
188 typename Mesh::Point p = mesh_.point(_vh);
189 typename Mesh::FaceHandle fh = bsp_->nearest(p).handle;
190
191 if ( ! fh.is_valid() )
192 return;
193
194 typename Mesh::CFVIter fv_it = refmesh_->cfv_iter(fh);
195
196
197 const typename Mesh::Point& p0 = refmesh_->point(*fv_it);
198 typename Mesh::Normal n0 = refmesh_->normal(*fv_it);
199 const typename Mesh::Point& p1 = refmesh_->point(*(++fv_it));
200 typename Mesh::Normal n1 = refmesh_->normal(*fv_it);
201 const typename Mesh::Point& p2 = refmesh_->point(*(++fv_it));
202 typename Mesh::Normal n2 = refmesh_->normal(*fv_it);
203
204
205 // project
206 ACG::Geometry::distPointTriangle(p, p0, p1, p2, p);
207 mesh_.set_point(_vh, p);
208
209
210 // get barycentric coordinates
211 if (!ACG::Geometry::baryCoord(p, p0, p1, p2, p))
212 p[0] = p[1] = p[2] = 1.0/3.0;
213
214
215 // interpolate normal
216 typename Mesh::Normal n;
217 n = (n0 *= p[0]);
218 n += (n1 *= p[1]);
219 n += (n2 *= p[2]);
220 n.normalize();
221 mesh_.set_normal(_vh, n);
222}
223
224
225//-----------------------------------------------------------------------------
226
227
228template <class Mesh>
229void
230BaseRemesherT<Mesh>::
231remesh(unsigned int _iters,
232 unsigned int _area_iters,
233 bool _use_projection,
234 Selection _selection) {
235
236 try
237 {
238 if (_selection == VERTEX_SELECTION)
239 prepare_vertex_selection();
240 else if (_selection == FACE_SELECTION)
241 prepare_face_selection();
242 remeshh(_iters, _area_iters, _use_projection);
243 }
244 catch (std::bad_alloc&)
245 {
246 mesh_.clear();
247 omerr() << "Remeshig: Out of memory\n";
248 }
249
250 cleanup();
251}
252
253//-----------------------------------------------------------------------------
254
255template <class Mesh>
256void
257BaseRemesherT<Mesh>::
258prepare_vertex_selection()
259{
260 typename Mesh::EIter e_it, e_end;
261 typename Mesh::VIter v_it, v_end;
262 typename Mesh::FIter f_it, f_end;
263 typename Mesh::CFVIter fv_it;
264 typename Mesh::CVOHIter vh_it;
265 typename Mesh::VHandle v0, v1;
266 typename Mesh::FHandle f0, f1;
267
268
269 // need vertex and edge status
270 mesh_.request_vertex_status();
271 mesh_.request_edge_status();
272 mesh_.request_face_status();
273
274
275
276 // if nothing selected -> select all
277 nothing_selected_ = true;
278 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
279 v_it!=v_end; ++v_it)
280 if (mesh_.status(*v_it).selected())
281 { nothing_selected_ = false; break; }
282
283 if (nothing_selected_)
284 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
285 v_it!=v_end; ++v_it)
286 mesh_.status(*v_it).set_selected(true);
287
288
289
290 // lock un-selected vertices & edges
291 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
292 v_it!=v_end; ++v_it)
293 mesh_.status(*v_it).set_locked(!mesh_.status(*v_it).selected());
294
295 for (e_it=mesh_.edges_begin(), e_end=mesh_.edges_end();
296 e_it!=e_end; ++e_it)
297 {
298 v0 = mesh_.to_vertex_handle(mesh_.halfedge_handle(*e_it, 0));
299 v1 = mesh_.to_vertex_handle(mesh_.halfedge_handle(*e_it, 1));
300 mesh_.status(*e_it).set_locked(mesh_.status(v0).locked() ||
301 mesh_.status(v1).locked());
302 }
303
304
305
306 // handle feature corners:
307 // lock corner vertices (>2 feature edges) and endpoints
308 // of feature lines (1 feature edge)
309 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
310 v_it!=v_end; ++v_it)
311 {
312 if (mesh_.status(*v_it).feature())
313 {
314 int c=0;
315 for (vh_it=mesh_.cvoh_iter(*v_it); vh_it.is_valid(); ++vh_it)
316 if (mesh_.status(mesh_.edge_handle(*vh_it)).feature())
317 ++c;
318 if (c!=2) mesh_.status(*v_it).set_locked(true);
319 }
320 }
321
322
323
324 // build reference mesh
325 init_reference();
326// if (emit_progress_) Progress().step(5);
327
328
329 // add properties
330 mesh_.add_property(valences_);
331 mesh_.add_property(update_);
332 mesh_.add_property(area_);
333}
334
335//-----------------------------------------------------------------------------
336
337
338template <class Mesh>
339void
340BaseRemesherT<Mesh>::
341prepare_face_selection()
342{
343 typename Mesh::EIter e_it, e_end;
344 typename Mesh::VIter v_it, v_end;
345 typename Mesh::FIter f_it, f_end;
346 typename Mesh::CFVIter fv_it;
347 typename Mesh::CVOHIter vh_it;
348 typename Mesh::VHandle v0, v1;
349 typename Mesh::FHandle f0, f1;
350
351
352 // need vertex and edge status
353 mesh_.request_vertex_status();
354 mesh_.request_edge_status();
355 mesh_.request_face_status();
356
357 // if nothing selected -> select all
358 nothing_selected_ = true;
359 for (f_it = mesh_.faces_begin(), f_end = mesh_.faces_end(); f_it != f_end;
360 ++f_it)
361 {
362 if (mesh_.status(*f_it).selected())
363 {
364 nothing_selected_ = false;
365 break;
366 }
367 }
368
369 if (nothing_selected_)
370 MeshSelection::selectAllFaces(&mesh_);
371
372
373
374 // lock un-selected vertices & edges
375 for (v_it = mesh_.vertices_begin(), v_end = mesh_.vertices_end();
376 v_it != v_end; ++v_it)
377 {
378 bool all_faces_selected = true;
379
380 for (typename Mesh::ConstVertexFaceIter vf_it = mesh_.cvf_iter(*v_it);
381 vf_it.is_valid(); ++vf_it)
382 {
383 if (!mesh_.status(*vf_it).selected())
384 {
385 all_faces_selected = false;
386 break;
387 }
388 }
389 mesh_.status(*v_it).set_locked(!all_faces_selected);
390 }
391
392 for (e_it = mesh_.edges_begin(), e_end = mesh_.edges_end();
393 e_it != e_end; ++e_it)
394 {
395 if (mesh_.is_boundary(*e_it))
396 {
397 mesh_.status(*e_it).set_locked(true);
398 }
399 else
400 {
401 f0 = mesh_.face_handle(mesh_.halfedge_handle(*e_it, 0));
402 f1 = mesh_.face_handle(mesh_.halfedge_handle(*e_it, 1));
403
404 mesh_.status(*e_it).set_locked(!(mesh_.status(f0).selected() && mesh_.status(f1).selected()));
405 }
406 }
407
408 MeshSelection::clearFaceSelection(&mesh_);
409
410 // handle feature corners:
411 // lock corner vertices (>2 feature edges) and endpoints
412 // of feature lines (1 feature edge)
413 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
414 v_it!=v_end; ++v_it)
415 {
416 if (mesh_.status(*v_it).feature())
417 {
418 int c=0;
419 for (vh_it=mesh_.cvoh_iter(*v_it); vh_it.is_valid(); ++vh_it)
420 if (mesh_.status(mesh_.edge_handle(*vh_it)).feature())
421 ++c;
422 if (c!=2) mesh_.status(*v_it).set_locked(true);
423 }
424 }
425
426
427
428 // build reference mesh
429 init_reference();
430// if (emit_progress_) Progress().step(5);
431
432
433 // add properties
434 mesh_.add_property(valences_);
435 mesh_.add_property(update_);
436 mesh_.add_property(area_);
437}
438
439
440//-----------------------------------------------------------------------------
441
442
443template <class Mesh>
444void
445BaseRemesherT<Mesh>::
446remeshh(unsigned int _iters,
447 unsigned int _area_iters, bool _use_projection) {
448
449 double progress_step = 100.0 / (((double)_iters/4.0 + (double)_area_iters));
450 double total_progress = 0.0;
451
452 for (unsigned int i = 0; i < _iters; ++i) {
453
454 split_long_edges();
455 if(progress_ != NULL) {
456 total_progress += progress_step;
457 progress_->sendProgressSignal(total_progress);
458 }
459
460 collapse_short_edges();
461 if(progress_ != NULL) {
462 total_progress += progress_step;
463 progress_->sendProgressSignal(total_progress);
464 }
465
466 flip_edges();
467 if(progress_ != NULL) {
468 total_progress += progress_step;
469 progress_->sendProgressSignal(total_progress);
470 }
471
472 tangential_smoothing(_use_projection);
473 if(progress_ != NULL) {
474 total_progress += progress_step;
475 progress_->sendProgressSignal(total_progress);
476 }
477 }
478
479 if (_area_iters) {
480 balanace_area(_area_iters, _use_projection);
481 if(progress_ != NULL) {
482 total_progress += progress_step;
483 progress_->sendProgressSignal(total_progress);
484 }
485 }
486
487 // feature edges block flips -> might lead to caps
488 remove_caps();
489}
490
491
492//-----------------------------------------------------------------------------
493
494
495template <class Mesh>
496void
497BaseRemesherT<Mesh>::
498cleanup()
499{
500 typename Mesh::EIter e_it, e_end;
501 typename Mesh::VIter v_it, v_end;
502 typename Mesh::FIter f_it, f_end;
503 typename Mesh::CFVIter fv_it;
504 typename Mesh::CVOHIter vh_it;
505 typename Mesh::VHandle v0, v1;
506 typename Mesh::FHandle f0, f1;
507
508
509 // remove properties
510 mesh_.remove_property(valences_);
511 mesh_.remove_property(update_);
512 mesh_.remove_property(area_);
513
514
515 // remove reference again
516 delete_reference();
517
518
519 // unlock all vertices & edges
520 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
521 v_it!=v_end; ++v_it)
522 {
523 mesh_.status(*v_it).set_locked(false);
524 }
525 for (e_it=mesh_.edges_begin(), e_end=mesh_.edges_end();
526 e_it!=e_end; ++e_it)
527 {
528 mesh_.status(*e_it).set_locked(false);
529 }
530
531
532 // de-select if nothing was selected before
533 if (nothing_selected_)
534 for (v_it=mesh_.vertices_begin(), v_end=mesh_.vertices_end();
535 v_it!=v_end; ++v_it)
536 mesh_.status(*v_it).set_selected(false);
537
538
539
540 // free vertex and edge status
541 mesh_.release_vertex_status();
542 mesh_.release_edge_status();
543 mesh_.release_face_status();
544}
545
546
547//-----------------------------------------------------------------------------
548
549
550template <class Mesh>
551void
552BaseRemesherT<Mesh>::
553split_long_edges()
554{
555 typename Mesh::VIter v_it, v_end;
556 typename Mesh::EIter e_it, e_end;
557 typename Mesh::VHandle v0, v1, vh;
558 typename Mesh::EHandle eh, e0, e1;
559 typename Mesh::FHandle f0, f1, f2, f3;
560
561 bool ok, is_feature, is_boundary;
562 int i;
563
564 // un-tagg
565 for (v_it = mesh_.vertices_begin(), v_end = mesh_.vertices_end(); v_it != v_end; ++v_it)
566 mesh_.status(*v_it).set_tagged(false);
567
568 // handle Nastran PIDs during refinement
569 OpenMesh::FPropHandleT<unsigned int> pids;
570 mesh_.get_property_handle(pids, "Nastran PIDs");
571
572 // split long edges
573 for (ok = false, i = 0; !ok && i < 100; ++i) {
574 ok = true;
575
576 for (e_it = mesh_.edges_begin(), e_end = mesh_.edges_end(); e_it != e_end; ++e_it) {
577 v0 = mesh_.to_vertex_handle(mesh_.halfedge_handle(*e_it, 0));
578 v1 = mesh_.to_vertex_handle(mesh_.halfedge_handle(*e_it, 1));
579
580 if (!mesh_.status(*e_it).locked() && is_too_long(v0, v1)) {
581 const typename Mesh::Point& p0 = mesh_.point(v0);
582 const typename Mesh::Point& p1 = mesh_.point(v1);
583
584 is_feature = mesh_.status(*e_it).feature();
585 is_boundary = mesh_.is_boundary(*e_it);
586
587 vh = mesh_.add_vertex((p0 + p1) * 0.5);
588 mesh_.split(*e_it, vh);
589
590 mesh_.status(vh).set_selected(true);
591
592 if (is_feature) {
593 eh = (is_boundary ? mesh_.edge_handle(mesh_.n_edges() - 2) : mesh_.edge_handle(mesh_.n_edges() - 3));
594
595 mesh_.status(eh).set_feature(true);
596 } else {
597 project_to_reference(vh);
598 }
599
600 if (pids.is_valid()) {
601 e0 = *e_it;
602 e1 = (is_boundary ? mesh_.edge_handle(mesh_.n_edges() - 2) : mesh_.edge_handle(mesh_.n_edges() - 3));
603
604 f0 = mesh_.face_handle(mesh_.halfedge_handle(e0, 0));
605 f1 = mesh_.face_handle(mesh_.halfedge_handle(e0, 1));
606 f2 = mesh_.face_handle(mesh_.halfedge_handle(e1, 0));
607 f3 = mesh_.face_handle(mesh_.halfedge_handle(e1, 1));
608
609 if (f0.is_valid() && f3.is_valid()) mesh_.property(pids, f3) = mesh_.property(pids, f0);
610
611 if (f1.is_valid() && f2.is_valid()) mesh_.property(pids, f2) = mesh_.property(pids, f1);
612 }
613
614 ok = false;
615 }
616 }
617 }
618
619 if (i == 100) omlog() << "split break\n";
620}
621
622
623//-----------------------------------------------------------------------------
624
625
626template <class Mesh>
627void
628BaseRemesherT<Mesh>::
629collapse_short_edges()
630{
631 typename Mesh::EIter e_it, e_end;
632 typename Mesh::CVVIter vv_it;
633 typename Mesh::VHandle v0, v1;
634 typename Mesh::HHandle h0, h1, h01, h10;
635 typename Mesh::Point p;
636 bool ok, skip, b0, b1, l0, l1, f0, f1;
637 int i;
638 bool hcol01, hcol10;
639
640 for (ok = false, i = 0; !ok && i < 100; ++i) {
641 ok = true;
642
643 for (e_it = mesh_.edges_begin(), e_end = mesh_.edges_end(); e_it != e_end; ++e_it) {
644 if (!mesh_.status(*e_it).deleted() && !mesh_.status(*e_it).locked()) {
645 h10 = mesh_.halfedge_handle(*e_it, 0);
646 h01 = mesh_.halfedge_handle(*e_it, 1);
647 v0 = mesh_.to_vertex_handle(h10);
648 v1 = mesh_.to_vertex_handle(h01);
649
650 if (is_too_short(v0, v1)) {
651 // get status
652 b0 = mesh_.is_boundary(v0);
653 b1 = mesh_.is_boundary(v1);
654 l0 = mesh_.status(v0).locked();
655 l1 = mesh_.status(v1).locked();
656 f0 = mesh_.status(v0).feature();
657 f1 = mesh_.status(v1).feature();
658 hcol01 = hcol10 = true;
659
660 // boundary rules
661 if (b0 && b1) {
662 if (!mesh_.is_boundary(*e_it))
663 continue;
664 } else if (b0)
665 hcol01 = false;
666 else if (b1)
667 hcol10 = false;
668
669 // locked rules
670 if (l0 && l1)
671 continue;
672 else if (l0)
673 hcol01 = false;
674 else if (l1)
675 hcol10 = false;
676
677 // feature vertex rules
678 if (f0 && f1)
679 continue;
680 else if (f0)
681 hcol01 = false;
682 else if (f1)
683 hcol10 = false;
684
685 // feature edge rules
686 // if vertex is incident to feature edge but collapse edge is not feature dont collapse
687
688 // hcol01
689 {
690 int feature_valence_0 = 0;
691 for (auto eh : mesh_.ve_range(v0))
692 if (mesh_.status(eh).feature())
693 ++feature_valence_0;
694 if (feature_valence_0 == 2)
695 {
696 if (!mesh_.status(*e_it).feature())
697 hcol01 = false;
698 }
699 else if (feature_valence_0 != 0)
700 hcol01 = false;
701 }
702 // hcol10
703 {
704 int feature_valence_1 = 0;
705 for (auto eh : mesh_.ve_range(v1))
706 if (mesh_.status(eh).feature())
707 ++feature_valence_1;
708 if (feature_valence_1 == 2)
709 {
710 if (!mesh_.status(*e_it).feature())
711 hcol10 = false;
712 }
713 else if (feature_valence_1 != 0)
714 hcol10 = false;
715 }
716
717 // the other two edges removed by collapse must not be features
718 h0 = mesh_.prev_halfedge_handle(h01);
719 h1 = mesh_.next_halfedge_handle(h10);
720 if (mesh_.status(mesh_.edge_handle(h0)).feature() || mesh_.status(mesh_.edge_handle(h1)).feature())
721 hcol01 = false;
722 h0 = mesh_.prev_halfedge_handle(h10);
723 h1 = mesh_.next_halfedge_handle(h01);
724 if (mesh_.status(mesh_.edge_handle(h0)).feature() || mesh_.status(mesh_.edge_handle(h1)).feature())
725 hcol10 = false;
726
727 // topological rules
728 if (hcol01)
729 hcol01 = mesh_.is_collapse_ok(h01);
730 if (hcol10)
731 hcol10 = mesh_.is_collapse_ok(h10);
732
733 // both collapses possible: collapse into vertex w/ higher valence
734 if (hcol01 && hcol10) {
735 if (mesh_.valence(v0) < mesh_.valence(v1))
736 hcol10 = false;
737 else
738 hcol01 = false;
739 }
740
741 // check if collapse flips triangles
742 if (collapse_flips_normal(h01))
743 hcol01 = false;
744 if (collapse_flips_normal(h10))
745 hcol10 = false;
746
747 // try v1 -> v0
748 if (hcol10) {
749 // don't create too long edges
750 skip = false;
751 for (vv_it = mesh_.cvv_iter(v1); vv_it.is_valid() && !skip; ++vv_it)
752 if (is_too_long(v0, *vv_it))
753 skip = true;
754
755 if (!skip) {
756 mesh_.collapse(h10);
757 ok = false;
758 }
759 }
760
761 // try v0 -> v1
762 else if (hcol01) {
763 // don't create too long edges
764 skip = false;
765 for (vv_it = mesh_.cvv_iter(v0); vv_it.is_valid() && !skip; ++vv_it)
766 if (is_too_long(v1, *vv_it))
767 skip = true;
768
769 if (!skip) {
770 mesh_.collapse(h01);
771 ok = false;
772 }
773 }
774 }
775 }
776 }
777 }
778
779 mesh_.garbage_collection();
780
781 if (i == 100)
782 omlog() << "collapse break\n";
783}
784
785
786//-----------------------------------------------------------------------------
787
788
789template <class Mesh>
790void
791BaseRemesherT<Mesh>::
792flip_edges()
793{
794 typename Mesh::EIter e_it, e_end;
795 typename Mesh::VIter v_it, v_end;
796 typename Mesh::VHandle v0, v1, v2, v3, vh;
797 typename Mesh::HHandle hh;
798 typename Mesh::Point p;
799 typename Mesh::FHandle fh;
800
801 int val0, val1, val2, val3;
802 int val_opt0, val_opt1, val_opt2, val_opt3;
803 int ve0, ve1, ve2, ve3, ve_before, ve_after;
804 bool ok;
805 int i;
806
807 // compute vertex valences
808 for (v_it = mesh_.vertices_begin(), v_end = mesh_.vertices_end(); v_it != v_end; ++v_it)
809 mesh_.property(valences_, *v_it) = mesh_.valence(*v_it);
810
811 // flip all edges
812 for (ok = false, i = 0; !ok && i < 100; ++i) {
813 ok = true;
814
815 for (e_it = mesh_.edges_begin(), e_end = mesh_.edges_end(); e_it != e_end; ++e_it) {
816 if (!mesh_.status(*e_it).locked() && !mesh_.status(*e_it).feature()) {
817 hh = mesh_.halfedge_handle(*e_it, 0);
818 v0 = mesh_.to_vertex_handle(hh);
819 v2 = mesh_.to_vertex_handle(mesh_.next_halfedge_handle(hh));
820 if ( !mesh_.next_halfedge_handle(hh).is_valid() ) {
821 std::cerr << "Error v2" << std::endl;
822 continue;
823 }
824 hh = mesh_.halfedge_handle(*e_it, 1);
825 v1 = mesh_.to_vertex_handle(hh);
826 v3 = mesh_.to_vertex_handle(mesh_.next_halfedge_handle(hh));
827 if ( !mesh_.next_halfedge_handle(hh).is_valid() ) {
828 std::cerr << "Error v3" << std::endl;
829 continue;
830 }
831
832 if ( !v2.is_valid())
833 continue;
834
835 if (!mesh_.status(v0).locked() && !mesh_.status(v1).locked() && !mesh_.status(v2).locked()
836 && !mesh_.status(v3).locked()) {
837 val0 = mesh_.property(valences_, v0);
838 val1 = mesh_.property(valences_, v1);
839 val2 = mesh_.property(valences_, v2);
840 val3 = mesh_.property(valences_, v3);
841
842 val_opt0 = (mesh_.is_boundary(v0) ? 4 : 6);
843 val_opt1 = (mesh_.is_boundary(v1) ? 4 : 6);
844 val_opt2 = (mesh_.is_boundary(v2) ? 4 : 6);
845 val_opt3 = (mesh_.is_boundary(v3) ? 4 : 6);
846
847 ve0 = (val0 - val_opt0);
848 ve0 *= ve0;
849 ve1 = (val1 - val_opt1);
850 ve1 *= ve1;
851 ve2 = (val2 - val_opt2);
852 ve2 *= ve2;
853 ve3 = (val3 - val_opt3);
854 ve3 *= ve3;
855
856 ve_before = ve0 + ve1 + ve2 + ve3;
857
858 --val0;
859 --val1;
860 ++val2;
861 ++val3;
862
863 ve0 = (val0 - val_opt0);
864 ve0 *= ve0;
865 ve1 = (val1 - val_opt1);
866 ve1 *= ve1;
867 ve2 = (val2 - val_opt2);
868 ve2 *= ve2;
869 ve3 = (val3 - val_opt3);
870 ve3 *= ve3;
871
872 ve_after = ve0 + ve1 + ve2 + ve3;
873
874 if (ve_before > ve_after && mesh_.is_flip_ok(*e_it) && !edge_flip_flips_normal(*e_it)) {
875 mesh_.flip(*e_it);
876 --mesh_.property(valences_, v0);
877 --mesh_.property(valences_, v1);
878 ++mesh_.property(valences_, v2);
879 ++mesh_.property(valences_, v3);
880 ok = false;
881 }
882 }
883 }
884 }
885 }
886
887 if (i == 100)
888 omlog() << "flip break\n";
889}
890
891
892//-----------------------------------------------------------------------------
893
894
895template <class Mesh>
896void
897BaseRemesherT<Mesh>::
898tangential_smoothing(bool _use_projection)
899{
900 typename Mesh::VIter v_it, v_end(mesh_.vertices_end());
901 typename Mesh::CVVIter vv_it;
902 typename Mesh::Point u, n;
903
904
905 // tag active vertices
906 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
907 mesh_.status(*v_it).set_tagged( !mesh_.status(*v_it).locked() &&
908 !mesh_.status(*v_it).feature() &&
909 !mesh_.is_boundary(*v_it) );
910
911 mesh_.update_face_normals();
912
913 double damping = 1.0;
914
915 // smooth
916 for (int iters=0; iters<10; ++iters)
917 {
918 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
919 {
920 if (mesh_.status(*v_it).tagged())
921 {
922 u.vectorize(0.0);
923 int valence = 0;
924 int feature_valence = 0;
925
926 for (auto heh : mesh_.voh_range(*v_it))
927 {
928 u += mesh_.point(mesh_.to_vertex_handle(heh));
929 if (mesh_.status(mesh_.edge_handle(heh)).feature())
930 ++feature_valence;
931 ++valence;
932 }
933
934 if (valence)
935 {
936 u *= (1.0/valence);
937 u -= mesh_.point(*v_it);
938 n = mesh_.normal(*v_it);
939 n *= (u | n);
940 u -= n;
941 }
942
943 if (feature_valence == 2)
944 {
945 // project update onto feature directions
946 typename Mesh::Point feature_dir(0.0,0.0,0.0);
947 for (auto heh : mesh_.voh_range(*v_it))
948 {
949 if (mesh_.status(mesh_.edge_handle(heh)).feature())
950 {
951 auto dir = mesh_.point(mesh_.to_vertex_handle(heh)) - mesh_.point(mesh_.from_vertex_handle(heh));
952 if ((dir | feature_dir) > 0)
953 feature_dir += dir;
954 else
955 feature_dir -= dir;
956 }
957 }
958 if (feature_dir.sqrnorm() > 0)
959 feature_dir.normalize();
960 u = (feature_dir | u) * feature_dir;
961 }
962 else if (feature_valence != 0)
963 {
964 // more than two or exactly one incident feature edge means vertex should be preserved
965 u.vectorize(0.0);
966 }
967
968 mesh_.property(update_, *v_it) = u;
969 }
970 }
971
972 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
973 if (mesh_.status(*v_it).tagged())
974 mesh_.point(*v_it) += damping*mesh_.property(update_, *v_it);
975
976 // check if normals changed
977 bool normals_changed = false;
978 for (auto fh : mesh_.faces())
979 {
980 if ((mesh_.calc_face_normal(fh) | mesh_.normal(fh)) < 0)
981 {
982 normals_changed = true;
983 break;
984 }
985 }
986
987 if (normals_changed)
988 {
989 // revert update and try again with smaller step
990 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
991 if (mesh_.status(*v_it).tagged())
992 mesh_.point(*v_it) -= damping*mesh_.property(update_, *v_it);
993 damping *= 0.5;
994 }
995 }
996
997
998 // reset tagged status
999 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1000 mesh_.status(*v_it).set_tagged(false);
1001
1002
1003 // project
1004 if (_use_projection)
1005 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1006 if (!mesh_.status(*v_it).locked() && !mesh_.status(*v_it).feature())
1007 project_to_reference(*v_it);
1008}
1009
1010
1011//-----------------------------------------------------------------------------
1012
1013
1014template <class Mesh>
1015void
1016BaseRemesherT<Mesh>::
1017balanace_area(unsigned int _iters, bool _use_projection)
1018{
1019 typename Mesh::VIter v_it, v_end(mesh_.vertices_end());
1020 typename Mesh::CVVIter vv_it;
1021 typename Mesh::Scalar w, ww;
1022 typename Mesh::Point u, n;
1023
1024
1025 DiffGeoT<Mesh> diffgeo(mesh_);
1026
1027
1028 // tag active vertices
1029 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1030 {
1031 bool active = ( !mesh_.status(*v_it).locked() &&
1032 !mesh_.status(*v_it).feature() &&
1033 !mesh_.is_boundary(*v_it) );
1034
1035 // don't move neighbors of boundary vertices
1036 for (vv_it=mesh_.cvv_iter(*v_it); active && vv_it.is_valid(); ++vv_it)
1037 if (mesh_.is_boundary(*vv_it))
1038 active = false;
1039
1040 mesh_.status(*v_it).set_tagged( active );
1041 }
1042
1043
1044 // tag2 vertices for which area has to be computed
1045 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1046 {
1047 mesh_.status(*v_it).set_tagged2(false);
1048 for (vv_it=mesh_.cvv_iter(*v_it); vv_it.is_valid(); ++vv_it)
1049 {
1050 if (mesh_.status(*vv_it).tagged())
1051 {
1052 mesh_.status(*v_it).set_tagged2(true);
1053 break;
1054 }
1055 }
1056 }
1057
1058
1059
1060 for (unsigned int bla=0; bla<_iters; ++bla)
1061 {
1062 // compute vertex areas
1063 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1064 if (mesh_.status(*v_it).tagged2())
1065 mesh_.property(area_, *v_it) = pow(diffgeo.compute_area(*v_it), 2);
1066
1067
1068
1069 // smooth
1070 for (int iters=0; iters<10; ++iters)
1071 {
1072 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1073 {
1074 if (mesh_.status(*v_it).tagged())
1075 {
1076 u.vectorize(0.0);
1077 ww = 0;
1078
1079 for (vv_it=mesh_.cvv_iter(*v_it); vv_it.is_valid(); ++vv_it)
1080 {
1081 w = mesh_.property(area_, *vv_it);
1082 u += mesh_.point(*vv_it) * w;
1083 ww += w;
1084 }
1085
1086 if (ww)
1087 {
1088 u *= (1.0/ww);
1089 u -= mesh_.point(*v_it);
1090 n = mesh_.normal(*v_it);
1091 n *= (u | n);
1092 u -= n;
1093
1094 u *= 0.1;
1095 }
1096
1097 mesh_.property(update_, *v_it) = u;
1098 }
1099 }
1100
1101 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1102 if (!mesh_.status(*v_it).locked() &&
1103 !mesh_.status(*v_it).feature() &&
1104 !mesh_.is_boundary(*v_it))
1105 mesh_.point(*v_it) += mesh_.property(update_, *v_it);
1106 }
1107
1108
1109 // project
1110 if (_use_projection)
1111 for (v_it=mesh_.vertices_begin(); v_it!=v_end; ++v_it)
1112 if (!mesh_.status(*v_it).locked() && !mesh_.status(*v_it).feature())
1113 project_to_reference(*v_it);
1114
1115
1116// if (emit_progress_)
1117// if (!Progress().step(3))
1118// break;
1119 }
1120}
1121
1122
1123//-----------------------------------------------------------------------------
1124
1125
1126template <class Mesh>
1127void
1128BaseRemesherT<Mesh>::
1129remove_caps()
1130{
1131 typename Mesh::EdgeIter e_it, e_end(mesh_.edges_end());
1132 typename Mesh::HalfedgeHandle h;
1133 typename Mesh::Scalar a0, a1, amin, aa(cos(170.0 * M_PI / 180.0));
1134 typename Mesh::VHandle v, vb, vd;
1135 typename Mesh::FHandle fb, fd;
1136 typename Mesh::Point a, b, c, d;
1137
1138
1139 // handle Nastran PIDs for edge flips
1140 OpenMesh::FPropHandleT<unsigned int> pids;
1141 mesh_.get_property_handle(pids, "Nastran PIDs");
1142
1143
1144 for (e_it=mesh_.edges_begin(); e_it!=e_end; ++e_it)
1145 {
1146 if (!mesh_.status(*e_it).locked() &&
1147 mesh_.is_flip_ok(*e_it))
1148 {
1149 h = mesh_.halfedge_handle(*e_it, 0);
1150 a = mesh_.point(mesh_.to_vertex_handle(h));
1151
1152 h = mesh_.next_halfedge_handle(h);
1153 b = mesh_.point(vb=mesh_.to_vertex_handle(h));
1154
1155 h = mesh_.halfedge_handle(*e_it, 1);
1156 c = mesh_.point(mesh_.to_vertex_handle(h));
1157
1158 h = mesh_.next_halfedge_handle(h);
1159 d = mesh_.point(vd=mesh_.to_vertex_handle(h));
1160
1161 a0 = ( (a-b).normalize() | (c-b).normalize() );
1162 a1 = ( (a-d).normalize() | (c-d).normalize() );
1163
1164 if (a0 < a1) { amin = a0; v = vb; }
1165 else { amin = a1; v = vd; }
1166
1167 // is it a cap?
1168 if (amin < aa)
1169 {
1170 // feature edge and feature vertex -> seems to be intended
1171 if (mesh_.status(*e_it).feature() && mesh_.status(v).feature())
1172 continue;
1173
1174 // handle PIDs: flip = split + collapse
1175 if (pids.is_valid())
1176 {
1177 fb = mesh_.face_handle(mesh_.halfedge_handle(*e_it, 0));
1178 fd = mesh_.face_handle(mesh_.halfedge_handle(*e_it, 1));
1179
1180 if (v == vb)
1181 mesh_.property(pids, fb) = mesh_.property(pids, fd);
1182 else
1183 mesh_.property(pids, fd) = mesh_.property(pids, fb);
1184 }
1185
1186 // project v onto feature edge
1187 if (mesh_.status(*e_it).feature())
1188 mesh_.set_point(v, (a+c)*0.5);
1189
1190 // flip
1191 mesh_.flip(*e_it);
1192 }
1193 }
1194 }
1195}
1196
1197
1198template<class Mesh>
1199bool
1200BaseRemesherT<Mesh>::
1201edge_flip_flips_normal(EdgeHandle _eh)
1202{
1203 if (mesh_.is_boundary(_eh))
1204 return true;
1205
1206 auto heh = mesh_.halfedge_handle(_eh, 0);
1207
1208 // get the four points of the two incident faces in ccw order
1209 auto p0 = mesh_.point(mesh_.to_vertex_handle(heh));
1210 auto p1 = mesh_.point(mesh_.to_vertex_handle(mesh_.next_halfedge_handle(heh)));
1211 auto p2 = mesh_.point(mesh_.from_vertex_handle(heh));
1212 auto p3 = mesh_.point(mesh_.to_vertex_handle(mesh_.next_halfedge_handle(mesh_.opposite_halfedge_handle(heh))));
1213
1214 // compute normals before flip
1215 auto n0_before = (p1-p0) % (p2-p0);
1216 auto n1_before = (p2-p0) % (p3-p0);
1217
1218 // compute normals after flip
1219 auto n0_after = (p1-p0) % (p3-p0);
1220 auto n1_after = (p3-p2) % (p1-p2);
1221
1222 // compare all pairs of before and after normals
1223 if ((n0_before | n0_after) < 0 || (n1_before | n1_after) < 0 ||
1224 (n0_before | n1_after) < 0 || (n1_before | n0_after) < 0)
1225 return true;
1226
1227 return false;
1228}
1229
1230
1231template<class Mesh>
1232bool
1233BaseRemesherT<Mesh>::
1234collapse_flips_normal(HalfedgeHandle _heh)
1235{
1236 if (!mesh_.is_collapse_ok(_heh))
1237 return true;
1238
1239 // get faces that are removed by the collapse
1240 auto fh0 = mesh_.face_handle(_heh);
1241 auto fh1 = mesh_.face_handle(mesh_.opposite_halfedge_handle(_heh));
1242
1243 auto collapsing_vertex = mesh_.from_vertex_handle(_heh);
1244 auto point_before = mesh_.point(collapsing_vertex);
1245
1246 // compute normals before collapse
1247 std::vector<decltype (mesh_.calc_face_normal(OpenMesh::FaceHandle(0)))> normals_before;
1248 for (auto fh : mesh_.vf_range(collapsing_vertex))
1249 normals_before.push_back(mesh_.calc_face_normal(fh));
1250
1251 // move point
1252 mesh_.point(collapsing_vertex) = mesh_.point(mesh_.to_vertex_handle(_heh));
1253
1254 bool collapse_ok = true;
1255
1256 int i = 0;
1257 for (auto fh : mesh_.vf_range(collapsing_vertex))
1258 {
1259 if (fh != fh0 && fh != fh1) // we don't care about faces that are removec by the collapse
1260 {
1261 auto normal_after = mesh_.calc_face_normal(fh);
1262 if ((normal_after | normals_before[i]) <= 0)
1263 collapse_ok = false;
1264 }
1265 ++i;
1266 }
1267
1268 // move point back
1269 mesh_.point(collapsing_vertex) = point_before;
1270
1271 return !collapse_ok;
1272}
1273
1274
1275//=============================================================================
1276} // namespace Remeshing
1277//=============================================================================
MeshSelectionT.hh
Functions for selection on a mesh.
MeshSelection::growVertexSelection
void growVertexSelection(MeshT *_mesh)
Grow vertex selection.
Definition MeshSelectionT_impl.hh:181
MeshSelection::clearVertexSelection
void clearVertexSelection(MeshT *_mesh)
Set all vertices to unselected.
Definition MeshSelectionT_impl.hh:114
OpenMesh::BaseHandle::is_valid
bool is_valid() const
The handle is valid iff the index is not negative.
Definition Handles.hh:72
OpenMesh::PolyMeshT::VertexHandle
Kernel::VertexHandle VertexHandle
Handle for referencing the corresponding item.
Definition PolyMeshT.hh:136
OpenMesh::PolyMeshT::Scalar
Kernel::Scalar Scalar
Scalar type.
Definition PolyMeshT.hh:110
OpenMesh::PolyMeshT::Normal
Kernel::Normal Normal
Normal type.
Definition PolyMeshT.hh:114
OpenMesh::PolyMeshT::update_face_normals
void update_face_normals()
Update normal vectors for all faces.
Definition PolyMeshT_impl.hh:335
OpenMesh::PolyMeshT::update_normals
void update_normals()
Compute normals for all primitives.
Definition PolyMeshT_impl.hh:317
OpenMesh::PolyMeshT::ConstVertexFaceIter
Kernel::ConstVertexFaceIter ConstVertexFaceIter
Circulator.
Definition PolyMeshT.hh:176
OpenMesh::PolyMeshT::add_vertex
SmartVertexHandle add_vertex(const Point _p)
Definition PolyMeshT.hh:255
OpenMesh::PolyMeshT::FaceHandle
Kernel::FaceHandle FaceHandle
Scalar type.
Definition PolyMeshT.hh:139
OpenMesh::PolyMeshT::HalfedgeHandle
Kernel::HalfedgeHandle HalfedgeHandle
Scalar type.
Definition PolyMeshT.hh:137
OpenMesh::PolyMeshT::EdgeIter
Kernel::EdgeIter EdgeIter
Scalar type.
Definition PolyMeshT.hh:145
OpenMesh::PolyMeshT::Point
Kernel::Point Point
Coordinate type.
Definition PolyMeshT.hh:112
OpenMesh::TriMeshT::calc_face_normal
Normal calc_face_normal(FaceHandle _fh) const
Definition TriMeshT_impl.hh:74
OpenMesh::TriMeshT::split
SmartVertexHandle split(EdgeHandle _eh, const Point &_p)
Edge split (= 2-to-4 split)
Definition TriMeshT.hh:275
OpenMesh::TriMesh_ArrayKernelT< CustomTraits >
ACG::Geometry::baryCoord
bool baryCoord(const VectorT< Scalar, 3 > &_p, const VectorT< Scalar, 3 > &_u, const VectorT< Scalar, 3 > &_v, const VectorT< Scalar, 3 > &_w, VectorT< Scalar, 3 > &_result)
Definition Algorithms.cc:83
ACG::Geometry::distPointTriangle
Vec::value_type distPointTriangle(const Vec &_p, const Vec &_v0, const Vec &_v1, const Vec &_v2, Vec &_nearestPoint)
distance from point _p to triangle (_v0, _v1, _v2)
Definition Algorithms.hh:326
OpenMesh::FaceHandle
Handle for a face entity.
Definition Handles.hh:142