Coverage for pygeodesy/utily.py: 91%

335 statements  

« prev     ^ index     » next       coverage.py v7.2.2, created at 2023-12-29 12:35 -0500

1 

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

3 

4u'''Various utility functions. 

5 

6After I{(C) Chris Veness 2011-2015} published under the same MIT Licence**, see 

7U{Latitude/Longitude<https://www.Movable-Type.co.UK/scripts/latlong.html>} and 

8U{Vector-based geodesy<https://www.Movable-Type.co.UK/scripts/latlong-vectors.html>}. 

9''' 

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

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

12 

13from pygeodesy.basics import _copysign, isint, isstr, neg, _passargs 

14from pygeodesy.constants import EPS, EPS0, INF, NAN, NEG0, PI, PI2, PI_2, R_M, \ 

15 _float as _F, _isfinite, isnan, isnear0, _over, \ 

16 _umod_360, _umod_PI2, _M_KM, _M_NM, _M_SM, _0_0, \ 

17 _1__90, _0_5, _1_0, _N_1_0, _2__PI, _10_0, _90_0, \ 

18 _N_90_0, _180_0, _N_180_0, _360_0, _400_0 

19from pygeodesy.errors import _ValueError, _xkwds, _xkwds_get, _ALL_LAZY, _MODS 

20from pygeodesy.interns import _edge_, _radians_, _semi_circular_, _SPACE_ 

21# from pygeodesy.lazily import _ALL_LAZY, _ALL_MODS as _MODS # from .errors 

22from pygeodesy.units import Degrees, Feet, Float, Lam, Lam_, Meter, Meter2, Radians 

23 

24from math import acos, asin, atan2, cos, degrees, fabs, radians, sin, tan # pow 

25 

26__all__ = _ALL_LAZY.utily 

27__version__ = '23.12.10' 

28 

29# read constant name "_M_Unit" as "meter per Unit" 

30_M_CHAIN = _F( 20.1168) # yard2m(1) * 22 

31_M_FATHOM = _F( 1.8288) # yard2m(1) * 2 or _M_NM * 1e-3 

32_M_FOOT = _F( 0.3048) # Int'l (1 / 3.2808398950131 = 10_000 / (254 * 12)) 

33_M_FOOT_GE = _F( 0.31608) # German Fuss (1 / 3.1637560111364) 

34_M_FOOT_FR = _F( 0.3248406) # French Pied-du-Roi or pied (1 / 3.0784329298739) 

35_M_FOOT_USVY = _F( 0.3048006096012192) # US Survey (1200 / 3937) 

36_M_FURLONG = _F( 201.168) # 220 * yard2m(1) == 10 * m2chain(1) 

37# _M_KM = _F(1000.0) # kilo meter 

38# _M_NM = _F(1852.0) # nautical mile 

39# _M_SM = _F(1609.344) # statute mile 

40_M_TOISE = _F( 1.9490436) # French toise, 6 pieds (6 / 3.0784329298739) 

41_M_YARD_UK = _F( 0.9144) # 254 * 12 * 3 / 10_000 == 3 * ft2m(1) Int'l 

42 

43 

44def _abs1nan(x): 

45 '''(INTERNAL) Bracket C{x}. 

46 ''' 

47 return _N_1_0 < x < _1_0 or isnan(x) 

48 

49 

50def acos1(x): 

51 '''Return C{math.acos(max(-1, min(1, B{x})))}. 

52 ''' 

53 return acos(x) if _abs1nan(x) else (PI if x < 0 else _0_0) 

54 

55 

56def acre2ha(acres): 

57 '''Convert acres to hectare. 

58 

59 @arg acres: Value in acres (C{scalar}). 

60 

61 @return: Value in C{hectare} (C{float}). 

62 

63 @raise ValueError: Invalid B{C{acres}}. 

64 ''' 

65 # 0.40468564224 == acre2m2(1) / 10_000 

66 return Float(ha=Float(acres) * 0.40468564224) 

67 

68 

69def acre2m2(acres): 

70 '''Convert acres to I{square} meter. 

71 

72 @arg acres: Value in acres (C{scalar}). 

73 

74 @return: Value in C{meter^2} (C{float}). 

75 

76 @raise ValueError: Invalid B{C{acres}}. 

77 ''' 

78 # 4046.8564224 == chain2m(1) * furlong2m(1) 

79 return Meter2(Float(acres) * 4046.8564224) 

80 

81 

82def asin1(x): 

83 '''Return C{math.asin(max(-1, min(1, B{x})))}. 

84 ''' 

85 return asin(x) if _abs1nan(x) else _copysign(PI_2, x) 

86 

87 

88def atan1(y, x=_1_0): 

89 '''Return C{atan(B{y} / B{x})} angle in C{radians} M{[-PI/2..+PI/2]} 

90 using C{atan2} for consistency and to avoid C{ZeroDivisionError}. 

91 ''' 

92 return atan2(-y, -x) if x < 0 else atan2(y, x or _0_0) # -0. to 0. 

93 

94 

95def atan1d(y, x=_1_0): 

96 '''Return C{atan(B{y} / B{x})} angle in C{degrees} M{[-90..+90]} 

97 using C{atan2d} for consistency and to avoid C{ZeroDivisionError}. 

98 

99 @see: Function L{pygeodesy.atan2d}. 

100 ''' 

101 return atan2d(-y, -x) if x < 0 else atan2d(y, x or _0_0) # -0. to 0. 

102 

103 

104def atan2b(y, x): 

105 '''Return C{atan2(B{y}, B{x})} in degrees M{[0..+360]}. 

106 

107 @see: Function L{pygeodesy.atan2d}. 

108 ''' 

109 b = atan2d(y, x) 

110 if b < 0: 

111 b += _360_0 

112 return b 

113 

114 

115def atan2d(y, x, reverse=False): 

116 '''Return C{atan2(B{y}, B{x})} in degrees M{[-180..+180]}, 

117 optionally I{reversed} (by 180 degrees for C{azimuth}s). 

118 

119 @see: I{Karney}'s C++ function U{Math.atan2d 

120 <https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Math.html>}. 

121 ''' 

122 if fabs(y) > fabs(x) > 0: 

123 if y < 0: # q = 3 

124 d = degrees(atan2(x, -y)) - _90_0 

125 else: # q = 2 

126 d = _90_0 - degrees(atan2(x, y)) 

127 elif isnan(x) or isnan(y): 

128 return NAN 

129 elif y: 

130 if x > 0: # q = 0 

131 d = degrees(atan2(y, x)) 

132 elif x < 0: # q = 1 

133 d = _copysign(_180_0, y) - degrees(atan2(y, -x)) 

134 else: # x == 0 

135 d = _copysign(_90_0, y) 

136 else: 

137 d = _180_0 if x < 0 else _0_0 

138 return _azireversed(d) if reverse else d 

139 

140 

141def _azireversed(azimuth): # in .rhumbBase 

142 '''(INTERNAL) Return the I{reverse} B{C{azimuth}} in degrees M{[-180..+180]}. 

143 ''' 

144 return azimuth + (_N_180_0 if azimuth > 0 else _180_0) 

145 

146 

147def chain2m(chains): 

148 '''Convert I{UK} chains to meter. 

149 

150 @arg chains: Value in chains (C{scalar}). 

151 

152 @return: Value in C{meter} (C{float}). 

153 

154 @raise ValueError: Invalid B{C{chains}}. 

155 ''' 

156 return Meter(Float(chains=chains) * _M_CHAIN) 

157 

158 

159def circle4(earth, lat): 

160 '''Get the equatorial or a parallel I{circle of latitude}. 

161 

162 @arg earth: The earth radius, ellipsoid or datum 

163 (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, 

164 L{Datum} or L{a_f2Tuple}). 

165 @arg lat: Geodetic latitude (C{degrees90}, C{str}). 

166 

167 @return: A L{Circle4Tuple}C{(radius, height, lat, beta)} 

168 instance. 

169 

170 @raise RangeError: Latitude B{C{lat}} outside valid range and 

171 L{pygeodesy.rangerrors} set to C{True}. 

172 

173 @raise TypeError: Invalid B{C{earth}}. 

174 

175 @raise ValueError: B{C{earth}} or B{C{lat}}. 

176 ''' 

177 E = _MODS.datums._earth_ellipsoid(earth) 

178 return E.circle4(lat) 

179 

180 

181def cot(rad, **error_kwds): 

182 '''Return the C{cotangent} of an angle in C{radians}. 

183 

184 @arg rad: Angle (C{radians}). 

185 @kwarg error_kwds: Error to raise (C{ValueError}). 

186 

187 @return: C{cot(B{rad})}. 

188 

189 @raise ValueError: L{pygeodesy.isnear0}C{(sin(B{rad})}. 

190 ''' 

191 s, c = sincos2(rad) 

192 if c: 

193 if isnear0(s): 

194 raise _valueError(cot, rad, **error_kwds) 

195 c = c / s # /= chokes PyChecker 

196 return c 

197 

198 

199def cot_(*rads, **error_kwds): 

200 '''Return the C{cotangent} of angle(s) in C{radiansresection}. 

201 

202 @arg rads: One or more angles (C{radians}). 

203 @kwarg error_kwds: Error to raise (C{ValueError}). 

204 

205 @return: Yield the C{cot(B{rad})} for each angle. 

206 

207 @raise ValueError: See L{pygeodesy.cot}. 

208 ''' 

209 try: 

210 for r in rads: 

211 yield cot(r) 

212 except ValueError: 

213 raise _valueError(cot_, r, **error_kwds) 

214 

215 

216def cotd(deg, **error_kwds): 

217 '''Return the C{cotangent} of an angle in C{degrees}. 

218 

219 @arg deg: Angle (C{degrees}). 

220 @kwarg error_kwds: Error to raise (C{ValueError}). 

221 

222 @return: C{cot(B{deg})}. 

223 

224 @raise ValueError: L{pygeodesy.isnear0}C{(sin(B{deg})}. 

225 ''' 

226 s, c = sincos2d(deg) 

227 if c: 

228 if isnear0(s): 

229 raise _valueError(cotd, deg, **error_kwds) 

230 c = c / s # /= chokes PyChecker 

231 elif s < 0: 

232 c = -c # negate-0 

233 return c 

234 

235 

236def cotd_(*degs, **error_kwds): 

237 '''Return the C{cotangent} of angle(s) in C{degrees}. 

238 

239 @arg degs: One or more angles (C{degrees}). 

240 @kwarg error_kwds: Error to raise (C{ValueError}). 

241 

242 @return: Yield the C{cot(B{deg})} for each angle. 

243 

244 @raise ValueError: See L{pygeodesy.cotd}. 

245 ''' 

246 try: 

247 for d in degs: 

248 yield cotd(d) 

249 except ValueError: 

250 raise _valueError(cotd_, d, **error_kwds) 

251 

252 

253def degrees90(rad): 

254 '''Convert radians to degrees and wrap M{[-90..+90)}. 

255 

256 @arg rad: Angle (C{radians}). 

257 

258 @return: Angle, wrapped (C{degrees90}). 

259 ''' 

260 return wrap90(degrees(rad)) 

261 

262 

263def degrees180(rad): 

264 '''Convert radians to degrees and wrap M{[-180..+180)}. 

265 

266 @arg rad: Angle (C{radians}). 

267 

268 @return: Angle, wrapped (C{degrees180}). 

269 ''' 

270 return wrap180(degrees(rad)) 

271 

272 

273def degrees360(rad): 

274 '''Convert radians to degrees and wrap M{[0..+360)}. 

275 

276 @arg rad: Angle (C{radians}). 

277 

278 @return: Angle, wrapped (C{degrees360}). 

279 ''' 

280 return _umod_360(degrees(rad)) 

281 

282 

283def degrees2grades(deg): 

284 '''Convert degrees to I{grades} (aka I{gons} or I{gradians}). 

285 

286 @arg deg: Angle (C{degrees}). 

287 

288 @return: Angle (C{grades}). 

289 ''' 

290 return Degrees(deg) * _400_0 / _360_0 

291 

292 

293def degrees2m(deg, radius=R_M, lat=0): 

294 '''Convert an angle to a distance along the equator or 

295 along the parallel at an other (geodetic) latitude. 

296 

297 @arg deg: The angle (C{degrees}). 

298 @kwarg radius: Mean earth radius, ellipsoid or datum 

299 (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, 

300 L{Datum} or L{a_f2Tuple}). 

301 @kwarg lat: Parallel latitude (C{degrees90}, C{str}). 

302 

303 @return: Distance (C{meter}, same units as B{C{radius}} 

304 or equatorial and polar radii) or C{0.0} for 

305 near-polar B{C{lat}}. 

306 

307 @raise RangeError: Latitude B{C{lat}} outside valid range and 

308 L{pygeodesy.rangerrors} set to C{True}. 

309 

310 @raise TypeError: Invalid B{C{radius}}. 

311 

312 @raise ValueError: Invalid B{C{deg}}, B{C{radius}} or 

313 B{C{lat}}. 

314 

315 @see: Function L{radians2m} and L{m2degrees}. 

316 ''' 

317 return _Radians2m(Lam_(deg=deg, clip=0), radius, lat) 

318 

319 

320def fathom2m(fathoms): 

321 '''Convert I{Imperial} fathom to meter. 

322 

323 @arg fathoms: Value in fathoms (C{scalar}). 

324 

325 @return: Value in C{meter} (C{float}). 

326 

327 @raise ValueError: Invalid B{C{fathoms}}. 

328 

329 @see: Function L{toise2m}, U{Fathom<https://WikiPedia.org/wiki/Fathom>} 

330 and U{Klafter<https://WikiPedia.org/wiki/Klafter>}. 

331 ''' 

332 return Meter(Float(fathoms=fathoms) * _M_FATHOM) 

333 

334 

335def ft2m(feet, usurvey=False, pied=False, fuss=False): 

336 '''Convert I{International}, I{US Survey}, I{French} or I{German} 

337 B{C{feet}} to C{meter}. 

338 

339 @arg feet: Value in feet (C{scalar}). 

340 @kwarg usurvey: If C{True}, convert I{US Survey} foot else ... 

341 @kwarg pied: If C{True}, convert French I{pied-du-Roi} else ... 

342 @kwarg fuss: If C{True}, convert German I{Fuss}, otherwise 

343 I{International} foot to C{meter}. 

344 

345 @return: Value in C{meter} (C{float}). 

346 

347 @raise ValueError: Invalid B{C{feet}}. 

348 ''' 

349 return Meter(Feet(feet) * (_M_FOOT_USVY if usurvey else 

350 (_M_FOOT_FR if pied else 

351 (_M_FOOT_GE if fuss else _M_FOOT)))) 

352 

353 

354def furlong2m(furlongs): 

355 '''Convert a furlong to meter. 

356 

357 @arg furlongs: Value in furlongs (C{scalar}). 

358 

359 @return: Value in C{meter} (C{float}). 

360 

361 @raise ValueError: Invalid B{C{furlongs}}. 

362 ''' 

363 return Meter(Float(furlongs=furlongs) * _M_FURLONG) 

364 

365 

366def grades(rad): 

367 '''Convert radians to I{grades} (aka I{gons} or I{gradians}). 

368 

369 @arg rad: Angle (C{radians}). 

370 

371 @return: Angle (C{grades}). 

372 ''' 

373 return Float(grades=Float(rad=rad) * _400_0 / PI2) 

374 

375 

376def grades400(rad): 

377 '''Convert radians to I{grades} (aka I{gons} or I{gradians}) and wrap M{[0..+400)}. 

378 

379 @arg rad: Angle (C{radians}). 

380 

381 @return: Angle, wrapped (C{grades}). 

382 ''' 

383 return Float(grades400=(grades(rad) % _400_0) or _0_0) # _umod_400 

384 

385 

386def grades2degrees(gon): 

387 '''Convert I{grades} (aka I{gons} or I{gradians}) to C{degrees}. 

388 

389 @arg gon: Angle (C{grades}). 

390 

391 @return: Angle (C{degrees}). 

392 ''' 

393 return Degrees(Float(gon=gon) * _360_0 / _400_0) 

394 

395 

396def grades2radians(gon): 

397 '''Convert I{grades} (aka I{gons} or I{gradians}) to C{radians}. 

398 

399 @arg gon: Angle (C{grades}). 

400 

401 @return: Angle (C{radians}). 

402 ''' 

403 return Radians(Float(gon=gon) * PI2 / _400_0) 

404 

405 

406def km2m(km): 

407 '''Convert kilo meter to meter (m). 

408 

409 @arg km: Value in kilo meter (C{scalar}). 

410 

411 @return: Value in meter (C{float}). 

412 

413 @raise ValueError: Invalid B{C{km}}. 

414 ''' 

415 return Meter(Float(km=km) * _M_KM) 

416 

417 

418def _loneg(lon): 

419 '''(INTERNAL) "Complement" of C{lon}. 

420 ''' 

421 return _180_0 - lon 

422 

423 

424def m2chain(meter): 

425 '''Convert meter to I{UK} chains. 

426 

427 @arg meter: Value in meter (C{scalar}). 

428 

429 @return: Value in C{chains} (C{float}). 

430 

431 @raise ValueError: Invalid B{C{meter}}. 

432 ''' 

433 return Float(chain=Meter(meter) / _M_CHAIN) # * 0.049709695378986715 

434 

435 

436def m2degrees(distance, radius=R_M, lat=0): 

437 '''Convert a distance to an angle along the equator or 

438 along the parallel at an other (geodetic) latitude. 

439 

440 @arg distance: Distance (C{meter}, same units as B{C{radius}}). 

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

442 an L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or 

443 L{a_f2Tuple}). 

444 @kwarg lat: Parallel latitude (C{degrees90}, C{str}). 

445 

446 @return: Angle (C{degrees}) or C{INF} for near-polar B{C{lat}}. 

447 

448 @raise RangeError: Latitude B{C{lat}} outside valid range and 

449 L{pygeodesy.rangerrors} set to C{True}. 

450 

451 @raise TypeError: Invalid B{C{radius}}. 

452 

453 @raise ValueError: Invalid B{C{distance}}, B{C{radius}} 

454 or B{C{lat}}. 

455 

456 @see: Function L{m2radians} and L{degrees2m}. 

457 ''' 

458 return degrees(m2radians(distance, radius=radius, lat=lat)) 

459 

460 

461def m2fathom(meter): 

462 '''Convert meter to I{Imperial} fathoms. 

463 

464 @arg meter: Value in meter (C{scalar}). 

465 

466 @return: Value in C{fathoms} (C{float}). 

467 

468 @raise ValueError: Invalid B{C{meter}}. 

469 

470 @see: Function L{m2toise}, U{Fathom<https://WikiPedia.org/wiki/Fathom>} 

471 and U{Klafter<https://WikiPedia.org/wiki/Klafter>}. 

472 ''' 

473 return Float(fathom=Meter(meter) / _M_FATHOM) # * 0.546806649 

474 

475 

476def m2ft(meter, usurvey=False, pied=False, fuss=False): 

477 '''Convert meter to I{International}, I{US Survey}, I{French} or 

478 or I{German} feet (C{ft}). 

479 

480 @arg meter: Value in meter (C{scalar}). 

481 @kwarg usurvey: If C{True}, convert to I{US Survey} foot else ... 

482 @kwarg pied: If C{True}, convert to French I{pied-du-Roi} else ... 

483 @kwarg fuss: If C{True}, convert to German I{Fuss}, otherwise to 

484 I{International} foot. 

485 

486 @return: Value in C{feet} (C{float}). 

487 

488 @raise ValueError: Invalid B{C{meter}}. 

489 ''' 

490 # * 3.2808333333333333, US Survey 3937 / 1200 

491 # * 3.2808398950131235, Int'l 10_000 / (254 * 12) 

492 return Float(feet=Meter(meter) / (_M_FOOT_USVY if usurvey else 

493 (_M_FOOT_FR if pied else 

494 (_M_FOOT_GE if fuss else _M_FOOT)))) 

495 

496 

497def m2furlong(meter): 

498 '''Convert meter to furlongs. 

499 

500 @arg meter: Value in meter (C{scalar}). 

501 

502 @return: Value in C{furlongs} (C{float}). 

503 

504 @raise ValueError: Invalid B{C{meter}}. 

505 ''' 

506 return Float(furlong=Meter(meter) / _M_FURLONG) # * 0.00497096954 

507 

508 

509def m2km(meter): 

510 '''Convert meter to kilo meter (Km). 

511 

512 @arg meter: Value in meter (C{scalar}). 

513 

514 @return: Value in Km (C{float}). 

515 

516 @raise ValueError: Invalid B{C{meter}}. 

517 ''' 

518 return Float(km=Meter(meter) / _M_KM) 

519 

520 

521def m2NM(meter): 

522 '''Convert meter to nautical miles (NM). 

523 

524 @arg meter: Value in meter (C{scalar}). 

525 

526 @return: Value in C{NM} (C{float}). 

527 

528 @raise ValueError: Invalid B{C{meter}}. 

529 ''' 

530 return Float(NM=Meter(meter) / _M_NM) # * 5.39956804e-4 

531 

532 

533def m2radians(distance, radius=R_M, lat=0): 

534 '''Convert a distance to an angle along the equator or along the 

535 parallel at an other (geodetic) latitude. 

536 

537 @arg distance: Distance (C{meter}, same units as B{C{radius}}). 

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

539 an L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or 

540 L{a_f2Tuple}). 

541 @kwarg lat: Parallel latitude (C{degrees90}, C{str}). 

542 

543 @return: Angle (C{radians}) or C{INF} for near-polar B{C{lat}}. 

544 

545 @raise RangeError: Latitude B{C{lat}} outside valid range and 

546 L{pygeodesy.rangerrors} set to C{True}. 

547 

548 @raise TypeError: Invalid B{C{radius}}. 

549 

550 @raise ValueError: Invalid B{C{distance}}, B{C{radius}} 

551 or B{C{lat}}. 

552 

553 @see: Function L{m2degrees} and L{radians2m}. 

554 ''' 

555 m = circle4(radius, lat).radius 

556 return INF if m < EPS0 else Radians(Float(distance=distance) / m) 

557 

558 

559def m2SM(meter): 

560 '''Convert meter to statute miles (SM). 

561 

562 @arg meter: Value in meter (C{scalar}). 

563 

564 @return: Value in C{SM} (C{float}). 

565 

566 @raise ValueError: Invalid B{C{meter}}. 

567 ''' 

568 return Float(SM=Meter(meter) / _M_SM) # * 6.21369949e-4 == 1 / 1609.344 

569 

570 

571def m2toise(meter): 

572 '''Convert meter to French U{toises<https://WikiPedia.org/wiki/Toise>}. 

573 

574 @arg meter: Value in meter (C{scalar}). 

575 

576 @return: Value in C{toises} (C{float}). 

577 

578 @raise ValueError: Invalid B{C{meter}}. 

579 

580 @see: Function L{m2fathom}. 

581 ''' 

582 return Float(toise=Meter(meter) / _M_TOISE) # * 0.513083632632119 

583 

584 

585def m2yard(meter): 

586 '''Convert meter to I{UK} yards. 

587 

588 @arg meter: Value in meter (C{scalar}). 

589 

590 @return: Value in C{yards} (C{float}). 

591 

592 @raise ValueError: Invalid B{C{meter}}. 

593 ''' 

594 return Float(yard=Meter(meter) / _M_YARD_UK) # * 1.0936132983377078 

595 

596 

597def NM2m(nm): 

598 '''Convert nautical miles to meter (m). 

599 

600 @arg nm: Value in nautical miles (C{scalar}). 

601 

602 @return: Value in meter (C{float}). 

603 

604 @raise ValueError: Invalid B{C{nm}}. 

605 ''' 

606 return Meter(Float(nm=nm) * _M_NM) 

607 

608 

609def radians2m(rad, radius=R_M, lat=0): 

610 '''Convert an angle to a distance along the equator or 

611 along the parallel at an other (geodetic) latitude. 

612 

613 @arg rad: The angle (C{radians}). 

614 @kwarg radius: Mean earth radius, ellipsoid or datum 

615 (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, 

616 L{Datum} or L{a_f2Tuple}). 

617 @kwarg lat: Parallel latitude (C{degrees90}, C{str}). 

618 

619 @return: Distance (C{meter}, same units as B{C{radius}} 

620 or equatorial and polar radii) or C{0.0} for 

621 near-polar B{C{lat}}. 

622 

623 @raise RangeError: Latitude B{C{lat}} outside valid range and 

624 L{pygeodesy.rangerrors} set to C{True}. 

625 

626 @raise TypeError: Invalid B{C{radius}}. 

627 

628 @raise ValueError: Invalid B{C{rad}}, B{C{radius}} or 

629 B{C{lat}}. 

630 

631 @see: Function L{degrees2m} and L{m2radians}. 

632 ''' 

633 return _Radians2m(Lam(rad=rad, clip=0), radius, lat) 

634 

635 

636def _Radians2m(rad, radius, lat): 

637 '''(INTERNAL) Helper for C{degrees2m} and C{radians2m}. 

638 ''' 

639 m = circle4(radius, lat).radius 

640 return _0_0 if m < EPS0 else (rad * m) 

641 

642 

643def radiansPI(deg): 

644 '''Convert and wrap degrees to radians M{[-PI..+PI]}. 

645 

646 @arg deg: Angle (C{degrees}). 

647 

648 @return: Radians, wrapped (C{radiansPI}) 

649 ''' 

650 return wrapPI(radians(deg)) 

651 

652 

653def radiansPI2(deg): 

654 '''Convert and wrap degrees to radians M{[0..+2PI)}. 

655 

656 @arg deg: Angle (C{degrees}). 

657 

658 @return: Radians, wrapped (C{radiansPI2}) 

659 ''' 

660 return _umod_PI2(radians(deg)) 

661 

662 

663def radiansPI_2(deg): 

664 '''Convert and wrap degrees to radians M{[-3PI/2..+PI/2]}. 

665 

666 @arg deg: Angle (C{degrees}). 

667 

668 @return: Radians, wrapped (C{radiansPI_2}) 

669 ''' 

670 return wrapPI_2(radians(deg)) 

671 

672 

673def _sin0cos2(q, r, sign): 

674 '''(INTERNAL) 2-tuple (C{sin(r), cos(r)}) in quadrant C{0 <= B{q} <= 3} 

675 and C{sin} zero I{signed} with B{C{sign}}. 

676 ''' 

677 if r < PI_2: 

678 s, c = sin(r), cos(r) 

679 t = s, c, -s, -c, s 

680 else: # r == PI_2 

681 t = _1_0, _0_0, _N_1_0, _0_0, _1_0 

682# else: # r == 0, testUtility failures 

683# t = _0_0, _1_0, _0_0, _N_1_0, _0_0 

684# q &= 3 

685 s = t[q] or (NEG0 if sign < 0 else _0_0) 

686 c = t[q + 1] or _0_0 

687 return s, c 

688 

689 

690def sincos2(rad): 

691 '''Return the C{sine} and C{cosine} of an angle in C{radians}. 

692 

693 @arg rad: Angle (C{radians}). 

694 

695 @return: 2-Tuple (C{sin(B{rad})}, C{cos(B{rad})}). 

696 

697 @see: U{GeographicLib<https://GeographicLib.SourceForge.io/C++/doc/ 

698 classGeographicLib_1_1Math.html#sincosd>} function U{sincosd 

699 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/ 

700 python/geographiclib/geomath.py#l155>} and C++ U{sincosd 

701 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/ 

702 include/GeographicLib/Math.hpp#l558>}. 

703 ''' 

704 if _isfinite(rad): 

705 q = int(rad * _2__PI) # int(math.floor) 

706 if q < 0: 

707 q -= 1 

708 t = _sin0cos2(q & 3, rad - q * PI_2, rad) 

709 else: 

710 t = NAN, NAN 

711 return t 

712 

713 

714def sincos2_(*rads): 

715 '''Return the C{sine} and C{cosine} of angle(s) in C{radians}. 

716 

717 @arg rads: One or more angles (C{radians}). 

718 

719 @return: Yield the C{sin(B{rad})} and C{cos(B{rad})} for each angle. 

720 

721 @see: function L{sincos2}. 

722 ''' 

723 for r in rads: 

724 s, c = sincos2(r) 

725 yield s 

726 yield c 

727 

728 

729def sincos2d(deg, **adeg): 

730 '''Return the C{sine} and C{cosine} of an angle in C{degrees}. 

731 

732 @arg deg: Angle (C{degrees}). 

733 @kwarg adeg: Optional correction (C{degrees}). 

734 

735 @return: 2-Tuple (C{sin(B{deg_})}, C{cos(B{deg_})}, C{B{deg_} = 

736 B{deg} + B{adeg}}). 

737 

738 @see: U{GeographicLib<https://GeographicLib.SourceForge.io/C++/doc/ 

739 classGeographicLib_1_1Math.html#sincosd>} function U{sincosd 

740 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/ 

741 python/geographiclib/geomath.py#l155>} and C++ U{sincosd 

742 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/ 

743 include/GeographicLib/Math.hpp#l558>}. 

744 ''' 

745 if _isfinite(deg): 

746 q = int(deg * _1__90) # int(math.floor) 

747 if q < 0: 

748 q -= 1 

749 d = deg - q * _90_0 

750 if adeg: 

751 t = _xkwds_get(adeg, adeg=_0_0) 

752 d = _MODS.karney._around(d + t) 

753 t = _sin0cos2(q & 3, radians(d), deg) 

754 else: 

755 t = NAN, NAN 

756 return t 

757 

758 

759def SinCos2(x): 

760 '''Get C{sin} and C{cos} of I{typed} angle. 

761 

762 @arg x: Angle (L{Degrees}, L{Radians} or C{radians}). 

763 

764 @return: 2-Tuple (C{sin(B{x})}, C{cos(B{x})}). 

765 ''' 

766 return sincos2d(x) if isinstance(x, Degrees) else ( 

767 sincos2(x) if isinstance(x, Radians) else 

768 sincos2(float(x))) # assume C{radians} 

769 

770 

771def sincos2d_(*degs): 

772 '''Return the C{sine} and C{cosine} of angle(s) in C{degrees}. 

773 

774 @arg degs: One or more angles (C{degrees}). 

775 

776 @return: Yield the C{sin(B{deg})} and C{cos(B{deg})} for each angle. 

777 

778 @see: Function L{sincos2d}. 

779 ''' 

780 for d in degs: 

781 s, c = sincos2d(d) 

782 yield s 

783 yield c 

784 

785 

786def sincostan3(rad): 

787 '''Return the C{sine}, C{cosine} and C{tangent} of an angle in C{radians}. 

788 

789 @arg rad: Angle (C{radians}). 

790 

791 @return: 3-Tuple (C{sin(B{rad})}, C{cos(B{rad})}, C{tan(B{rad})}). 

792 

793 @see: Function L{sincos2}. 

794 ''' 

795 s, c = sincos2(float(rad)) 

796 t = NAN if s is NAN else (_over(s, c) if s else neg(s, neg0=c < 0)) 

797 return s, c, t 

798 

799 

800def SM2m(sm): 

801 '''Convert statute miles to meter (m). 

802 

803 @arg sm: Value in statute miles (C{scalar}). 

804 

805 @return: Value in meter (C{float}). 

806 

807 @raise ValueError: Invalid B{C{sm}}. 

808 ''' 

809 return Meter(Float(sm=sm) * _M_SM) 

810 

811 

812def tan_2(rad, **semi): # edge=1 

813 '''Compute the tangent of half angle. 

814 

815 @arg rad: Angle (C{radians}). 

816 @kwarg semi: Angle or edge name and index 

817 for semi-circular error. 

818 

819 @return: M{tan(rad / 2)} (C{float}). 

820 

821 @raise ValueError: If B{C{rad}} is semi-circular 

822 and B{C{semi}} is given. 

823 ''' 

824 # .formy.excessKarney_, .sphericalTrigonometry.areaOf 

825 if semi and isnear0(fabs(rad) - PI): 

826 for n, v in semi.items(): 

827 break 

828 n = _SPACE_(n, _radians_) if not isint(v) else \ 

829 _SPACE_(_MODS.streprs.Fmt.SQUARE(**semi), _edge_) 

830 raise _ValueError(n, rad, txt=_semi_circular_) 

831 

832 return tan(rad * _0_5) if _isfinite(rad) else NAN 

833 

834 

835def tand(deg, **error_kwds): 

836 '''Return the C{tangent} of an angle in C{degrees}. 

837 

838 @arg deg: Angle (C{degrees}). 

839 @kwarg error_kwds: Error to raise (C{ValueError}). 

840 

841 @return: C{tan(B{deg})}. 

842 

843 @raise ValueError: If L{pygeodesy.isnear0}C{(cos(B{deg})}. 

844 ''' 

845 s, c = sincos2d(deg) 

846 if s: 

847 if isnear0(c): 

848 raise _valueError(tand, deg, **error_kwds) 

849 s = s / c # /= chokes PyChecker 

850 elif c < 0: 

851 s = -s # negate-0 

852 return s 

853 

854 

855def tand_(*degs, **error_kwds): 

856 '''Return the C{tangent} of angle(s) in C{degrees}. 

857 

858 @arg degs: One or more angles (C{degrees}). 

859 @kwarg error_kwds: Error to raise (C{ValueError}). 

860 

861 @return: Yield the C{tan(B{deg})} for each angle. 

862 

863 @raise ValueError: See L{pygeodesy.tand}. 

864 ''' 

865 for d in degs: 

866 yield tand(d, **error_kwds) 

867 

868 

869def tanPI_2_2(rad): 

870 '''Compute the tangent of half angle, 90 degrees rotated. 

871 

872 @arg rad: Angle (C{radians}). 

873 

874 @return: M{tan((rad + PI/2) / 2)} (C{float}). 

875 ''' 

876 return tan((rad + PI_2) * _0_5) if _isfinite(rad) else ( 

877 NAN if isnan(rad) else (_N_90_0 if rad < 0 else _90_0)) 

878 

879 

880def toise2m(toises): 

881 '''Convert French U{toises<https://WikiPedia.org/wiki/Toise>} to meter. 

882 

883 @arg toises: Value in toises (C{scalar}). 

884 

885 @return: Value in C{meter} (C{float}). 

886 

887 @raise ValueError: Invalid B{C{toises}}. 

888 

889 @see: Function L{fathom2m}. 

890 ''' 

891 return Meter(Float(toises=toises) * _M_TOISE) 

892 

893 

894def truncate(x, ndigits=None): 

895 '''Truncate to the given number of digits. 

896 

897 @arg x: Value to truncate (C{scalar}). 

898 @kwarg ndigits: Number of digits (C{int}), 

899 aka I{precision}. 

900 

901 @return: Truncated B{C{x}} (C{float}). 

902 

903 @see: Python function C{round}. 

904 ''' 

905 if isint(ndigits): 

906 p = _10_0**ndigits 

907 x = int(x * p) / p 

908 return x 

909 

910 

911def unroll180(lon1, lon2, wrap=True): 

912 '''Unroll longitudinal delta and wrap longitude in degrees. 

913 

914 @arg lon1: Start longitude (C{degrees}). 

915 @arg lon2: End longitude (C{degrees}). 

916 @kwarg wrap: If C{True}, wrap and unroll to the M{(-180..+180]} 

917 range (C{bool}). 

918 

919 @return: 2-Tuple C{(B{lon2}-B{lon1}, B{lon2})} unrolled (C{degrees}, 

920 C{degrees}). 

921 

922 @see: Capability C{LONG_UNROLL} in U{GeographicLib 

923 <https://GeographicLib.SourceForge.io/html/python/interface.html#outmask>}. 

924 ''' 

925 d = lon2 - lon1 

926 if wrap: 

927 u = wrap180(d) 

928 if u != d: 

929 return u, (lon1 + u) 

930 return d, lon2 

931 

932 

933def _unrollon(p1, p2, wrap=False): # unroll180 == .karney._unroll2 

934 '''(INTERNAL) Wrap/normalize, unroll and replace longitude. 

935 ''' 

936 lat, lon = p2.lat, p2.lon 

937 if wrap and _Wrap.normal: 

938 lat, lon = _Wrap.latlon(lat, lon) 

939 _, lon = unroll180(p1.lon, lon, wrap=True) 

940 if lat != p2.lat or fabs(lon - p2.lon) > EPS: 

941 p2 = p2.dup(lat=lat, lon=wrap180(lon)) 

942 # p2 = p2.copy(); p2.latlon = lat, wrap180(lon) 

943 return p2 

944 

945 

946def _unrollon3(p1, p2, p3, wrap=False): 

947 '''(INTERNAL) Wrap/normalize, unroll 2 points. 

948 ''' 

949 w = wrap 

950 if w: 

951 w = _Wrap.normal 

952 p2 = _unrollon(p1, p2, wrap=w) 

953 p3 = _unrollon(p1, p3, wrap=w) 

954 p2 = _unrollon(p2, p3) 

955 return p2, p3, w # was wrapped? 

956 

957 

958def unrollPI(rad1, rad2, wrap=True): 

959 '''Unroll longitudinal delta and wrap longitude in radians. 

960 

961 @arg rad1: Start longitude (C{radians}). 

962 @arg rad2: End longitude (C{radians}). 

963 @kwarg wrap: If C{True}, wrap and unroll to the M{(-PI..+PI]} 

964 range (C{bool}). 

965 

966 @return: 2-Tuple C{(B{rad2}-B{rad1}, B{rad2})} unrolled 

967 (C{radians}, C{radians}). 

968 

969 @see: Capability C{LONG_UNROLL} in U{GeographicLib 

970 <https://GeographicLib.SourceForge.io/html/python/interface.html#outmask>}. 

971 ''' 

972 r = rad2 - rad1 

973 if wrap: 

974 u = wrapPI(r) 

975 if u != r: 

976 return u, (rad1 + u) 

977 return r, rad2 

978 

979 

980def _valueError(where, x, **kwds): 

981 '''(INTERNAL) Return a C{_ValueError}. 

982 ''' 

983 t = _MODS.streprs.Fmt.PAREN(where.__name__, x) 

984 return _ValueError(t, **kwds) 

985 

986 

987class _Wrap(object): 

988 

989 _normal = False # default 

990 

991 @property 

992 def normal(self): 

993 '''Get the current L{normal} setting (C{True}, 

994 C{False} or C{None}). 

995 ''' 

996 return self._normal 

997 

998 @normal.setter # PYCHOK setter! 

999 def normal(self, setting): 

1000 '''Set L{normal} to C{True}, C{False} or C{None}. 

1001 ''' 

1002 t = {True: (_MODS.formy.normal, _MODS.formy.normal_), 

1003 False: (self.wraplatlon, self.wraphilam), 

1004 None: (_passargs, _passargs)}.get(setting, ()) 

1005 if t: 

1006 self.latlon, self.philam = t 

1007 self._normal = setting 

1008 

1009 def latlonDMS2(self, lat, lon, **DMS2_kwds): 

1010 if isstr(lat) or isstr(lon): 

1011 kwds = _xkwds(DMS2_kwds, clipLon=0, clipLat=0) 

1012 lat, lon = _MODS.dms.parseDMS2(lat, lon, **kwds) 

1013 return self.latlon(lat, lon) 

1014 

1015# def normalatlon(self, *latlon): 

1016# return _MODS.formy.normal(*latlon) 

1017 

1018# def normalamphi(self, *philam): 

1019# return _MODS.formy.normal_(*philam) 

1020 

1021 def wraplatlon(self, lat, lon): 

1022 return wrap90(lat), wrap180(lon) 

1023 

1024 latlon = wraplatlon # default 

1025 

1026 def latlon3(self, lon1, lat2, lon2, wrap): 

1027 if wrap: 

1028 lat2, lon2 = self.latlon(lat2, lon2) 

1029 lon21, lon2 = unroll180(lon1, lon2) 

1030 else: 

1031 lon21 = lon2 - lon1 

1032 return lon21, lat2, lon2 

1033 

1034 def _latlonop(self, wrap): 

1035 if wrap and self._normal is not None: 

1036 return self.latlon 

1037 else: 

1038 return _passargs 

1039 

1040 def wraphilam(self, phi, lam): 

1041 return wrapPI_2(phi), wrapPI(lam) 

1042 

1043 philam = wraphilam # default 

1044 

1045 def philam3(self, lam1, phi2, lam2, wrap): 

1046 if wrap: 

1047 phi2, lam2 = self.philam(phi2, lam2) 

1048 lam21, lam2 = unrollPI(lam1, lam2) 

1049 else: 

1050 lam21 = lam2 - lam1 

1051 return lam21, phi2, lam2 

1052 

1053 def _philamop(self, wrap): 

1054 if wrap and self._normal is not None: 

1055 return self.philam 

1056 else: 

1057 return _passargs 

1058 

1059 def point(self, ll, wrap=True): # in .points._fractional, -.PointsIter.iterate, ... 

1060 '''Return C{ll} or a copy, I{normalized} or I{wrap}'d. 

1061 ''' 

1062 if wrap and self._normal is not None: 

1063 lat, lon = ll.latlon 

1064 if fabs(lon) > 180 or fabs(lat) > 90: 

1065 _n = self.latlon 

1066 ll = ll.copy(name=_n.__name__) 

1067 ll.latlon = _n(lat, lon) 

1068 return ll 

1069 

1070_Wrap = _Wrap() # PYCHOK singleton 

1071 

1072 

1073# def _wrap(angle, wrap, modulo): 

1074# '''(INTERNAL) Angle wrapper M{((wrap-modulo)..+wrap]}. 

1075# 

1076# @arg angle: Angle (C{degrees}, C{radians} or C{grades}). 

1077# @arg wrap: Range (C{degrees}, C{radians} or C{grades}). 

1078# @arg modulo: Upper limit (360 C{degrees}, PI2 C{radians} or 400 C{grades}). 

1079# 

1080# @return: The B{C{angle}}, wrapped (C{degrees}, C{radians} or C{grades}). 

1081# ''' 

1082# a = float(angle) 

1083# if not (wrap - modulo) <= a < wrap: 

1084# # math.fmod(-1.5, 3.14) == -1.5, but -1.5 % 3.14 == 1.64 

1085# # math.fmod(-1.5, 360) == -1.5, but -1.5 % 360 == 358.5 

1086# a %= modulo 

1087# if a > wrap: 

1088# a -= modulo 

1089# return a 

1090 

1091 

1092def wrap90(deg): 

1093 '''Wrap degrees to M{[-90..+90]}. 

1094 

1095 @arg deg: Angle (C{degrees}). 

1096 

1097 @return: Degrees, wrapped (C{degrees90}). 

1098 ''' 

1099 w = wrap180(deg) 

1100 return (w - _180_0) if w > 90 else ((w + _180_0) if w < -90 else w) 

1101 

1102 

1103def wrap180(deg): 

1104 '''Wrap degrees to M{[-180..+180]}. 

1105 

1106 @arg deg: Angle (C{degrees}). 

1107 

1108 @return: Degrees, wrapped (C{degrees180}). 

1109 ''' 

1110 d = float(deg) 

1111 w = _umod_360(d) 

1112 if w > 180: 

1113 w -= _360_0 

1114 elif d < 0 and w == 180: 

1115 w = -w 

1116 return w 

1117 

1118 

1119def wrap360(deg): # see .streprs._umod_360 

1120 '''Wrap degrees to M{[0..+360)}. 

1121 

1122 @arg deg: Angle (C{degrees}). 

1123 

1124 @return: Degrees, wrapped (C{degrees360}). 

1125 ''' 

1126 return _umod_360(float(deg)) 

1127 

1128 

1129def wrapPI(rad): 

1130 '''Wrap radians to M{[-PI..+PI]}. 

1131 

1132 @arg rad: Angle (C{radians}). 

1133 

1134 @return: Radians, wrapped (C{radiansPI}). 

1135 ''' 

1136 r = float(rad) 

1137 w = _umod_PI2(r) 

1138 if w > PI: 

1139 w -= PI2 

1140 elif r < 0 and w == PI: 

1141 w = -PI 

1142 return w 

1143 

1144 

1145def wrapPI2(rad): 

1146 '''Wrap radians to M{[0..+2PI)}. 

1147 

1148 @arg rad: Angle (C{radians}). 

1149 

1150 @return: Radians, wrapped (C{radiansPI2}). 

1151 ''' 

1152 return _umod_PI2(float(rad)) 

1153 

1154 

1155def wrapPI_2(rad): 

1156 '''Wrap radians to M{[-PI/2..+PI/2]}. 

1157 

1158 @arg rad: Angle (C{radians}). 

1159 

1160 @return: Radians, wrapped (C{radiansPI_2}). 

1161 ''' 

1162 w = wrapPI(rad) 

1163 return (w - PI) if w > PI_2 else ((w + PI) if w < (-PI_2) else w) 

1164 

1165 

1166# def wraplatlon(lat, lon): 

1167# '''Both C{wrap90(B{lat})} and C{wrap180(B{lon})}. 

1168# ''' 

1169# return wrap90(lat), wrap180(lon) 

1170 

1171 

1172def wrap_normal(*normal): 

1173 '''Define the operation for the keyword argument C{B{wrap}=True}, 

1174 across L{pygeodesy}: I{wrap}, I{normalize} or I{no-op}. For 

1175 backward compatibility, the default is I{wrap}. 

1176 

1177 @arg normal: If C{True}, I{normalize} lat- and longitude using 

1178 L{normal} or L{normal_}, if C{False}, I{wrap} the 

1179 lat- and longitude individually by L{wrap90} or 

1180 L{wrapPI_2} respectively L{wrap180}, L{wrapPI} or 

1181 if C{None}, leave lat- and longitude I{unchanged}. 

1182 Do not supply any value to get the current setting. 

1183 

1184 @return: The previous L{wrap_normal} setting (C{bool} or C{None}). 

1185 ''' 

1186 t = _Wrap.normal 

1187 if normal: 

1188 _Wrap.normal = normal[0] 

1189 return t 

1190 

1191 

1192# def wraphilam(phi, lam,): 

1193# '''Both C{wrapPI_2(B{phi})} and C{wrapPI(B{lam})}. 

1194# ''' 

1195# return wrapPI_2(phi), wrapPI(lam) 

1196 

1197 

1198def yard2m(yards): 

1199 '''Convert I{UK} yards to meter. 

1200 

1201 @arg yards: Value in yards (C{scalar}). 

1202 

1203 @return: Value in C{meter} (C{float}). 

1204 

1205 @raise ValueError: Invalid B{C{yards}}. 

1206 ''' 

1207 return Float(yards=yards) * _M_YARD_UK 

1208 

1209# **) MIT License 

1210# 

1211# Copyright (C) 2016-2024 -- mrJean1 at Gmail -- All Rights Reserved. 

1212# 

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

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

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

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

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

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

1219# 

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

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

1222# 

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

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

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

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

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

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

1229# OTHER DEALINGS IN THE SOFTWARE.