1

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

3

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

5

6 I{Local cartesian} and I{local tangent plane} classes L{LocalCartesian}, approximations L{ChLVa}

7 and L{ChLVe} and L{Ltp}, L{ChLV}, L{LocalError}, L{Attitude} and L{Frustum}.

8

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

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

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

12'''

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

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

15

16 from pygeodesy.basics import isscalar, issubclassof, map1, _xargs_names

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

18 _2_0, _60_0, _90_0, _100_0, _180_0, _3600_0, \

19 _N_1_0 # PYCHOK used!

20 from pygeodesy.datums import _WGS84, _xinstanceof

21 from pygeodesy.ecef import _EcefBase, EcefKarney, _llhn4, _xyzn4

22 from pygeodesy.errors import _NotImplementedError, _TypesError, _ValueError, _xkwds

23 from pygeodesy.fmath import fabs, fdot, Fhorner

24 from pygeodesy.fsums import _floor, Fsum, fsumf_, fsum1f_

25 from pygeodesy.interns import NN, _0_, _COMMASPACE_, _DOT_, _ecef_, _height_, \

26 _invalid_, _lat0_, _lon0_, _ltp_, _M_, _name_, _too_

27# from pygeodesy.lazily import _ALL_LAZY # from vector3d

28 from pygeodesy.ltpTuples import Attitude4Tuple, ChLVEN2Tuple, ChLV9Tuple, \

29 ChLVYX2Tuple, Footprint5Tuple, Local9Tuple, \

30 ChLVyx2Tuple, _XyzLocals4, _XyzLocals5, Xyz4Tuple

31 from pygeodesy.named import _NamedBase, notOverloaded

32 from pygeodesy.namedTuples import LatLon3Tuple, LatLon4Tuple, Vector3Tuple

33 from pygeodesy.props import Property, Property_RO, property_doc_, property_RO, \

34 _update_all

35 from pygeodesy.streprs import Fmt, strs, unstr

36 from pygeodesy.units import Bearing, Degrees, Meter

37 from pygeodesy.utily import cotd, sincos2d, sincos2d_, tand, tand_, wrap180, wrap360

38 from pygeodesy.vector3d import _ALL_LAZY, Vector3d

39

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

41

42__all__ = _ALL_LAZY.ltp

43__version__ = '23.05.26'

44

45_height0_ = _height_ + _0_

46_narrow_ = 'narrow'

47_wide_ = 'wide'

48_Xyz_ = 'Xyz'

49

50

51 def _fov_2(**fov):

52 # Half a field-of-view angle in C{degrees}.

53 f = Degrees(Error=LocalError, **fov) * _0_5

54 if EPS < f < _90_0:

55 return f

56 t = _invalid_ if f < 0 else _too_(_wide_ if f > EPS else _narrow_)

57 raise LocalError(txt=t, **fov)

58

59

60 class Attitude(_NamedBase):

61 '''The orientation of a plane or camera in space.

62 '''

63 _alt = Meter( alt =_0_0)

64 _roll = Degrees(roll=_0_0)

65 _tilt = Degrees(tilt=_0_0)

66 _yaw = Bearing(yaw =_0_0)

67

68 def __init__(self, alt_attitude=INT0, tilt=INT0, yaw=INT0, roll=INT0, name=NN):

69 '''New L{Attitude}.

70

71 @kwarg alt_attitude: An altitude (C{meter}) above earth or an attitude

72 (L{Attitude} or L{Attitude4Tuple}) with the

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

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

75 (clockwise rotation along and around the x- or East axis).

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

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

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

79 (clockwise rotation along and around the y- or North axis).

80 @kwarg name: Optional name C{str}).

81

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

83 B{C{roll}}.

84

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

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

87 '''

88 if isscalar(alt_attitude):

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

90 else:

91 try:

92 t = alt_attitude.atyr

93 except AttributeError:

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

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

96 if v:

97 setattr(self, n, v)

98 n = name or t.name

99 if n:

100 self.name = n

101

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

103 def alt(self):

104 return self._alt

105

106 @alt.setter # PYCHOK setter!

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

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

109 if self._alt != a:

110 _update_all(self)

111 self._alt = a

112

113 altitude = alt

114

115 @Property_RO

116 def atyr(self):

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

118 '''

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

120

121 @Property_RO

122 def matrix(self):

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

124

125 @see: The matrix M of case 10 in U{Appendix A

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

127 '''

128 def _5to3(x, y, _y, z, _z):

129 return x, fsum1f_(y, _y), fsum1f_(z, _z)

130

131 r0, r1, r2 = self._rows3

132 return _5to3(*r0), _5to3(*r1), r2

133

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

135 def roll(self):

136 return self._roll

137

138 @roll.setter # PYCHOK setter!

139 def roll(self, roll):

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

141 if self._roll != r:

142 _update_all(self)

143 self._roll = r

144

145 bank = roll

146

147 @Property_RO

148 def _rows3(self):

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

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

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

152 sa, ca, sb, cb, sg, cg = sincos2d_(-self.yaw, self.roll, self.tilt)

153 return ((ca * cb, ca * sb * sg, -sa * cg, ca * sb * cg, sa * sg),

154 (sa * cb, sa * sb * sg, ca * cg, sa * sb * cg, -ca * sg),

155 ( -sb, cb * sg, cb * cg))

156

157 def rotate(self, x_xyz, y=None, z=None, Vector=None, **Vector_kwds):

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

159

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

161 (C{Cartesian}, L{Vector3d} or L{Vector3Tuple}).

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

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

164 @kwarg Vector: Class to return transformed point (C{Cartesian},

165 L{Vector3d} or C{Vector3Tuple}) or C{None}.

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

167 ignored if C{B{Vector} is None}.

168

169 @return: A B{C{Vector}} instance or a L{Vector3Tuple}C{(x, y, z)} if

170 C{B{Vector}=None}.

171

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

173 '''

174 try:

175 x, y, z = map( float, x_xyz.xyz)

176 except AttributeError:

177 x, y, z = map1(float, x_xyz, y, z)

178

179 r0, r1, r2 = self._rows3

180 X = fdot(r0, x, y, y, z, z)

181 Y = fdot(r1, x, y, y, z, z)

182 Z = fdot(r2, x, y, z)

183 return Vector3Tuple(X, Y, Z, name=self.name) if Vector is None else \

184 Vector(X, Y, Z, **_xkwds(Vector_kwds, name=self.name))

185

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

187 def tilt(self):

188 return self._tilt

189

190 @tilt.setter # PYCHOK setter!

191 def tilt(self, tilt):

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

193 if self._tilt != t:

194 _update_all(self)

195 self._tilt = t

196

197 elevation = pitch = tilt

198

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

200 '''Format this attitude as string.

201

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

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

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

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

206

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

208 '''

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

210

211 @Property_RO

212 def tyr3d(self):

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

214

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

216 '''

217 def _r2d(r):

218 return fsumf_(_N_1_0, *r)

219

220 return Vector3d(*map(_r2d, self._rows3), name=tyr3d.__name__)

221

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

223 def yaw(self):

224 return self._yaw

225

226 @yaw.setter # PYCHOK setter!

227 def yaw(self, yaw):

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

229 if self._yaw != y:

230 _update_all(self)

231 self._yaw = y

232

233 bearing = heading = yaw

234

235

236 class AttitudeError(_ValueError):

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

238 '''

239 pass

240

241

242 class Frustum(_NamedBase):

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

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

245

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

247 '''

248 _h_2 = _0_0 # half hfov in degrees

249 _ltp = None # local tangent plane

250 _tan_h_2 = _0_0 # tan(_h_2)

251 _v_2 = _0_0 # half vfov in degrees

252

253 def __init__(self, hfov, vfov, ltp=None):

254 '''New L{Frustum}.

255

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

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

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

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, fov_2=h)

266

267 if ltp:

268 self._ltp = _xLtp(ltp)

269

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

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

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

273

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

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

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

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

278 (clockwise rotation along and around the x- or East axis).

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

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

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

282 (clockwise rotation along and around the y- or North axis).

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

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

285 frustum's C{ltp}.

286

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

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

289

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

291

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

293

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

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

296 over the horizon.

297

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

299 '''

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

301 # left and right corners, or swapped

302 if r < EPS: # no roll

303 r = a * tan_h_2

304 l = -r # PYCHOK l is ell

305 else: # roll

306 r, l = tand_(r - h_2, r + h_2, roll_hfov=r) # PYCHOK l is ell

307 r *= -a # negate right positive

308 l *= -a # PYCHOK l is ell

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

310 return (l, y), (r, y)

311

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

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

314 s, c = sincos2d(b)

315 for x, y in xy5:

316 yield Xyz4Tuple(fsum1f_(x * c, y * s),

317 fsum1f_(y * c, -x * s), z, ltp)

318

319 try:

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

321 except AttributeError:

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

323

324 a = Meter(altitude=a)

325 if a < EPS: # too low

326 raise _ValueError(altitude=a)

327 if z: # PYCHOK no cover

328 z = Meter(z=z)

329 a -= z

330 if a < EPS: # z above a

331 raise _ValueError(altitude_z=a)

332 else:

333 z = _0_0

334

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

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

337 if not EPS < e < _180_0:

338 raise _ValueError(tilt=t)

339 if e > _90_0:

340 e = _180_0 - e

341 b = _umod_360(b + _180_0)

342

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

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

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

346 if fabs(y) < EPS:

347 y = _0_0

348

349 # center and corners, clockwise from upperleft, rolled

350 xy5 = ((x, y),) + _xy2(a, e - self._v_2, self._h_2, self._tan_h_2, r) \

351 + _xy2(a, e + self._v_2, -self._h_2, -self._tan_h_2, r) # swapped

352 # turn center and corners by yaw, clockwise

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

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

355

356 @Property_RO

357 def hfov(self):

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

359 '''

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

361

362 @Property_RO

363 def ltp(self):

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

365 '''

366 return self._ltp

367

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

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

370

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

372 @kwarg fmt: Optional, C{float} format (C{str}).

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

374

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

376 '''

377 t = self.hfov, self.vfov

378 if self.ltp:

379 t += self.ltp,

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

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

382

383 @Property_RO

384 def vfov(self):

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

386 '''

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

388

389

390 class LocalError(_ValueError):

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

392 '''

393 pass

394

395

396 class LocalCartesian(_NamedBase):

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

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

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

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

401 classGeographicLib_1_1LocalCartesian.html>}.

402

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

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

405

406 The conversions all take place via geocentric coordinates using a

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

408

409 @see: Class L{Ltp}.

410 '''

411 _ecef = EcefKarney(_WGS84)

412 _t0 = None # origin (..., lat0, lon0, height0, ...) L{Ecef9Tuple}

413 _9Tuple = Local9Tuple

414

415 def __init__(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, name=NN):

416 '''New L{LocalCartesian} converter.

417

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

419 L{Ltp} or L{Ecef9Tuple}) or the C{scalar} latitude

420 of the (goedetic) origin (C{degrees}).

421 @kwarg lon0: Optional longitude of the (goedetic) origin for

422 C{scalar} B{C{latlonh0}} and B{C{height0}} (C{degrees}).

423 @kwarg height0: Optional origin height (C{meter}), vertically

424 above (or below) the surface of the ellipsoid.

425 @kwarg ecef: An ECEF converter (L{EcefKarney} I{only}).

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

427

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

429 L{Ltp} or L{Ecef9Tuple} or B{C{latlonh0}}, B{C{lon0}}

430 or B{C{height0}} invalid, non-C{scalar}.

431

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

433

434 @note: If BC{latlonh0} is an L{Ltp}, only the lat-, longitude and

435 height are duplicated, I{not} the ECEF converter.

436 '''

437 if isinstance(latlonh0, LocalCartesian):

438 self._ecef = latlonh0.ecef

439 self._t0 = latlonh0._t0

440 self.name = name or latlonh0.name

441 else:

442 self.reset(latlonh0, lon0, height0, name=name)

443 if ecef: # PYCHOK no cover

444 _xinstanceof(EcefKarney, ecef=ecef)

445 self._ecef = ecef

446

447 def __eq__(self, other):

448 '''Compare this and an other instance.

449

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

451

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

453 '''

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

455 other.ecef == self.ecef and

456 other._t0 == self._t0)

457

458 @Property_RO

459 def datum(self):

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

461 '''

462 return self.ecef.datum

463

464 @Property_RO

465 def ecef(self):

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

467 '''

468 return self._ecef

469

470 def _ecef2local(self, ecef, Xyz, Xyz_kwds):

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

472

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

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

475 @arg Xyz_kwds: B{C{Xyz}} keyword arguments, ignored if C{B{Xyz} is None}.

476

477 @return: An C{B{Xyz}(x, y, z, ltp, **B{Xyz_kwds}} instance or if

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

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

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

481 always.

482 '''

483 ltp = self

484 if ecef.datum != ltp.datum:

485 ecef = ecef.toDatum(ltp.datum)

486 x, y, z = self.M.rotate(ecef.xyz, *ltp._xyz0)

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

488 ltp, ecef, None, name=ecef.name)

489 if Xyz:

490 if not issubclassof(Xyz, *_XyzLocals4): # Vector3d

491 raise _TypesError(_Xyz_, Xyz, *_XyzLocals4)

492 r = r.toXyz(Xyz=Xyz, **Xyz_kwds)

493 return r

494

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

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

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

498

499 @arg latlonh: Either a C{LatLon}, an L{Ltp}, an L{Ecef9Tuple} or

500 C{scalar} (geodetic) latitude (C{degrees}).

501 @kwarg lon: Optional C{scalar} (geodetic) longitude for C{scalar}

502 B{C{latlonh}} (C{degrees}).

503 @kwarg height: Optional height (C{meter}), vertically above (or below)

504 the surface of the ellipsoid.

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

506 iff available (C{bool}).

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

508

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

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

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

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

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

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

515

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

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

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

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

520 '''

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

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

523 x, y, z = self.M.rotate(t.xyz, *self._xyz0)

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

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

526

527 @Property_RO

528 def height0(self):

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

530 '''

531 return self._t0.height

532

533 @Property_RO

534 def lat0(self):

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

536 '''

537 return self._t0.lat

538

539 @Property_RO

540 def latlonheight0(self):

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

542 '''

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

544

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

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

547

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

549 @kwarg nine: Return 3- or 9-tuple (C{bool}).

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

551

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

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

554 optionally including rotation matrix C{M} or C{None}.

555 '''

556 t = self.M.unrotate(local.xyz, *self._xyz0)

557 if nine:

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

559 return t

560

561 @Property_RO

562 def lon0(self):

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

564 '''

565 return self._t0.lon

566

567 @Property_RO

568 def M(self):

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

570 '''

571 return self._t0.M

572

573 def reset(self, latlonh0=INT0, lon0=INT0, height0=INT0, name=NN):

574 '''Reset the (geodetic) origin.

575

576 @kwarg latlonh0: Either a C{LatLon}, an L{Ecef9Tuple} or C{scalar}

577 latitude of the origin (C{degrees}).

578 @kwarg lon0: Optional C{scalar} longitude of the origin for

579 C{scalar} B{C{latlonh0}} (C{degrees}).

580 @kwarg height0: Optional origin height (C{meter}), vertically

581 above (or below) the surface of the ellipsoid.

582 @kwarg name: Optional, new name (C{str}).

583

584 @raise LocalError: If B{C{latlonh0}} not C{LatLon}, L{Ecef9Tuple},

585 C{scalar} or invalid or if B{C{lon0}} not

586 C{scalar} for C{scalar} B{C{latlonh0}} or

587 invalid or if B{C{height0}} invalid.

588 '''

589 _update_all(self) # force reset

590

591 lat0, lon0, height0, n = _llhn4(latlonh0, lon0, height0,

592 suffix=_0_, Error=LocalError, name=name)

593 if n:

594 self.rename(n)

595 else:

596 n = self.name

597 self._t0 = self.ecef._forward(lat0, lon0, height0, n, M=True)

598

599 def reverse(self, xyz, y=None, z=None, M=False, name=NN):

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

601

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

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

604 @kwarg y: Local C{y} coordinate for C{scalar} B{C{xyz}} and B{C{z}} (C{meter}).

605 @kwarg z: Local C{z} coordinate for C{scalar} B{C{xyz}} and B{C{y}} (C{meter}).

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

607 available (C{bool}).

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

609

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

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

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

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

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

615

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

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

618 '''

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

620 c = self.M.unrotate((x, y, z), *self._xyz0)

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

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

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

624

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

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

627

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

629

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

631 '''

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

633

634 @Property_RO

635 def _xyz0(self):

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

637 '''

638 return self._t0.xyz

639

640

641 class Ltp(LocalCartesian):

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

643 (re-)configurable ECEF converter.

644 '''

645 def __init__(self, latlonh0=0, lon0=0, height0=0, ecef=None, name=NN):

646 '''New C{Ltp}.

647

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

649 L{Ltp} or L{Ecef9Tuple}) or the C{scalar} latitude

650 of the origin (C{degrees}).

651 @kwarg lon0: Optional longitude of the (goedetic) origin for

652 C{scalar} B{C{latlonh0}} and B{C{height0}} (C{degrees}).

653 @kwarg height0: Optional origin height (C{meter}), vertically

654 above (or below) the surface of the ellipsoid.

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

656 L{EcefFarrell22}, L{EcefSudano}, L{EcefVeness} or

657 L{EcefYou} I{instance}), overriding the default

658 L{EcefKarney}C{(datum=Datums.WGS84)}.

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

660

661 @return: New instance (C{Ltp}).

662

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

664 L{Ltp} or L{Ecef9Tuple} or B{C{latlonh0}}, B{C{lon0}}

665 or B{C{height0}} invalid, non-C{scalar}.

666

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

668

669 @note: If BC{latlonh0} is an L{Ltp}, only the lat-, longitude and

670 height are duplicated, I{not} the ECEF converter.

671 '''

672 LocalCartesian.__init__(self, latlonh0, lon0=lon0, height0=height0, name=name)

673 if ecef:

674 self.ecef = ecef

675

676 @Property

677 def ecef(self):

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

679 '''

680 return self._ecef

681

682 @ecef.setter # PYCHOK setter!

683 def ecef(self, ecef):

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

685

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

687 '''

688 _xinstanceof(_EcefBase, ecef=ecef)

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

690 self.reset(self._t0)

691 self._ecef = ecef

692

693

694 class _ChLV(object):

695 '''(INTERNAL) Base class for C{ChLV*} classes.

696 '''

697 _03_falsing = ChLVyx2Tuple(0.6e6, 0.2e6)

698# _92_falsing = ChLVYX2Tuple(2.0e6, 1.0e6) # _95_ - _03_

699 _95_falsing = ChLVEN2Tuple(2.6e6, 1.2e6)

700

701 def _ChLV9Tuple(self, fw, M, name, *Y_X_h_lat_lon_h):

702 '''(INTERNAL) Helper for C{ChLVa/e.forward} and C{.reverse}.

703 '''

704 if bool(M): # PYCHOK no cover

705 m = self.forward if fw else self.reverse # PYCHOK attr

706 n = _DOT_(self.__class__.__name__, m.__name__)

707 raise _NotImplementedError(unstr(n, M=M), txt=None)

708 t = Y_X_h_lat_lon_h + (self, self._t0, None) # PYCHOK _t0

709 return ChLV9Tuple(t, name=name)

710

711 @property_RO

712 def _enh_n_h(self):

713 '''(INTERNAL) Get C{ChLV*.reverse} args[1:4] names, I{once}.

714 '''

715 t = _xargs_names(_ChLV.reverse)[1:4]

716 _ChLV._enh_n_h = t # overwrite this property_RO

717 # assert _xargs_names( ChLV.reverse)[1:4] == t

718 # assert _xargs_names(ChLVa.reverse)[1:4] == t

719 # assert _xargs_names(ChLVe.reverse)[1:4] == t

720 return t

721

722 def forward(self, latlonh, lon=None, height=0, M=None, name=NN):

723 '''Convert WGS84 geodetic to I{Swiss} projection coordinates.

724

725 @arg latlonh: Either a C{LatLon}, L{Ltp} or C{scalar} (geodetic) latitude (C{degrees}).

726 @kwarg lon: Optional, C{scalar} (geodetic) longitude for C{scalar} B{C{latlonh}} (C{degrees}).

727 @kwarg height: Optional, height, vertically above (or below) the surface of the ellipsoid

728 (C{meter}) for C{scalar} B{C{latlonh}} and B{C{lon}}.

729 @kwarg M: If C{True}, return the I{concatenated} rotation L{EcefMatrix} iff available

730 for C{ChLV} only, C{None} otherwise (C{bool}).

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

732

733 @return: A L{ChLV9Tuple}C{(Y, X, h_, lat, lon, height, ltp, ecef, M)} with the unfalsed

734 I{Swiss Y, X} coordinates, I{Swiss h_} height, the given I{geodetic} C{lat},

735 C{lon} and C{height}, this C{ChLV*} instance and C{ecef} (L{Ecef9Tuple}) at

736 I{Bern, Ch} and rotation matrix C{M}. The returned C{ltp} is this C{ChLV},

737 C{ChLVa} or C{ChLVe} instance.

738

739 @raise LocalError: Invalid or non-C{scalar} B{C{latlonh}}, B{C{lon}} or B{C{height}}.

740 '''

741 notOverloaded(self, latlonh, lon=lon, height=height, M=M, name=name)

742

743 def reverse(self, enh_, n=None, h_=0, M=None, name=NN):

744 '''Convert I{Swiss} projection to WGS84 geodetic coordinates.

745

746 @arg enh_: A Swiss projection (L{ChLV9Tuple}) or the C{scalar}, falsed I{Swiss E_LV95}

747 or I{y_LV03} easting (C{meter}).

748 @kwarg n: Falsed I{Swiss N_LV85} or I{x_LV03} northing for C{scalar} B{C{enh_}} and

749 B{C{h_}} (C{meter}).

750 @kwarg h_: I{Swiss h'} height for C{scalar} B{C{enh_}} and B{C{n}} (C{meter}).

751 @kwarg M: If C{True}, return the I{concatenated} rotation L{EcefMatrix} iff available

752 for C{ChLV} only, C{None} otherwise (C{bool}).

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

754

755 @return: A L{ChLV9Tuple}C{(Y, X, h_, lat, lon, height, ltp, ecef, M)} with the unfalsed

756 I{Swiss Y, X} coordinates, I{Swiss h_} height, the given I{geodetic} C{lat},

757 C{lon} and C{height}, this C{ChLV*} instance and C{ecef} (L{Ecef9Tuple}) at

758 I{Bern, Ch} and rotation matrix C{M}. The returned C{ltp} is this C{ChLV},

759 C{ChLVa} or C{ChLVe} instance.

760

761 @raise LocalError: Invalid or non-C{scalar} B{C{enh_}}, B{C{n}} or B{C{h_}}.

762 '''

763 notOverloaded(self, enh_, n=n, h_=h_, M=M, name=name)

764

765 @staticmethod

766 def _falsing2(LV95):

767 '''(INTERNAL) Get the C{LV95} or C{LV03} falsing.

768 '''

769 return _ChLV._95_falsing if LV95 in (True, 95) else (

770 _ChLV._03_falsing if LV95 in (False, 3) else ChLVYX2Tuple(0, 0))

771

772 @staticmethod

773 def _llh2abh_3(lat, lon, h):

774 '''(INTERNAL) Helper for C{ChLVa/e.forward}.

775 '''

776 def _deg2ab(deg, sLL):

777 # convert degrees to arc-seconds

778 def _dms(ds, p, q, swap):

779 d = _floor(ds)

780 t = (ds - d) * p

781 m = _floor(t)

782 s = (t - m) * p

783 if swap:

784 d, s = s, d

785 return d + (m + s * q) * q

786

787 s = _dms(deg, _60_0, _0_01, False) # deg2sexag

788 s = _dms( s, _100_0, _60_0, True) # sexag2asec

789 return (s - sLL) / ChLV._s_ab

790

791 a = _deg2ab(lat, ChLV._sLat) # phi', lat_aux

792 b = _deg2ab(lon, ChLV._sLon) # lam', lng_aux

793 h_ = fsumf_(h, -ChLV.Bern.height, 2.73 * b, 6.94 * a)

794 return a, b, h_

795

796 @staticmethod

797 def _YXh_2abh3(Y, X, h_):

798 '''(INTERNAL) Helper for C{ChLVa/e.reverse}.

799 '''

800 def _YX2ab(YX):

801 return YX * ChLV._ab_m

802

803 a, b = map1(_YX2ab, Y, X)

804 h = fsumf_(h_, ChLV.Bern.height, -12.6 * a, -22.64 * b)

805 return a, b, h

806

807 def _YXh_n4(self, enh_, n, h_, name):

808 '''(INTERNAL) Helper for C{ChLV*.reverse}.

809 '''

810 Y, X, h_, name = _xyzn4(enh_, n, h_, ChLV9Tuple, name=name,

811 _xyz_y_z_names=self._enh_n_h)

812 if isinstance(enh_, ChLV9Tuple):

813 Y, X = enh_.Y, enh_.X

814 else: # isscalar(enh_)

815 Y, X = ChLV.unfalse2(Y, X) # PYCHOK ChLVYX2Tuple

816 return Y, X, h_, name

817

818

819 class ChLV(_ChLV, Ltp):

820 '''Conversion between I{WGS84 geodetic} and I{Swiss} projection coordinates using

821 L{pygeodesy.EcefKarney}'s Earth-Centered, Earth-Fixed (ECEF) methods.

822

823 @see: U{Swiss projection formulas<https://www.SwissTopo.admin.CH/en/maps-data-online/

824 calculation-services.html>}, page 7ff, U{NAVREF<https://www.SwissTopo.admin.CH/en/

825 maps-data-online/calculation-services/navref.html>}, U{REFRAME<https://www.SwissTopo.admin.CH/

826 en/maps-data-online/calculation-services/reframe.html>} and U{SwissTopo Scripts GPS WGS84

827 <-> LV03<https://GitHub.com/ValentinMinder/Swisstopo-WGS84-LV03>}.

828 '''

829 _9Tuple = ChLV9Tuple

830

831 _ab_d = 0.36 # a, b units per degree, ...

832 _ab_m = 1.0e-6 # ... per meter and ...

833 _ab_M = _1_0 # ... per 1000 kilometer

834 _s_d = _3600_0 # arc-seconds per degree ...

835 _s_ab = _s_d / _ab_d # ... and per a, b unit

836 _sLat = 169028.66 # Bern, Ch in ...

837 _sLon = 26782.5 # ... arc-seconds ...

838 # lat, lon, height == 46°57'08.66", 7°26'22.50", 49.55m ("new" 46°57'07.89", 7°26'22.335")

839 Bern = LatLon4Tuple(_sLat / _s_d, _sLon / _s_d, 49.55, _WGS84, name='Bern')

840

841 def __init__(self, latlonh0=Bern, **other_Ltp_kwds):

842 '''New ECEF-based I{WGS84-Swiss} L{ChLV} converter, centered at I{Bern, Ch}.

843

844 @kwarg latlonh0: The I{geodetic} origin and height, overriding C{Bern, Ch}.

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

846

847 @see: L{Ltp.__init__} for more information.

848 '''

849 Ltp.__init__(self, latlonh0, **_xkwds(other_Ltp_kwds, ecef=None, name=ChLV.Bern.name))

850

851 def forward(self, latlonh, lon=None, height=0, M=None, name=NN): # PYCHOK unused M

852 # overloaded for the _ChLV.forward.__doc__

853 return Ltp.forward(self, latlonh, lon=lon, height=height, M=M, name=name)

854

855 def reverse(self, enh_, n=None, h_=0, M=None, name=NN):

856 # overloaded for the _ChLV.reverse.__doc__

857 Y, X, h_, name = self._YXh_n4(enh_, n, h_, name=name)

858 return Ltp.reverse(self, Y, X, h_, M=M, name=name)

859

860 @staticmethod

861 def false2(Y, X, LV95=True, name=NN):

862 '''Add the I{Swiss LV95} or I{LV03} falsing.

863

864 @arg Y: Unfalsed I{Swiss Y} easting (C{meter}).

865 @arg X: Unfalsed I{Swiss X} northing (C{meter}).

866 @kwarg LV95: If C{True} add C{LV95} falsing, if C{False} add

867 C{LV03} falsing, otherwise leave unfalsed.

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

869

870 @return: A L{ChLVEN2Tuple}C{(E_LV95, N_LV95)} or a

871 L{ChLVyx2Tuple}C{(y_LV03, x_LV03)} with falsed B{C{Y}}

872 and B{C{X}}, otherwise a L{ChLVYX2Tuple}C{(Y, X)}

873 with B{C{Y}} and B{C{X}} as-is.

874 '''

875 e, n = t = _ChLV._falsing2(LV95)

876 return t.classof(e + Y, n + X, name=name)

877

878 @staticmethod

879 def isLV03(e, n):

880 '''Is C{(B{e}, B{n})} a valid I{Swiss LV03} projection?

881

882 @arg e: Falsed (or unfalsed) I{Swiss} easting (C{meter}).

883 @arg n: Falsed (or unfalsed) I{Swiss} northing (C{meter}).

884

885 @return: C{True} if C{(B{e}, B{n})} is a valid, falsed I{Swiss

886 LV03}, projection C{False} otherwise.

887 '''

888 # @see: U{Map<https://www.SwissTopo.admin.CH/en/knowledge-facts/

889 # surveying-geodesy/reference-frames/local/lv95.html>}

890 return 400.e3 < e < 900.e3 and 40.e3 < n < 400.e3

891

892 @staticmethod

893 def isLV95(e, n, raiser=True):

894 '''Is C{(B{e}, B{n})} a valid I{Swiss LV95} or I{LV03} projection?

895

896 @arg e: Falsed (or unfalsed) I{Swiss} easting (C{meter}).

897 @arg n: Falsed (or unfalsed) I{Swiss} northing (C{meter}).

898 @kwarg raiser: If C{True}, throw a L{LocalError} if B{C{e}} and

899 B{C{n}} are invalid I{Swiss LV95} nor I{LV03}.

900

901 @return: C{True} or C{False} if C{(B{e}, B{n})} is a valid I{Swiss

902 LV95} respectively I{LV03} projection, C{None} otherwise.

903 '''

904 if ChLV.isLV03(e, n):

905 return False

906 elif ChLV.isLV03(e - 2.e6, n - 1.e6): # _92_falsing = _95_ - _03_

907 return True

908 elif raiser: # PYCHOK no cover

909 raise LocalError(unstr(ChLV.isLV95, e=e, n=n))

910 return None

911

912 @staticmethod

913 def unfalse2(e, n, LV95=None, name=NN):

914 '''Remove the I{Swiss LV95} or I{LV03} falsing.

915

916 @arg e: Falsed I{Swiss E_LV95} or I{y_LV03} easting (C{meter}).

917 @arg n: Falsed I{Swiss N_LV95} or I{x_LV03} northing (C{meter}).

918 @kwarg LV95: If C{True} remove I{LV95} falsing, if C{False} remove

919 I{LV03} falsing, otherwise use method C{isLV95(B{e},

920 B{n})}.

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

922

923 @return: A L{ChLVYX2Tuple}C{(Y, X)} with the unfalsed B{C{e}}

924 respectively B{C{n}}.

925 '''

926 Y, X = _ChLV._falsing2(ChLV.isLV95(e, n) if LV95 is None else LV95)

927 return ChLVYX2Tuple(e - Y, n - X, name=name)

928

929

930 class ChLVa(_ChLV, LocalCartesian):

931 '''Conversion between I{WGS84 geodetic} and I{Swiss} projection coordinates

932 using the U{Approximate<https://www.SwissTopo.admin.CH/en/maps-data-online/

933 calculation-services.html>} formulas, page 13.

934

935 @see: Older U{references<https://GitHub.com/alphasldiallo/Swisstopo-WGS84-LV03>}.

936 '''

937 def __init__(self, name=ChLV.Bern.name):

938 '''New I{Approximate WGS84-Swiss} L{ChLVa} converter, centered at I{Bern, Ch}.

939

940 @kwarg name: Optional name (C{str}), overriding C{Bern.name}.

941 '''

942 LocalCartesian.__init__(self, latlonh0=ChLV.Bern, name=name)

943

944 def forward(self, latlonh, lon=None, height=0, M=None, name=NN):

945 # overloaded for the _ChLV.forward.__doc__

946 lat, lon, h, name = _llhn4(latlonh, lon, height, name=name)

947 a, b, h_ = _ChLV._llh2abh_3(lat, lon, h)

948 a2, b2 = a**2, b**2

949

950 Y = fsumf_( 72.37, 211455.93 * b,

951 -10938.51 * b * a,

952 -0.36 * b * a2,

953 -44.54 * b * b2) # + 600_000

954 X = fsumf_(147.07, 308807.95 * a,

955 3745.25 * b2,

956 76.63 * a2,

957 -194.56 * b2 * a,

958 119.79 * a2 * a) # + 200_000

959 return self._ChLV9Tuple(True, M, name, Y, X, h_, lat, lon, h)

960

961 def reverse(self, enh_, n=None, h_=0, M=None, name=NN):

962 # overloaded for the _ChLV.reverse.__doc__

963 Y, X, h_, name = self._YXh_n4(enh_, n, h_, name=name)

964 a, b, h = _ChLV._YXh_2abh3(Y, X, h_)

965 a2, b2 = a**2, b**2

966

967 lon = Fsum( 2.6779094, 4.728982 * a,

968 0.791484 * a * b,

969 0.1306 * a * b2,

970 -0.0436 * a * a2).fover(ChLV._ab_d)

971 lat = Fsum(16.9023892, 3.238272 * b,

972 -0.270978 * a2,

973 -0.002528 * b2,

974 -0.0447 * a2 * b,

975 -0.014 * b2 * b).fover(ChLV._ab_d)

976 return self._ChLV9Tuple(False, M, name, Y, X, h_, lat, lon, h)

977

978

979 class ChLVe(_ChLV, LocalCartesian):

980 '''Conversion between I{WGS84 geodetic} and I{Swiss} projection coordinates

981 using the U{Ellipsoidal approximate<https://www.SwissTopo.admin.CH/en/

982 maps-data-online/calculation-services.html>} formulas, pp 10-11 and U{Bolliger,

983 J.<https://eMuseum.GGGS.CH/literatur-lv/liste-Dateien/1967_Bolliger_a.pdf>}

984 pp 148-151 (also U{GGGS<https://eMuseum.GGGS.CH/literatur-lv/liste.htm>}).

985

986 @note: Methods L{ChLVe.forward} and L{ChLVe.reverse} have an additional keyword

987 argument C{B{gamma}=False} to approximate the I{meridian convergence}.

988 If C{B{gamma}=True} a 2-tuple C{(t, gamma)} is returned with C{t} the

989 usual result (C{ChLV9Tuple}) and C{gamma}, the I{meridian convergence}

990 (decimal C{degrees}). To convert C{gamma} to C{grades} or C{gons},

991 use function L{pygeodesy.degrees2grades}.

992

993 @see: Older U{references<https://GitHub.com/alphasldiallo/Swisstopo-WGS84-LV03>}.

994 '''

995 def __init__(self, name=ChLV.Bern.name):

996 '''New I{Approximate WGS84-Swiss} L{ChLVe} converter, centered at I{Bern, Ch}.

997

998 @kwarg name: Optional name (C{str}), overriding C{Bern.name}.

999 '''

1000 LocalCartesian.__init__(self, latlonh0=ChLV.Bern, name=name)

1001

1002 def forward(self, latlonh, lon=None, height=0, M=None, name=NN, gamma=False): # PYCHOK gamma

1003 # overloaded for the _ChLV.forward.__doc__

1004 lat, lon, h, name = _llhn4(latlonh, lon, height, name=name)

1005 a, b, h_ = _ChLV._llh2abh_3(lat, lon, h)

1006 F = Fhorner

1007

1008 B1 = F(a, 211428.533991, -10939.608605, -2.658213, -8.539078, -0.00345, -0.007992)

1009 B3 = F(a, -44.232717, 4.291740, -0.309883, 0.013924)

1010 B5 = F(a, 0.019784, -0.004277)

1011 Y = F(b, 0, B1, 0, B3, 0, B5).fover(ChLV._ab_M) # 1000 Km!

1012

1013 B0 = F(a, 0, 308770.746371, 75.028131, 120.435227, 0.009488, 0.070332, -0.00001)

1014 B2 = F(a, 3745.408911, -193.792705, 4.340858, -0.376174, 0.004053)

1015 B4 = F(a, -0.734684, 0.144466, -0.011842)

1016 B6 = 0.000488

1017 X = F(b, B0, 0, B2, 0, B4, 0, B6).fover(ChLV._ab_M) # 1000 Km!

1018

1019 t = self._ChLV9Tuple(True, M, name, Y, X, h_, lat, lon, h)

1020 if gamma:

1021 U1 = F(a, 2255515.207166, 2642.456961, 1.284180, 2.577486, 0.001165)

1022 U3 = F(a, -412.991934, 64.106344, -2.679566, 0.123833)

1023 U5 = F(a, 0.204129, -0.037725)

1024 g = F(b, 0, U1, 0, U3, 0, U5).fover(ChLV._ab_m) # * ChLV._ab_d degrees?

1025 t = t, g

1026 return t

1027

1028 def reverse(self, enh_, n=None, h_=0, M=None, name=NN, gamma=False): # PYCHOK gamma

1029 # overloaded for the _ChLV.reverse.__doc__

1030 Y, X, h_, name = self._YXh_n4(enh_, n, h_, name=name)

1031 a, b, h = _ChLV._YXh_2abh3(Y, X, h_)

1032 F = Fhorner

1033

1034 A1 = F(b, 47297.3056722, 7925.714783, 1328.129667, 255.02202, 48.17474, 9.0243)

1035 A3 = F(b, -442.709889, -255.02202, -96.34947, -30.0808)

1036 A5 = F(b, 9.63495, 9.0243)

1037 lon = F(a, ChLV._sLon, A1, 0, A3, 0, A5).fover(ChLV._s_d)

1038 # == (ChLV._sLon + a * (A1 + a**2 * (A3 + a**2 * A5))) / ChLV._s_d

1039

1040 A0 = F(b, ChLV._sLat, 32386.4877666, -25.486822, -132.457771, 0.48747, 0.81305, -0.0069)

1041 A2 = F(b, -2713.537919, -450.442705, -75.53194, -14.63049, -2.7604)

1042 A4 = F(b, 24.42786, 13.20703, 4.7476)

1043 A6 = -0.4249

1044 lat = F(a, A0, 0, A2, 0, A4, 0, A6).fover(ChLV._s_d)

1045

1046 t = self._ChLV9Tuple(False, M, name, Y, X, h_, lat, lon, h)

1047 if gamma:

1048 U1 = F(b, 106679.792202, 17876.57022, 4306.5241, 794.87772, 148.1545, 27.8725)

1049 U3 = F(b, -1435.508, -794.8777, -296.309, -92.908)

1050 U5 = F(b, 29.631, 27.873)

1051 g = F(a, 0, U1, 0, U3, 0, U5).fover(ChLV._s_ab) # degrees

1052 t = t, g

1053 return t

1054

1055

1056 def tyr3d(tilt=INT0, yaw=INT0, roll=INT0, Vector=Vector3d, **Vector_kwds):

1057 '''Convert an attitude oriention into a (3-D) direction vector.

1058

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

1060 (clockwise rotation along and around the x-axis).

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

1062 (counter-clockwise rotation along and around the z-axis).

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

1064 (clockwise rotation along and around the y-axis).

1065

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

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

1068

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

1070 and function L{pygeodesy.hartzell} argument C{los}.

1071 '''

1072 d = Attitude4Tuple(_0_0, tilt, yaw, roll).tyr3d

1073 return d if Vector is type(d) else (

1074 Vector3Tuple(d.x, d.y, d.z, name=d.name) if Vector is None else

1075 Vector(d.x, d.y, d.z, **_xkwds(Vector_kwds, name=d.name))) # PYCHOK indent

1076

1077

1078 def _xLtp(ltp, *dflt):

1079 '''(INTERNAL) Validate B{C{ltp}}.

1080 '''

1081 if dflt and ltp is None:

1082 ltp = dflt[0]

1083 if isinstance(ltp, (LocalCartesian, Ltp)):

1084 return ltp

1085 raise _TypesError(_ltp_, ltp, Ltp, LocalCartesian)

1086

1087# **) MIT License

1088#

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

1090#

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

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

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

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

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

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

1097#

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

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

1100#

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

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

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

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

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

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

1107# OTHER DEALINGS IN THE SOFTWARE.