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, isint, isstr, neg

14 from 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

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

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

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

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

23

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

25

26__all__ = _ALL_LAZY.utily

27__version__ = '23.10.24'

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

44 def _abs1nan(x):

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

46 '''

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

48

49

50 def 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

56 def 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

69 def 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

82 def 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

88 def 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

95 def 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

104 def 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

115 def 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

141 def _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

147 def 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

159 def 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

181 def 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

199 def 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

216 def 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

236 def 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

253 def 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

263 def 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

273 def 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

283 def 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

293 def 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

320 def 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

335 def 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

354 def 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

366 def 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

376 def 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

386 def 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

396 def 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

406 def 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

418 def _loneg(lon):

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

420 '''

421 return _180_0 - lon

422

423

424 def 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

436 def 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

461 def 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

476 def 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

497 def 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

509 def 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

521 def 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

533 def 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

559 def 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

571 def 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

585 def 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

597 def 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

609 def _passarg(arg): # in .auxilats.auxLat, .formy

610 '''(INTERNAL) Helper, no-op.

611 '''

612 return arg

613

614

615 def _passargs(*args): # in .formy

616 '''(INTERNAL) Helper, no-op.

617 '''

618 return args

619

620

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

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

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

624

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

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

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

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

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

630

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

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

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

634

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

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

637

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

639

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

641 B{C{lat}}.

642

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

644 '''

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

646

647

648 def _Radians2m(rad, radius, lat):

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

650 '''

651 m = circle4(radius, lat).radius

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

653

654

655 def radiansPI(deg):

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

657

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

659

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

661 '''

662 return wrapPI(radians(deg))

663

664

665 def radiansPI2(deg):

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

667

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

669

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

671 '''

672 return _umod_PI2(radians(deg))

673

674

675 def radiansPI_2(deg):

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

677

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

679

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

681 '''

682 return wrapPI_2(radians(deg))

683

684

685 def _sin0cos2(q, r, sign):

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

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

688 '''

689 if r < PI_2:

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

691 t = s, c, -s, -c, s

692 else: # r == PI_2

693 t = _1_0, _0_0, _N_1_0, _0_0, _1_0

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

695# t = _0_0, _1_0, _0_0, _N_1_0, _0_0

696# q &= 3

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

698 c = t[q + 1] or _0_0

699 return s, c

700

701

702 def sincos2(rad):

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

704

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

706

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

708

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

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

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

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

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

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

715 '''

716 if _isfinite(rad):

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

718 if q < 0:

719 q -= 1

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

721 else:

722 t = NAN, NAN

723 return t

724

725

726 def sincos2_(*rads):

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

728

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

730

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

732

733 @see: function L{sincos2}.

734 '''

735 for r in rads:

736 s, c = sincos2(r)

737 yield s

738 yield c

739

740

741 def sincos2d(deg, **adeg):

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

743

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

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

746

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

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

749

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

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

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

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

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

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

756 '''

757 if _isfinite(deg):

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

759 if q < 0:

760 q -= 1

761 d = deg - q * _90_0

762 if adeg:

763 t = _xkwds_get(adeg, adeg=_0_0)

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

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

766 else:

767 t = NAN, NAN

768 return t

769

770

771 def SinCos2(x):

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

773

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

775

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

777 '''

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

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

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

781

782

783 def sincos2d_(*degs):

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

785

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

787

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

789

790 @see: Function L{sincos2d}.

791 '''

792 for d in degs:

793 s, c = sincos2d(d)

794 yield s

795 yield c

796

797

798 def sincostan3(rad):

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

800

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

802

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

804

805 @see: Function L{sincos2}.

806 '''

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

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

809 return s, c, t

810

811

812 def SM2m(sm):

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

814

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

816

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

818

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

820 '''

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

822

823

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

825 '''Compute the tangent of half angle.

826

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

828 @kwarg semi: Angle or edge name and index

829 for semi-circular error.

830

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

832

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

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

835 '''

836 # .formy.excessKarney_, .sphericalTrigonometry.areaOf

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

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

839 break

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

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

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

843

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

845

846

847 def tand(deg, **error_kwds):

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

849

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

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

852

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

854

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

856 '''

857 s, c = sincos2d(deg)

858 if s:

859 if isnear0(c):

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

861 s = s / c # /= chokes PyChecker

862 elif c < 0:

863 s = -s # negate-0

864 return s

865

866

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

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

869

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

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

872

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

874

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

876 '''

877 for d in degs:

878 yield tand(d, **error_kwds)

879

880

881 def tanPI_2_2(rad):

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

883

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

885

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

887 '''

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

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

890

891

892 def toise2m(toises):

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

894

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

896

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

898

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

900

901 @see: Function L{fathom2m}.

902 '''

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

904

905

906 def truncate(x, ndigits=None):

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

908

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

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

911 aka I{precision}.

912

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

914

915 @see: Python function C{round}.

916 '''

917 if isint(ndigits):

918 p = _10_0**ndigits

919 x = int(x * p) / p

920 return x

921

922

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

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

925

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

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

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

929 range (C{bool}).

930

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

932 C{degrees}).

933

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

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

936 '''

937 d = lon2 - lon1

938 if wrap:

939 u = wrap180(d)

940 if u != d:

941 return u, (lon1 + u)

942 return d, lon2

943

944

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

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

947 '''

948 lat, lon = p2.lat, p2.lon

949 if wrap and _Wrap.normal:

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

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

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

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

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

955 return p2

956

957

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

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

960 '''

961 w = wrap

962 if w:

963 w = _Wrap.normal

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

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

966 p2 = _unrollon(p2, p3)

967 return p2, p3, w # was wrapped?

968

969

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

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

972

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

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

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

976 range (C{bool}).

977

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

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

980

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

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

983 '''

984 r = rad2 - rad1

985 if wrap:

986 u = wrapPI(r)

987 if u != r:

988 return u, (rad1 + u)

989 return r, rad2

990

991

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

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

994 '''

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

996 return _ValueError(t, **kwds)

997

998

999 class _Wrap(object):

1000

1001 _normal = False # default

1002

1003 @property

1004 def normal(self):

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

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

1007 '''

1008 return self._normal

1009

1010 @normal.setter # PYCHOK setter!

1011 def normal(self, setting):

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

1013 '''

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

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

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

1017 if t:

1018 self.latlon, self.philam = t

1019 self._normal = setting

1020

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

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

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

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

1025 return self.latlon(lat, lon)

1026

1027# def normalatlon(self, *latlon):

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

1029

1030# def normalamphi(self, *philam):

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

1032

1033 def wraplatlon(self, lat, lon):

1034 return wrap90(lat), wrap180(lon)

1035

1036 latlon = wraplatlon # default

1037

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

1039 if wrap:

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

1041 lon21, lon2 = unroll180(lon1, lon2)

1042 else:

1043 lon21 = lon2 - lon1

1044 return lon21, lat2, lon2

1045

1046 def _latlonop(self, wrap):

1047 if wrap and self._normal is not None:

1048 return self.latlon

1049 else:

1050 return _passargs

1051

1052 def wraphilam(self, phi, lam):

1053 return wrapPI_2(phi), wrapPI(lam)

1054

1055 philam = wraphilam # default

1056

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

1058 if wrap:

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

1060 lam21, lam2 = unrollPI(lam1, lam2)

1061 else:

1062 lam21 = lam2 - lam1

1063 return lam21, phi2, lam2

1064

1065 def _philamop(self, wrap):

1066 if wrap and self._normal is not None:

1067 return self.philam

1068 else:

1069 return _passargs

1070

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

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

1073 '''

1074 if wrap and self._normal is not None:

1075 lat, lon = ll.latlon

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

1077 _n = self.latlon

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

1079 ll.latlon = _n(lat, lon)

1080 return ll

1081

1082_Wrap = _Wrap() # PYCHOK singleton

1083

1084

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

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

1087#

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

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

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

1091#

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

1093# '''

1094# a = float(angle)

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

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

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

1098# a %= modulo

1099# if a > wrap:

1100# a -= modulo

1101# return a

1102

1103

1104 def wrap90(deg):

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

1106

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

1108

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

1110 '''

1111 w = wrap180(deg)

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

1113

1114

1115 def wrap180(deg):

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

1117

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

1119

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

1121 '''

1122 d = float(deg)

1123 w = _umod_360(d)

1124 if w > 180:

1125 w -= _360_0

1126 elif d < 0 and w == 180:

1127 w = -w

1128 return w

1129

1130

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

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

1133

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

1135

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

1137 '''

1138 return _umod_360(float(deg))

1139

1140

1141 def wrapPI(rad):

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

1143

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

1145

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

1147 '''

1148 r = float(rad)

1149 w = _umod_PI2(r)

1150 if w > PI:

1151 w -= PI2

1152 elif r < 0 and w == PI:

1153 w = -PI

1154 return w

1155

1156

1157 def wrapPI2(rad):

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

1159

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

1161

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

1163 '''

1164 return _umod_PI2(float(rad))

1165

1166

1167 def wrapPI_2(rad):

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

1169

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

1171

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

1173 '''

1174 w = wrapPI(rad)

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

1176

1177

1178# def wraplatlon(lat, lon):

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

1180# '''

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

1182

1183

1184 def wrap_normal(*normal):

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

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

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

1188

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

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

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

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

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

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

1195

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

1197 '''

1198 t = _Wrap.normal

1199 if normal:

1200 _Wrap.normal = normal[0]

1201 return t

1202

1203

1204# def wraphilam(phi, lam,):

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

1206# '''

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

1208

1209

1210 def yard2m(yards):

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

1212

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

1214

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

1216

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

1218 '''

1219 return Float(yards=yards) * _M_YARD_UK

1220

1221# **) MIT License

1222#

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

1224#

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

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

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

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

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

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

1231#

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

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

1234#

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

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

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

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

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

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

1241# OTHER DEALINGS IN THE SOFTWARE.