Math.cpp
Go to the documentation of this file.
1 
10 #include <GeographicLib/Math.hpp>
11 
12 #if defined(_MSC_VER)
13 // Squelch warnings about constant conditional expressions
14 # pragma warning (disable: 4127)
15 #endif
16 
17 namespace GeographicLib {
18 
19  using namespace std;
20 
21  template<typename T> T Math::eatanhe(T x, T es) {
22  return es > T(0) ? es * atanh(es * x) : -es * atan(es * x);
23  }
24 
25  template<typename T> T Math::taupf(T tau, T es) {
26  T tau1 = hypot(T(1), tau),
27  sig = sinh( eatanhe(tau / tau1, es ) );
28  return hypot(T(1), sig) * tau - sig * tau1;
29  }
30 
31  template<typename T> T Math::tauf(T taup, T es) {
32  static const int numit = 5;
33  static const T tol = sqrt(numeric_limits<T>::epsilon()) / T(10);
34  T e2m = T(1) - sq(es),
35  // To lowest order in e^2, taup = (1 - e^2) * tau = _e2m * tau; so use
36  // tau = taup/_e2m as a starting guess. (This starting guess is the
37  // geocentric latitude which, to first order in the flattening, is equal
38  // to the conformal latitude.) Only 1 iteration is needed for |lat| <
39  // 3.35 deg, otherwise 2 iterations are needed. If, instead, tau = taup
40  // is used the mean number of iterations increases to 1.99 (2 iterations
41  // are needed except near tau = 0).
42  tau = taup/e2m,
43  stol = tol * max(T(1), abs(taup));
44  // min iterations = 1, max iterations = 2; mean = 1.94
45  for (int i = 0; i < numit || GEOGRAPHICLIB_PANIC; ++i) {
46  T taupa = taupf(tau, es),
47  dtau = (taup - taupa) * (1 + e2m * sq(tau)) /
48  ( e2m * hypot(T(1), tau) * hypot(T(1), taupa) );
49  tau += dtau;
50  if (!(abs(dtau) >= stol))
51  break;
52  }
53  return tau;
54  }
55 
57  // Instantiate
58  template Math::real Math::eatanhe<Math::real>(Math::real, Math::real);
59  template Math::real Math::taupf<Math::real>(Math::real, Math::real);
60  template Math::real Math::tauf<Math::real>(Math::real, Math::real);
62 
63 } // namespace GeographicLib
es
EigenSolver< MatrixXf > es
Definition: EigenSolver_compute.cpp:1
atan
const EIGEN_DEVICE_FUNC AtanReturnType atan() const
Definition: ArrayCwiseUnaryOps.h:283
GeographicLib
Namespace for GeographicLib.
Definition: JacobiConformal.hpp:15
x
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
Definition: gnuplot_common_settings.hh:12
GeographicLib::Math::tauf
static T tauf(T taup, T es)
Definition: Math.cpp:31
T
Eigen::Triplet< double > T
Definition: Tutorial_sparse_example.cpp:6
epsilon
static double epsilon
Definition: testRot3.cpp:37
GeographicLib::Math::real
double real
Definition: Math.hpp:129
GEOGRAPHICLIB_PANIC
#define GEOGRAPHICLIB_PANIC
Definition: Math.hpp:87
Eigen::Triplet< double >
Math.hpp
Header for GeographicLib::Math class.
std
Definition: BFloat16.h:88
GeographicLib::Math::taupf
static T taupf(T tau, T es)
Definition: Math.cpp:25
gtsam::tol
const G double tol
Definition: Group.h:79
abs
#define abs(x)
Definition: datatypes.h:17
GeographicLib::Math::eatanhe
static T eatanhe(T x, T es)
Definition: Math.cpp:21
sinh
const EIGEN_DEVICE_FUNC SinhReturnType sinh() const
Definition: ArrayCwiseUnaryOps.h:339
max
#define max(a, b)
Definition: datatypes.h:20
ceres::sqrt
Jet< T, N > sqrt(const Jet< T, N > &f)
Definition: jet.h:418
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9


gtsam
Author(s):
autogenerated on Sat Nov 16 2024 04:03:00