Coverage for pygeodesy/ltp.py: 95%

1 

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

3 

4u'''I{Local Tangent Plane} (LTP) and I{local} cartesian coordinates aka U{topocentric 

5<https://ssd.JPL.NASA.gov/glossary/topocentric.html>} coordinates. 

6 

7I{Local cartesian} class L{LocalCartesian}, I{local tangent plane} classes L{Ltp} and 

8L{LqRD} and classes L{Attitude}, L{Frustum} and L{LocalError}. 

9 

10@see: U{Local tangent plane coordinates<https://WikiPedia.org/wiki/Local_tangent_plane_coordinates>} 

11 and class L{LocalCartesian}, transcoded from I{Charles Karney}'s C++ classU{LocalCartesian 

12 <https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1LocalCartesian.html>}. 

13''' 

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

15from __future__ import division as _; del _ # noqa: E702 ; 

16 

17from pygeodesy.basics import map1, map2, _xinstanceof, _xsubclassof 

18from pygeodesy.constants import EPS, INT0, _umod_360, _0_0, _0_01, _0_5, _1_0, \ 

19 _2_0, _8_0, _60_0, _90_0, _100_0, _180_0, _3600_0, \ 

20 _N_1_0 # PYCHOK used! 

21from pygeodesy.datums import Datums, _WGS84 

22# from pygeodesy.dms import parseDMS # from .units 

23from pygeodesy.ecef import _EcefBase, EcefKarney, Ecef9Tuple, _llhn4, _xyzn4 

24from pygeodesy.errors import _ValueError, _xattr, _xkwds, _xkwds_get, _xkwds_pop2 

25from pygeodesy.fmath import fabs, fdot, fdot_ 

26from pygeodesy.fsums import fsumf_ 

27# from pygeodesy.internals import typename # from .basics 

28from pygeodesy.interns import _0_, _COMMASPACE_, _ecef_, _height_, _invalid_, \ 

29 _lat0_, _lon0_, _M_, _name_, _too_ 

30# from pygeodesy.lazily import _ALL_LAZY # from vector3d 

31from pygeodesy.ltpTuples import Attitude4Tuple, Footprint5Tuple, Local9Tuple, \ 

32 _XyzLocals4, _XyzLocals5, Xyz4Tuple 

33from pygeodesy.named import _name__, _name2__, _NamedBase 

34from pygeodesy.namedTuples import LatLon3Tuple, LatLon4Tuple, Vector3Tuple, \ 

35 Bounds4Tuple, RD4Tuple # PYCHOK used! 

36from pygeodesy.props import deprecated_property_RO, Property, Property_RO, \ 

37 property_doc_, property_ROver, _update_all 

38from pygeodesy.streprs import Fmt, strs 

39from pygeodesy.units import Bearing, Degrees, _isHeight, Meter, parseDMS 

40from pygeodesy.utily import cotd, _loneg, sincos2d, sincos2d_, tand, tand_, \ 

41 wrap180, wrap360 

42from pygeodesy.vector3d import _ALL_LAZY, Vector3d 

43 

44# from math import fabs, floor as _floor # from .fmath, .fsums 

45 

46__all__ = _ALL_LAZY.ltp 

47__version__ = '26.10.08' 

48 

49_GRS80 = Datums.GRS80 

50_height0_ = _height_ + _0_ 

51_narrow_ = 'narrow' 

52_wide_ = 'wide' 

53# del Datums 

54 

55 

56class Attitude(_NamedBase): 

57 '''The pose of a plane or camera in space. 

58 ''' 

59 _alt = Meter( alt =_0_0) 

60 _roll = Degrees(roll=_0_0) 

61 _tilt = Degrees(tilt=_0_0) 

62 _yaw = Bearing(yaw =_0_0) 

63 

64 def __init__(self, alt_attitude=INT0, tilt=INT0, yaw=INT0, roll=INT0, **name): 

65 '''New L{Attitude}. 

66 

67 @kwarg alt_attitude: Altitude (C{meter}) above earth or a previous attitude 

68 (L{Attitude} or L{Attitude4Tuple}) with the C{B{alt}itude}, 

69 B{C{tilt}}, B{C{yaw}} and B{C{roll}}. 

70 @kwarg tilt: Pitch, elevation from horizontal (C{degrees180}), negative down 

71 (clockwise rotation along and around the x- or East axis), iff 

72 B{C{alt_attitude}} is C{meter}, ignored otherwise. 

73 @kwarg yaw: Bearing, heading (compass C{degrees360}), clockwise from North 

74 (counter-clockwise rotation along and around the z- or Up axis) 

75 iff B{C{alt_attitude}} is C{meter}, ignored otherwise. 

76 @kwarg roll: Roll, bank (C{degrees180}), positive to the right and down 

77 (clockwise rotation along and around the y- or North axis), iff 

78 B{C{alt_attitude}} is C{meter}, ignored otherwise. 

79 @kwarg name: Optional C{B{name}=NN} C{str}). 

80 

81 @raise AttitudeError: Invalid B{C{alt_attitude}}, B{C{tilt}}, B{C{yaw}} or 

82 B{C{roll}}. 

83 

84 @see: U{Principal axes<https://WikiPedia.org/wiki/Aircraft_principal_axes>} and 

85 U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}. 

86 ''' 

87 if _isHeight(alt_attitude): 

88 t = Attitude4Tuple(alt_attitude, tilt, yaw, roll) 

89 else: 

90 try: 

91 t = alt_attitude.atyr 

92 except AttributeError: 

93 raise AttitudeError(alt=alt_attitude, tilt=tilt, yaw=yaw, rol=roll) 

94 for n, v in t.items(): 

95 if v: 

96 setattr(self, n, v) 

97 n = _name__(name, _or_nameof=t) 

98 if n: 

99 self.name = n 

100 

101 @property_doc_(' altitude above earth in C{meter}.') 

102 def alt(self): 

103 return self._alt 

104 

105 @alt.setter # PYCHOK setter! 

106 def alt(self, alt): # PYCHOK no cover 

107 a = Meter(alt=alt, Error=AttitudeError) 

108 if self._alt != a: 

109 _update_all(self) 

110 self._alt = a 

111 

112 altitude = alt 

113 

114 @Property_RO 

115 def atyr(self): 

116 '''Return this attitude's alt[itude], tilt, yaw and roll as an L{Attitude4Tuple}. 

117 ''' 

118 return Attitude4Tuple(self.alt, self.tilt, self.yaw, self.roll, name=self.name) 

119 

120 @Property_RO 

121 def matrix(self): 

122 '''Get the 3x3 rotation matrix C{R(yaw)·R(tilt)·R(roll)}, aka I{ZYX} (C{float}, row-order). 

123 

124 @see: Matrix M of case 10 in U{Appendix A 

125 <https://ntrs.NASA.gov/api/citations/19770019231/downloads/19770019231.pdf>}. 

126 ''' 

127 # to follow the definitions of rotation angles alpha, beta and gamma: 

128 # negate yaw since yaw is counter-clockwise around the z-axis, swap 

129 # tilt and roll since tilt is around the x- and roll around the y-axis 

130 sa, ca, sb, cb, sg, cg = sincos2d_(-self.yaw, self.roll, self.tilt) # PYLINT E1130 

131 return ((ca * cb, fdot_(ca, sb * sg, -sa, cg), fdot_(ca, sb * cg, sa, sg)), # PYLINT E1130 

132 (sa * cb, fdot_(sa, sb * sg, ca, cg), fdot_(sa, sb * cg, -ca, sg)), # PYLINT E1130 

133 ( -sb, cb * sg, cb * cg)) # PYLINT E1130 

134 

135 @property_doc_(' roll/bank in C{degrees180}, positive to the right and down.') 

136 def roll(self): 

137 return self._roll 

138 

139 @roll.setter # PYCHOK setter! 

140 def roll(self, roll): 

141 r = Degrees(roll=roll, wrap=wrap180, Error=AttitudeError) 

142 if self._roll != r: 

143 _update_all(self) 

144 self._roll = r 

145 

146 bank = roll 

147 

148 def rotate(self, x_xyz, y=None, z=None, Vector=None, **name_Vector_kwds): 

149 '''Transform a (local) cartesian by this attitude's matrix. 

150 

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

152 L{Vector3d} or L{Vector3Tuple}). 

153 @kwarg y: Y component of vector (C{scalar}), same units as C{scalar} B{C{x}}, 

154 ignored otherwise. 

155 @kwarg z: Z component of vector (C{scalar}), same units as C{sclar} B{C{x}}, 

156 ignored otherwise. 

157 @kwarg Vector: Class to return transformed point (C{Cartesian}, L{Vector3d} 

158 or C{Vector3Tuple}) or C{None}. 

159 @kwarg name_Vector_kwds: Optional C{B{name}=NN} (C{str}) and optionally, 

160 additional B{C{Vector}} keyword arguments, ignored if C{B{Vector} 

161 is None}. 

162 

163 @return: A named B{C{Vector}} instance or if C{B{Vector} is None}, 

164 a named L{Vector3Tuple}C{(x, y, z)}. 

165 

166 @raise AttitudeError: Invalid B{C{x_xyz}}, B{C{y}} or B{C{z}}. 

167 

168 @raise TypeError: Invalid B{C{Vector}} or B{C{name_Vector_kwds}} item. 

169 

170 @see: U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}. 

171 ''' 

172 try: 

173 try: 

174 xyz = map2(float, x_xyz.xyz3) 

175 except AttributeError: 

176 xyz = map1(float, x_xyz, y, z) 

177 except (TypeError, ValueError) as x: 

178 raise AttitudeError(x_xyz=x_xyz, y=y, z=z, cause=x) 

179 

180 x, y, z = (fdot(r, *xyz) for r in self.matrix) 

181 n, kwds = _name2__(name_Vector_kwds, _or_nameof=self) 

182 return Vector3Tuple(x, y, z, name=n) if Vector is None else \ 

183 Vector(x, y, z, name=n, **kwds) 

184 

185 @property_doc_(' tilt/pitch/elevation from horizontal in C{degrees180}, negative down.') 

186 def tilt(self): 

187 return self._tilt 

188 

189 @tilt.setter # PYCHOK setter! 

190 def tilt(self, tilt): 

191 t = Degrees(tilt=tilt, wrap=wrap180, Error=AttitudeError) 

192 if self._tilt != t: 

193 _update_all(self) 

194 self._tilt = t 

195 

196 elevation = pitch = tilt 

197 

198 def toStr(self, prec=6, sep=_COMMASPACE_, **unused): # PYCHOK signature 

199 '''Format this attitude as string. 

200 

201 @kwarg prec: The C{float} precision, number of decimal digits (0..9). 

202 Trailing zero decimals are stripped for B{C{prec}} values 

203 of 1 and above, but kept for negative B{C{prec}} values. 

204 @kwarg sep: Separator to join (C{str}). 

205 

206 @return: This attitude (C{str}). 

207 ''' 

208 return self.atyr.toStr(prec=prec, sep=sep) 

209 

210 @Property_RO 

211 def tyr3d(self): 

212 '''Get this attitude's (3-D) directional vector (L{Vector3d}). 

213 

214 @see: U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}. 

215 ''' 

216 def _r2d(r): 

217 return fsumf_(_N_1_0, *r) 

218 

219 return Vector3d(*map(_r2d, self.matrix), name__=tyr3d) 

220 

221 @property_doc_(' yaw/bearing/heading in compass C{degrees360}, clockwise from North.') 

222 def yaw(self): 

223 return self._yaw 

224 

225 @yaw.setter # PYCHOK setter! 

226 def yaw(self, yaw): 

227 y = Bearing(yaw=yaw, Error=AttitudeError) 

228 if self._yaw != y: 

229 _update_all(self) 

230 self._yaw = y 

231 

232 bearing = heading = yaw # azimuth 

233 

234 

235class AttitudeError(_ValueError): 

236 '''An L{Attitude} or L{Attitude4Tuple} issue. 

237 ''' 

238 pass # PYLINT W0107 

239 

240 

241class Frustum(_NamedBase): 

242 '''A rectangular pyramid, typically representing a camera's I{field-of-view} 

243 (fov) and the intersection with (or projection to) a I{local tangent plane}. 

244 

245 @see: U{Viewing frustum<https://WikiPedia.org/wiki/Viewing_frustum>}. 

246 ''' 

247 _h_2 = _0_0 # half hfov in degrees 

248 _ltp = None # local tangent plane 

249 _tan_h_2 = _0_0 # tan(_h_2) 

250 _v_2 = _0_0 # half vfov in degrees 

251 

252 def __init__(self, hfov, vfov, ltp=None, **name): 

253 '''New L{Frustum}. 

254 

255 @arg hfov: Horizontal field-of-view (C{degrees180}). 

256 @arg vfov: Vertical field-of-view (C{degrees180}). 

257 @kwarg ltp: Optional I{local tangent plane} (L{Ltp}). 

258 @kwarg name: Optional C{B{name}=NN} (C{str}). 

259 

260 @raise LocalError: Invalid B{C{hfov}} or B{C{vfov}}. 

261 ''' 

262 self._h_2 = h = _fov_2(hfov=hfov) 

263 self._v_2 = _fov_2(vfov=vfov) 

264 

265 self._tan_h_2 = tand(h, hfov_2=h) 

266 

267 if ltp: 

268 self._ltp = _xLtp(ltp) 

269 if name: 

270 self.name # PYCHOK effect 

271 

272 def footprint5(self, alt_attitude, tilt=0, yaw=0, roll=0, z=_0_0, ltp=None, **name): # MCCABE 15 

273 '''Compute the center and corners of the intersection with (or projection 

274 to) the I{local tangent plane} (LTP). 

275 

276 @arg alt_attitude: An altitude (C{meter}) above I{local tangent plane} or 

277 an attitude (L{Attitude} or L{Attitude4Tuple}) with the 

278 C{B{alt}itude}, B{C{tilt}}, B{C{yaw}} and B{C{roll}}. 

279 @kwarg tilt: Pitch, elevation from horizontal (C{degrees}), negative down 

280 (clockwise rotation along and around the x- or East axis) iff 

281 B{C{alt_attitude}} is C{meter}, ignored otherwise. 

282 @kwarg yaw: Bearing, heading (compass C{degrees}), clockwise from North 

283 (counter-clockwise rotation along and around the z- or Up axis) 

284 iff B{C{alt_attitude}} is C{meter}, ignored otherwise. 

285 @kwarg roll: Roll, bank (C{degrees}), positive to the right and down 

286 (clockwise rotation along and around the y- or North axis) iff 

287 B{C{alt_attitude}} is C{meter}, ignored otherwise. 

288 @kwarg z: Optional height of the footprint (C{meter}) above I{local tangent plane}. 

289 @kwarg ltp: The I{local tangent plane} (L{Ltp}), overriding this 

290 frustum's C{ltp}. 

291 @kwarg name: Optional C{B{name}=NN} (C{str}). 

292 

293 @return: A L{Footprint5Tuple}C{(center, upperleft, upperight, loweright, 

294 lowerleft)} with the C{center} and 4 corners, each an L{Xyz4Tuple}. 

295 

296 @raise TypeError: Invalid B{C{ltp}}. 

297 

298 @raise UnitError: Invalid B{C{altitude}}, B{C{tilt}}, B{C{roll}} or B{C{z}}. 

299 

300 @raise ValueError: If B{C{altitude}} too low, B{C{z}} too high or B{C{tilt}} 

301 or B{C{roll}} -including B{C{vfov}} respectively B{C{hfov}}- 

302 over the horizon. 

303 

304 @see: U{Principal axes<https://WikiPedia.org/wiki/Aircraft_principal_axes>}. 

305 ''' 

306 def _xy2(a, e, h_2, tan_h_2, r): 

307 # left and right corners, or swapped 

308 if r < EPS: # no roll 

309 r = a * tan_h_2 

310 l = -r # noqa: E741 l is ell 

311 else: # roll 

312 r, l = tand_(r - h_2, r + h_2, roll_hfov=r) # noqa: E741 l is ell 

313 r *= -a # negate right positive 

314 l *= -a # noqa: E741 l is ell 

315 y = a * cotd(e, tilt_vfov=e) 

316 return (l, y), (r, y) 

317 

318 def _xyz5(b, xy5, z, ltp): 

319 # rotate (x, y)'s by bearing, clockwise 

320 sc = sincos2d(b) 

321 for x, y in xy5: 

322 yield Xyz4Tuple(fdot(sc, x, y), 

323 fdot(sc, -x, y), z, ltp) 

324 

325 try: 

326 a, t, y, r = alt_attitude.atyr 

327 except AttributeError: 

328 a, t, y, r = alt_attitude, tilt, yaw, roll 

329 

330 a = Meter(altitude=a) 

331 if a < EPS: # too low 

332 raise _ValueError(altitude=a) 

333 if z: # PYCHOK no cover 

334 z = Meter(z=z) 

335 a -= z 

336 if a < EPS: # z above a 

337 raise _ValueError(altitude_z=a) 

338 else: 

339 z = _0_0 

340 

341 b = Degrees(yaw=y, wrap=wrap360) # bearing 

342 e = -Degrees(tilt=t, wrap=wrap180) # elevation, pitch 

343 if not EPS < e < _180_0: 

344 raise _ValueError(tilt=t) 

345 if e > _90_0: 

346 e = _loneg(e) 

347 b = _umod_360(b + _180_0) 

348 

349 r = Degrees(roll=r, wrap=wrap180) # roll center 

350 x = (-a * tand(r, roll=r)) if r else _0_0 

351 y = a * cotd(e, tilt=t) # ground range 

352 if fabs(y) < EPS: 

353 y = _0_0 

354 

355 v, h, t = self._v_2, self._h_2, self._tan_h_2 

356 # center and corners, clockwise from upperleft, rolled 

357 xy5 = ((x, y),) + _xy2(a, e - v, h, t, r) \ 

358 + _xy2(a, e + v, -h, -t, r) # swapped 

359 # turn center and corners by yaw, clockwise 

360 p = self.ltp if ltp is None else ltp # None OK 

361 return Footprint5Tuple(_xyz5(b, xy5, z, p), **name) # *_xyz5 

362 

363 @Property_RO 

364 def hfov(self): 

365 '''Get the horizontal C{fov} (C{degrees}). 

366 ''' 

367 return Degrees(hfov=self._h_2 * _2_0) 

368 

369 @Property_RO 

370 def ltp(self): 

371 '''Get the I{local tangent plane} (L{Ltp}) or C{None}. 

372 ''' 

373 return self._ltp 

374 

375 def toStr(self, prec=3, fmt=Fmt.F, sep=_COMMASPACE_): # PYCHOK signature 

376 '''Convert this frustum to a "hfov, vfov, ltp" string. 

377 

378 @kwarg prec: Number of (decimal) digits, unstripped (0..8 or C{None}). 

379 @kwarg fmt: Optional, C{float} format (C{letter}). 

380 @kwarg sep: Separator to join (C{str}). 

381 

382 @return: Frustum in the specified form (C{str}). 

383 ''' 

384 t = self.hfov, self.vfov 

385 if self.ltp: 

386 t += self.ltp, 

387 t = strs(t, prec=prec, fmt=fmt) 

388 return sep.join(t) if sep else t 

389 

390 @Property_RO 

391 def vfov(self): 

392 '''Get the vertical C{fov} (C{degrees}). 

393 ''' 

394 return Degrees(vfov=self._v_2 * _2_0) 

395 

396 

397class LocalError(_ValueError): 

398 '''A L{LocalCartesian} or L{Ltp} related issue. 

399 ''' 

400 pass # PYLINT W0107 

401 

402 

403class LocalCartesian(_NamedBase): 

404 '''Conversion between geodetic C{(lat, lon, height)} and I{local 

405 cartesian} C{(x, y, z)} coordinates with I{geodetic} origin 

406 C{(lat0, lon0, height0)}, transcoded from I{Karney}'s C++ class 

407 U{LocalCartesian<https://GeographicLib.SourceForge.io/C++/doc/ 

408 classGeographicLib_1_1LocalCartesian.html>}. 

409 

410 The C{z} axis is normal to the ellipsoid, the C{y} axis points due 

411 North. The plane C{z = -height0} is tangent to the ellipsoid. 

412 

413 The conversions all take place via geocentric coordinates using a 

414 geocentric L{EcefKarney}, by default the WGS84 datum/ellipsoid. 

415 

416 @see: Class L{Ltp}. 

417 ''' 

418 _Ecef = EcefKarney 

419 _ecef = EcefKarney(_WGS84) 

420 _lon00 = INT0 # self.lon0 

421 _9t0 = None # origin (..., lat0, lon0, height0, ...) L{Ecef9Tuple} 

422 _9Tuple = Local9Tuple 

423 

424 def __init__(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, **lon00_name): 

425 '''New L{LocalCartesian} converter. 

426 

427 @kwarg latlonh0: The (geodetic) origin (C{LatLon}, L{LatLon4Tuple}, L{Ltp} 

428 L{LocalCartesian} or L{Ecef9Tuple}) or the C{scalar} 

429 latitude of the (goedetic) origin (C{degrees}). 

430 @kwarg lon0: Longitude of the (goedetic) origin (C{degrees}), required if 

431 B{C{latlonh0}} is C{scalar}, ignored otherwise. 

432 @kwarg height0: Optional height (C{meter}, conventionally) at the (goedetic) 

433 origin perpendicular to and above (or below) the ellipsoid's 

434 surface, like B{C{lon0}}. 

435 @kwarg ecef: An ECEF converter (L{EcefKarney} I{only}), like B{C{lon0}}. 

436 @kwarg lon00_name: Optional C{B{name}=NN} (C{str}) and keyword argument 

437 C{B{lon00}=B{lon0}} for the arbitrary I{polar} longitude 

438 (C{degrees}), see method C{reverse} and property C{lon00} 

439 for further details. 

440 

441 @raise LocalError: If B{C{latlonh0}} not C{LatLon}, L{LatLon4Tuple}, L{Ltp}, 

442 L{LocalCartesian} or L{Ecef9Tuple} or B{C{latlonh0}}, 

443 B{C{lon0}}, B{C{height0}} or B{C{lon00}} invalid. 

444 

445 @raise TypeError: Invalid B{C{ecef}} or not L{EcefKarney}. 

446 

447 @note: If BC{latlonh0} is an L{Ltp} or L{LocalCartesian}, only C{lat0}, C{lon0}, 

448 C{height0} and I{polar} C{lon00} are copied, I{not} the ECEF converter. 

449 ''' 

450 self.reset(latlonh0, lon0=lon0, height0=height0, ecef=ecef, **lon00_name) 

451 

452 def __eq__(self, other): 

453 '''Compare this and an other instance. 

454 

455 @arg other: The other ellipsoid (L{LocalCartesian} or L{Ltp}). 

456 

457 @return: C{True} if equal, C{False} otherwise. 

458 ''' 

459 return other is self or (isinstance(other, self.__class__) and 

460 other.ecef == self.ecef and 

461 other._9t0 == self._9t0) 

462 

463 @Property_RO 

464 def datum(self): 

465 '''Get the ECEF converter's datum (L{Datum}). 

466 ''' 

467 return self.ecef.datum 

468 

469 @Property_RO 

470 def ecef(self): 

471 '''Get the ECEF converter (L{EcefKarney}). 

472 ''' 

473 return self._ecef 

474 

475 def _ecef2local(self, ecef, Xyz, name_Xyz_kwds): # in _EcefLocal._Ltp_ecef2local 

476 '''(INTERNAL) Convert geocentric/geodetic to local, like I{forward}. 

477 

478 @arg ecef: Geocentric (and geodetic) (L{Ecef9Tuple}). 

479 @arg Xyz: An L{XyzLocal}, L{Aer}, L{Enu} or L{Ned} I{class} or C{None}. 

480 @arg name_Xyz_kwds: Optional C{B{name}=NN} (C{str}) and optionally, 

481 additional B{C{Xyz}} keyword arguments, ignored if C{B{Xyz} 

482 is None}. 

483 

484 @return: An C{B{Xyz}(x, y, z, ltp, **B{name_Xyz_kwds}} instance or 

485 if C{B{Xyz} is None}, a L{Local9Tuple}C{(x, y, z, lat, lon, 

486 height, ltp, ecef, M)} with this C{ltp}, B{C{ecef}} 

487 (L{Ecef9Tuple}) converted to this C{datum} and C{M=None}, 

488 always. 

489 

490 @raise TypeError: Invalid B{C{Xyz}} or B{C{name_Xyz_kwds}} item. 

491 ''' 

492 _xinstanceof(Ecef9Tuple, ecef=ecef) 

493 if ecef.datum != self.datum: 

494 ecef = ecef.toDatum(self.datum) 

495 n, kwds = _name2__(name_Xyz_kwds, _or_nameof=ecef) 

496 x, y, z = self.M.rotate(ecef.xyz, *self._9t0_xyz) 

497 r = Local9Tuple(x, y, z, ecef.lat, ecef.lon, ecef.height, 

498 self, ecef, None, name=n) 

499 if Xyz: 

500 _xsubclassof(*_XyzLocals4, Xyz=Xyz) # Vector3d 

501 r = r.toXyz(Xyz=Xyz, name=n, **kwds) 

502 return r 

503 

504 @Property_RO 

505 def ellipsoid(self): 

506 '''Get the ECEF converter's ellipsoid (L{Ellipsoid}). 

507 ''' 

508 return self.ecef.datum.ellipsoid 

509 

510 def forward(self, latlonh, lon=None, height=0, M=False, **name): 

511 '''Convert I{geodetic} C{(lat, lon, height)} to I{local} cartesian 

512 C{(x, y, z)}. 

513 

514 @arg latlonh: Either a C{LatLon}, L{Ltp}, L{Ecef9Tuple} or C{scalar} 

515 (geodetic) latitude (C{degrees}). 

516 @kwarg lon: Optional C{scalar} (geodetic) longitude (C{degrees}) iff 

517 B{C{latlonh}} is C{scalar}, ignored otherwise. 

518 @kwarg height: Optional height (C{meter}, conventionally) perpendicular 

519 to and above (or below) the ellipsoid's surface, iff 

520 B{C{latlonh}} is C{scalar}, ignored otherwise. 

521 @kwarg M: Optionally, return the I{concatenated} rotation L{EcefMatrix}, 

522 iff available (C{bool}). 

523 @kwarg name: Optional C{B{name}=NN} (C{str}). 

524 

525 @return: A L{Local9Tuple}C{(x, y, z, lat, lon, height, ltp, ecef, M)} 

526 with I{local} C{x}, C{y}, C{z}, I{geodetic} C{(lat}, C{lon}, 

527 C{height}, this C{ltp}, C{ecef} (L{Ecef9Tuple}) with 

528 I{geocentric} C{x}, C{y}, C{z} (and I{geodetic} C{lat}, 

529 C{lon}, C{height}) and the I{concatenated} rotation matrix 

530 C{M} (L{EcefMatrix}) if requested. 

531 

532 @raise LocalError: If B{C{latlonh}} not C{scalar}, C{LatLon}, L{Ltp}, 

533 L{Ecef9Tuple} or invalid or if B{C{lon}} not 

534 C{scalar} for C{scalar} B{C{latlonh}} or invalid 

535 or if B{C{height}} invalid. 

536 ''' 

537 lat, lon, h, n = _llhn4(latlonh, lon, height, Error=LocalError, **name) 

538 t = self.ecef._forward(lat, lon, h, n, M=M) 

539 x, y, z = self.M.rotate(t.xyz, *self._9t0_xyz) 

540 m = self.M.multiply(t.M) if M else None 

541 return self._9Tuple(x, y, z, lat, lon, h, self, t, m, name=n or self.name) 

542 

543 @Property_RO 

544 def height0(self): 

545 '''Get the origin's height (C{meter}). 

546 ''' 

547 return self._9t0.height 

548 

549 @Property_RO 

550 def lat0(self): 

551 '''Get the origin's latitude (C{degrees}). 

552 ''' 

553 return self._9t0.lat 

554 

555 @Property_RO 

556 def latlonheight0(self): 

557 '''Get the origin's lat-, longitude and height (L{LatLon3Tuple}C{(lat, lon, height)}). 

558 ''' 

559 return LatLon3Tuple(self.lat0, self.lon0, self.height0, name=self.name) 

560 

561 def _local2ecef(self, local, nine=False, M=False): 

562 '''(INTERNAL) Convert I{local} to geocentric/geodetic, like I{.reverse}. 

563 

564 @arg local: Local (L{XyzLocal}, L{Enu}, L{Ned}, L{Aer} or L{Local9Tuple}). 

565 @kwarg nine: If C{True}, return a 9-, otherwise a 3-tuple (C{bool}). 

566 @kwarg M: Include the rotation matrix (C{bool}). 

567 

568 @return: A I{geocentric} 3-tuple C{(x, y, z)} or if C{B{nine}=True}, an 

569 L{Ecef9Tuple}C{(x, y, z, lat, lon, height, C, M, datum)} with 

570 rotation matrix C{M} (L{EcefMatrix}) if requested. 

571 ''' 

572 _xinstanceof(*_XyzLocals5, local=local) 

573 t = self.M.unrotate(local.xyz, *self._9t0_xyz) 

574 if nine: 

575 t = self.ecef.reverse(*t, M=M) 

576 return t 

577 

578 @Property_RO 

579 def lon0(self): 

580 '''Get the origin's longitude (C{degrees}). 

581 ''' 

582 return self._9t0.lon 

583 

584 @Property 

585 def lon00(self): 

586 '''Get the arbitrary, I{polar} longitude (C{degrees}). 

587 ''' 

588 return self._lon00 

589 

590 @lon00.setter # PYCHOK setter! 

591 def lon00(self, lon00): 

592 '''Set the arbitrary, I{polar} longitude (C{degrees}). 

593 ''' 

594 # lon00 <https://GitHub.com/mrJean1/PyGeodesy/issues/77> 

595 self._lon00 = Degrees(lon00=lon00) 

596 

597 @Property_RO 

598 def M(self): 

599 '''Get the rotation matrix (C{EcefMatrix}). 

600 ''' 

601 return self._9t0.M 

602 

603 def reset(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, **lon00_name): 

604 '''Reset this converter, see L{LocalCartesian.__init__} for further details. 

605 ''' 

606 _, name = _xkwds_pop2(lon00_name, lon00=None) # PYCHOK get **name 

607 if isinstance(latlonh0, LocalCartesian): 

608 if self._9t0: 

609 _update_all(self) 

610 self._ecef = latlonh0.ecef 

611 self._lon00 = latlonh0.lon00 

612 self._9t0 = latlonh0._9t0 

613 n = _name__(name, _or_nameof=latlonh0) 

614 else: 

615 n = _name__(name, _or_nameof=self) 

616 lat0, lon0, height0, n = _llhn4(latlonh0, lon0, height0, suffix=_0_, 

617 Error=LocalError, name=n) 

618 if ecef: # PYCHOK no cover 

619 _xinstanceof(self._Ecef, ecef=ecef) 

620 _update_all(self) 

621 self._ecef = ecef 

622 elif self._9t0: 

623 _update_all(self) 

624 self._9t0 = self.ecef._forward(lat0, lon0, height0, n, M=True) 

625 self.lon00 = _xattr(latlonh0, lon00=_xkwds_get(lon00_name, lon00=lon0)) 

626 if n: 

627 self.rename(n) 

628 

629 def reverse(self, xyz, y=None, z=None, M=False, **lon00_name): 

630 '''Convert I{local} C{(x, y, z)} to I{geodetic} C{(lat, lon, height)}. 

631 

632 @arg xyz: A I{local} (L{XyzLocal}, L{Enu}, L{Ned}, L{Aer}, L{Local9Tuple}) or 

633 local C{x} coordinate (C{scalar}). 

634 @kwarg y: Local C{y} coordinate (C{meter}), iff B{C{xyz}} is C{scalar}, 

635 ignored otherwise. 

636 @kwarg z: Local C{z} coordinate (C{meter}), iff B{C{xyz}} is C{scalar}, 

637 ignored otherwise. 

638 @kwarg M: Optionally, return the I{concatenated} rotation L{EcefMatrix}, iff 

639 available (C{bool}). 

640 @kwarg lon00_name: Optional C{B{name}=NN} (C{str}) and keyword argument 

641 C{B{lon00}=B{lon0}} for the arbitrary I{polar} longitude 

642 (C{degrees}), overriding see the property C{B{lon00}=B{lon0}} 

643 value. The I{polar} longitude (C{degrees}) is returned with 

644 I{polar} latitudes C{abs(B{lat0}) == 90} for local C{B{x}=0} 

645 and C{B{y}=0} locations. 

646 

647 @return: A L{Local9Tuple}C{(x, y, z, lat, lon, height, ltp, ecef, M)} with 

648 I{local} C{x}, C{y}, C{z}, I{geodetic} C{lat}, C{lon}, C{height}, 

649 this C{ltp}, an C{ecef} (L{Ecef9Tuple}) with the I{geocentric} C{x}, 

650 C{y}, C{z} (and I{geodetic} C{lat}, C{lon}, C{height}) and the 

651 I{concatenated} rotation matrix C{M} (L{EcefMatrix}) if requested. 

652 

653 @raise LocalError: Invalid B{C{xyz}} or C{scalar} C{x} or B{C{y}} and/or B{C{z}} 

654 not C{scalar} for C{scalar} B{C{xyz}}. 

655 ''' 

656 lon00, name =_xkwds_pop2(lon00_name, lon00=self.lon00) 

657 x, y, z, n = _xyzn4(xyz, y, z, _XyzLocals5, Error=LocalError, name=name) 

658 c = self.M.unrotate((x, y, z), *self._9t0_xyz) 

659 t = self.ecef.reverse(*c, M=M, lon00=lon00) 

660 m = self.M.multiply(t.M) if M else None 

661 return self._9Tuple(x, y, z, t.lat, t.lon, t.height, self, t, m, name=n or self.name) 

662 

663 @Property_RO 

664 def _9t0_xyz(self): 

665 '''(INTERNAL) Get C{(x0, y0, z0)} as L{Vector3Tuple}. 

666 ''' 

667 return self._9t0.xyz 

668 

669 def toStr(self, prec=9, **unused): # PYCHOK signature 

670 '''Return this L{LocalCartesian} as a string. 

671 

672 @kwarg prec: Precision, number of (decimal) digits (0..9). 

673 

674 @return: This L{LocalCartesian} representation (C{str}). 

675 ''' 

676 return self.attrs(_lat0_, _lon0_, _height0_, _M_, _ecef_, _name_, prec=prec) 

677 

678 

679class Ltp(LocalCartesian): 

680 '''A I{local tangent plan} (LTP), a sub-class of C{LocalCartesian} with 

681 (re-)configurable ECEF converter. 

682 ''' 

683 _Ecef = _EcefBase 

684 

685 def __init__(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, **lon00_name): # PYLINT W0231 

686 '''New C{Ltp}, see L{LocalCartesian.__init__} for more details. 

687 

688 @kwarg ecef: Optional ECEF converter (L{EcefKarney}, L{EcefFarrell21}, 

689 L{EcefFarrell22}, L{EcefSudano}, L{EcefVeness} or L{EcefYou} 

690 I{instance}), overriding the default L{EcefKarney}C{(datum=Datums.WGS84)} 

691 for C{scalar} B{C{latlonh0}}. 

692 

693 @see: Class L{LocalCartesian<LocalCartesian.__init__>} for further details. 

694 

695 @raise TypeError: Invalid B{C{ecef}}. 

696 ''' 

697 LocalCartesian.reset(self, latlonh0, lon0=lon0, height0=height0, 

698 ecef=ecef, **lon00_name) 

699 

700 @Property 

701 def ecef(self): 

702 '''Get this LTP's ECEF converter (C{Ecef...} I{instance}). 

703 ''' 

704 return self._ecef 

705 

706 @ecef.setter # PYCHOK setter! 

707 def ecef(self, ecef): 

708 '''Set this LTP's ECEF converter (C{Ecef...} I{instance}). 

709 

710 @raise TypeError: Invalid B{C{ecef}}. 

711 ''' 

712 _xinstanceof(_EcefBase, ecef=ecef) 

713 if self._ecef != ecef: # PYCHOK no cover 

714 self.reset(self._9t0) 

715 self._ecef = ecef 

716 

717 

718class LqRD(Ltp): 

719 '''A I{local tangent plan} (LTP) for conversion between I{GRS80 (ETRS89) geodetic} and 

720 I{local Netherlands}' C{quasi-B{R}ijksB{D}riehoeksmeting (RD)} coordinates. 

721 

722 This C{quasi-RD} transformer B{does not} implement any U{RD NAP<https://www.NSGI.NL/ 

723 coordinatenstelsels-en-transformaties/coordinatentransformaties/rdnap-etrs89-rdnaptrans>} 

724 specification and B{does not} provide I{Netherlands}' C{B{N}ormaal B{A}msterdams B{P}eil 

725 (NAP)} quasi-geodetic-height. 

726 

727 The L{LqRD.forward} C{x} and C{y} results differ 3 meter near the center up to 600 meter 

728 at the corners of the L{RD region<LqRD.region4>} with C{RDx} and C{RDy} values from 

729 formal C{RD NAP 2018} implementations like U{pyrdnap<https://PyPI.org/project/pyrdnap>}. 

730 

731 The L{LqRD.forward} C{z} values represent perpendicular distances to this local tangent 

732 plane (LTP). Other heights in L{LqRD} are I{GRS80 (ETRS89) heights} above (or below) 

733 the ellipsoid's surface. B{None} are C{NAP} quasi-geodetic-heights. 

734 ''' 

735 Amersfoort = LatLon4Tuple(parseDMS('52 9 22.178N'), # height=0.0, not .height0_ETRS! 

736 parseDMS(' 5 23 15.5E'), _0_0, _GRS80, name='Amersfoort') 

737# _Ecef = _EcefBase 

738 _ecef = EcefKarney(_GRS80) 

739 _x0 = Meter(x0=155029.8) # 155000.0 see pyrdnap -v1 -forward Amersfoort.latlon 

740 _y0 = Meter(y0=463109.9) # 463000.0 see pyrdnap -v1 -forward Amersfoort.latlon 

741 

742 def __init__(self, latlonh0=Amersfoort, **other_Ltp_kwds): 

743 '''New ECEF-based I{GRS80 (ETRS89)} L{LqRD} converter, centered at I{Amersfoort, NL}. 

744 

745 @kwarg latlonh0: The I{geodetic} origin and height, overriding C{Amersfoort}. 

746 @kwarg other_Ltp_kwds: Optional, other L{Ltp.__init__} keyword arguments. 

747 

748 @see: Class L{Ltp<Ltp.__init__>} for more information. 

749 ''' 

750 Ltp.__init__(self, latlonh0, **_xkwds(other_Ltp_kwds, ecef=None, name=LqRD.Amersfoort.name)) 

751 

752 def _asRD(self, b4): # bounds of C{r} as C{RD4Tuple} 

753 _xinstanceof(Bounds4Tuple, b4=b4) 

754 S, W, N, E = b4 

755 b = self.forward(S, W) 

756 t = self.forward(N, E) 

757 # assert b.x < t.x and b.y < t.y 

758 return RD4Tuple(b.x, b.y, t.x, t.y, name=b4.name) 

759 

760 def bounds4(self, asRD=False): 

761 '''Get the South, West, North and East bounds of the Netherlands' U{EEZ 

762 <http://MarineRegions.org/mrgid/5668>} and U{EPSG:28992<https://EPSG.io/28992>}. 

763 

764 @kwarg asRd: Use C{B{asRD}=True} for the bounds in C{meter}, otherwise in 

765 C{degrees} (C{bool}). 

766 

767 @return: A L{Bounds4Tuple}C{(latS, lonW, latN, lonE)} with lat- and longitudes 

768 in C{degrees} or an L{RD4Tuple}C{(minRDx, minRDy, maxRDx, maxRDy)} with 

769 the C{quasi-RD} bounds in C{meter}. 

770 

771 @see: U{EEZ<https://NL.WikiPedia.org/wiki/Nederlandse_Exclusieve_Economische_Zone>} 

772 ''' 

773 b = self._bounds4 

774 return self._asRD(b) if asRD else b 

775 

776 @property_ROver 

777 def _bounds4(self): # as RD-Bessel L{Bounds4Tuple} 

778 b = Bounds4Tuple('51 19 48.6', '2 32 21.6', '55 45 54', '7 12 37') # EEZ 

779 b = b.toUnits().union(50.75, 3.2, 53.7, 7.22) # EPSG:28992 

780 return b.toUnits(name='RD bounds ') 

781 

782 def forward(self, lat_latlonh, lon=None, height=0, **M_name): # PYCHOK signature 

783 '''Convert I{geodetic} C{(lat, lon, height)} to I{local} C{quasi-RD (x, y, z)}. 

784 

785 @return: A L{Local9Tuple}C{(x, y, z, lat, lon, height, ltp, ecef, M)}. 

786 

787 @see: Method L{LocalCartesian.forward} for more information. 

788 ''' 

789 r = Ltp.forward(self, lat_latlonh, lon, height, **M_name) 

790 return r.dup(x=r.x + self.x0, y=r.y + self.y0, name=r.name) 

791 

792 @property_ROver 

793 def height0_ETRS(self): 

794 '''Get C{Amersfoort}'s I{GRS80 (ETRS89) height} (C{Meter}). 

795 ''' 

796 return Meter(height0_ETRS=43.0) # see pyrdnap h0_ETRS 

797 

798 @deprecated_property_RO 

799 def region(self): # PYCHOK no cover 

800 '''DEPRECATED on 2026.06.12, use method L{pygeodesy.LqRD.region4()}.''' 

801 return self._region4 

802 

803 def region4(self, asRD=False): # in pyrdnap.rd0._RD 

804 '''Get the South, West, North and East bounds of the C{RD} region. 

805 

806 @kwarg asRd: Use C{B{asRD}=True} for the bounds in C{meter}, 

807 otherwise in C{degrees} (C{bool}). 

808 

809 @return: A L{Bounds4Tuple}C{(latS, lonW, latN, lonE)} with lat- and 

810 longitudes in C{degrees} or an L{RD4Tuple}C{(minRDx, minRDy, 

811 maxRDx, maxRDy)} with the C{quasi-RD} bounds in C{meter}. 

812 ''' 

813 r = self._region4 

814 return self._asRD(r) if asRD else r 

815 

816# @property_ROver 

817# def _region4ETRS(self): # as ETRS (ETRS89) L{Bounds4Tuple} 

818# return Bounds4Tuple('49 59 57.39', '2 0 0.12', 

819# '55 59 54.82', '7 59 56.97').toUnits(name='ETRS region ') 

820 

821 @property_ROver 

822 def _region4(self): # as RD-Bessel L{Bounds4Tuple} 

823 return Bounds4Tuple(50.0, _2_0, 56.0, _8_0).toUnits(name='RD region ') 

824 

825 def reverse(self, x_xyz, y=None, z=None, **M_name): # PYCHOK signature 

826 '''Convert I{local} C{quasi-RD (x, y, z)} to I{geodetic} C{(lat, lon, height)}. 

827 

828 @return: A L{Local9Tuple}C{(x, y, z, lat, lon, height, ltp, ecef, M)}. 

829 

830 @see: Method L{LocalCartesian.reverse} for more information. 

831 ''' 

832 x, y, z = x_xyz.xyz if y is z is None else map1(Meter, x_xyz, y, z) 

833 r = Ltp.reverse(self, x - self.x0, y - self.y0, z, **M_name) 

834 return r.dup(x=x, y=y, name=r.name) 

835 

836 @property_doc_(' the C{quasi-RD} false Easting (C{meter}).') 

837 def x0(self): 

838 return self._x0 

839 

840 @x0.setter # PYCHOK setter! 

841 def x0(self, meter): 

842 self._x0 = Meter(x0=meter) 

843 

844 @property_doc_(' the C{quasi-RD} false Northing (C{meter}).') 

845 def y0(self): 

846 return self._y0 

847 

848 @y0.setter # PYCHOK setter! 

849 def y0(self, meter): 

850 self._y0 = Meter(y0=meter) 

851