1

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

3

4 u'''Various utility functions.

5

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

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

8 U{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

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

12

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

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

15 _M_KM, _M_NM, _M_SM, _0_0, _1__90, _0_5, _1_0, \

16 _N_1_0, _2__PI, _10_0, _90_0, _N_90_0, _180_0, \

17 _N_180_0, _360_0, _400_0, _copysign_0_0, \

18 _float as _F, _isfinite, isnan, isnear0, \

19 _over, _umod_360, _umod_PI2

20 from pygeodesy.errors import _ValueError, _xkwds, _xkwds_get1, _ALL_LAZY, _MODS

21 from pygeodesy.internals import _passargs # , _MODS?

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

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

24 from pygeodesy.units import Degrees, Degrees_, Feet, Float, Lam, Lamd, \

25 Meter, Meter2, Radians, Radians_

26

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

28

29__all__ = _ALL_LAZY.utily

30__version__ = '24.06.15'

31

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

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

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

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

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

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

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

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

40# _M_KM = _F(1000.0) # kilo meter

41# _M_NM = _F(1852.0) # nautical mile

42# _M_SM = _F(1609.344) # statute mile

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

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

45

46

47 def _abs1nan(x):

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

49 '''

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

51

52

53 def acos1(x):

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

55 '''

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

57

58

59 def acre2ha(acres):

60 '''Convert acres to hectare.

61

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

63

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

65

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

67 '''

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

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

70

71

72 def acre2m2(acres):

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

74

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

76

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

78

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

80 '''

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

82 return Meter2(Float(acres) * 4046.8564224)

83

84

85 def asin1(x):

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

87 '''

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

89

90

91 def atan1(y, x=_1_0):

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

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

94 '''

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

96

97

98 def atan1d(y, x=_1_0):

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

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

101

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

103 '''

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

105

106

107 def atan2b(y, x):

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

109

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

111 '''

112 b = atan2d(y, x)

113 if b < 0:

114 b += _360_0

115 return b

116

117

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

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

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

121

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

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

124 '''

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

126 if y < 0: # q = 3

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

128 else: # q = 2

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

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

131 return NAN

132 elif y:

133 if x > 0: # q = 0

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

135 elif x < 0: # q = 1

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

137 else: # x == 0

138 d = _copysign(_90_0, y)

139 else:

140 d = _180_0 if x < 0 else _0_0

141 return _azireversed(d) if reverse else d

142

143

144 def _azireversed(azimuth): # in .rhumbBase

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

146 '''

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

148

149

150 def chain2m(chains):

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

152

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

154

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

156

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

158 '''

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

160

161

162 def circle4(earth, lat):

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

164

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

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

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

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

169

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

171 instance.

172

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

174 L{rangerrors<pygeodesy.rangerrors>} is C{True}.

175

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

177

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

179 '''

180 E = _MODS.datums._earth_ellipsoid(earth)

181 return E.circle4(lat)

182

183

184 def cot(rad, **error_kwds):

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

186

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

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

189

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

191

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

193 '''

194 s, c = sincos2(rad)

195 if c:

196 if isnear0(s):

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

198 c = c / s # /= chokes PyChecker

199 return c

200

201

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

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

204

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

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

207

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

209

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

211 '''

212 try:

213 for r in rads:

214 yield cot(r)

215 except ValueError:

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

217

218

219 def cotd(deg, **error_kwds):

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

221

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

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

224

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

226

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

228 '''

229 s, c = sincos2d(deg)

230 if c:

231 if isnear0(s):

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

233 c = c / s # /= chokes PyChecker

234 elif s < 0:

235 c = -c # negate-0

236 return c

237

238

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

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

241

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

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

244

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

246

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

248 '''

249 try:

250 for d in degs:

251 yield cotd(d)

252 except ValueError:

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

254

255

256 def degrees90(rad):

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

258

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

260

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

262 '''

263 return wrap90(degrees(rad))

264

265

266 def degrees180(rad):

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

268

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

270

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

272 '''

273 return wrap180(degrees(rad))

274

275

276 def degrees360(rad):

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

278

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

280

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

282 '''

283 return _umod_360(degrees(rad))

284

285

286 def degrees2grades(deg):

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

288

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

290

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

292 '''

293 return Degrees(deg) * _400_0 / _360_0

294

295

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

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

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

299

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

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

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

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

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

305

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

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

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

309

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

311 L{rangerrors<pygeodesy.rangerrors>} is C{True}.

312

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

314

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

316 B{C{lat}}.

317

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

319 '''

320 return _Radians2m(Lamd(deg=deg, clip=0), radius, lat)

321

322

323 def fathom2m(fathoms):

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

325

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

327

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

329

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

331

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

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

334 '''

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

336

337

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

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

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

341

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

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

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

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

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

347

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

349

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

351 '''

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

353 (_M_FOOT_FR if pied else

354 (_M_FOOT_GE if fuss else _M_FOOT))))

355

356

357 def furlong2m(furlongs):

358 '''Convert a furlong to meter.

359

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

361

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

363

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

365 '''

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

367

368

369 def grades(rad):

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

371

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

373

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

375 '''

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

377

378

379 def grades400(rad):

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

381

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

383

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

385 '''

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

387

388

389 def grades2degrees(gon):

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

391

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

393

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

395 '''

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

397

398

399 def grades2radians(gon):

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

401

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

403

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

405 '''

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

407

408

409 def km2m(km):

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

411

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

413

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

415

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

417 '''

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

419

420

421 def _loneg(lon):

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

423 '''

424 return _180_0 - lon

425

426

427 def m2chain(meter):

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

429

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

431

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

433

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

435 '''

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

437

438

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

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

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

442

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

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

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

446 L{a_f2Tuple}).

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

448

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

450

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

452 L{rangerrors<pygeodesy.rangerrors>} is C{True}.

453

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

455

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

457 or B{C{lat}}.

458

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

460 '''

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

462

463

464 def m2fathom(meter):

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

466

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

468

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

470

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

472

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

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

475 '''

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

477

478

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

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

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

482

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

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

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

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

487 I{International} foot.

488

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

490

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

492 '''

493 # * 3.2808333333333333, US Survey 3937 / 1200

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

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

496 (_M_FOOT_FR if pied else

497 (_M_FOOT_GE if fuss else _M_FOOT))))

498

499

500 def m2furlong(meter):

501 '''Convert meter to furlongs.

502

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

504

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

506

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

508 '''

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

510

511

512 def m2km(meter):

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

514

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

516

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

518

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

520 '''

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

522

523

524 def m2NM(meter):

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

526

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

528

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

530

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

532 '''

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

534

535

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

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

538 parallel at an other (geodetic) latitude.

539

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

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

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

543 L{a_f2Tuple}).

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

545

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

547

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

549 L{rangerrors<pygeodesy.rangerrors>} is C{True}.

550

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

552

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

554 or B{C{lat}}.

555

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

557 '''

558 m = circle4(radius, lat).radius

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

560

561

562 def m2SM(meter):

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

564

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

566

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

568

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

570 '''

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

572

573

574 def m2toise(meter):

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

576

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

578

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

580

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

582

583 @see: Function L{m2fathom}.

584 '''

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

586

587

588 def m2yard(meter):

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

590

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

592

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

594

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

596 '''

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

598

599

600 def NM2m(nm):

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

602

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

604

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

606

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

608 '''

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

610

611

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

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

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

615

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

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

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

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

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

621

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

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

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

625

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

627 L{rangerrors<pygeodesy.rangerrors>} is C{True}.

628

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

630

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

632 B{C{lat}}.

633

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

635 '''

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

637

638

639 def _Radians2m(rad, radius, lat):

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

641 '''

642 m = circle4(radius, lat).radius

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

644

645

646 def radiansPI(deg):

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

648

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

650

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

652 '''

653 return wrapPI(radians(deg))

654

655

656 def radiansPI2(deg):

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

658

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

660

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

662 '''

663 return _umod_PI2(radians(deg))

664

665

666 def radiansPI_2(deg):

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

668

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

670

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

672 '''

673 return wrapPI_2(radians(deg))

674

675

676 def _sin0cos2(q, r, sign):

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

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

679 '''

680 if r < PI_2:

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

682 t = s, c, -s, -c, s

683 else: # r == PI_2

684 t = _1_0, _0_0, _N_1_0, _0_0, _1_0

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

686# t = _0_0, _1_0, _0_0, _N_1_0, _0_0

687# q &= 3

688 s = t[q] or _copysign_0_0(sign)

689 c = t[q + 1] or _0_0

690 return s, c

691

692

693 def SinCos2(x):

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

695

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

697

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

699 '''

700 return sincos2d(x) if isinstanceof(x, Degrees, Degrees_) else (

701 sincos2(x) if isinstanceof(x, Radians, Radians_) else

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

703

704

705 def sincos2(rad):

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

707

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

709

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

711

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

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

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

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

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

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

718 '''

719 if _isfinite(rad):

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

721 if q < 0:

722 q -= 1

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

724 else:

725 t = NAN, NAN

726 return t

727

728

729 def sincos2_(*rads):

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

731

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

733

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

735

736 @see: function L{sincos2}.

737 '''

738 for r in rads:

739 s, c = sincos2(r)

740 yield s

741 yield c

742

743

744 def sincos2d(deg, **adeg):

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

746

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

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

749

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

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

752

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

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

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

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

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

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

759 '''

760 if _isfinite(deg):

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

762 if q < 0:

763 q -= 1

764 d = deg - q * _90_0

765 if adeg:

766 t = _xkwds_get1(adeg, adeg=_0_0)

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

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

769 else:

770 t = NAN, NAN

771 return t

772

773

774 def sincos2d_(*degs):

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

776

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

778

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

780

781 @see: Function L{sincos2d}.

782 '''

783 for d in degs:

784 s, c = sincos2d(d)

785 yield s

786 yield c

787

788

789 def sincostan3(rad):

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

791

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

793

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

795

796 @see: Function L{sincos2}.

797 '''

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

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

800 return s, c, t

801

802

803 def SM2m(sm):

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

805

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

807

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

809

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

811 '''

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

813

814

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

816 '''Compute the tangent of half angle.

817

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

819 @kwarg semi: Angle or edge name and index

820 for semi-circular error.

821

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

823

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

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

826 '''

827 # .formy.excessKarney_, .sphericalTrigonometry.areaOf

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

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

830 break

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

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

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

834

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

836

837

838 def tand(deg, **error_kwds):

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

840

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

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

843

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

845

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

847 '''

848 s, c = sincos2d(deg)

849 if s:

850 if isnear0(c):

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

852 s = s / c # /= chokes PyChecker

853 elif c < 0:

854 s = -s # negate-0

855 return s

856

857

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

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

860

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

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

863

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

865

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

867 '''

868 for d in degs:

869 yield tand(d, **error_kwds)

870

871

872 def tanPI_2_2(rad):

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

874

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

876

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

878 '''

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

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

881

882

883 def toise2m(toises):

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

885

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

887

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

889

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

891

892 @see: Function L{fathom2m}.

893 '''

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

895

896

897 def truncate(x, ndigits=None):

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

899

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

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

902 aka I{precision}.

903

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

905

906 @see: Python function C{round}.

907 '''

908 if isint(ndigits):

909 p = _10_0**ndigits

910 x = int(x * p) / p

911 return x

912

913

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

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

916

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

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

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

920 range (C{bool}).

921

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

923 C{degrees}).

924

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

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

927 '''

928 d = lon2 - lon1

929 if wrap:

930 u = wrap180(d)

931 if u != d:

932 return u, (lon1 + u)

933 return d, lon2

934

935

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

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

938 '''

939 lat, lon = p2.lat, p2.lon

940 if wrap and _Wrap.normal:

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

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

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

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

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

946 return p2

947

948

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

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

951 '''

952 w = wrap

953 if w:

954 w = _Wrap.normal

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

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

957 p2 = _unrollon(p2, p3)

958 return p2, p3, w # was wrapped?

959

960

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

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

963

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

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

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

967 range (C{bool}).

968

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

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

971

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

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

974 '''

975 r = rad2 - rad1

976 if wrap:

977 u = wrapPI(r)

978 if u != r:

979 return u, (rad1 + u)

980 return r, rad2

981

982

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

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

985 '''

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

987 return _ValueError(t, **kwds)

988

989

990 class _Wrap(object):

991

992 _normal = False # default

993

994 @property

995 def normal(self):

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

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

998 '''

999 return self._normal

1000

1001 @normal.setter # PYCHOK setter!

1002 def normal(self, setting):

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

1004 '''

1005 m = _MODS.formy

1006 t = {True: (m.normal, m.normal_),

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

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

1009 if t:

1010 self.latlon, self.philam = t

1011 self._normal = setting

1012

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

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

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

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

1017 return self.latlon(lat, lon)

1018

1019# def normalatlon(self, *latlon):

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

1021

1022# def normalamphi(self, *philam):

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

1024

1025 def wraplatlon(self, lat, lon):

1026 return wrap90(lat), wrap180(lon)

1027

1028 latlon = wraplatlon # default

1029

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

1031 if wrap:

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

1033 lon21, lon2 = unroll180(lon1, lon2)

1034 else:

1035 lon21 = lon2 - lon1

1036 return lon21, lat2, lon2

1037

1038 def _latlonop(self, wrap):

1039 if wrap and self._normal is not None:

1040 return self.latlon

1041 else:

1042 return _passargs

1043

1044 def wraphilam(self, phi, lam):

1045 return wrapPI_2(phi), wrapPI(lam)

1046

1047 philam = wraphilam # default

1048

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

1050 if wrap:

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

1052 lam21, lam2 = unrollPI(lam1, lam2)

1053 else:

1054 lam21 = lam2 - lam1

1055 return lam21, phi2, lam2

1056

1057 def _philamop(self, wrap):

1058 if wrap and self._normal is not None:

1059 return self.philam

1060 else:

1061 return _passargs

1062

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

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

1065 '''

1066 if wrap and self._normal is not None:

1067 lat, lon = ll.latlon

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

1069 _n = self.latlon

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

1071 ll.latlon = _n(lat, lon)

1072 return ll

1073

1074_Wrap = _Wrap() # PYCHOK singleton

1075

1076

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

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

1079#

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

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

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

1083#

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

1085# '''

1086# a = float(angle)

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

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

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

1090# a %= modulo

1091# if a > wrap:

1092# a -= modulo

1093# return a

1094

1095

1096 def wrap90(deg):

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

1098

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

1100

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

1102 '''

1103 w = wrap180(deg)

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

1105

1106

1107 def wrap180(deg):

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

1109

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

1111

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

1113 '''

1114 d = float(deg)

1115 w = _umod_360(d)

1116 if w > 180:

1117 w -= _360_0

1118 elif d < 0 and w == 180:

1119 w = -w

1120 return w

1121

1122

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

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

1125

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

1127

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

1129 '''

1130 return _umod_360(float(deg))

1131

1132

1133 def wrapPI(rad):

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

1135

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

1137

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

1139 '''

1140 r = float(rad)

1141 w = _umod_PI2(r)

1142 if w > PI:

1143 w -= PI2

1144 elif r < 0 and w == PI:

1145 w = -PI

1146 return w

1147

1148

1149 def wrapPI2(rad):

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

1151

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

1153

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

1155 '''

1156 return _umod_PI2(float(rad))

1157

1158

1159 def wrapPI_2(rad):

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

1161

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

1163

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

1165 '''

1166 w = wrapPI(rad)

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

1168

1169

1170# def wraplatlon(lat, lon):

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

1172# '''

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

1174

1175

1176 def wrap_normal(*normal):

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

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

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

1180

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

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

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

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

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

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

1187

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

1189 '''

1190 t = _Wrap.normal

1191 if normal:

1192 _Wrap.normal = normal[0]

1193 return t

1194

1195

1196# def wraphilam(phi, lam,):

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

1198# '''

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

1200

1201

1202 def yard2m(yards):

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

1204

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

1206

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

1208

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

1210 '''

1211 return Float(yards=yards) * _M_YARD_UK

1212

1213# **) MIT License

1214#

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

1216#

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

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

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

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

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

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

1223#

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

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

1226#

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

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

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

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

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

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

1233# OTHER DEALINGS IN THE SOFTWARE.