Coverage for pygeodesy/frechet.py: 96%
205 statements
« prev ^ index » next coverage.py v7.2.2, created at 2023-10-04 14:05 -0400
« prev ^ index » next coverage.py v7.2.2, created at 2023-10-04 14:05 -0400
2# -*- coding: utf-8 -*-
4u'''Fréchet distances.
6Classes L{Frechet}, L{FrechetDegrees}, L{FrechetRadians},
7L{FrechetCosineAndoyerLambert}, L{FrechetCosineForsytheAndoyerLambert},
8L{FrechetCosineLaw}, L{FrechetDistanceTo}< L{FrechetEquirectangular},
9L{FrechetEuclidean}, L{FrechetExact}, L{FrechetFlatLocal}, L{FrechetFlatPolar},
10L{FrechetHaversine}, L{FrechetHubeny}, L{FrechetKarney}, L{FrechetThomas}
11and L{FrechetVincentys} to compute I{discrete} U{Fréchet
12<https://WikiPedia.org/wiki/Frechet_distance>} distances between two sets
13of C{LatLon}, C{NumPy}, C{tuples} or other types of points.
15Only L{FrechetDistanceTo} -iff used with L{ellipsoidalKarney.LatLon}
16points- and L{FrechetKarney} requires installation of I{Charles Karney}'s
17U{geographiclib<https://PyPI.org/project/geographiclib>}.
19Typical usage is as follows. First, create a C{Frechet} calculator
20from one set of C{LatLon} points.
22C{f = FrechetXyz(point1s, ...)}
24Get the I{discrete} Fréchet distance to another set of C{LatLon} points
25by
27C{t6 = f.discrete(point2s)}
29Or, use function C{frechet_} with a proper C{distance} function passed
30as keyword arguments as follows
32C{t6 = frechet_(point1s, point2s, ..., distance=...)}.
34In both cases, the returned result C{t6} is a L{Frechet6Tuple}.
36For C{(lat, lon, ...)} points in a C{NumPy} array or plain C{tuples},
37wrap the points in a L{Numpy2LatLon} respectively L{Tuple2LatLon}
38instance, more details in the documentation thereof.
40For other points, create a L{Frechet} sub-class with the appropriate
41C{distance} method overloading L{Frechet.distance} as in this example.
43 >>> from pygeodesy import Frechet, hypot_
44 >>>
45 >>> class F3D(Frechet):
46 >>> """Custom Frechet example.
47 >>> """
48 >>> def distance(self, p1, p2):
49 >>> return hypot_(p1.x - p2.x, p1.y - p2.y, p1.z - p2.z)
50 >>>
51 >>> f3D = F3D(xyz1, ..., units="...")
52 >>> t6 = f3D.discrete(xyz2)
54Transcribed from the original U{Computing Discrete Fréchet Distance
55<https://www.kr.TUWien.ac.AT/staff/eiter/et-archive/cdtr9464.pdf>} by
56Eiter, T. and Mannila, H., 1994, April 25, Technical Report CD-TR 94/64,
57Information Systems Department/Christian Doppler Laboratory for Expert
58Systems, Technical University Vienna, Austria.
60This L{Frechet.discrete} implementation optionally generates intermediate
61points for each point set separately. For example, using keyword argument
62C{fraction=0.5} adds one additional point halfway between each pair of
63points. Or using C{fraction=0.1} interpolates nine additional points
64between each points pair.
66The L{Frechet6Tuple} attributes C{fi1} and/or C{fi2} will be I{fractional}
67indices of type C{float} if keyword argument C{fraction} is used. Otherwise,
68C{fi1} and/or C{fi2} are simply type C{int} indices into the respective
69points set.
71For example, C{fractional} index value 2.5 means an intermediate point
72halfway between points[2] and points[3]. Use function L{fractional}
73to obtain the intermediate point for a I{fractional} index in the
74corresponding set of points.
76The C{Fréchet} distance was introduced in 1906 by U{Maurice Fréchet
77<https://WikiPedia.org/wiki/Maurice_Rene_Frechet>}, see U{reference
78[6]<https://www.kr.TUWien.ac.AT/staff/eiter/et-archive/cdtr9464.pdf>}.
79It is a measure of similarity between curves that takes into account the
80location and ordering of the points. Therefore, it is often a better metric
81than the well-known C{Hausdorff} distance, see the L{hausdorff} module.
82'''
84# from pygeodesy.basics import isscalar # from .points
85from pygeodesy.constants import EPS, EPS1, INF, NINF
86from pygeodesy.datums import _ellipsoidal_datum, _WGS84
87from pygeodesy.errors import _IsnotError, PointsError, _xattr, \
88 _xkwds, _xkwds_get
89import pygeodesy.formy as _formy
90from pygeodesy.interns import NN, _DOT_, _n_, _units_
91# from pygeodesy.iters import points2 as _points2 # from .points
92from pygeodesy.lazily import _ALL_LAZY, _FOR_DOCS
93from pygeodesy.named import _Named, _NamedTuple, notOverloaded, _Pass
94# from pygeodesy.namedTuples import PhiLam2Tuple # from .points
95from pygeodesy.points import _distanceTo, _fractional, isscalar, \
96 PhiLam2Tuple, points2 as _points2, radians
97from pygeodesy.props import property_doc_, property_RO
98from pygeodesy.units import FIx, Float, Number_, _xUnit, _xUnits
99from pygeodesy.unitsBase import _Str_degrees, _Str_meter, _Str_NN, \
100 _Str_radians, _Str_radians2
102from collections import defaultdict as _defaultdict
103# from math import radians # from .points
105__all__ = _ALL_LAZY.frechet
106__version__ = '23.09.22'
109def _fraction(fraction, n):
110 f = 1 # int, no fractional indices
111 if fraction in (None, 1):
112 pass
113 elif not (isscalar(fraction) and EPS < fraction < EPS1
114 and (float(n) - fraction) < n):
115 raise FrechetError(fraction=fraction)
116 elif fraction < EPS1:
117 f = float(fraction)
118 return f
121class FrechetError(PointsError):
122 '''Fréchet issue.
123 '''
124 pass
127class Frechet(_Named):
128 '''Frechet base class, requires method L{Frechet.distance} to
129 be overloaded.
130 '''
131 _datum = _WGS84
132 _func = None # formy function
133 _f1 = 1
134 _kwds = {} # func_ options
135 _n1 = 0
136 _ps1 = None
137 _units = _Str_NN # XXX Str to _Pass and for backward compatibility
139 def __init__(self, point1s, fraction=None, name=NN, units=NN, **kwds):
140 '''New C{Frechet...} calculator/interpolator.
142 @arg point1s: First set of points (C{LatLon}[], L{Numpy2LatLon}[],
143 L{Tuple2LatLon}[] or C{other}[]).
144 @kwarg fraction: Index fraction (C{float} in L{EPS}..L{EPS1}) to
145 interpolate intermediate B{C{point1s}} or use
146 C{None}, C{0} or C{1} for no intermediate
147 B{C{point1s}} and no I{fractional} indices.
148 @kwarg name: Optional calculator/interpolator name (C{str}).
149 @kwarg units: Optional distance units (C{Unit} or C{str}).
150 @kwarg kwds: Optional keyword argument for distance function,
151 retrievable with property C{kwds}.
153 @raise FrechetError: Insufficient number of B{C{point1s}} or
154 an invalid B{C{point1}}, B{C{fraction}}
155 or B{C{units}}.
156 '''
157 self._n1, self._ps1 = self._points2(point1s)
158 if fraction:
159 self.fraction = fraction
160 if name:
161 self.name = name
162 if units: # and not self.units:
163 self.units = units
164 if kwds:
165 self._kwds = kwds
167 @property_RO
168 def adjust(self):
169 '''Get the C{adjust} setting (C{bool} or C{None}).
170 '''
171 return _xkwds_get(self._kwds, adjust=None)
173 @property_RO
174 def datum(self):
175 '''Get the datum (L{Datum} or C{None} if not applicable).
176 '''
177 return self._datum
179 def _datum_setter(self, datum):
180 '''(INTERNAL) Set the datum.
181 '''
182 d = datum or _xattr(self._ps1[0], datum=None)
183 if d and d is not self._datum: # PYCHOK no cover
184 self._datum = _ellipsoidal_datum(d, name=self.name)
186 def discrete(self, point2s, fraction=None):
187 '''Compute the C{forward, discrete Fréchet} distance.
189 @arg point2s: Second set of points (C{LatLon}[], L{Numpy2LatLon}[],
190 L{Tuple2LatLon}[] or C{other}[]).
191 @kwarg fraction: Index fraction (C{float} in L{EPS}..L{EPS1}) to
192 interpolate intermediate B{C{point2s}} or use
193 C{None}, C{0} or C{1} for no intermediate
194 B{C{point2s}} and no I{fractional} indices.
196 @return: A L{Frechet6Tuple}C{(fd, fi1, fi2, r, n, units)}.
198 @raise FrechetError: Insufficient number of B{C{point2s}} or
199 an invalid B{C{point2}} or B{C{fraction}}.
201 @raise RecursionError: Recursion depth exceeded, see U{sys.getrecursionlimit
202 <https://docs.Python.org/3/library/sys.html#sys.getrecursionlimit>}.
203 '''
204 return self._discrete(point2s, fraction, self.distance)
206 def _discrete(self, point2s, fraction, distance):
207 '''(INTERNAL) Detailed C{discrete} with C{disance}.
208 '''
209 n2, ps2 = self._points2(point2s)
211 f2 = _fraction(fraction, n2)
212 p2 = self.points_fraction if f2 < EPS1 else self.points_ # PYCHOK expected
214 f1 = self.fraction
215 p1 = self.points_fraction if f1 < EPS1 else self.points_ # PYCHOK expected
217 def _dF(fi1, fi2):
218 return distance(p1(self._ps1, fi1), p2(ps2, fi2))
220 try:
221 return _frechet_(self._n1, f1, n2, f2, _dF, self.units)
222 except TypeError as x:
223 t = _DOT_(self.classname, self.discrete.__name__)
224 raise FrechetError(t, cause=x)
226 def distance(self, point1, point2):
227 '''Return the distance between B{C{point1}} and B{C{point2s}},
228 subject to the supplied optional keyword arguments, see
229 property C{kwds}.
230 '''
231 return self._func(point1.lat, point1.lon,
232 point2.lat, point2.lon, **self._kwds)
234 @property_doc_(''' the index fraction (C{float}).''')
235 def fraction(self):
236 '''Get the index fraction (C{float} or C{1}).
237 '''
238 return self._f1
240 @fraction.setter # PYCHOK setter!
241 def fraction(self, fraction):
242 '''Set the index fraction (C{float} in C{EPS}..L{EPS1}) to interpolate
243 intermediate B{C{point1s}} or use C{None}, C{0} or C{1} for no
244 intermediate B{C{point1s}} and no I{fractional} indices.
246 @raise FrechetError: Invalid B{C{fraction}}.
247 '''
248 self._f1 = _fraction(fraction, self._n1)
250 def _func(self, *args, **kwds): # PYCHOK no cover
251 '''(INTERNAL) I{Must be overloaded}.'''
252 notOverloaded(self, *args, **kwds)
254 @property_RO
255 def kwds(self):
256 '''Get the supplied, optional keyword arguments (C{dict}).
257 '''
258 return self._kwds
260 def point(self, point):
261 '''Convert a point for the C{.distance} method.
263 @arg point: The point to convert ((C{LatLon}, L{Numpy2LatLon},
264 L{Tuple2LatLon} or C{other}).
266 @return: The converted B{C{point}}.
267 '''
268 return point # pass thru
270 def points_(self, points, i):
271 '''Get and convert a point for the C{.distance} method.
273 @arg points: The orignal B{C{points}} to convert ((C{LatLon}[],
274 L{Numpy2LatLon}[], L{Tuple2LatLon}[] or C{other}[]).
275 @arg i: The B{C{points}} index (C{int}).
277 @return: The converted B{C{points[i]}}.
278 '''
279 return self.point(points[i])
281 def points_fraction(self, points, fi):
282 '''Get and convert a I{fractional} point for the C{.distance} method.
284 @arg points: The orignal B{C{points}} to convert ((C{LatLon}[],
285 L{Numpy2LatLon}[], L{Tuple2LatLon}[] or C{other}[]).
286 @arg fi: The I{fractional} index in B{C{points}} (C{float} or C{int}).
288 @return: The interpolated, converted, intermediate B{C{points[fi]}}.
289 '''
290 return self.point(_fractional(points, fi, None, wrap=None)) # was=self.wrap
292 def _points2(self, points):
293 '''(INTERNAL) Check a set of points, overloaded in L{FrechetDistanceTo}.
294 '''
295 return _points2(points, closed=False, Error=FrechetError)
297 @property_doc_(''' the distance units (C{Unit} or C{str}).''')
298 def units(self):
299 '''Get the distance units (C{Unit} or C{str}).
300 '''
301 return self._units
303 @units.setter # PYCHOK setter!
304 def units(self, units):
305 '''Set the distance units (C{Unit} or C{str}).
307 @raise TypeError: Invalid B{C{units}}.
308 '''
309 self._units = _xUnits(units, Base=Float)
311 @property_RO
312 def wrap(self):
313 '''Get the C{wrap} setting (C{bool} or C{None}).
314 '''
315 return _xkwds_get(self._kwds, wrap=None)
318class FrechetDegrees(Frechet):
319 '''DEPRECATED, use an other C{Frechet*} class.
320 '''
321 _units = _Str_degrees
323 if _FOR_DOCS:
324 __init__ = Frechet.__init__
325 discrete = Frechet.discrete
327 def distance(self, point1, point2, *args, **kwds): # PYCHOK no cover
328 '''Return the distance in C{degrees} between B{C{point1}} and B{C{point2}}.
329 I{Must be overloaded}.'''
330 notOverloaded(self, point1, point2, *args, **kwds)
333class FrechetRadians(Frechet):
334 '''DEPRECATED, use an other C{Frechet*} class.
335 '''
336 _units = _Str_radians
338 if _FOR_DOCS:
339 __init__ = Frechet.__init__
340 discrete = Frechet.discrete
342 def distance(self, point1, point2, *args, **kwds): # PYCHOK no cover
343 '''Return the distance in C{radians} between B{C{point1}} and B{C{point2}}.
344 I{Must be overloaded}.'''
345 notOverloaded(self, point1, point2, *args, **kwds)
347 def point(self, point):
348 '''Return B{C{point}} as L{PhiLam2Tuple} to maintain
349 I{backward compatibility} of L{FrechetRadians}.
351 @return: A L{PhiLam2Tuple}C{(phi, lam)}.
352 '''
353 try:
354 return point.philam
355 except AttributeError:
356 return PhiLam2Tuple(radians(point.lat), radians(point.lon))
359class _FrechetMeterRadians(Frechet):
360 '''(INTERNAL) Returning C{meter} or C{radians} depending on
361 the optional keyword arguments supplied at instantiation
362 of the C{Frechet*} sub-class.
363 '''
364 _units = _Str_meter
365 _units_ = _Str_radians
367 def discrete(self, point2s, fraction=None):
368 '''Overloaded method L{Frechet.discrete} to determine
369 the distance function and units from the optional
370 keyword arguments given at this instantiation, see
371 property C{kwds}.
373 @see: L{Frechet.discrete} for other details.
374 '''
375 return self._discrete(point2s, fraction, _formy._radistance(self))
377 def _func_(self, *args, **kwds): # PYCHOK no cover
378 '''(INTERNAL) I{Must be overloaded}.'''
379 notOverloaded(self, *args, **kwds)
382class FrechetCosineAndoyerLambert(_FrechetMeterRadians):
383 '''Compute the C{Frechet} distance based on the I{angular} distance
384 in C{radians} from function L{pygeodesy.cosineAndoyerLambert}.
385 '''
386 def __init__(self, point1s, fraction=None, name=NN, **datum_wrap):
387 '''New L{FrechetCosineAndoyerLambert} calculator/interpolator.
389 @kwarg datum_wrap: Optional keyword arguments for function
390 L{pygeodesy.cosineAndoyerLambert}.
392 @see: L{Frechet.__init__} for details about B{C{point1s}},
393 B{C{fraction}}, B{C{name}} and other exceptions.
394 '''
395 Frechet.__init__(self, point1s, fraction=fraction, name=name,
396 **datum_wrap)
397 self._func = _formy.cosineAndoyerLambert
398 self._func_ = _formy.cosineAndoyerLambert_
400 if _FOR_DOCS:
401 discrete = Frechet.discrete
404class FrechetCosineForsytheAndoyerLambert(_FrechetMeterRadians):
405 '''Compute the C{Frechet} distance based on the I{angular} distance
406 in C{radians} from function L{pygeodesy.cosineForsytheAndoyerLambert}.
407 '''
408 def __init__(self, point1s, fraction=None, name=NN, **datum_wrap):
409 '''New L{FrechetCosineForsytheAndoyerLambert} calculator/interpolator.
411 @kwarg datum_wrap: Optional keyword arguments for function
412 L{pygeodesy.cosineAndoyerLambert}.
414 @see: L{Frechet.__init__} for details about B{C{point1s}},
415 B{C{fraction}}, B{C{name}} and other exceptions.
416 '''
417 Frechet.__init__(self, point1s, fraction=fraction, name=name,
418 **datum_wrap)
419 self._func = _formy.cosineForsytheAndoyerLambert
420 self._func_ = _formy.cosineForsytheAndoyerLambert_
422 if _FOR_DOCS:
423 discrete = Frechet.discrete
426class FrechetCosineLaw(_FrechetMeterRadians):
427 '''Compute the C{Frechet} distance based on the I{angular} distance
428 in C{radians} from function L{pygeodesy.cosineLaw}.
430 @note: See note at function L{pygeodesy.vincentys_}.
431 '''
432 def __init__(self, point1s, fraction=None, name=NN, **radius_wrap):
433 '''New L{FrechetCosineLaw} calculator/interpolator.
435 @kwarg radius_wrap: Optional keyword arguments for function
436 L{pygeodesy.cosineLaw}.
438 @see: L{Frechet.__init__} for details about B{C{point1s}},
439 B{C{fraction}}, B{C{name}} and other exceptions.
440 '''
441 Frechet.__init__(self, point1s, fraction=fraction, name=name,
442 **radius_wrap)
443 self._func = _formy.cosineLaw
444 self._func_ = _formy.cosineLaw_
446 if _FOR_DOCS:
447 discrete = Frechet.discrete
450class FrechetDistanceTo(Frechet): # FrechetMeter
451 '''Compute the C{Frechet} distance based on the distance from the
452 point1s' C{LatLon.distanceTo} method, conventionally in C{meter}.
453 '''
454 _units = _Str_meter
456 def __init__(self, point1s, fraction=None, name=NN, **distanceTo_kwds):
457 '''New L{FrechetDistanceTo} calculator/interpolator.
459 @kwarg distanceTo_kwds: Optional keyword arguments for each
460 B{C{point1s}}' C{LatLon.distanceTo}
461 method.
463 @see: L{Frechet.__init__} for details about B{C{point1s}},
464 B{C{fraction}}, B{C{name}} and other exceptions.
466 @note: All B{C{point1s}} I{must} be instances of the same
467 ellipsoidal or spherical C{LatLon} class.
468 '''
469 Frechet.__init__(self, point1s, fraction=fraction, name=name,
470 **distanceTo_kwds)
472 if _FOR_DOCS:
473 discrete = Frechet.discrete
475 def distance(self, p1, p2):
476 '''Return the distance in C{meter}.
477 '''
478 return p1.distanceTo(p2, **self._kwds)
480 def _points2(self, points):
481 '''(INTERNAL) Check a set of points.
482 '''
483 np, ps = Frechet._points2(self, points)
484 return np, _distanceTo(FrechetError, points=ps)
487class FrechetEquirectangular(Frechet):
488 '''Compute the C{Frechet} distance based on the I{equirectangular}
489 distance in C{radians squared} like function L{pygeodesy.equirectangular}.
490 '''
491 _units = _Str_radians2
493 def __init__(self, point1s, fraction=None, name=NN, **adjust_limit_wrap):
494 '''New L{FrechetEquirectangular} calculator/interpolator.
496 @kwarg adjust_limit_wrap: Optional keyword arguments for function
497 L{pygeodesy.equirectangular_} I{with default}
498 C{B{limit}=0} for I{backward compatibility}.
500 @see: L{Frechet.__init__} for details about B{C{point1s}},
501 B{C{fraction}}, B{C{name}} and other exceptions.
502 '''
503 adjust_limit_wrap = _xkwds(adjust_limit_wrap, limit=0)
504 Frechet.__init__(self, point1s, fraction=fraction, name=name,
505 **adjust_limit_wrap)
506 self._func = _formy._equirectangular # helper
508 if _FOR_DOCS:
509 discrete = Frechet.discrete
512class FrechetEuclidean(_FrechetMeterRadians):
513 '''Compute the C{Frechet} distance based on the I{Euclidean}
514 distance in C{radians} from function L{pygeodesy.euclidean}.
515 '''
516 def __init__(self, point1s, fraction=None, name=NN, **adjust_radius_wrap): # was=True
517 '''New L{FrechetEuclidean} calculator/interpolator.
519 @kwarg adjust_radius_wrap: Optional keyword arguments for
520 function L{pygeodesy.euclidean}.
522 @see: L{Frechet.__init__} for details about B{C{point1s}},
523 B{C{fraction}}, B{C{name}} and other exceptions.
524 '''
525 Frechet.__init__(self, point1s, fraction=fraction, name=name,
526 **adjust_radius_wrap)
527 self._func = _formy.euclidean
528 self._func_ = _formy.euclidean_
530 if _FOR_DOCS:
531 discrete = Frechet.discrete
534class FrechetExact(Frechet):
535 '''Compute the C{Frechet} distance based on the I{angular} distance
536 in C{degrees} from method L{GeodesicExact}C{.Inverse}.
537 '''
538 _units = _Str_degrees
540 def __init__(self, point1s, fraction=None, name=NN, datum=None, **wrap):
541 '''New L{FrechetExact} calculator/interpolator.
543 @kwarg datum: Datum to override the default C{Datums.WGS84} and
544 first B{C{point1s}}' datum (L{Datum}, L{Ellipsoid},
545 L{Ellipsoid2} or L{a_f2Tuple}).
546 @kwarg wrap: Optional keyword argument for method C{Inverse1}
547 of class L{geodesicx.GeodesicExact}.
549 @raise TypeError: Invalid B{C{datum}}.
551 @see: L{Frechet.__init__} for details about B{C{point1s}},
552 B{C{fraction}}, B{C{name}} and other exceptions.
553 '''
554 Frechet.__init__(self, point1s, fraction=fraction, name=name,
555 **wrap)
556 self._datum_setter(datum)
557 self._func = self.datum.ellipsoid.geodesicx.Inverse1 # note -x
559 if _FOR_DOCS:
560 discrete = Frechet.discrete
563class FrechetFlatLocal(_FrechetMeterRadians):
564 '''Compute the C{Frechet} distance based on the I{angular} distance in
565 C{radians squared} like function L{pygeodesy.flatLocal_}/L{pygeodesy.hubeny}.
566 '''
567 _units_ = _Str_radians2 # see L{flatLocal_}
569 def __init__(self, point1s, fraction=None, name=NN, **datum_scaled_wrap):
570 '''New L{FrechetFlatLocal}/L{FrechetHubeny} calculator/interpolator.
572 @kwarg datum_scaled_wrap: Optional keyword arguments for
573 function L{pygeodesy.flatLocal}.
575 @see: L{Frechet.__init__} for details about B{C{point1s}},
576 B{C{fraction}}, B{C{name}} and other exceptions.
578 @note: The distance C{units} are C{radians squared}, not C{radians}.
579 '''
580 Frechet.__init__(self, point1s, fraction=fraction, name=name,
581 **datum_scaled_wrap)
582 self._func = _formy.flatLocal
583 self._func_ = self.datum.ellipsoid._hubeny_2
585 if _FOR_DOCS:
586 discrete = Frechet.discrete
589class FrechetFlatPolar(_FrechetMeterRadians):
590 '''Compute the C{Frechet} distance based on the I{angular} distance
591 in C{radians} from function L{flatPolar_}.
592 '''
593 def __init__(self, point1s, fraction=None, name=NN, **radius_wrap):
594 '''New L{FrechetFlatPolar} calculator/interpolator.
596 @kwarg radius_wrap: Optional keyword arguments for function
597 L{pygeodesy.flatPolar}.
599 @see: L{Frechet.__init__} for details about B{C{point1s}},
600 B{C{fraction}}, B{C{name}} and other exceptions.
601 '''
602 Frechet.__init__(self, point1s, fraction=fraction, name=name,
603 **radius_wrap)
604 self._func = _formy.flatPolar
605 self._func_ = _formy.flatPolar_
607 if _FOR_DOCS:
608 discrete = Frechet.discrete
611class FrechetHaversine(_FrechetMeterRadians):
612 '''Compute the C{Frechet} distance based on the I{angular}
613 distance in C{radians} from function L{pygeodesy.haversine_}.
615 @note: See note at function L{pygeodesy.vincentys_}.
616 '''
617 def __init__(self, point1s, fraction=None, name=NN, **radius_wrap):
618 '''New L{FrechetHaversine} calculator/interpolator.
620 @kwarg radius_wrap: Optional keyword arguments for function
621 L{pygeodesy.haversine}.
623 @see: L{Frechet.__init__} for details about B{C{point1s}},
624 B{C{fraction}}, B{C{name}} and other exceptions.
625 '''
626 Frechet.__init__(self, point1s, fraction=fraction, name=name,
627 **radius_wrap)
628 self._func = _formy.haversine
629 self._func_ = _formy.haversine_
631 if _FOR_DOCS:
632 discrete = Frechet.discrete
635class FrechetHubeny(FrechetFlatLocal): # for Karl Hubeny
636 if _FOR_DOCS:
637 __doc__ = FrechetFlatLocal.__doc__
638 __init__ = FrechetFlatLocal.__init__
639 discrete = FrechetFlatLocal.discrete
640 distance = FrechetFlatLocal.discrete
643class FrechetKarney(Frechet):
644 '''Compute the C{Frechet} distance based on the I{angular}
645 distance in C{degrees} from I{Karney}'s U{geographiclib
646 <https://PyPI.org/project/geographiclib>} U{geodesic.Geodesic
647 <https://GeographicLib.SourceForge.io/Python/doc/code.html>}
648 C{Inverse} method.
649 '''
650 _units = _Str_degrees
652 def __init__(self, point1s, fraction=None, name=NN, datum=None, **wrap):
653 '''New L{FrechetKarney} calculator/interpolator.
655 @kwarg datum: Datum to override the default C{Datums.WGS84} and
656 first B{C{knots}}' datum (L{Datum}, L{Ellipsoid},
657 L{Ellipsoid2} or L{a_f2Tuple}).
658 @kwarg wrap: Optional keyword argument for method C{Inverse1}
659 of class L{geodesicw.Geodesic}.
661 @raise ImportError: Package U{geographiclib
662 <https://PyPI.org/project/geographiclib>} missing.
664 @raise TypeError: Invalid B{C{datum}}.
666 @see: L{Frechet.__init__} for details about B{C{point1s}},
667 B{C{fraction}}, B{C{name}} and other exceptions.
668 '''
669 Frechet.__init__(self, point1s, fraction=fraction, name=name,
670 **wrap)
671 self._datum_setter(datum)
672 self._func = self.datum.ellipsoid.geodesic.Inverse1
674 if _FOR_DOCS:
675 discrete = Frechet.discrete
678class FrechetThomas(_FrechetMeterRadians):
679 '''Compute the C{Frechet} distance based on the I{angular} distance
680 in C{radians} from function L{pygeodesy.thomas_}.
681 '''
682 def __init__(self, point1s, fraction=None, name=NN, **datum_wrap):
683 '''New L{FrechetThomas} calculator/interpolator.
685 @kwarg datum_wrap: Optional keyword argument for function
686 L{pygeodesy.thomas}.
688 @see: L{Frechet.__init__} for details about B{C{point1s}},
689 B{C{fraction}}, B{C{name}} and other exceptions.
690 '''
691 Frechet.__init__(self, point1s, fraction=fraction, name=name,
692 **datum_wrap)
693 self._func = _formy.thomas
694 self._func_ = _formy.thomas_
696 if _FOR_DOCS:
697 discrete = Frechet.discrete
700class FrechetVincentys(_FrechetMeterRadians):
701 '''Compute the C{Frechet} distance based on the I{angular}
702 distance in C{radians} from function L{pygeodesy.vincentys_}.
704 @note: See note at function L{pygeodesy.vincentys_}.
705 '''
706 def __init__(self, point1s, fraction=None, name=NN, **radius_wrap):
707 '''New L{FrechetVincentys} calculator/interpolator.
709 @kwarg radius_wrap: Optional keyword arguments for function
710 L{pygeodesy.vincentys}.
712 @see: L{Frechet.__init__} for details about B{C{point1s}},
713 B{C{fraction}}, B{C{name}} and other exceptions.
714 '''
715 Frechet.__init__(self, point1s, fraction=fraction, name=name,
716 **radius_wrap)
717 self._func = _formy.vincentys
718 self._func_ = _formy.vincentys_
720 if _FOR_DOCS:
721 discrete = Frechet.discrete
724def _frechet_(ni, fi, nj, fj, dF, units): # MCCABE 14
725 '''(INTERNAL) Recursive core of function L{frechet_}
726 and method C{discrete} of C{Frechet...} classes.
727 '''
728 iFs = {}
730 def iF(i): # cache index, depth ints and floats
731 return iFs.setdefault(i, i)
733 cF = _defaultdict(dict)
735 def _rF(i, j, r): # recursive Fréchet
736 i = iF(i)
737 j = iF(j)
738 try:
739 t = cF[i][j]
740 except KeyError:
741 r = iF(r + 1)
742 try:
743 if i > 0:
744 if j > 0:
745 t = min(_rF(i - fi, j, r),
746 _rF(i - fi, j - fj, r),
747 _rF(i, j - fj, r))
748 elif j < 0:
749 raise IndexError
750 else: # j == 0
751 t = _rF(i - fi, 0, r)
752 elif i < 0:
753 raise IndexError
755 elif j > 0: # i == 0
756 t = _rF(0, j - fj, r)
757 elif j < 0: # i == 0
758 raise IndexError
759 else: # i == j == 0
760 t = (NINF, i, j, r)
762 d = dF(i, j)
763 if d > t[0]:
764 t = (d, i, j, r)
765 except IndexError:
766 t = (INF, i, j, r)
767 cF[i][j] = t
768 return t
770 t = _rF(ni - 1, nj - 1, 0)
771 t += (sum(map(len, cF.values())), units)
772# del cF, iFs
773 return Frechet6Tuple(t) # *t
776def frechet_(point1s, point2s, distance=None, units=NN):
777 '''Compute the I{discrete} U{Fréchet<https://WikiPedia.org/wiki/Frechet_distance>}
778 distance between two paths, each given as a set of points.
780 @arg point1s: First set of points (C{LatLon}[], L{Numpy2LatLon}[],
781 L{Tuple2LatLon}[] or C{other}[]).
782 @arg point2s: Second set of points (C{LatLon}[], L{Numpy2LatLon}[],
783 L{Tuple2LatLon}[] or C{other}[]).
784 @kwarg distance: Callable returning the distance between a B{C{point1s}}
785 and a B{C{point2s}} point (signature C{(point1, point2)}).
786 @kwarg units: Optional, the distance units (C{Unit} or C{str}).
788 @return: A L{Frechet6Tuple}C{(fd, fi1, fi2, r, n, units)} where C{fi1}
789 and C{fi2} are type C{int} indices into B{C{point1s}} respectively
790 B{C{point2s}}.
792 @raise FrechetError: Insufficient number of B{C{point1s}} or B{C{point2s}}.
794 @raise RecursionError: Recursion depth exceeded, see U{sys.getrecursionlimit()
795 <https://docs.Python.org/3/library/sys.html#sys.getrecursionlimit>}.
797 @raise TypeError: If B{C{distance}} is not a callable.
799 @note: Function L{frechet_} does I{not} support I{fractional} indices
800 for intermediate B{C{point1s}} and B{C{point2s}}.
801 '''
802 if not callable(distance):
803 raise _IsnotError(callable.__name__, distance=distance)
805 n1, ps1 = _points2(point1s, closed=False, Error=FrechetError)
806 n2, ps2 = _points2(point2s, closed=False, Error=FrechetError)
808 def _dF(i1, i2):
809 return distance(ps1[i1], ps2[i2])
811 return _frechet_(n1, 1, n2, 1, _dF, units)
814class Frechet6Tuple(_NamedTuple):
815 '''6-Tuple C{(fd, fi1, fi2, r, n, units)} with the I{discrete}
816 U{Fréchet<https://WikiPedia.org/wiki/Frechet_distance>} distance
817 C{fd}, I{fractional} indices C{fi1} and C{fi2} as C{FIx}, the
818 recursion depth C{r}, the number of distances computed C{n} and
819 the L{units} class or class or name of the distance C{units}.
821 If I{fractional} indices C{fi1} and C{fi2} are C{int}, the
822 returned C{fd} is the distance between C{point1s[fi1]} and
823 C{point2s[fi2]}. For C{float} indices, the distance is
824 between an intermediate point along C{point1s[int(fi1)]} and
825 C{point1s[int(fi1) + 1]} respectively an intermediate point
826 along C{point2s[int(fi2)]} and C{point2s[int(fi2) + 1]}.
828 Use function L{fractional} to compute the point at a
829 I{fractional} index.
830 '''
831 _Names_ = ('fd', 'fi1', 'fi2', 'r', _n_, _units_)
832 _Units_ = (_Pass, FIx, FIx, Number_, Number_, _Pass)
834 def toUnits(self, **Error): # PYCHOK expected
835 '''Overloaded C{_NamedTuple.toUnits} for C{fd} units.
836 '''
837 U = _xUnit(self.units, Float) # PYCHOK expected
838 self._Units_ = (U,) + Frechet6Tuple._Units_[1:]
839 return _NamedTuple.toUnits(self, **Error)
841# def __gt__(self, other):
842# _xinstanceof(Frechet6Tuple, other=other)
843# return self if self.fd > other.fd else other # PYCHOK .fd=[0]
844#
845# def __lt__(self, other):
846# _xinstanceof(Frechet6Tuple, other=other)
847# return self if self.fd < other.fd else other # PYCHOK .fd=[0]
849# **) MIT License
850#
851# Copyright (C) 2016-2023 -- mrJean1 at Gmail -- All Rights Reserved.
852#
853# Permission is hereby granted, free of charge, to any person obtaining a
854# copy of this software and associated documentation files (the "Software"),
855# to deal in the Software without restriction, including without limitation
856# the rights to use, copy, modify, merge, publish, distribute, sublicense,
857# and/or sell copies of the Software, and to permit persons to whom the
858# Software is furnished to do so, subject to the following conditions:
859#
860# The above copyright notice and this permission notice shall be included
861# in all copies or substantial portions of the Software.
862#
863# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
864# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
865# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
866# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
867# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
868# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
869# OTHER DEALINGS IN THE SOFTWARE.