GeographicLib  2.6
Intersect.hpp
Go to the documentation of this file.
1 /**
2  * \file Intersect.hpp
3  * \brief Header for GeographicLib::Intersect class
4  *
5  * Copyright (c) Charles Karney (2023-2024) <karney@alum.mit.edu> and licensed
6  * under the MIT/X11 License. For more information, see
7  * https://geographiclib.sourceforge.io/
8  **********************************************************************/
9 
10 #if !defined(GEOGRAPHICLIB_INTERSECT_HPP)
11 #define GEOGRAPHICLIB_INTERSECT_HPP 1
12 
13 #include <vector>
14 #include <set>
15 #include <GeographicLib/Math.hpp>
16 #include <GeographicLib/Geodesic.hpp>
17 #include <GeographicLib/GeodesicLine.hpp>
18 
19 namespace GeographicLib {
20 
21  /**
22  * \brief %Geodesic intersections
23  *
24  * Find the intersections of two geodesics \e X and \e Y. Four calling
25  * sequences are supported.
26  * - The geodesics are defined by a position (latitude and longitude) and an
27  * azimuth. In this case the \e closest intersection is found.
28  * - The geodesics are defined by two endpoints. The intersection of the two
29  * segments is found. If they don't intersect, then the closest
30  * intersection is returned.
31  * - The geodesics are defined as an intersection point, a single position
32  * and two azimuths. In this case, the next closest intersection is found.
33  * - The geodesics are defined as in the first case and all intersection
34  * within a specified distance are returned.
35  * .
36  * In all cases the position of the intersection is given by the signed
37  * displacements \e x and \e y along the geodesics from the starting point
38  * (the first point in the case of a geodesic segment). The closest
39  * itersection is defined as the one that minimizes the L1 distance,
40  * Intersect::Dist([<i>x</i>, <i>y</i>) = |<i>x</i>| + |<i>y</i>|.
41  *
42  * The routines also optionally return a coincidence indicator \e c. This is
43  * typically 0. However if the geodesics lie on top of one another at the
44  * point of intersection, then \e c is set to +1, if they are parallel, and
45  * &minus;1, if they are antiparallel.
46  *
47  * Example of use:
48  * \include example-Intersect.cpp
49  *
50  * <a href="IntersectTool.1.html">IntersectTool</a> is a command-line utility
51  * providing access to the functionality of this class.
52  *
53  * This solution for intersections is described in
54  * - C. F. F. Karney,<br>
55  * <a href="https://doi.org/10.1061/JSUED2.SUENG-1483">
56  * Geodesic intersections</a>,
57  * J. Surveying Eng. <b>150</b>(3), 04024005:1--9 (2024);
58  * preprint
59  * <a href="https://arxiv.org/abs/2308.00495">arxiv:2308.00495</a>.
60  * .
61  * It is based on the work of
62  * - S. Baseldga and J. C. Martinez-Llario,
63  * <a href="https://doi.org/10.1007/s11200-017-1020-z">
64  * Intersection and point-to-line solutions for geodesics
65  * on the ellipsoid</a>,
66  * Stud. Geophys. Geod. <b>62</b>, 353--363 (2018);
67  * DOI: <a href="https://doi.org/10.1007/s11200-017-1020-z">
68  * 10.1007/s11200-017-1020-z</a>.
69  **********************************************************************/
70 
71  class GEOGRAPHICLIB_EXPORT Intersect {
72  private:
73  typedef Math::real real;
74  public:
75  /**
76  * The type used to hold the two displacement along the geodesics. This is
77  * just a std::pair with \e x = \e first and \e y = \e second.
78  **********************************************************************/
79  typedef std::pair<Math::real, Math::real> Point;
80  /**
81  * The minimum capabilities for GeodesicLine objects which are passed to
82  * this class.
83  **********************************************************************/
84  static const unsigned LineCaps = Geodesic::LATITUDE | Geodesic::LONGITUDE |
85  Geodesic::AZIMUTH | Geodesic::REDUCEDLENGTH | Geodesic::GEODESICSCALE |
86  Geodesic::DISTANCE_IN;
87  private:
88  static const int numit_ = 100;
89  const Geodesic _geod;
90  real _a, _f, // equatorial radius, flattening
91  _rR, // authalic radius
92  _d, // pi*_rR
93  _eps, // criterion for intersection + coincidence
94  _tol, // convergence for Newton in Solve1
95  _delta, // for equality tests, safety margin for tiling
96  _t1, // min distance between intersections
97  _t2, // furthest dist to closest intersection
98  _t3, // 1/2 furthest min dist to next intersection
99  _t4, // capture radius for spherical sol in Solve0
100  _t5, // longest shortest geodesic
101  _d1, // tile spacing for Closest
102  _d2, // tile spacing for Next
103  _d3; // tile spacing for All
104  // The L1 distance
105  static Math::real d1(Math::real x, Math::real y)
106  { using std::fabs; return fabs(x) + fabs(y); }
107  // An internal version of Point with a little more functionality
108  class XPoint {
109  public:
110  real x, y;
111  int c;
112  XPoint(Math::real x, Math::real y, int c = 0)
113  : x(x), y(y), c(c)
114  {}
115  XPoint()
116  : x(Math::NaN()), y(Math::NaN()), c(0)
117  {}
118  XPoint(const Point& p)
119  : x(p.first), y(p.second), c(0)
120  {}
121  XPoint& operator+=(const XPoint& p) {
122  x += p.x; y += p.y;
123  if (p.c) c = p.c; // pass along a nonzero c from either operand
124  return *this;
125  }
126  XPoint operator+(const XPoint& p) const {
127  XPoint t = *this; t += p; return t;
128  }
129  Math::real Dist() const { return d1(x, y); }
130  Math::real Dist(const XPoint& p) const { return d1(x - p.x, y - p.y); }
131  Point data() const { return Point(x, y); }
132  };
133  // Comparing XPoints for insertions into sets, but ensure that close
134  // XPoints test equal.
135  class GEOGRAPHICLIB_EXPORT SetComp {
136  private:
137  const real _delta;
138  public:
139  SetComp(Math::real delta) : _delta(delta) {}
140  bool eq(const XPoint& p, const XPoint& q) const {
141  return d1(p.x - q.x, p.y - q.y) <= _delta;
142  }
143  bool operator()(const XPoint& p, const XPoint& q) const {
144  return !eq(p, q) && ( p.x != q.x ? p.x < q.x : p.y < q.y );
145  }
146  };
147  SetComp _comp;
148  // For ranking XPoints by closeness
149  class RankPoint {
150  private:
151  const real _x, _y;
152  public:
153  RankPoint(const Point& p0) : _x(p0.first), _y(p0.second) {}
154  RankPoint(const XPoint& p0) : _x(p0.x), _y(p0.y) {}
155  bool operator()(const XPoint& p, const XPoint& q) const {
156  real dp = d1(p.x - _x, p.y - _y),
157  dq = d1(q.x - _x, q.y - _y);
158  return dp != dq ? (dp < dq) :
159  (p.x != q.x ? (p.x < q.x) : (p.y < q.y));
160  }
161  };
162  // The spherical solution
163  XPoint Spherical(const GeodesicLine& lineX, const GeodesicLine& lineY,
164  const XPoint& p) const;
165  // The basic algorithm
166  XPoint Basic(const GeodesicLine& lineX, const GeodesicLine& lineY,
167  const XPoint& p0) const;
168  // The closest intersecton
169  XPoint ClosestInt(const GeodesicLine& lineX, const GeodesicLine& lineY,
170  const XPoint& p0) const;
171  // The next intersecton
172  XPoint NextInt(const GeodesicLine& lineX, const GeodesicLine& lineY) const;
173  // Segment intersecton
174  XPoint SegmentInt(const GeodesicLine& lineX, const GeodesicLine& lineY,
175  int& segmode) const;
176  // All intersectons
177  std::vector<XPoint>
178  AllInt0(const GeodesicLine& lineX, const GeodesicLine& lineY,
179  Math::real maxdist, const XPoint& p0) const;
180  std::vector<Point>
181  AllInternal(const GeodesicLine& lineX, const GeodesicLine& lineY,
182  Math::real maxdist, const Point& p0,
183  std::vector<int>& c, bool cp) const;
184  // Find {semi-,}conjugate point which is close to s3. Optional m12, M12,
185  // M21 use {semi-,}conjugacy relative to point 2
186  Math::real ConjugateDist(const GeodesicLine& line, Math::real s3, bool semi,
187  Math::real m12 = 0, Math::real M12 = 1,
188  Math::real M21 = 1) const;
189  Math::real distpolar(Math::real lat1, Math::real* lat2 = nullptr) const;
190  Math::real polarb(Math::real* lata = nullptr, Math::real* latb = nullptr)
191  const;
192  Math::real conjdist(Math::real azi, Math::real* ds = nullptr,
193  Math::real* sp = nullptr, Math::real* sm = nullptr)
194  const;
195  Math::real distoblique(Math::real* azi = nullptr, Math::real* sp = nullptr,
196  Math::real* sm = nullptr) const;
197  // p is intersection point on coincident lines orientation = c; p0 is
198  // origin point. Change p to center point wrt p0, i.e, abs((p-p0)_x) =
199  // abs((p-p0)_y)
200  static XPoint fixcoincident(const XPoint& p0, const XPoint& p);
201  static XPoint fixcoincident(const XPoint& p0, const XPoint& p, int c);
202  static XPoint fixsegment(Math::real sx, Math::real sy, const XPoint& p);
203  static int segmentmode(Math::real sx, Math::real sy, const XPoint& p) {
204  return (p.x < 0 ? -1 : p.x <= sx ? 0 : 1) * 3
205  + (p.y < 0 ? -1 : p.y <= sy ? 0 : 1);
206  }
207  mutable long long _cnt0, _cnt1, _cnt2, _cnt3, _cnt4;
208  public:
209  /** \name Constructor
210  **********************************************************************/
211  ///@{
212  /**
213  * Constructor for an ellipsoid with
214  *
215  * @param[in] geod a Geodesic object. This sets the parameters \e a and \e
216  * f for the ellipsoid.
217  * @exception GeographicErr if the eccentricity of the elliposdoid is too
218  * large.
219  *
220  * \note This class has been validated for -1/4 &le; \e f &le; 1/5. It may
221  * give satisfactory results slightly outside this range; however
222  * sufficient far outside the range, some internal checks will fail and an
223  * exception thrown.
224  *
225  * \note If |<i>f</i>| > 1/50, then the Geodesic object should be
226  * constructed with \e exact = true.
227  **********************************************************************/
228  Intersect(const Geodesic& geod);
229  ///@}
230 
231  /** \name Finding intersections
232  **********************************************************************/
233  ///@{
234  /**
235  * Find the closest intersection point, with each geodesic specified by
236  * position and azimuth.
237  *
238  * @param[in] latX latitude of starting point for geodesic \e X (degrees).
239  * @param[in] lonX longitude of starting point for geodesic \e X (degrees).
240  * @param[in] aziX azimuth at starting point for geodesic \e X (degrees).
241  * @param[in] latY latitude of starting point for geodesic \e Y (degrees).
242  * @param[in] lonY longitude of starting point for geodesic \e Y (degrees).
243  * @param[in] aziY azimuth at starting point for geodesic \e Y (degrees).
244  * @param[in] p0 an optional offset for the starting points (meters),
245  * default = [0,0].
246  * @param[out] c optional pointer to an integer coincidence indicator.
247  * @return \e p the intersection point closest to \e p0.
248  *
249  * The returned intersection minimizes Intersect::Dist(\e p, \e p0).
250  **********************************************************************/
251  Point Closest(Math::real latX, Math::real lonX, Math::real aziX,
252  Math::real latY, Math::real lonY, Math::real aziY,
253  const Point& p0 = Point(0, 0), int* c = nullptr) const;
254  /**
255  * Find the closest intersection point, with each geodesic given as a
256  * GeodesicLine.
257  *
258  * @param[in] lineX geodesic \e X.
259  * @param[in] lineY geodesic \e Y.
260  * @param[in] p0 an optional offset for the starting points (meters),
261  * default = [0,0].
262  * @param[out] c optional pointer to an integer coincidence indicator.
263  * @return \e p the intersection point closest to \e p0.
264  *
265  * The returned intersection minimizes Intersect::Dist(\e p, \e p0).
266  *
267  * \note \e lineX and \e lineY should be created with minimum capabilities
268  * Intersect::LineCaps. The methods for creating a GeodesicLine include
269  * all these capabilities by default.
270  **********************************************************************/
271  Point Closest(const GeodesicLine& lineX, const GeodesicLine& lineY,
272  const Point& p0 = Point(0, 0), int* c = nullptr) const;
273  /**
274  * Find the intersection of two geodesic segments defined by their
275  * endpoints.
276  *
277  * @param[in] latX1 latitude of first point for segment \e X (degrees).
278  * @param[in] lonX1 longitude of first point for segment \e X (degrees).
279  * @param[in] latX2 latitude of second point for segment \e X (degrees).
280  * @param[in] lonX2 longitude of second point for segment \e X (degrees).
281  * @param[in] latY1 latitude of first point for segment \e Y (degrees).
282  * @param[in] lonY1 longitude of first point for segment \e Y (degrees).
283  * @param[in] latY2 latitude of second point for segment \e Y (degrees).
284  * @param[in] lonY2 longitude of second point for segment \e Y (degrees).
285  * @param[out] segmode an indicator equal to zero if the segments
286  * intersect (see below).
287  * @param[out] c optional pointer to an integer coincidence indicator.
288  * @return \e p the intersection point if the segments intersect, otherwise
289  * the intersection point closest to the midpoints of the two
290  * segments.
291  *
292  * \warning The results are only well defined if there's a \e unique
293  * shortest geodesic between the endpoints of the two segments.
294  *
295  * \e segmode codes up information about the closest intersection in the
296  * case where the segments intersect. Let <i>x</i><sub>12</sub> be the
297  * length of the segment \e X and \e x = <i>p</i>.first, the position of
298  * the intersection on segment \e X. Define
299  * - \e k<sub><i>x</i></sub> = &minus;1, if \e x < 0,
300  * - \e k<sub><i>x</i></sub> = 0,
301  * if 0 &le; \e x &le; <i>x</i><sub>12</sub>,
302  * - \e k<sub><i>x</i></sub> = 1, if <i>x</i><sub>12</sub> < \e x.
303  * .
304  * and similarly for segment \e Y. Then
305  * \e segmode = 3 \e k<sub><i>x</i></sub> + \e k<sub><i>y</i></sub>.
306  **********************************************************************/
307  Point Segment(Math::real latX1, Math::real lonX1,
308  Math::real latX2, Math::real lonX2,
309  Math::real latY1, Math::real lonY1,
310  Math::real latY2, Math::real lonY2,
311  int& segmode, int* c = nullptr) const;
312  /**
313  * Find the intersection of two geodesic segments each defined by a
314  * GeodesicLine.
315  *
316  * @param[in] lineX segment \e X.
317  * @param[in] lineY segment \e Y.
318  * @param[out] segmode an indicator equal to zero if the segments
319  * intersect (see below).
320  * @param[out] c optional pointer to an integer coincidence indicator.
321  * @return \e p the intersection point if the segments intersect, otherwise
322  * the intersection point closest to the midpoints of the two
323  * segments.
324  *
325  * \warning \e lineX and \e lineY must represent shortest geodesics, e.g.,
326  * they can be created by Geodesic::InverseLine. The results are only well
327  * defined if there's a \e unique shortest geodesic between the endpoints
328  * of the two segments.
329  *
330  * \note \e lineX and \e lineY should be created with minimum capabilities
331  * Intersect::LineCaps. The methods for creating a GeodesicLine include
332  * all these capabilities by default.
333  *
334  * See previous definition of Intersect::Segment for more information on \e
335  * segmode.
336  **********************************************************************/
337  Point Segment(const GeodesicLine& lineX, const GeodesicLine& lineY,
338  int& segmode, int* c = nullptr) const;
339  /**
340  * Find the next closest intersection point to a given intersection,
341  * specified by position and two azimuths.
342  *
343  * @param[in] latX latitude of starting points for geodesics \e X and \e Y
344  * (degrees).
345  * @param[in] lonX longitude of starting points for geodesics \e X and \e Y
346  * (degrees).
347  * @param[in] aziX azimuth at starting point for geodesic \e X (degrees).
348  * @param[in] aziY azimuth at starting point for geodesic \e Y (degrees).
349  * @param[out] c optional pointer to an integer coincidence indicator.
350  * @return \e p the next closest intersection point.
351  *
352  * The returned intersection minimizes Intersect::Dist(\e p) (excluding \e
353  * p = [0,0]).
354  *
355  * \note Equidistant closest intersections are surprisingly common. If
356  * this may be a problem, use Intersect::All with a sufficiently large \e
357  * maxdist to capture close intersections.
358  **********************************************************************/
359  Point Next(Math::real latX, Math::real lonX,
360  Math::real aziX, Math::real aziY, int* c = nullptr) const;
361  /**
362  * Find the next closest intersection point to a given intersection,
363  * with each geodesic specified a GeodesicLine.
364  *
365  * @param[in] lineX geodesic \e X.
366  * @param[in] lineY geodesic \e Y.
367  * @param[out] c optional pointer to an integer coincidence indicator.
368  * @return \e p the next closest intersection point.
369  *
370  * \warning \e lineX and \e lineY must both have the same starting point,
371  * i.e., the distance between [<i>lineX</i>.Latitude(),
372  * <i>lineX</i>.Longitude()] and [<i>lineY</i>.Latitude(),
373  * <i>lineY</i>.Longitude()] must be zero.
374  *
375  * \note \e lineX and \e lineY should be created with minimum capabilities
376  * Intersect::LineCaps. The methods for creating a GeodesicLine include
377  * all these capabilities by default.
378  *
379  * \note Equidistant closest intersections are surprisingly common. If
380  * this may be a problem, use Intersect::All with a sufficiently large \e
381  * maxdist to capture close intersections.
382  **********************************************************************/
383  Point Next(const GeodesicLine& lineX, const GeodesicLine& lineY,
384  int* c = nullptr) const;
385  ///@}
386 
387  /** \name Finding all intersections
388  **********************************************************************/
389  ///@{
390  /**
391  * Find all intersections within a certain distance, with each geodesic
392  * specified by position and azimuth.
393  *
394  * @param[in] latX latitude of starting point for geodesic \e X (degrees).
395  * @param[in] lonX longitude of starting point for geodesic \e X (degrees).
396  * @param[in] aziX azimuth at starting point for geodesic \e X (degrees).
397  * @param[in] latY latitude of starting point for geodesic \e Y (degrees).
398  * @param[in] lonY longitude of starting point for geodesic \e Y (degrees).
399  * @param[in] aziY azimuth at starting point for geodesic \e Y (degrees).
400  * @param[in] maxdist the maximum distance for the returned intersections
401  * (meters).
402  * @param[out] c vector of coincidences.
403  * @param[in] p0 an optional offset for the starting points (meters),
404  * default = [0,0].
405  * @return \e plist a vector for the intersections closest to \e p0.
406  *
407  * Each intersection point satisfies Intersect::Dist(\e p, \e p0) &le; \e
408  * maxdist. The vector of returned intersections is sorted on the distance
409  * from \e p0.
410  **********************************************************************/
411  std::vector<Point> All(Math::real latX, Math::real lonX, Math::real aziX,
412  Math::real latY, Math::real lonY, Math::real aziY,
413  Math::real maxdist, std::vector<int>& c,
414  const Point& p0 = Point(0, 0))
415  const;
416  /**
417  * Find all intersections within a certain distance, with each geodesic
418  * specified by position and azimuth. Don't return vector of
419  * coincidences.
420  *
421  * @param[in] latX latitude of starting point for geodesic \e X (degrees).
422  * @param[in] lonX longitude of starting point for geodesic \e X (degrees).
423  * @param[in] aziX azimuth at starting point for geodesic \e X (degrees).
424  * @param[in] latY latitude of starting point for geodesic \e Y (degrees).
425  * @param[in] lonY longitude of starting point for geodesic \e Y (degrees).
426  * @param[in] aziY azimuth at starting point for geodesic \e Y (degrees).
427  * @param[in] maxdist the maximum distance for the returned intersections
428  * (meters).
429  * @param[in] p0 an optional offset for the starting points (meters),
430  * default = [0,0].
431  * @return \e plist a vector for the intersections closest to \e p0.
432  *
433  * Each intersection point satisfies Intersect::Dist(\e p, \e p0) &le; \e
434  * maxdist. The vector of returned intersections is sorted on the distance
435  * from \e p0.
436  **********************************************************************/
437  std::vector<Point> All(Math::real latX, Math::real lonX, Math::real aziX,
438  Math::real latY, Math::real lonY, Math::real aziY,
439  Math::real maxdist, const Point& p0 = Point(0, 0))
440  const;
441  /**
442  * Find all intersections within a certain distance, with each geodesic
443  * specified by a GeodesicLine.
444  *
445  * @param[in] lineX geodesic \e X.
446  * @param[in] lineY geodesic \e Y.
447  * @param[in] maxdist the maximum distance for the returned intersections
448  * (meters).
449  * @param[out] c vector of coincidences.
450  * @param[in] p0 an optional offset for the starting points (meters),
451  * default = [0,0].
452  * @return \e plist a vector for the intersections closest to \e p0.
453  *
454  * Each intersection point satisfies Intersect::Dist(\e p, \e p0) &le; \e
455  * maxdist. The vector of returned intersections is sorted on the distance
456  * from \e p0.
457  *
458  * \note \e lineX and \e lineY should be created with minimum capabilities
459  * Intersect::LineCaps. The methods for creating a GeodesicLine include
460  * all these capabilities by default.
461  **********************************************************************/
462  std::vector<Point> All(const GeodesicLine& lineX, const GeodesicLine& lineY,
463  Math::real maxdist, std::vector<int>& c,
464  const Point& p0 = Point(0, 0))
465  const;
466  /**
467  * Find all intersections within a certain distance, with each geodesic
468  * specified by a GeodesicLine. Don't return vector or coincidences.
469  *
470  * @param[in] lineX geodesic \e X.
471  * @param[in] lineY geodesic \e Y.
472  * @param[in] maxdist the maximum distance for the returned intersections
473  * (meters).
474  * @param[in] p0 an optional offset for the starting points (meters),
475  * default = [0,0].
476  * @return \e plist a vector for the intersections closest to \e p0.
477  *
478  * Each intersection point satisfies Intersect::Dist(\e p, \e p0) &le; \e
479  * maxdist. The vector of returned intersections is sorted on the distance
480  * from \e p0.
481  *
482  * \note \e lineX and \e lineY should be created with minimum capabilities
483  * Intersect::LineCaps. The methods for creating a GeodesicLine include
484  * all these capabilities by default.
485  **********************************************************************/
486  std::vector<Point> All(const GeodesicLine& lineX, const GeodesicLine& lineY,
487  Math::real maxdist, const Point& p0 = Point(0, 0))
488  const;
489  ///@}
490 
491  /** \name Diagnostic counters
492  **********************************************************************/
493  ///@{
494  /**
495  * @return the cumulative number of invocations of **h**.
496  *
497  * This is a count of the number of times the spherical triangle needs to
498  * be solved. Each involves a call to Geodesic::Inverse and this is a good
499  * metric for the overall cost. This counter is set to zero by the
500  * constructor.
501  *
502  * \warning The counter is a mutable variable and so is not thread safe.
503  **********************************************************************/
504  long long NumInverse() const { return _cnt0; }
505  /**
506  * @return the cumulative number of invocations of **b**.
507  *
508  * This is a count of the number of invocations of the basic algorithm,
509  * which is used by all the intersection methods. This counter is set to
510  * zero by the constructor.
511  *
512  * \warning The counter is a mutable variable and so is not thread safe.
513  **********************************************************************/
514  long long NumBasic() const { return _cnt1; }
515  /**
516  * @return the number of times intersection point was changed in
517  * Intersect::Closest and Intersect::Next.
518  *
519  * If this counter is incremented by just 1 in Intersect::Closest, then the
520  * initial result of the basic algorithm was eventually accepted. This
521  * counter is set to zero by the constructor.
522  *
523  * \note This counter is also incremented by Intersect::Segment, which
524  * calls Intersect::Closest.
525  *
526  * \warning The counter is a mutable variable and so is not thread safe.
527  **********************************************************************/
528  long long NumChange() const { return _cnt2; }
529  /**
530  * @return the number of times a corner point is checked in
531  * Intersect::Segment.
532  *
533  * This counter is set to zero by the constructor.
534  *
535  * \warning The counter is a mutable variable and so is not thread safe.
536  **********************************************************************/
537  long long NumCorner() const { return _cnt3; }
538  /**
539  * @return the number of times a corner point is returned by
540  * Intersect::Segment.
541  *
542  * This counter is set to zero by the constructor.
543  *
544  * \note A conjecture is that a corner point never results in an
545  * intersection that overrides the intersection closest to the midpoints of
546  * the segments; i.e., NumCorner() always returns 0.
547  *
548  * \warning The counter is a mutable variable and so is not thread safe.
549  **********************************************************************/
550  long long NumOverride() const { return _cnt4; }
551  ///@}
552 
553  /** \name Insepctor function
554  **********************************************************************/
555  ///@{
556  /**
557  * @return \e geod the Geodesic object used in the constructor.
558  *
559  * This can be used to query Geodesic::EquatorialRadius(),
560  * Geodesic::Flattening(), Geodesic::Exact(), and
561  * Geodesic::EllipsoidArea().
562  **********************************************************************/
563  const Geodesic& GeodesicObject() const { return _geod; }
564  ///@}
565 
566  /**
567  * The L1 distance.
568  *
569  * @param[in] p the position along geodesics \e X and \e Y.
570  * @param[in] p0 [optional] the reference position, default = [0, 0].
571 
572  * @return the L1 distance of \e p from \e p0, i.e.,
573  * |<i>p</i><sub><i>x</i></sub> &minus; <i>p0</i><sub><i>x</i></sub>| +
574  * |<i>p</i><sub><i>y</i></sub> &minus; <i>p0</i><sub><i>y</i></sub>|.
575  **********************************************************************/
576  static Math::real Dist(const Point& p, const Point& p0 = Point(0, 0)) {
577  using std::fabs;
578  return fabs(p.first - p0.first) + fabs(p.second - p0.second);
579  }
580  };
581 
582 } // namespace GeographicLib
583 
584 #endif // GEOGRAPHICLIB_INTERSECT_HPP
GeodesicLine.hpp
Header for GeographicLib::GeodesicLine class.
GeographicLib::Intersect::GeodesicObject
const Geodesic & GeodesicObject() const
Definition: Intersect.hpp:563
GEOGRAPHICLIB_EXPORT
#define GEOGRAPHICLIB_EXPORT
Definition: Constants.hpp:59
GeographicLib::Intersect
Geodesic intersections
Definition: Intersect.hpp:71
GeographicLib::Intersect::NumInverse
long long NumInverse() const
Definition: Intersect.hpp:504
GeographicLib::Intersect::Point
std::pair< Math::real, Math::real > Point
Definition: Intersect.hpp:79
GeographicLib::Intersect::Dist
static Math::real Dist(const Point &p, const Point &p0=Point(0, 0))
Definition: Intersect.hpp:576
Geodesic.hpp
Header for GeographicLib::Geodesic class.
Math.hpp
Header for GeographicLib::Math class.
GeographicLib::Intersect::NumOverride
long long NumOverride() const
Definition: Intersect.hpp:550
GeographicLib
Namespace for GeographicLib.
Definition: Accumulator.cpp:12
GeographicLib::Intersect::NumCorner
long long NumCorner() const
Definition: Intersect.hpp:537
real
GeographicLib::Math::real real
Definition: Geod3Solve.cpp:25
GeographicLib::Intersect::NumChange
long long NumChange() const
Definition: Intersect.hpp:528
GeographicLib::Intersect::NumBasic
long long NumBasic() const
Definition: Intersect.hpp:514
GeographicLib::Geodesic
Geodesic calculations
Definition: Geodesic.hpp:175