File: Linear_cell_complex_3_subdivision.cpp

package info (click to toggle)
cgal 6.1.1-2
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 144,952 kB
  • sloc: cpp: 811,597; ansic: 208,576; sh: 493; python: 411; makefile: 286; javascript: 174
file content (193 lines) | stat: -rw-r--r-- 5,548 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
// Copyright (c) 2011 CNRS and LIRIS' Establishments (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org)
//
// $URL$
// $Id$
// SPDX-License-Identifier: LGPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s)     : Guillaume Damiand <guillaume.damiand@liris.cnrs.fr>
//
#include "typedefs.h"

#define PI 3.1415926535897932

// Smoth a vertex depending on the vertices of its incident facet.
class Smooth_old_vertex
{
public:
  /**  Constructor.
   * @param alcc is the lcc to smooth
   * @param amark is a mark designing old darts (i.e. darts not created during
   *        the triangulation step)
   */
  Smooth_old_vertex (LCC & alcc, LCC::size_type /* TODO amark*/):mlcc (alcc)
  {
  }

  std::pair<Point_3, Dart_descriptor> operator  () (Vertex & v) const
  {
    Dart_descriptor d = v.dart ();

    int degree = 0;
    bool open = false;

    LCC::One_dart_per_incident_cell_range<1,0>::iterator it (mlcc, d),
      itend(mlcc.one_dart_per_incident_cell<1,0>(d).end());
    for (; it != itend; ++it)
    {
      ++degree;
      if (mlcc.is_free(it, 2)) open = true;
    }

    if (open)
    { return std::make_pair(v.point(), d); }

    LCC::FT alpha = (4.0f - 2.0f *
                     (LCC::FT) cos (2.0f * PI / (LCC::FT) degree)) / 9.0f;
    LCC::Vector vec =
      LCC::Traits::Construct_scaled_vector()
      ( LCC::Traits::Construct_vector() (CGAL::ORIGIN, v.point()), (1.0f - alpha));

    for (it.rewind (); it != itend; ++it)
    {
      CGAL_assertion (!mlcc.is_free(it,2));
      vec = vec + (mlcc.point(mlcc.beta(it,2)) - CGAL::ORIGIN)
        * alpha / degree;
    }

    std::pair<Point_3, Dart_descriptor> res=std::make_pair
      (LCC::Traits::Construct_translated_point() (CGAL::ORIGIN, vec), d);

    return res;
  }
private:
  LCC & mlcc;
};

// Flip an edge, work only in 2D and 3D
Dart_descriptor
flip_edge (LCC & m, Dart_descriptor d)
{
  CGAL_assertion ( !m.is_free(d,2) );
  CGAL_assertion ( !m.is_free(d,1) && !m.is_free(d,0) );
  CGAL_assertion ( !m.is_free(m.beta(d,2), 0) && !m.is_free(m.beta(d, 2), 1) );

  if (!m.is_removable<1>(d)) return LCC::null_descriptor;

  Dart_descriptor d1 = m.beta(d,1);
  Dart_descriptor d2 = m.beta(d,2,0);

  CGAL_assertion ( !m.is_free(d1,1) && !m.is_free(d2,0) );

  Dart_descriptor d3 = m.beta(d1,1);
  Dart_descriptor d4 = m.beta(d2, 0);

  // We isolated the edge
  m.basic_link_beta_1(m.beta(d,0), m.beta(d,2,1));
  m.basic_link_beta_0(m.beta(d,1), m.beta(d,2,0));
  if ( !m.is_free(d,3) )
  {
    m.basic_link_beta_0(m.beta(d,0,3), m.beta(d,2,1,3));
    m.basic_link_beta_1(m.beta(d,1,3), m.beta(d,2,0,3));
  }

  // Then we push the two extremities.
  m.basic_link_beta_0(d3, d);
  m.basic_link_beta_0(d2, m.beta(d,2));
  m.link_beta_1(d4, d);
  m.link_beta_1(d1, m.beta(d,2));

  if ( !m.is_free(d,3) )
  {
    m.basic_link_beta_0(m.beta(d4,3), m.beta(d,3));
    m.basic_link_beta_0(m.beta(d1,3), m.beta(d,2,3));
    m.link_beta_1(m.beta(d3,3), m.beta(d,3));
    m.link_beta_1(m.beta(d2,3), m.beta(d,2,3));
  }

  // CGAL::remove_cell<LCC,1>(m, d);
  // insert_cell_1_in_cell_2(m, d1, d1->beta(1)->beta(1));

  return d;
}

// Subdivide each facet of the lcc by using sqrt(3)-subdivision.
void
subdivide_lcc_3 (LCC & m)
{
  if (m.number_of_darts () == 0)
    return;

  LCC::size_type mark = m.get_new_mark ();
  LCC::size_type treated = m.get_new_mark ();
  m.negate_mark (mark);  // All the old darts are marked in O(1).

  // 1) We smoth the old vertices.
  std::vector <std::pair<Point_3, Dart_descriptor> > vertices;  // smooth the old vertices
  vertices.reserve (m.number_of_attributes<0> ());  // get intermediate space
  std::transform (m.vertex_attributes().begin (),
                  m.vertex_attributes().end (),
                  std::back_inserter (vertices),
                  Smooth_old_vertex (m, mark));

  // 2) We subdivide each facet.
  m.negate_mark (treated);  // All the darts are marked in O(1).

  for (LCC::Dart_range::iterator it (m.darts().begin ());
       m.number_of_marked_darts (treated) > 0; ++it)
  {
    if (m.is_marked (it, treated))
    {
      // We unmark the darts of the facet to process only once dart/facet.
      CGAL::unmark_cell < LCC, 2 > (m, it, treated);
      // We triangulate the facet.
      m.insert_barycenter_in_cell<2>(it);
    }
  }

  CGAL_assertion (m.is_whole_map_unmarked (treated));
  CGAL_assertion (m.is_valid ());
  m.free_mark (treated);

  // 3) We update the coordinates of old vertices.
  for (std::vector<std::pair<Point_3, Dart_descriptor> >::iterator
         vit=vertices.begin(); vit!=vertices.end(); ++vit)
  {
    m.point(vit->second)=vit->first;
  }

  // 4) We flip all the old edges.
  m.negate_mark (mark);  // Now only new darts are marked.
  Dart_descriptor d2 =LCC::null_descriptor;
  for (LCC::Dart_range::iterator it (m.darts().begin ());
       it != m.darts().end ();)
  {
    d2 = it++;
    if (!m.is_marked (d2, mark))  // This is an old dart.
    {
      // We process only the last dart of a same edge.
      if (!m.is_free(d2,2) && (m.beta(d2,2,3)==m.beta(d2,3,2)))
      {
        if (m.is_marked(m.beta(d2,2), mark) &&
            (m.is_free(d2,3) ||
             (m.is_marked(m.beta(d2,3), mark) &&
              m.is_marked(m.beta(d2,2,3), mark))))
        {
          flip_edge (m, d2);
          m.mark(d2, mark);
        }
        else
          m.mark (d2, mark);
      }
      else
        m.mark (d2, mark);
    }
  }

  CGAL_assertion (m.is_whole_map_marked (mark));
  m.free_mark (mark);

  CGAL_postcondition ( m.is_valid ());
}