14 # pragma warning (disable: 4127) 27 , _omega2(
Math::sq(_omega))
28 , _aomega2(
Math::sq(_omega * _a))
33 throw GeographicErr(
"Gravitational constants is not positive");
34 if (!(_omega == 0 && _f == 0 && _J2 == 0)) {
36 if (_J2 > 0 && Math::isfinite(_J2) && flatp)
38 if (!(_J2 > 0 && Math::isfinite(_J2)) && !flatp)
40 if (!(Math::isfinite(_omega) && _omega != 0))
48 _ep2 = _e2 / (1 - _e2);
53 _U0 = _GM / _E * atan(sqrt(_ep2)) + _aomega2 / 3;
56 _m = _aomega2 * _b / _GM;
58 Q = _m * sqrt(_ep2) * qpf(_ep2) / (3 * _q0),
60 _gammae = _GM / (_a * _b) * G;
61 _gammap = _GM / (_a * _a) * (1 + Q);
63 _k = (_m + 3 * Q / 2 - _e2 * (1 + Q)) / G;
65 _fstar = (_m + 3 * Q / 2 - _f * (1 + Q)) / G;
86 if (abs(x) >= real(0.5)) {
87 real y = sqrt(abs(x)), x2 =
Math::sq(x);
88 return ((x > 0 ? atan(y) :
Math::atanh(y))- y * (1 - x / 3 + x2 / 5)) /
92 for (
int n = 7; ; n += 2) {
105 {
return 1/real(5) + x * atan7(x); }
115 return sqrt(ep2) * ep2 * (3 * (3 + ep2) * atan5(ep2) - 1) / 6;
121 return sqrt(ep2) * (5 - 3 * (1 + ep2) * (1 + 5 * ep2 * atan7(ep2))) /
133 return ep2 * (1 - 3 * (1 + ep2) * atan5(ep2));
142 for (
int j = n; j--;)
145 -3 * e2n * ((1 - n) + 5 * n * _J2 / _e2) / ((2 * n + 1) * (2 * n + 3));
151 sphi2 = abs(lat) == 90 ? 1 :
Math::sq(sin(phi));
153 return _gammae * (1 + _k * sphi2) / sqrt(1 - _e2 * sphi2);
157 real& GammaX, real& GammaY, real& GammaZ)
170 u = sqrt((Q >= 0 ? (Q + disc) : t2 / (disc - Q)) / 2),
176 cbet = s ? cbet/s : 0;
177 sbet = s ? sbet/s : 1;
185 Vres = (_GM / _E * atan(_E / u)
186 + _aomega2 * q * (
Math::sq(sbet) - 1/real(3)) / 2),
189 + (_aomega2 * _E * qp
191 gamb = _aomega2 * q * sbet * cbet * invw / uE,
194 GammaX = t * cbet * gamu - invw * sbet * gamb;
195 GammaY = GammaX * slam;
197 GammaZ = invw * sbet * gamu + t * cbet * gamb;
210 real& gammaX, real& gammaY, real& gammaZ)
213 real Ures =
V0(X, Y, Z, gammaX, gammaY, gammaZ) +
Phi(X, Y, fX, fY);
220 real& gammay, real& gammaz)
223 real M[Geocentric::dim2_];
224 _earth.IntForward(lat, 0, h, X, Y, Z, M);
225 real gammaX, gammaY, gammaZ,
226 Ures =
U(X, Y, Z, gammaX, gammaY, gammaZ);
228 gammay = M[1] * gammaX + M[4] * gammaY + M[7] * gammaZ;
229 gammaz = M[2] * gammaX + M[5] * gammaY + M[8] * gammaZ;
234 real omega, real J2) {
236 K = 2 *
Math::sq(a * omega) * a / (15 * GM),
244 h = e2 * (1 - sqrt(e2) * K / q0) - 3 * J2,
245 dh = 1 - sqrt(e2) * K * (3 * q0 - 2 * e2 * dq0) / (2 *
Math::sq(q0)),
251 return e2 / (1 + sqrt(1 - e2));
255 real omega, real f) {
257 K = 2 *
Math::sq(a * omega) * a / (15 * GM),
259 q0 = qf(e2 / (1 - e2));
260 return e2 * (1 - K * sqrt(e2) / q0) / 3;
Math::real SurfaceGravity(real lat) const
Math::real Gravity(real lat, real h, real &gammay, real &gammaz) const
The normal gravity of the earth.
static bool isfinite(T x)
Mathematical functions needed by GeographicLib.
static Math::real FlatteningToJ2(real a, real GM, real omega, real f)
Namespace for GeographicLib.
static Math::real J2ToFlattening(real a, real GM, real omega, real J2)
static const NormalGravity & GRS80()
Header for GeographicLib::NormalGravity class.
Exception handling for GeographicLib.
Math::real U(real X, real Y, real Z, real &gammaX, real &gammaY, real &gammaZ) const
static const NormalGravity & WGS84()
Math::real Phi(real X, real Y, real &fX, real &fY) const
#define GEOGRAPHICLIB_PANIC
Math::real V0(real X, real Y, real Z, real &GammaX, real &GammaY, real &GammaZ) const