File: simplify_polyline.h

package info (click to toggle)
cgal 6.1-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 144,912 kB
  • sloc: cpp: 810,858; ansic: 208,477; sh: 493; python: 411; makefile: 286; javascript: 174
file content (189 lines) | stat: -rw-r--r-- 6,184 bytes parent folder | download
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
// Copyright (c) 2020 GeometryFactory (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL: https://github.com/CGAL/cgal/blob/v6.1/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/simplify_polyline.h $
// $Id: include/CGAL/Polygon_mesh_processing/internal/simplify_polyline.h b26b07a1242 $
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
//
// Author(s)     : Sébastien Loriot

#ifndef CGAL_POLYGON_MESH_PROCESSING_SIMPLIFY_POLYLINE_H
#define CGAL_POLYGON_MESH_PROCESSING_SIMPLIFY_POLYLINE_H

#include <CGAL/license/Polygon_mesh_processing/geometric_repair.h>

#include <CGAL/Named_function_parameters.h>
#include <CGAL/boost/graph/named_params_helper.h>

#include <type_traits>


namespace CGAL {
namespace Polygon_mesh_processing {
namespace experimental {

enum Polyline_simplification_algorithms { DOUGLAS_PEUCKER, ITERATIVE };

template <typename PointRangeIn, typename PointRangeOut,
          typename NamedParametersIn = parameters::Default_named_parameters,
          typename NamedParametersOut = parameters::Default_named_parameters>
void simplify_polyline(const PointRangeIn& input,
                             PointRangeOut& output,
                       const double max_squared_frechet_distance,
                       const NamedParametersIn& np_in = parameters::default_values(),
                       const NamedParametersOut& np_out = parameters::default_values())
{
  using parameters::choose_parameter;
  using parameters::get_parameter;

  static_assert(std::is_same<
    typename std::iterator_traits<typename PointRangeIn::const_iterator>::value_type,
    typename std::iterator_traits<typename PointRangeOut::const_iterator>::value_type >::value, "");

  typedef typename GetPointMap<PointRangeIn, NamedParametersIn>::type Point_map_in;
  typedef typename GetPointMap<PointRangeOut, NamedParametersOut>::type Point_map_out;
  typedef typename Point_set_processing_3_np_helper<PointRangeIn, NamedParametersIn>::Geom_traits Kernel;

  Point_map_in in_pm = choose_parameter<Point_map_in>(get_parameter(np_in, internal_np::point_map));
  Point_map_out out_pm = choose_parameter<Point_map_out>(get_parameter(np_out, internal_np::point_map));

  const Polyline_simplification_algorithms algorithm =
    choose_parameter(get_parameter(np_in, internal_np::algorithm), DOUGLAS_PEUCKER);

  switch(algorithm)
  {
    case ITERATIVE:
    {
      const bool is_closed = input.front()==input.back();

      std::size_t nb_points = is_closed ? input.size()-1 : input.size();

      // skip points in the input range that do not contains any information
      if (nb_points<=2)
      {
        output.reserve(input.size());
        for (const auto& p : input)
        {
          output.push_back(p);
          put(out_pm, output.back(), get(in_pm, p));
        }
        return;
      }

      auto is_valid_approx = [&input, &in_pm, max_squared_frechet_distance](
        std::size_t b, std::size_t e,
        const typename Kernel::Line_3& line)
      {
        typename Kernel::Compare_squared_distance_3 compare_squared_distance;
        for (std::size_t i=b+1; i<e; ++i)
        {

          if (compare_squared_distance(get(in_pm, input[i]), line, max_squared_frechet_distance) == LARGER)
            return false;
        }
        return true;
      };

      std::size_t bi=0;
      while(bi!=nb_points)
      {
        std::size_t ei=bi+2;

        output.push_back(input[bi]);
        put(out_pm, output.back(), get(in_pm, input[bi]));

        while(ei<nb_points)
        {
          typename Kernel::Line_3 sl(get(in_pm,input[bi]), get(in_pm,input[ei]));

          if (is_valid_approx(bi,ei,sl))
            ++ei; // we skip ei-1
          else
          {
            bi=ei-1; // ei-1 shall not be skipped
            break;
          }
        }
        if(ei>=nb_points) break;
      }
      output.push_back(input[nb_points-1]);
      put(out_pm, output.back(), get(in_pm, input[nb_points-1]));
      if (is_closed)
      {
        output.push_back(input.back());
        put(out_pm, output.back(), get(in_pm, input.back()));
      }
      return;
    }
    case DOUGLAS_PEUCKER:
    {
      const bool is_closed = input.front()==input.back();

      std::size_t nb_points = is_closed ? input.size()-1 : input.size();

      if (nb_points<=2)
      {
        output.reserve(input.size());
        for (const auto& p : input)
        {
          output.push_back(p);
          put(out_pm, output.back(), get(in_pm, p));
        }
        return;
      }

      std::vector< std::pair<std::size_t, std::size_t> > ranges;
      ranges.push_back(std::make_pair(0, nb_points-1));

      std::vector<bool> kept(input.size(), false);
      if (is_closed) kept[nb_points]=true;
      while( !ranges.empty() )
      {
        std::size_t rb, re;
        std::tie(rb, re) = ranges.back();
        ranges.pop_back();
        kept[rb]=true;
        kept[re]=true;
        if (rb+1==re) continue;

        typename Kernel::Line_3 line(get(in_pm, input[rb]), get(in_pm, input[re]));
        double max_d = max_squared_frechet_distance;
        std::size_t max_i = 0;
        for (std::size_t i=rb; i<re; ++i)
        {
          double d = squared_distance(line, get(in_pm, input[i]));
          if (d > max_d)
          {
            max_d = d;
            max_i = i;
          }
        }
        if (max_i != 0)
        {
          ranges.push_back( std::make_pair(max_i, re) );
          ranges.push_back( std::make_pair(rb, max_i) );
        }
      }
      std::size_t nb_kept=0;
      for (std::size_t i=0; i<input.size(); ++i)
        if (kept[i]) ++nb_kept;

      output.reserve(nb_kept);
      for (std::size_t i=0; i<input.size(); ++i)
        if (kept[i])
        {
          output.push_back(input[i]);
          put(out_pm, output.back(), get(in_pm, input[i]));
        }
      //TODO if is_closed-==true, shall we add en extra step to see if we can remove output.front() and output[output.size()-2] (initial endpoints)
    }
  }
}

} } } // end of CGAL::Polygon_mesh_processing::experimental namespace


#endif // CGAL_POLYGON_MESH_PROCESSING_SIMPLIFY_POLYLINE_H