# A.I.M.S algorithms

geodesic_algorithm_dijkstra_alternative.h
Go to the documentation of this file.
2 #ifndef GEODESIC_ALGORITHM_DIJKSTRA_ALTERNATIVE_010506
3 #define GEODESIC_ALGORITHM_DIJKSTRA_ALTERNATIVE_010506
4
7 #include <vector>
8 #include <set>
9 #include <assert.h>
10
11 namespace geodesic{
12
14 {
15  typedef DijkstraNode1* node_pointer;
16 public:
19
20  double& distance_from_source(){return m_distance;};
21  node_pointer& previous(){return m_previous;};
22  unsigned& source_index(){return m_source_index;};
23  vertex_pointer& vertex(){return m_vertex;};
24
25  void clear()
26  {
27  m_distance = GEODESIC_INF;
28  m_previous = NULL;
29  }
30
31  bool operator()(node_pointer const s1, node_pointer const s2) const
32  {
33  return s1->distance_from_source() != s2->distance_from_source() ?
35  s1->vertex()->id() < s2->vertex()->id();
36  };
37
38 private:
39  double m_distance; //distance to the closest source
40  unsigned m_source_index; //closest source index
41  node_pointer m_previous; //previous node in the geodesic path
42  vertex_pointer m_vertex; //correspoding vertex
43 };
44
46 {
47 public:
49  typedef Node* node_pointer;
50
53  m_nodes(mesh->vertices().size())
54  {
55  m_type = DIJKSTRA;
56  for(unsigned i=0; i<m_nodes.size(); ++i)
57  {
58  m_nodes[i].vertex() = &m_mesh->vertices()[i];
59  }
60  };
61
63
64  virtual void propagate(std::vector<SurfacePoint>& sources,
65  double max_propagation_distance = GEODESIC_INF, //propagation algorithm stops after reaching the certain distance from the source
66  std::vector<SurfacePoint>* stop_points = NULL); //or after ensuring that all the stop_points are covered
67
68  virtual void trace_back(SurfacePoint& destination, //trace back piecewise-linear path
69  std::vector<SurfacePoint>& path);
70
71  virtual unsigned best_source(SurfacePoint& point, //quickly find what source this point belongs to and what is the distance to this source
72  double& best_source_distance);
73 private:
74
75  void set_sources(std::vector<SurfacePoint>& sources)
76  {
77  m_sources = sources;
78  }
79
80  bool check_stop_conditions(unsigned& index);
81
82  std::vector<Node> m_nodes; //nodes of the subdivision graph located on the vertices
83
84  std::vector<SurfacePoint> m_sources; //for simplicity, we keep sources as they are
85
86  typedef std::set<node_pointer, Node> queue_type;
87  queue_type m_queue;
88 };
89
90 inline unsigned GeodesicAlgorithmDijkstraAlternative::best_source(SurfacePoint& point, //quickly find what source this point belongs to and what is the distance to this source
91  double& best_source_distance)
92 {
93  if(point.type() == VERTEX)
94  {
96  node_pointer node = &m_nodes[v->id()];
97  best_source_distance = node->distance_from_source();
98  return node->source_index();
99  }
100  else
101  {
102  std::vector<vertex_pointer> possible_vertices; //initialize vertices directly visible from sources
103  m_mesh->closest_vertices(&point, &possible_vertices);
104
105  best_source_distance = GEODESIC_INF;
106  vertex_pointer min_vertex = NULL;
107  for(unsigned i=0; i<possible_vertices.size(); ++i)
108  {
109  vertex_pointer v = possible_vertices[i];
110
111  double distance_from_source = m_nodes[v->id()].distance_from_source();
112  double distance_from_dest = v->distance(&point);
113  if(distance_from_source + distance_from_dest < best_source_distance)
114  {
115  best_source_distance = distance_from_source + distance_from_dest;
116  min_vertex = v;
117  }
118  }
119  assert(min_vertex);
120  node_pointer node = &m_nodes[min_vertex->id()];
121  return node->source_index();
122  }
123 }
124
125 inline void GeodesicAlgorithmDijkstraAlternative::trace_back(SurfacePoint& destination, //trace back piecewise-linear path
126  std::vector<SurfacePoint>& path)
127 {
128  path.clear();
129  path.push_back(destination);
130
131  if(destination.type() != VERTEX) //snap to the closest vertex
132  {
133  std::vector<vertex_pointer> possible_vertices; //initialize vertices directly visible from sources
134  m_mesh->closest_vertices(&destination, &possible_vertices);
135
136  double min_distance = GEODESIC_INF;
137  vertex_pointer min_vertex = NULL;
138  for(unsigned i=0; i<possible_vertices.size(); ++i)
139  {
140  vertex_pointer v = possible_vertices[i];
141
142  double distance_from_source = m_nodes[v->id()].distance_from_source();
143  double distance_from_dest = v->distance(&destination);
144  if(distance_from_source + distance_from_dest < min_distance)
145  {
146  min_distance = distance_from_source + distance_from_dest;
147  min_vertex = v;
148  }
149  }
150  assert(min_vertex);
151  path.push_back(SurfacePoint(min_vertex));
152  }
153
154  node_pointer node = &m_nodes[path.back().base_element()->id()];
156  {
157  node = node->previous();
158  path.push_back(SurfacePoint(node->vertex()));
159  }
160
161  SurfacePoint& source = m_sources[node->source_index()]; //add source to the path if it is not already there
162  if(source.type() != VERTEX ||
163  source.base_element()->id() != node->vertex()->id())
164  {
165  path.push_back(source);
166  }
167 }
168
169 inline void GeodesicAlgorithmDijkstraAlternative::propagate(std::vector<SurfacePoint>& sources,
170  double max_propagation_distance, //propagation algorithm stops after reaching the certain distance from the source
171  std::vector<SurfacePoint>* stop_points) //or after ensuring that all the stop_points are covered
172 {
173  set_stop_conditions(stop_points, max_propagation_distance);
174  set_sources(sources);
175
176  m_queue.clear(); //clear everything
177  for(unsigned i=0; i<m_nodes.size(); ++i)
178  {
179  m_nodes[i].clear();
180  }
181
182  clock_t start = clock();
183
184  std::vector<vertex_pointer> visible_vertices; //initialize vertices directly visible from sources
185  for(unsigned i=0; i<m_sources.size(); ++i)
186  {
187  SurfacePoint* s = &m_sources[i];
188  m_mesh->closest_vertices(s, &visible_vertices);
189  for(unsigned j=0; j<visible_vertices.size(); ++j)
190  {
191  vertex_pointer v = visible_vertices[j];
192  double distance = s->distance(v);
193
194  node_pointer node = &m_nodes[v->id()];
195  if(distance < node->distance_from_source())
196  {
197  node->distance_from_source() = distance;
198  node->source_index() = i;
199  node->previous() = NULL;
200  }
201  }
202  visible_vertices.clear();
203  }
204
205  for(unsigned i=0; i<m_nodes.size(); ++i) //initialize the queue
206  {
207  if(m_nodes[i].distance_from_source() < GEODESIC_INF)
208  {
209  m_queue.insert(&m_nodes[i]);
210  }
211  }
212
213  unsigned counter = 0;
214  unsigned satisfied_index = 0;
215
216  while(!m_queue.empty()) //main cycle
217  {
218  if(counter++ % 10 == 0) //check if we covered all required vertices
219  {
220  if (check_stop_conditions(satisfied_index))
221  {
222  break;
223  }
224  //std::cout << counter << " " << m_queue.size() << " " << (*m_queue.begin())->distance_from_source()<< std::endl;
225  }
226
227  node_pointer min_node = *m_queue.begin();
228  m_queue.erase(m_queue.begin());
229
230  vertex_pointer v = min_node->vertex();
231
233  {
235  vertex_pointer next_v = e->opposite_vertex(v);
236  node_pointer next_node = &m_nodes[next_v->id()];
237  //double current_distance =
238  if(next_node->distance_from_source() > min_node->distance_from_source() + e->length())
239  {
240  if(next_node->distance_from_source() < GEODESIC_INF) //remove it from the queue
241  {
242  queue_type::iterator iter = m_queue.find(next_node);
243  assert(iter != m_queue.end());
244  m_queue.erase(iter);
245  }
246  next_node->distance_from_source() = min_node->distance_from_source() + e->length();
247  next_node->source_index() = min_node->source_index();
248  next_node->previous() = min_node;
249  m_queue.insert(next_node);
250  }
251  }
252  }
253  //std::cout << std::endl;
254
255  clock_t finish = clock();
256  m_time_consumed = (static_cast<double>(finish)-static_cast<double>(start))/CLOCKS_PER_SEC;
257 }
258
259 inline bool GeodesicAlgorithmDijkstraAlternative::check_stop_conditions(unsigned& index)
260 {
261  double queue_min_distance = (*m_queue.begin())->distance_from_source();
262
263  if(queue_min_distance < m_max_propagation_distance)
264  {
265  return false;
266  }
267
268  while(index < m_stop_vertices.size())
269  {
270  vertex_pointer v = m_stop_vertices[index].first;
271  Node& node = m_nodes[v->id()];
272  if(queue_min_distance < node.distance_from_source() + m_stop_vertices[index].second)
273  {
274  return false;
275  }
276  ++index;
277  }
278  return true;
279 }
280
281 } //geodesic
282
283 #endif //GEODESIC_ALGORITHM_DIJKSTRA_ALTERNATIVE_010506
virtual unsigned best_source(SurfacePoint &point, double &best_source_distance)