Coverage for pygeodesy/ellipsoidalNvector.py: 97%

130 statements  

« prev     ^ index     » next       coverage.py v7.2.2, created at 2023-09-15 09:43 -0400

1 

2# -*- coding: utf-8 -*- 

3 

4u'''Ellipsoidal, C{N-vector}-based geodesy. 

5 

6Ellipsoidal classes geodetic (lat-/longitude) L{LatLon}, geocentric 

7(ECEF) L{Cartesian}, DEPRECATED L{Ned} and L{Nvector} and functions 

8L{meanOf}, L{sumOf} and DEPRECATED L{toNed}. 

9 

10Pure Python implementation of n-vector-based geodetic (lat-/longitude) 

11methods by I{(C) Chris Veness 2011-2016} published under the same MIT 

12Licence**, see U{Vector-based geodesy 

13<https://www.Movable-Type.co.UK/scripts/latlong-vectors.html>}. 

14 

15These classes and functions work with: (a) geodesic (polar) lat-/longitude 

16points on the earth's surface and (b) 3-D vectors used as n-vectors 

17representing points on the earth's surface or vectors normal to the plane 

18of a great circle. 

19 

20See also Kenneth Gade U{'A Non-singular Horizontal Position Representation' 

21<https://www.NavLab.net/Publications/A_Nonsingular_Horizontal_Position_Representation.pdf>}, 

22The Journal of Navigation (2010), vol 63, nr 3, pp 395-417. 

23''' 

24# make sure int/int division yields float quotient, see .basics 

25from __future__ import division as _; del _ # PYCHOK semicolon 

26 

27from pygeodesy.basics import issubclassof, map2, _xinstanceof 

28from pygeodesy.datums import _ellipsoidal_datum, _spherical_datum, _WGS84 

29# from pygeodesy.dms import toDMS # _MODS 

30from pygeodesy.ellipsoidalBase import CartesianEllipsoidalBase, \ 

31 _TOL_M, LatLonEllipsoidalBase, \ 

32 _nearestOn, _Wrap 

33from pygeodesy.errors import _IsnotError, _xkwds 

34# from pygeodesy.fmath import fdot # from .nvectorBase 

35from pygeodesy.interns import NN, _Nv00_, _COMMASPACE_ 

36from pygeodesy.interns import _down_, _east_, _north_, _pole_ # PYCHOK used! 

37from pygeodesy.lazily import _ALL_LAZY, _ALL_MODS as _MODS, _ALL_OTHER 

38# from pygeodesy.ltp import Ltp # _MODS 

39from pygeodesy.ltpTuples import Aer as _Aer, Ned as _Ned, Ned4Tuple, \ 

40 sincos2d_, _xnamed 

41# from pygeodesy.named import _xnamed # from .ltpTuples 

42from pygeodesy.nvectorBase import fabs, fdot, NorthPole, LatLonNvectorBase, \ 

43 NvectorBase, sumOf as _sumOf 

44from pygeodesy.props import deprecated_class, deprecated_function, \ 

45 deprecated_method, Property_RO 

46from pygeodesy.streprs import Fmt, fstr, _xzipairs 

47from pygeodesy.units import Bearing, Distance, Height, Scalar 

48# from pygeodesy.utily import sincos2d_, _Wrap # from .ltpTuples, .ellipsoidalBase 

49 

50# from math import fabs # from .nvectorBase 

51 

52__all__ = _ALL_LAZY.ellipsoidalNvector 

53__version__ = '23.05.04' 

54 

55 

56class Ned(_Ned): 

57 '''DEPRECATED, use class L{pygeodesy.Ned}.''' 

58 

59 def __init__(self, north, east, down, name=NN): 

60 deprecated_class(self.__class__) 

61 _Ned.__init__(self, north, east, down, name=name) 

62 

63 @deprecated_method # PYCHOK expected 

64 def toRepr(self, prec=None, fmt=Fmt.SQUARE, sep=_COMMASPACE_, **unused): 

65 '''DEPRECATED, use class L{pygeodesy.Ned}. 

66 

67 @kwarg prec: Number of (decimal) digits, unstripped (C{int}). 

68 @kwarg fmt: Enclosing backets format (C{str}). 

69 @kwarg sep: Separator between NEDs (C{str}). 

70 

71 @return: This Ned as "[L:f, B:degrees360, E:degrees90]" (C{str}) 

72 with length or slantrange C{L}, bearing or azimuth C{B} 

73 and elevation C{E}. 

74 ''' 

75 dms = _MODS.dms 

76 t = (fstr(self.slantrange, prec=3 if prec is None else prec), 

77 dms.toDMS(self.azimuth, form=dms.F_D, prec=prec, ddd=0), 

78 dms.toDMS(self.elevation, form=dms.F_D, prec=prec, ddd=0)) 

79 return _xzipairs('LBE', t, sep=sep, fmt=fmt) 

80 

81 

82class Cartesian(CartesianEllipsoidalBase): 

83 '''Extended to convert geocentric, L{Cartesian} points to 

84 L{Nvector} and n-vector-based, geodetic L{LatLon}. 

85 ''' 

86 @Property_RO 

87 def Ecef(self): 

88 '''Get the ECEF I{class} (L{EcefVeness}), I{lazily}. 

89 ''' 

90 return _MODS.ecef.EcefVeness 

91 

92 def toLatLon(self, **LatLon_and_kwds): # PYCHOK LatLon=LatLon, datum=None 

93 '''Convert this cartesian to an C{Nvector}-based geodetic point. 

94 

95 @kwarg LatLon_and_kwds: Optional L{LatLon}, B{C{datum}} and other 

96 keyword arguments. Use C{B{LatLon}=...} to 

97 override this L{LatLon} class or specify 

98 C{B{LatLon} is None}. 

99 

100 @return: The geodetic point (L{LatLon}) or if B{C{LatLon}} is set 

101 to C{None}, an L{Ecef9Tuple}C{(x, y, z, lat, lon, height, 

102 C, M, datum)} with C{C} and C{M} if available. 

103 

104 @raise TypeError: Invalid B{C{LatLon_and_kwds}}. 

105 ''' 

106 kwds = _xkwds(LatLon_and_kwds, LatLon=LatLon, datum=self.datum) 

107 return CartesianEllipsoidalBase.toLatLon(self, **kwds) 

108 

109 def toNvector(self, **Nvector_and_kwds): # PYCHOK Datums.WGS84 

110 '''Convert this cartesian to L{Nvector} components, I{including height}. 

111 

112 @kwarg Nvector_and_kwds: Optional L{Nvector}, B{C{datum}} and other 

113 keyword arguments. Use C{B{Nvector}=...} to 

114 override this L{Nvector} class or specify 

115 C{B{Nvector} is None}. 

116 

117 @return: The C{n-vector} components (L{Nvector}) or if B{C{Nvector}} 

118 is set to C{None}, a L{Vector4Tuple}C{(x, y, z, h)} 

119 

120 @raise TypeError: Invalid B{C{Nvector_and_kwds}}. 

121 

122 @example: 

123 

124 >>> from ellipsoidalNvector import LatLon 

125 >>> c = Cartesian(3980581, 97, 4966825) 

126 >>> n = c.toNvector() # (0.62282, 0.000002, 0.78237, +0.24) 

127 ''' 

128 kwds = _xkwds(Nvector_and_kwds, Nvector=Nvector, datum=self.datum) 

129 return CartesianEllipsoidalBase.toNvector(self, **kwds) 

130 

131 

132class LatLon(LatLonNvectorBase, LatLonEllipsoidalBase): 

133 '''An n-vector-based, ellipsoidal L{LatLon} point. 

134 

135 @example: 

136 

137 >>> from ellipsoidalNvector import LatLon 

138 >>> p = LatLon(52.205, 0.119) # height=0, datum=Datums.WGS84 

139 ''' 

140 _Nv = None # cached toNvector (L{Nvector}) 

141 

142 def _update(self, updated, *attrs, **setters): # PYCHOK args 

143 '''(INTERNAL) Zap cached attributes if updated. 

144 ''' 

145 if updated: 

146 LatLonNvectorBase._update(self, updated, _Nv=self._Nv) # special case 

147 LatLonEllipsoidalBase._update(self, updated, *attrs, **setters) 

148 

149# def crossTrackDistanceTo(self, start, end, radius=R_M): 

150# '''Return the (signed) distance from this point to the great 

151# circle defined by a start point and an end point or bearing. 

152# 

153# @arg start: Start point of great circle line (L{LatLon}). 

154# @arg end: End point of great circle line (L{LatLon}) or 

155# initial bearing (compass C{degrees360}) at the 

156# start point. 

157# @kwarg radius: Mean earth radius (C{meter}). 

158# 

159# @return: Distance to great circle, negative if to left or 

160# positive if to right of line (C{meter}, same units 

161# as B{C{radius}}). 

162# 

163# @raise TypeError: If B{C{start}} or B{C{end}} point is not L{LatLon}. 

164# 

165# @example: 

166# 

167# >>> p = LatLon(53.2611, -0.7972) 

168# 

169# >>> s = LatLon(53.3206, -1.7297) 

170# >>> b = 96.0 

171# >>> d = p.crossTrackDistanceTo(s, b) # -305.7 

172# 

173# >>> e = LatLon(53.1887, 0.1334) 

174# >>> d = p.crossTrackDistanceTo(s, e) # -307.5 

175# ''' 

176# self.others(start=start) 

177# 

178# if isscalar(end): # gc from point and bearing 

179# gc = start.greatCircle(end) 

180# else: # gc by two points 

181# gc = start.toNvector().cross(end.toNvector()) 

182# 

183# # (signed) angle between point and gc normal vector 

184# v = self.toNvector() 

185# a = gc.angleTo(v, vSign=v.cross(gc)) 

186# a = (-PI_2 - a) if a < 0 else (PI_2 - a) 

187# return a * float(radius) 

188 

189 def deltaTo(self, other, Ned=Ned, wrap=False): 

190 '''Calculate the local delta from this to an other point. 

191 

192 @note: This is a linear delta, I{unrelated} to a geodesic 

193 on the ellipsoid. 

194 

195 @arg other: The other point (L{LatLon}). 

196 @kwarg Ned: Class to use (L{pygeodesy.Ned} or 

197 L{pygeodesy.Ned4Tuple}), overriding the 

198 default DEPRECATED L{Ned}. 

199 @kwarg wrap: If C{True}, wrap or I{normalize} the 

200 B{C{other}} point (C{bool}). 

201 

202 @return: Delta from this to the other point (B{C{Ned}}). 

203 

204 @raise TypeError: The B{C{other}} point is not L{LatLon} or 

205 B{C{Ned}} is not L{pygeodesy.Ned} nor 

206 L{pygeodesy.Ned4Tuple} nor DEPRECATED L{Ned}. 

207 

208 @raise ValueError: If ellipsoids are incompatible. 

209 

210 @example: 

211 

212 >>> a = LatLon(49.66618, 3.45063) 

213 >>> b = LatLon(48.88667, 2.37472) 

214 >>> delta = a.deltaTo(b) # [N:-86126, E:-78900, D:1069] 

215 >>> d = delta.length # 116807.681 m 

216 >>> b = delta.bearing # 222.493° 

217 >>> e = delta.elevation # -0.5245° 

218 ''' 

219 self.ellipsoids(other) # throws TypeError and ValueError 

220 

221 p = self.others(other) 

222 if wrap: 

223 p = _Wrap.point(p) 

224 # get delta in cartesian frame 

225 dc = p.toCartesian().minus(self.toCartesian()) 

226 # rotate dc to get delta in n-vector reference 

227 # frame using the rotation matrix row vectors 

228 ned_ = map2(dc.dot, self._rotation3) 

229 if issubclassof(Ned, Ned4Tuple): 

230 ned_ += (_MODS.ltp.Ltp(self, ecef=self.Ecef(self.datum)),) 

231 elif not issubclassof(Ned, _Ned): 

232 raise _IsnotError(Fmt.sub_class(_Ned, Ned4Tuple), Ned=Ned) 

233 return Ned(*ned_, name=self.name) 

234 

235# def destination(self, distance, bearing, radius=R_M, height=None): 

236# '''Return the destination point after traveling from this 

237# point the given distance on the given initial bearing. 

238# 

239# @arg distance: Distance traveled (C{meter}, same units as 

240# given earth B{C{radius}}). 

241# @arg bearing: Initial bearing (compass C{degrees360}). 

242# @kwarg radius: Mean earth radius (C{meter}). 

243# @kwarg height: Optional height at destination point, 

244# overriding default (C{meter}, same units 

245# as B{C{radius}}). 

246# 

247# @return: Destination point (L{LatLon}). 

248# 

249# @example: 

250# 

251# >>> p = LatLon(51.4778, -0.0015) 

252# >>> q = p.destination(7794, 300.7) 

253# >>> q.toStr() # '51.5135°N, 000.0983°W' ? 

254# ''' 

255# r = _angular(distance, radius) # angular distance in radians 

256# # great circle by starting from this point on given bearing 

257# gc = self.greatCircle(bearing) 

258# 

259# v1 = self.toNvector() 

260# x = v1.times(cos(r)) # component of v2 parallel to v1 

261# y = gc.cross(v1).times(sin(r)) # component of v2 perpendicular to v1 

262# 

263# v2 = x.plus(y).unit() 

264# return v2.toLatLon(height=self.height if height is C{None} else height) 

265 

266 def destinationNed(self, delta): 

267 '''Calculate the destination point using the supplied NED delta 

268 from this point. 

269 

270 @arg delta: Delta from this to the other point in the local 

271 tangent plane (LTP) of this point (L{Ned}). 

272 

273 @return: Destination point (L{LatLon}). 

274 

275 @raise TypeError: If B{C{delta}} is not L{pygeodesy.Ned} or 

276 DEPRECATED L{Ned}. 

277 

278 @example: 

279 

280 >>> a = LatLon(49.66618, 3.45063) 

281 >>> delta = Ned(-86126, -78900, 1069) # from Aer(222.493, -0.5245, 116807.681) 

282 >>> b = a.destinationNed(delta) # 48.886669°N, 002.374721°E or 48°53′12.01″N, 002°22′29.0″E +0.20m 

283 ''' 

284 _xinstanceof(_Ned, delta=delta) 

285 

286 nv, ev, dv = self._rotation3 

287 # convert NED delta to standard coordinate frame of n-vector 

288 dn = delta.ned 

289 # rotate dn to get delta in cartesian (ECEF) coordinate 

290 # reference frame using the rotation matrix column vectors 

291 dc = Cartesian(fdot(dn, nv.x, ev.x, dv.x), 

292 fdot(dn, nv.y, ev.y, dv.y), 

293 fdot(dn, nv.z, ev.z, dv.z)) 

294 

295 # apply (cartesian) delta to this Cartesian to obtain destination as cartesian 

296 v = self.toCartesian().plus(dc) 

297 return v.toLatLon(datum=self.datum, LatLon=self.classof) # Cartesian(v.x, v.y, v.z).toLatLon(...) 

298 

299 def distanceTo(self, other, radius=None, wrap=False): 

300 '''I{Approximate} the distance from this to an other point. 

301 

302 @arg other: The other point (L{LatLon}). 

303 @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, 

304 L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or 

305 L{a_f2Tuple}), overriding the mean radius C{R1} 

306 of this point's datum.. 

307 @kwarg wrap: If C{True}, wrap or I{normalize} and unroll the 

308 B{C{other}} and angular distance (C{bool}). 

309 

310 @return: Distance (C{meter}, same units as B{C{radius}}). 

311 

312 @raise TypeError: The B{C{other}} point is not L{LatLon}. 

313 

314 @raise ValueError: Invalid B{C{radius}}. 

315 

316 @example: 

317 

318 >>> p = LatLon(52.205, 0.119) 

319 >>> q = LatLon(48.857, 2.351); 

320 >>> d = p.distanceTo(q) # 404300 

321 ''' 

322 p = self.others(other) 

323 if wrap: 

324 p = _Wrap.point(p) 

325 a = self._N_vector.angleTo(p._N_vector, wrap=wrap) 

326 d = self.datum if radius is None else _spherical_datum(radius) 

327 return fabs(a) * d.ellipsoid.R1 # see .utily.radians2m 

328 

329 @Property_RO 

330 def Ecef(self): 

331 '''Get the ECEF I{class} (L{EcefVeness}), I{lazily}. 

332 ''' 

333 return _MODS.ecef.EcefVeness 

334 

335 @deprecated_method 

336 def equals(self, other, eps=None): # PYCHOK no cover 

337 '''DEPRECATED, use method L{isequalTo}. 

338 ''' 

339 return self.isequalTo(other, eps=eps) 

340 

341 def isequalTo(self, other, eps=None): 

342 '''Compare this point with an other point. 

343 

344 @arg other: The other point (L{LatLon}). 

345 @kwarg eps: Optional margin (C{float}). 

346 

347 @return: C{True} if points are identical, including 

348 datum, I{ignoring height}, C{False} otherwise. 

349 

350 @raise TypeError: The B{C{other}} point is not L{LatLon}. 

351 

352 @raise ValueError: Invalid B{C{eps}}. 

353 

354 @see: Method C{isequalTo3} to include I{height}. 

355 

356 @example: 

357 

358 >>> p = LatLon(52.205, 0.119) 

359 >>> q = LatLon(52.205, 0.119) 

360 >>> e = p.isequalTo(q) # True 

361 ''' 

362 return self.datum == self.others(other).datum and \ 

363 _MODS.formy._isequalTo(self, other, eps=eps) 

364 

365# def greatCircle(self, bearing): 

366# '''Return the great circle heading on the given bearing 

367# from this point. 

368# 

369# Direction of vector is such that initial bearing vector 

370# b = c × p, where p is representing this point. 

371# 

372# @arg bearing: Bearing from this point (compass C{degrees360}). 

373# 

374# @return: N-vector representing great circle (L{Nvector}). 

375# 

376# @example: 

377# 

378# >>> p = LatLon(53.3206, -1.7297) 

379# >>> g = p.greatCircle(96.0) 

380# >>> g.toStr() # '(-0.794, 0.129, 0.594)' 

381# ''' 

382# a, b, _ = self.philamheight 

383# t = radians(bearing) 

384# 

385# sa, ca, sb, cb, st, ct = sincos2_(a, b, t) 

386# return self._xnamed(Nvector(sb * ct - sa * cb * st, 

387# -cb * ct - sa * sb * st, 

388# ca * st) 

389 

390# def initialBearingTo(self, other, wrap=False): 

391# '''Return the initial bearing (forward azimuth) from 

392# this to an other point. 

393# 

394# @arg other: The other point (L{LatLon}). 

395# @kwarg wrap: If C{True}, wrap or I{normalize} 

396# and unroll the B{C{other}} (C{bool}). 

397# 

398# @return: Initial bearing (compass C{degrees360}). 

399# 

400# @raise TypeError: The B{C{other}} point is not L{LatLon}. 

401# 

402# @example: 

403# 

404# >>> p1 = LatLon(52.205, 0.119) 

405# >>> p2 = LatLon(48.857, 2.351) 

406# >>> b = p1.bearingTo(p2) # 156.2 

407# ''' 

408# p = self.others(other) 

409# if wrap: 

410# p = _Wrap.point(p) 

411# v1 = self.toNvector() 

412# 

413# gc1 = v1.cross(p.toNvector()) # gc through v1 & v2 

414# gc2 = v1.cross(_NP3) # gc through v1 & North pole 

415# 

416# # bearing is (signed) angle between gc1 & gc2 

417# return degrees360(gc1.angleTo(gc2, vSign=v1)) 

418 

419 def intermediateTo(self, other, fraction, height=None, wrap=False): 

420 '''Return the point at given fraction between this and 

421 an other point. 

422 

423 @arg other: The other point (L{LatLon}). 

424 @arg fraction: Fraction between both points (C{scalar}, 

425 0.0 at this to 1.0 at the other point. 

426 @kwarg height: Optional height, overriding the fractional 

427 height (C{meter}). 

428 @kwarg wrap: If C{True}, wrap or I{normalize} the 

429 B{C{other}} point (C{bool}). 

430 

431 @return: Intermediate point (L{LatLon}). 

432 

433 @raise TypeError: The B{C{other}} point is not L{LatLon}. 

434 

435 @example: 

436 

437 >>> p = LatLon(52.205, 0.119) 

438 >>> q = LatLon(48.857, 2.351) 

439 >>> p = p.intermediateTo(q, 0.25) # 51.3721°N, 000.7073°E 

440 ''' 

441 p = self.others(other) 

442 if wrap: 

443 p = _Wrap.point(p) 

444 f = Scalar(fraction=fraction) 

445 h = self._havg(other, f=f, h=height) 

446 i = self.toNvector().intermediateTo(p.toNvector(), f) 

447 return i.toLatLon(height=h, LatLon=self.classof) # Nvector(i.x, i.y, i.z).toLatLon(...) 

448 

449 @Property_RO 

450 def _rotation3(self): 

451 '''(INTERNAL) Get the rotation matrix from n-vector coordinate frame axes. 

452 ''' 

453 nv = self.toNvector() # local (n-vector) coordinate frame 

454 

455 dv = nv.negate() # down, opposite to n-vector 

456 ev = NorthPole.cross(nv, raiser=_pole_).unit() # east, pointing perpendicular to the plane 

457 nv = ev.cross(dv) # north, by right hand rule 

458 return nv, ev, dv # matrix rows 

459 

460 def toCartesian(self, **Cartesian_and_kwds): # PYCHOK Cartesian=Cartesian, datum=None 

461 '''Convert this point to an C{Nvector}-based geodetic point. 

462 

463 @kwarg Cartesian_and_kwds: Optional L{Cartesian}, B{C{datum}} and other 

464 keyword arguments. Use C{B{Cartesian}=...} 

465 to override this L{Cartesian} class or specify 

466 C{B{Cartesian}=None}. 

467 

468 @return: The geodetic point (L{Cartesian}) or if B{C{Cartesian}} is set 

469 to C{None}, an L{Ecef9Tuple}C{(x, y, z, lat, lon, height, C, M, 

470 datum)} with C{C} and C{M} if available. 

471 

472 @raise TypeError: Invalid B{C{Cartesian}} or other B{C{Cartesian_and_kwds}}. 

473 ''' 

474 kwds = _xkwds(Cartesian_and_kwds, Cartesian=Cartesian, datum=self.datum) 

475 return LatLonEllipsoidalBase.toCartesian(self, **kwds) 

476 

477 def toNvector(self, **Nvector_and_kwds): # PYCHOK signature 

478 '''Convert this point to L{Nvector} components, I{including height}. 

479 

480 @kwarg Nvector_and_kwds: Optional L{Nvector}, B{C{datum}} and other 

481 keyword arguments. Use C{B{Nvector}=...} 

482 to override this L{Nvector} class or specify 

483 C{B{Nvector}=None}. 

484 

485 @return: The C{n-vector} components (L{Nvector}) or if B{C{Nvector}} 

486 is set to C{None}, a L{Vector4Tuple}C{(x, y, z, h)}. 

487 

488 @raise TypeError: Invalid B{C{Nvector}} or other B{C{Nvector_and_kwds}}. 

489 

490 @example: 

491 

492 >>> p = LatLon(45, 45) 

493 >>> n = p.toNvector() 

494 >>> n.toStr() # [0.50, 0.50, 0.70710] 

495 ''' 

496 kwds = _xkwds(Nvector_and_kwds, Nvector=Nvector, datum=self.datum) 

497 return LatLonNvectorBase.toNvector(self, **kwds) 

498 

499 

500_Nvll = LatLon(0, 0, name=_Nv00_) # reference instance (L{LatLon}) 

501 

502 

503class Nvector(NvectorBase): 

504 '''An n-vector is a position representation using a (unit) vector 

505 normal to the earth ellipsoid. Unlike lat-/longitude points, 

506 n-vectors have no singularities or discontinuities. 

507 

508 For many applications, n-vectors are more convenient to work 

509 with than other position representations like lat-/longitude, 

510 earth-centred earth-fixed (ECEF) vectors, UTM coordinates, etc. 

511 

512 Note commonality with L{sphericalNvector.Nvector}. 

513 ''' 

514 _datum = _WGS84 # default datum (L{Datum}) 

515 

516 def __init__(self, x_xyz, y=None, z=None, h=0, datum=None, ll=None, name=NN): 

517 '''New n-vector normal to the earth's surface. 

518 

519 @arg x_xyz: X component of vector (C{scalar}) or (3-D) vector 

520 (C{Nvector}, L{Vector3d}, L{Vector3Tuple} or 

521 L{Vector4Tuple}). 

522 @kwarg y: Y component of vector (C{scalar}), ignored if B{C{x_xyz}} 

523 is not C{scalar}, otherwise same units as B{C{x_xyz}}. 

524 @kwarg z: Z component of vector (C{scalar}), ignored if B{C{x_xyz}} 

525 is not C{scalar}, otherwise same units as B{C{x_xyz}}. 

526 @kwarg h: Optional height above model surface (C{meter}). 

527 @kwarg datum: Optional datum this n-vector is defined in 

528 (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or 

529 L{a_f2Tuple}). 

530 @kwarg ll: Optional, original latlon (C{LatLon}). 

531 @kwarg name: Optional name (C{str}). 

532 

533 @raise TypeError: If B{C{datum}} is not a L{Datum}. 

534 

535 @example: 

536 

537 >>> from ellipsoidalNvector import Nvector 

538 >>> v = Nvector(0.5, 0.5, 0.7071, 1) 

539 >>> v.toLatLon() # 45.0°N, 045.0°E, +1.00m 

540 ''' 

541 NvectorBase.__init__(self, x_xyz, y=y, z=z, h=h, ll=ll, name=name) 

542 if datum not in (None, self._datum): 

543 self._datum = _ellipsoidal_datum(datum, name=name) 

544 

545 @Property_RO 

546 def datum(self): 

547 '''Get this n-vector's datum (L{Datum}). 

548 ''' 

549 return self._datum 

550 

551 def toCartesian(self, **Cartesian_and_kwds): # PYCHOK Cartesian=Cartesian 

552 '''Convert this n-vector to C{Nvector}-based cartesian (ECEF) coordinates. 

553 

554 @kwarg Cartesian_and_kwds: Optional L{Cartesian}, B{C{h}}, B{C{datum}} and 

555 other keyword arguments. Use C{B{Cartesian}=...} 

556 to override this L{Cartesian} class or specify 

557 C{B{Cartesian} is None}. 

558 

559 @return: The cartesian point (L{Cartesian}) or if B{C{Cartesian}} is set 

560 to C{None}, an L{Ecef9Tuple}C{(x, y, z, lat, lon, height, C, M, 

561 datum)} with C{C} and C{M} if available. 

562 

563 @raise TypeError: Invalid B{C{Cartesian_and_kwds}}. 

564 

565 @example: 

566 

567 >>> v = Nvector(0.5, 0.5, 0.7071) 

568 >>> c = v.toCartesian() # [3194434, 3194434, 4487327] 

569 >>> p = c.toLatLon() # 45.0°N, 45.0°E 

570 ''' 

571 kwds = _xkwds(Cartesian_and_kwds, h=self.h, Cartesian=Cartesian, 

572 datum=self.datum) 

573 return NvectorBase.toCartesian(self, **kwds) # class or .classof 

574 

575 def toLatLon(self, **LatLon_and_kwds): # PYCHOK height=None, LatLon=LatLon 

576 '''Convert this n-vector to an C{Nvector}-based geodetic point. 

577 

578 @kwarg LatLon_and_kwds: Optional L{LatLon}, B{C{height}}, B{C{datum}} 

579 and other keyword arguments. Use C{B{LatLon}=...} 

580 to override this L{LatLon} class or specify 

581 C{B{LatLon} is None}. 

582 

583 @return: The geodetic point (L{LatLon}) or if B{C{LatLon}} is set 

584 to C{None}, an L{Ecef9Tuple}C{(x, y, z, lat, lon, height, 

585 C, M, datum)} with C{C} and C{M} if available. 

586 

587 @raise TypeError: Invalid B{C{LatLon_and_kwds}}. 

588 

589 @example: 

590 

591 >>> v = Nvector(0.5, 0.5, 0.7071) 

592 >>> p = v.toLatLon() # 45.0°N, 45.0°E 

593 ''' 

594 kwds = _xkwds(LatLon_and_kwds, height=self.h, datum=self.datum, LatLon=LatLon) 

595 return NvectorBase.toLatLon(self, **kwds) # class or .classof 

596 

597 def unit(self, ll=None): 

598 '''Normalize this vector to unit length. 

599 

600 @kwarg ll: Optional, original latlon (C{LatLon}). 

601 

602 @return: Normalised vector (L{Nvector}). 

603 ''' 

604 u = NvectorBase.unit(self, ll=ll) 

605 if u.datum != self.datum: 

606 u._update(False, datum=self.datum) 

607 return u 

608 

609 

610def meanOf(points, datum=_WGS84, height=None, wrap=False, 

611 **LatLon_and_kwds): 

612 '''Compute the geographic mean of several points. 

613 

614 @arg points: Points to be averaged (L{LatLon}[]). 

615 @kwarg datum: Optional datum to use (L{Datum}). 

616 @kwarg height: Optional height at mean point, overriding 

617 the mean height (C{meter}). 

618 @kwarg wrap: If C{True}, wrap or I{normalize} B{C{points}} 

619 (C{bool}). 

620 @kwarg LatLon_and_kwds: Optional B{C{LatLon}} class to return 

621 the mean points and overriding this L{LatLon} 

622 (or C{None}) and additional B{C{LatLon}} 

623 keyword arguments, ignored if C{B{LatLon} 

624 is None}. 

625 

626 @return: Geographic mean point and mean height (B{C{LatLon}}) 

627 or if B{C{LatLon}} is C{None}, an L{Ecef9Tuple}C{(x, 

628 y, z, lat, lon, height, C, M, datum)} with C{C} and 

629 C{M} if available. 

630 

631 @raise ValueError: Insufficient number of B{C{points}}. 

632 ''' 

633 Ps = _Nvll.PointsIter(points, wrap=wrap) 

634 # geographic mean 

635 m = sumOf(p._N_vector for p in Ps.iterate(closed=False)) 

636 kwds = _xkwds(LatLon_and_kwds, height=height, datum=datum, 

637 LatLon=LatLon, name=meanOf.__name__) 

638 return m.toLatLon(**kwds) 

639 

640 

641def nearestOn(point, point1, point2, within=True, height=None, wrap=False, 

642 equidistant=None, tol=_TOL_M, LatLon=LatLon, **LatLon_kwds): 

643 '''I{Iteratively} locate the closest point on the geodesic between 

644 two other (ellipsoidal) points. 

645 

646 @arg point: Reference point (C{LatLon}). 

647 @arg point1: Start point of the geodesic (C{LatLon}). 

648 @arg point2: End point of the geodesic (C{LatLon}). 

649 @kwarg within: If C{True} return the closest point I{between} 

650 B{C{point1}} and B{C{point2}}, otherwise the 

651 closest point elsewhere on the geodesic (C{bool}). 

652 @kwarg height: Optional height for the closest point (C{meter}, 

653 conventionally) or C{None} or C{False} for the 

654 interpolated height. If C{False}, the closest 

655 takes the heights of the points into account. 

656 @kwarg wrap: If C{True}, wrap or I{normalize} and unroll I{only} 

657 B{C{point1}} and B{C{point2}} (C{bool}). 

658 @kwarg equidistant: An azimuthal equidistant projection (I{class} 

659 or function L{pygeodesy.equidistant}) or C{None} 

660 for the preferred C{B{point}.Equidistant}. 

661 @kwarg tol: Convergence tolerance (C{meter}). 

662 @kwarg LatLon: Optional class to return the closest point 

663 (L{LatLon}) or C{None}. 

664 @kwarg LatLon_kwds: Optional, additional B{C{LatLon}} keyword 

665 arguments, ignored if C{B{LatLon} is None}. 

666 

667 @return: Closest point, a B{C{LatLon}} instance or if C{B{LatLon} 

668 is None}, a L{LatLon4Tuple}C{(lat, lon, height, datum)}. 

669 

670 @raise ImportError: Package U{geographiclib 

671 <https://PyPI.org/project/geographiclib>} 

672 not installed or not found. 

673 

674 @raise TypeError: Invalid or non-ellipsoidal B{C{point}}, B{C{point1}} 

675 or B{C{point2}} or invalid B{C{equidistant}}. 

676 

677 @raise ValueError: No convergence for the B{C{tol}}. 

678 

679 @see: U{The B{ellipsoidal} case<https://GIS.StackExchange.com/questions/48937/ 

680 calculating-intersection-of-two-circles>} and U{Karney's paper 

681 <https://ArXiv.org/pdf/1102.1215.pdf>}, pp 20-21, section B{14. MARITIME 

682 BOUNDARIES} for more details about the iteration algorithm. 

683 ''' 

684 return _nearestOn(point, point1, point2, within=within, height=height, wrap=wrap, 

685 equidistant=equidistant, tol=tol, LatLon=LatLon, **LatLon_kwds) 

686 

687 

688def sumOf(nvectors, Vector=Nvector, h=None, **Vector_kwds): 

689 '''Return the vectorial sum of two or more n-vectors. 

690 

691 @arg nvectors: Vectors to be added (L{Nvector}[]). 

692 @kwarg Vector: Optional class for the vectorial sum (L{Nvector}). 

693 @kwarg h: Optional height, overriding the mean height (C{meter}). 

694 @kwarg Vector_kwds: Optional, additional B{C{Vector}} keyword 

695 arguments, ignored if C{B{Vector} is None}. 

696 

697 @return: Vectorial sum (B{C{Vector}}). 

698 

699 @raise VectorError: No B{C{nvectors}}. 

700 ''' 

701 return _sumOf(nvectors, Vector=Vector, h=h, **Vector_kwds) 

702 

703 

704@deprecated_function 

705def toNed(distance, bearing, elevation, Ned=Ned, name=NN): 

706 '''DEPRECATED, use L{pygeodesy.Aer}C{(bearing, elevation, 

707 distance).xyzLocal.toNed(B{Ned}, name=B{name})} or 

708 L{XyzLocal}C{(pygeodesy.Aer(bearing, elevation, 

709 distance)).toNed(B{Ned}, name=B{name})}. 

710 

711 Create an NED vector from distance, bearing and elevation 

712 (in local coordinate system). 

713 

714 @arg distance: NED vector length (C{meter}). 

715 @arg bearing: NED vector bearing (compass C{degrees360}). 

716 @arg elevation: NED vector elevation from local coordinate 

717 frame horizontal (C{degrees}). 

718 @kwarg Ned: Optional class to return the NED (C{Ned}) or 

719 C{None}. 

720 @kwarg name: Optional name (C{str}). 

721 

722 @return: An NED vector equivalent to this B{C{distance}}, 

723 B{C{bearing}} and B{C{elevation}} (DEPRECATED L{Ned}) 

724 or a DEPRECATED L{Ned3Tuple}C{(north, east, down)} 

725 if C{B{Ned} is None}. 

726 

727 @raise ValueError: Invalid B{C{distance}}, B{C{bearing}} 

728 or B{C{elevation}}. 

729 ''' 

730 if True: # use new Aer class 

731 n, e, d, _ = _Aer(bearing, elevation, distance).xyz4 

732 else: # DEPRECATED 

733 d = Distance(distance) 

734 

735 sb, cb, se, ce = sincos2d_(Bearing(bearing), 

736 Height(elevation=elevation)) 

737 n = cb * d * ce 

738 e = sb * d * ce 

739 d *= se 

740 

741 r = _MODS.deprecated.Ned3Tuple(n, e, -d) if Ned is None else \ 

742 Ned(n, e, -d) 

743 return _xnamed(r, name) 

744 

745 

746__all__ += _ALL_OTHER(Cartesian, LatLon, Ned, Nvector, # classes 

747 meanOf, sumOf, toNed) 

748 

749# **) MIT License 

750# 

751# Copyright (C) 2016-2023 -- mrJean1 at Gmail -- All Rights Reserved. 

752# 

753# Permission is hereby granted, free of charge, to any person obtaining a 

754# copy of this software and associated documentation files (the "Software"), 

755# to deal in the Software without restriction, including without limitation 

756# the rights to use, copy, modify, merge, publish, distribute, sublicense, 

757# and/or sell copies of the Software, and to permit persons to whom the 

758# Software is furnished to do so, subject to the following conditions: 

759# 

760# The above copyright notice and this permission notice shall be included 

761# in all copies or substantial portions of the Software. 

762# 

763# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS 

764# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 

765# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL 

766# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR 

767# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, 

768# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR 

769# OTHER DEALINGS IN THE SOFTWARE.