$darkmode
VCG Library
geodesic.h
1 /****************************************************************************
2 * VCGLib o o *
3 * Visual and Computer Graphics Library o o *
4 * _ O _ *
5 * Copyright(C) 2004-2016 \/)\/ *
6 * Visual Computing Lab /\/| *
7 * ISTI - Italian National Research Council | *
8 * \ *
9 * All rights reserved. *
10 * *
11 * This program is free software; you can redistribute it and/or modify *
12 * it under the terms of the GNU General Public License as published by *
13 * the Free Software Foundation; either version 2 of the License, or *
14 * (at your option) any later version. *
15 * *
16 * This program is distributed in the hope that it will be useful, *
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
19 * GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
20 * for more details. *
21 * *
22 ****************************************************************************/
23 #include <vcg/simplex/face/pos.h>
24 #include <vcg/simplex/face/topology.h>
25 #include <vcg/complex/algorithms/update/quality.h>
26 #include <deque>
27 #include <functional>
28 #ifndef __VCGLIB_GEODESIC
29 #define __VCGLIB_GEODESIC
30 
31 namespace vcg{
32 namespace tri{
33 
34 template <class MeshType>
36  typedef typename MeshType::VertexType VertexType;
37  typedef typename MeshType::ScalarType ScalarType;
38  typedef typename MeshType::FacePointer FacePointer;
39 
41 
42  ScalarType operator()(const VertexType * v0, const VertexType * v1) const
43  {return vcg::Distance(v0->cP(),v1->cP());}
44 
45  ScalarType operator()(const FacePointer f0, const FacePointer f1) const
46  {return vcg::Distance(Barycenter(*f0),Barycenter(*f1));}
47 };
48 
49 
50 template <class MeshType>
52 private:
53  // The only member of this class. The attribute handle used to
54  typename MeshType::template PerVertexAttributeHandle<float> wH;
55 public:
56  typedef typename MeshType::VertexType VertexType;
57  typedef typename MeshType::ScalarType ScalarType;
58  typedef typename MeshType::FacePointer FacePointer;
59 
65 
66  IsotropicDistance(MeshType &m, float variance)
67  {
68  // the wH attribute store the scaling factor to be applied to the distance
69  wH = tri::Allocator<MeshType>:: template GetPerVertexAttribute<float> (m,"distW");
70  float qmin = 1.0f/variance;
71  float qmax = variance;
72  float qrange = qmax-qmin;
73  std::pair<float,float> minmax = Stat<MeshType>::ComputePerVertexQualityMinMax(m);
74  float range = minmax.second-minmax.first;
75  for(size_t i=0;i<m.vert.size();++i)
76  wH[i]=qmin+((m.vert[i].Q()-minmax.first)/range)*qrange;
77 
78 // qDebug("Range %f %f %f",minmax.first,minmax.second,range);
79  }
80 
81  ScalarType operator()( VertexType * v0, VertexType * v1)
82  {
83  float scale = (wH[v0]+wH[v1])/2.0f;
84  return (1.0f/scale)*vcg::Distance(v0->cP(),v1->cP());
85  }
86 };
87 
88 
89 template <class MeshType>
91 {
92  BasicCrossFunctor(MeshType &m) { tri::RequirePerVertexCurvatureDir(m); }
93  typedef typename MeshType::VertexType VertexType;
94 
95  typename MeshType::CoordType D1(VertexType &v) { return v.PD1(); }
96  typename MeshType::CoordType D2(VertexType &v) { return v.PD2(); }
97 };
98 
107 template <class MeshType>
109  typedef typename MeshType::VertexType VertexType;
110  typedef typename MeshType::ScalarType ScalarType;
111  typedef typename MeshType::CoordType CoordType;
112  typedef typename MeshType::VertexIterator VertexIterator;
113 
114  typename MeshType::template PerVertexAttributeHandle<CoordType> wxH,wyH;
115 public:
116  template <class CrossFunctor > AnisotropicDistance(MeshType &m, CrossFunctor &cf)
117  {
118  wxH = tri::Allocator<MeshType>:: template GetPerVertexAttribute<CoordType> (m,"distDirX");
119  wyH = tri::Allocator<MeshType>:: template GetPerVertexAttribute<CoordType> (m,"distDirY");
120 
121  for(VertexIterator vi=m.vert.begin();vi!=m.vert.end();++vi)
122  {
123  wxH[vi]=cf.D1(*vi);
124  wyH[vi]=cf.D2(*vi);
125  }
126  }
127 
128  ScalarType operator()( VertexType * v0, VertexType * v1)
129  {
130  CoordType dd = v0->cP()-v1->cP();
131  float x = (fabs(dd * wxH[v0])+fabs(dd *wxH[v1]))/2.0f;
132  float y = (fabs(dd * wyH[v0])+fabs(dd *wyH[v1]))/2.0f;
133 
134  return sqrt(x*x+y*y);
135  }
136 };
137 
138 
139 
140 
147 template <class MeshType>
148 class Geodesic{
149 
150 public:
151 
152  typedef typename MeshType::VertexType VertexType;
153  typedef typename MeshType::VertexIterator VertexIterator;
154  typedef typename MeshType::VertexPointer VertexPointer;
155  typedef typename MeshType::FacePointer FacePointer;
156  typedef typename MeshType::FaceType FaceType;
157  typedef typename MeshType::CoordType CoordType;
158  typedef typename MeshType::ScalarType ScalarType;
159 
160 
161 
162 /* Auxiliary class for keeping the heap of vertices to visit and their estimated distance */
163  struct VertDist{
164  VertDist(){}
165  VertDist(VertexPointer _v, ScalarType _d):v(_v),d(_d){}
166 
167  VertexPointer v;
168  ScalarType d;
169  };
170 
171 
172  struct DIJKDist{
173  DIJKDist(VertexPointer _v):v(_v), q(_v->Q()){}
174  VertexPointer v;
175  ScalarType q;
176 
177  bool operator < (const DIJKDist &o) const
178  {
179  if( q != o.q)
180  return q > o.q;
181  return v<o.v;
182  }
183  };
184 
185  /* Auxiliary class for keeping the heap of vertices to visit and their estimated distance */
186  struct FaceDist{
187  FaceDist(FacePointer _f):f(_f){}
188  FacePointer f;
189  bool operator < (const FaceDist &o) const
190  {
191  if( f->Q() != o.f->Q())
192  return f->Q() > o.f->Q();
193  return f<o.f;
194  }
195  };
196 
197  /* Temporary data to associate to all the vertices: estimated distance and boolean flag */
198  struct TempData{
199  TempData(){}
200  TempData(const ScalarType & _d):d(_d),source(0),parent(0){}
201 
202  ScalarType d;
203  VertexPointer source;//closest source
204  VertexPointer parent;
205  };
206 
207  typedef SimpleTempData<std::vector<VertexType>, TempData > TempDataType;
208 
209 
210  struct pred {
211  pred() {};
212  bool operator()(const VertDist& v0, const VertDist& v1) const {
213  return (v0.d > v1.d);
214  }
215  };
216 
217  /*
218  *
219  curr: vertex for which distance should be estimated
220  d_pw1: distance of pw1 from the source
221  d_curr: distance of curr from the source
222 
223 The function estimates the distance of pw from the source
224 in the assumption the mesh is developable (and without holes)
225 along the path, so that (source,pw1,curr) from a triangle.
226 All the math is to comput the angles at pw1 and curr with the Erone formula.
227 
228 The if cases take care of the cases where the angles are obtuse.
229 
230  curr
231  d_pw1 +
232  | +pw
233 source+ |
234  d_curr +
235  pw1
236 
237  */
238  template <class DistanceFunctor>
239  static ScalarType Distance(DistanceFunctor &distFunc,
240  const VertexPointer &pw,
241  const VertexPointer &pw1,
242  const VertexPointer &curr,
243  const ScalarType &d_pw1,
244  const ScalarType &d_curr)
245  {
246  ScalarType curr_d=0;
247 
248  ScalarType ew_c = distFunc(pw,curr);
249  ScalarType ew_w1 = distFunc(pw,pw1);
250  ScalarType ec_w1 = distFunc(pw1,curr);
251  CoordType w_c = (pw->cP()-curr->cP()).Normalize() * ew_c;
252  CoordType w_w1 = (pw->cP() - pw1->cP()).Normalize() * ew_w1;
253  CoordType w1_c = (pw1->cP() - curr->cP()).Normalize() * ec_w1;
254 
255  ScalarType alpha,alpha_, beta,beta_,theta,h,delta,s,a,b;
256 
257  alpha = acos((w_c.dot(w1_c))/(ew_c*ec_w1));
258  s = (d_curr + d_pw1+ec_w1)/2;
259  a = s/ec_w1;
260  b = a*s;
261  alpha_ = 2*acos ( std::min<ScalarType>(1.0,sqrt( (b- a* d_pw1)/d_curr)));
262 
263  if ( alpha+alpha_ > M_PI){
264  curr_d = d_curr + ew_c;
265  }else
266  {
267  beta_ = 2*acos ( std::min<ScalarType>(1.0,sqrt( (b- a* d_curr)/d_pw1)));
268  beta = acos((w_w1).dot(-w1_c)/(ew_w1*ec_w1));
269 
270  if ( beta+beta_ > M_PI)
271  curr_d = d_pw1 + ew_w1;
272  else
273  {
274  theta = ScalarType(M_PI)-alpha-alpha_;
275  delta = cos(theta)* ew_c;
276  h = sin(theta)* ew_c;
277  curr_d = sqrt( pow(h,2)+ pow(d_curr + delta,2));
278  }
279  }
280  return (curr_d);
281  }
282 
283 
284 
285 
286 /*
287 This is the low level version of the geodesic computation framework.
288 Starting from the seeds, it assign a distance value to each vertex. The distance of a vertex is its
289 approximated geodesic distance to the closest seeds.
290 This is function is not meant to be called (although is not prevented). Instead, it is invoked by
291 wrapping function.
292 */
293 
294  template <class DistanceFunctor>
295  static VertexPointer Visit(
296  MeshType & m,
297  std::vector<VertDist> & seedVec, // the set of seeds to start from
298  DistanceFunctor &distFunc,
299  ScalarType distance_threshold = std::numeric_limits<ScalarType>::max(), // cut off distance (do no compute anything farther than this value)
300  typename MeshType::template PerVertexAttributeHandle<VertexPointer> * vertSource = NULL, // if present we put in this attribute the closest source for each vertex
301  typename MeshType::template PerVertexAttributeHandle<VertexPointer> * vertParent = NULL, // if present we put in this attribute the parent in the path that goes from the vertex to the closest source
302  std::vector<VertexPointer> *InInterval=NULL)
303  {
304  VertexPointer farthest=0;
305 // int t0=clock();
306  //Requirements
307  tri::RequireVFAdjacency(m);
308  tri::RequirePerVertexQuality(m);
309 
310  assert(!seedVec.empty());
311 
312  TempDataType TD(m.vert, std::numeric_limits<ScalarType>::max());
313 
314  // initialize Heap
315  std::vector<VertDist> frontierHeap;
316  typename std::vector <VertDist >::iterator ifr;
317  for(ifr = seedVec.begin(); ifr != seedVec.end(); ++ifr){
318  TD[(*ifr).v].d = (*ifr).d;
319  TD[(*ifr).v].source = (*ifr).v;
320  TD[(*ifr).v].parent = (*ifr).v;
321  frontierHeap.push_back(*ifr);
322  }
323  make_heap(frontierHeap.begin(),frontierHeap.end(),pred());
324 
325  ScalarType curr_d,d_curr = 0.0,d_heap;
326  ScalarType max_distance=0.0;
327 // int t1=clock();
328  while(!frontierHeap.empty() && max_distance < distance_threshold)
329  {
330  pop_heap(frontierHeap.begin(),frontierHeap.end(),pred());
331  VertexPointer curr = (frontierHeap.back()).v;
332  if (InInterval!=NULL) InInterval->push_back(curr);
333 
334  if(vertSource!=NULL) (*vertSource)[curr] = TD[curr].source;
335  if(vertParent!=NULL) (*vertParent)[curr] = TD[curr].parent;
336 
337  d_heap = (frontierHeap.back()).d;
338  frontierHeap.pop_back();
339 
340  assert(TD[curr].d <= d_heap);
341  if(TD[curr].d < d_heap ) // a vertex whose distance has been improved after it was inserted in the queue
342  continue;
343  assert(TD[curr].d == d_heap);
344 
345  d_curr = TD[curr].d;
346 
347  for(face::VFIterator<FaceType> vfi(curr) ; vfi.f!=0; ++vfi )
348  {
349  for(int k=0;k<2;++k)
350  {
351  VertexPointer pw,pw1;
352  if(k==0) {
353  pw = vfi.f->V1(vfi.z);
354  pw1= vfi.f->V2(vfi.z);
355  }
356  else {
357  pw = vfi.f->V2(vfi.z);
358  pw1= vfi.f->V1(vfi.z);
359  }
360 
361  const ScalarType & d_pw1 = TD[pw1].d;
362  {
363  const ScalarType inter = distFunc(curr,pw1);//(curr->P() - pw1->P()).Norm();
364  const ScalarType tol = (inter + d_curr + d_pw1)*.0001f;
365 
366  if ( (TD[pw1].source != TD[curr].source)||// not the same source
367  (inter + d_curr < d_pw1 +tol ) ||
368  (inter + d_pw1 < d_curr +tol ) ||
369  (d_curr + d_pw1 < inter +tol ) // triangular inequality
370  )
371  curr_d = d_curr + distFunc(pw,curr);//(pw->P()-curr->P()).Norm();
372  else
373  curr_d = Distance(distFunc,pw,pw1,curr,d_pw1,d_curr);
374  }
375 
376  if(TD[pw].d > curr_d){
377  TD[pw].d = curr_d;
378  TD[pw].source = TD[curr].source;
379  TD[pw].parent = curr;
380  frontierHeap.push_back(VertDist(pw,curr_d));
381  push_heap(frontierHeap.begin(),frontierHeap.end(),pred());
382  }
383 // if(isLeaf){
384  if(d_curr > max_distance){
385  max_distance = d_curr;
386  farthest = curr;
387  }
388 // }
389  }
390  } // end for VFIterator
391  }// end while
392 // int t2=clock();
393 
394  // Copy found distance onto the Quality (\todo parametric!)
395  if (InInterval==NULL)
396  {
397  for(VertexIterator vi = m.vert.begin(); vi != m.vert.end(); ++vi) if(!(*vi).IsD())
398  (*vi).Q() = TD[&(*vi)].d;
399  }
400  else
401  {
402  assert(InInterval->size()>0);
403  for(size_t i=0;i<InInterval->size();i++)
404  (*InInterval)[i]->Q() = TD[(*InInterval)[i]].d;
405  }
406 // int t3=clock();
407 // printf("Init %6.3f\nVisit %6.3f\nFinal %6.3f\n",float(t1-t0)/CLOCKS_PER_SEC,float(t2-t1)/CLOCKS_PER_SEC,float(t3-t2)/CLOCKS_PER_SEC);
408  return farthest;
409  }
410 
411 public:
443  static bool Compute( MeshType & m,
444  const std::vector<VertexPointer> & seedVec)
445  {
447  return Compute(m,seedVec,dd);
448  }
449 
450  template <class DistanceFunctor>
451  static bool Compute( MeshType & m,
452  const std::vector<VertexPointer> & seedVec,
453  DistanceFunctor &distFunc,
454  ScalarType maxDistanceThr = std::numeric_limits<ScalarType>::max(),
455  std::vector<VertexPointer> *withinDistanceVec=NULL,
456  typename MeshType::template PerVertexAttributeHandle<VertexPointer> * sourceSeed = NULL,
457  typename MeshType::template PerVertexAttributeHandle<VertexPointer> * parentSeed = NULL
458  )
459  {
460  if(seedVec.empty()) return false;
461  std::vector<VertDist> vdSeedVec;
462  typename std::vector<VertexPointer>::const_iterator fi;
463  for( fi = seedVec.begin(); fi != seedVec.end() ; ++fi)
464  vdSeedVec.push_back(VertDist(*fi,0.0));
465  Visit(m, vdSeedVec, distFunc, maxDistanceThr, sourceSeed, parentSeed, withinDistanceVec);
466  return true;
467  }
468 
469  /* \brief Assigns to each vertex of the mesh its distance to the closest vertex on the boundary
470 
471 It is just a simple wrapper of the basic Compute()
472 
473  Note: update the field Q() of the vertices
474  Note: it needs the border bit set.
475  */
476  static bool DistanceFromBorder( MeshType & m, typename MeshType::template PerVertexAttributeHandle<VertexPointer> * sources = NULL)
477  {
478  std::vector<VertexPointer> fro;
479  for(VertexIterator vi = m.vert.begin(); vi != m.vert.end(); ++vi)
480  if( (*vi).IsB())
481  fro.push_back(&(*vi));
482  if(fro.empty()) return false;
483  EuclideanDistance<MeshType> dd;
485  return Compute(m,fro,dd,std::numeric_limits<ScalarType>::max(),0,sources);
486  }
487 
488 
489  static bool ConvertPerVertexSeedToPerFaceSeed(MeshType &m, const std::vector<VertexPointer> &vertexSeedVec,
490  std::vector<FacePointer> &faceSeedVec)
491  {
492  tri::RequireVFAdjacency(m);
493  tri::RequirePerFaceMark(m);
494 
495  faceSeedVec.clear();
496  tri::UnMarkAll(m);
497  for(size_t i=0;i<vertexSeedVec.size();++i)
498  {
499  for(face::VFIterator<FaceType> vfi(vertexSeedVec[i]);!vfi.End();++vfi)
500  {
501  if(tri::IsMarked(m,vfi.F())) return false;
502  faceSeedVec.push_back(vfi.F());
503  tri::Mark(m,vfi.F());
504  }
505  }
506  return true;
507  }
508 
509  static inline std::string sourcesAttributeName(void) { return "sources"; }
510  static inline std::string parentsAttributeName(void) { return "parent"; }
511 
512  template <class DistanceFunctor>
513  static void PerFaceDijkstraCompute(MeshType &m, const std::vector<FacePointer> &seedVec,
514  DistanceFunctor &distFunc,
515  ScalarType maxDistanceThr = std::numeric_limits<ScalarType>::max(),
516  std::vector<FacePointer> *InInterval=NULL,
517  FacePointer FaceTarget=NULL,
518  bool avoid_selected=false)
519  {
520  tri::RequireFFAdjacency(m);
521  tri::RequirePerFaceMark(m);
522  tri::RequirePerFaceQuality(m);
523 
524  typename MeshType::template PerFaceAttributeHandle<FacePointer> sourceHandle
525  = tri::Allocator<MeshType>::template GetPerFaceAttribute<FacePointer>(m, sourcesAttributeName());
526 
527  typename MeshType::template PerFaceAttributeHandle<FacePointer> parentHandle
528  = tri::Allocator<MeshType>::template GetPerFaceAttribute<FacePointer>(m, parentsAttributeName());
529 
530  std::vector<FaceDist> Heap;
531  tri::UnMarkAll(m);
532  for(size_t i=0;i<seedVec.size();++i)
533  {
534  tri::Mark(m,seedVec[i]);
535  seedVec[i]->Q()=0;
536  sourceHandle[seedVec[i]]=seedVec[i];
537  parentHandle[seedVec[i]]=seedVec[i];
538  Heap.push_back(FaceDist(seedVec[i]));
539  if (InInterval!=NULL) InInterval->push_back(seedVec[i]);
540  }
541 
542  std::make_heap(Heap.begin(),Heap.end());
543  while(!Heap.empty())
544  {
545  pop_heap(Heap.begin(),Heap.end());
546  FacePointer curr = (Heap.back()).f;
547  if ((FaceTarget!=NULL)&&(curr==FaceTarget))return;
548  Heap.pop_back();
549 
550  for(int i=0;i<3;++i)
551  {
552  if(!face::IsBorder(*curr,i) )
553  {
554  FacePointer nextF = curr->FFp(i);
555  ScalarType nextDist = curr->Q() + distFunc(curr,nextF);
556  if( (nextDist < maxDistanceThr) &&
557  (!tri::IsMarked(m,nextF) || nextDist < nextF->Q()) )
558  {
559  nextF->Q() = nextDist;
560  if ((avoid_selected)&&(nextF->IsS()))continue;
561  tri::Mark(m,nextF);
562  Heap.push_back(FaceDist(nextF));
563  push_heap(Heap.begin(),Heap.end());
564  if (InInterval!=NULL) InInterval->push_back(nextF);
565  sourceHandle[nextF] = sourceHandle[curr];
566  parentHandle[nextF] = curr;
567 // printf("Heapsize %i nextDist = %f curr face %i next face %i \n",Heap.size(), nextDist, tri::Index(m,curr), tri::Index(m,nextF));
568  }
569  }
570  }
571  }
572  }
573 
574 
575 
576 
577  template <class DistanceFunctor>
578  static void PerVertexDijkstraCompute(MeshType &m, const std::vector<VertexPointer> &seedVec,
579  DistanceFunctor &distFunc,
580  ScalarType maxDistanceThr = std::numeric_limits<ScalarType>::max(),
581  std::vector<VertexPointer> *InInterval=NULL,
582  typename MeshType::template PerVertexAttributeHandle<VertexPointer> * sourceHandle= NULL,
583  typename MeshType::template PerVertexAttributeHandle<VertexPointer> * parentHandle=NULL,
584  bool avoid_selected=false,
585  VertexPointer target=NULL)
586  {
587  tri::RequireVFAdjacency(m);
588  tri::RequirePerVertexMark(m);
589  tri::RequirePerVertexQuality(m);
590 
591  std::vector<DIJKDist> Heap;
592  tri::UnMarkAll(m);
593 
594  for(size_t i=0;i<seedVec.size();++i)
595  {
596  assert(!tri::IsMarked(m,seedVec[i]));
597  tri::Mark(m,seedVec[i]);
598  seedVec[i]->Q()=0;
599  if (sourceHandle!=NULL)
600  (*sourceHandle)[seedVec[i]]=seedVec[i];
601  if (parentHandle!=NULL)
602  (*parentHandle)[seedVec[i]]=seedVec[i];
603  Heap.push_back(DIJKDist(seedVec[i]));
604  if (InInterval!=NULL) InInterval->push_back(seedVec[i]);
605  }
606 
607  std::make_heap(Heap.begin(),Heap.end());
608  while(!Heap.empty())
609  {
610  pop_heap(Heap.begin(),Heap.end());
611  VertexPointer curr = (Heap.back()).v;
612  if ((target!=NULL)&&(target==curr))return;
613  Heap.pop_back();
614  std::vector<VertexPointer> vertVec;
615  face::VVStarVF<FaceType>(curr,vertVec);
616  for(size_t i=0;i<vertVec.size();++i)
617  {
618  VertexPointer nextV = vertVec[i];
619  if ((avoid_selected)&&(nextV->IsS()))continue;
620  ScalarType nextDist = curr->Q() + distFunc(curr,nextV);
621  if( (nextDist < maxDistanceThr) &&
622  (!tri::IsMarked(m,nextV) || nextDist < nextV->Q()) )
623  {
624  nextV->Q() = nextDist;
625  tri::Mark(m,nextV);
626  Heap.push_back(DIJKDist(nextV));
627  push_heap(Heap.begin(),Heap.end());
628  if (InInterval!=NULL) InInterval->push_back(nextV);
629  if (sourceHandle!=NULL)
630  (*sourceHandle)[nextV] = (*sourceHandle)[curr];
631  if (parentHandle!=NULL)
632  (*parentHandle)[nextV] = curr;
633 // printf("Heapsize %i nextDist = %f curr vert %i next vert %i \n",Heap.size(), nextDist, tri::Index(m,curr), tri::Index(m,nextV));
634  }
635  }
636  }
637  }
638 
639 
640 };// end class
641 }// end namespace tri
642 }// end namespace vcg
643 #endif
Class to safely add and delete elements in a mesh.
Definition: allocate.h:97
Definition: geodesic.h:108
class for computing approximate geodesic distances on a mesh
Definition: geodesic.h:148
static bool Compute(MeshType &m, const std::vector< VertexPointer > &seedVec)
Given a set of source vertices compute the approximate geodesic distance to all the other vertices.
Definition: geodesic.h:443
Definition: geodesic.h:51
IsotropicDistance(MeshType &m, float variance)
Definition: geodesic.h:66
static void VertexConstant(MeshType &m, VertexQualityType q)
Definition: quality.h:67
bool IsBorder(FaceType const &f, const int j)
Definition: topology.h:55
Definition: namespaces.dox:6
Definition: geodesic.h:91
Definition: geodesic.h:35
Definition: geodesic.h:172
Definition: geodesic.h:186
Definition: geodesic.h:198
Definition: geodesic.h:163
Definition: geodesic.h:210