GeographicLib  2.1.2
AuxLatitude.cpp
Go to the documentation of this file.
1 /**
2  * \file AuxLatitude.cpp
3  * \brief Implementation for the GeographicLib::AuxLatitude and
4  * GeographicLib::AuxAngle classes.
5  *
6  * \note This is just sample code. It is not part of GeographicLib itself.
7  *
8  * This file is an implementation of the methods described in
9  * - C. F. F. Karney,
10  * On auxiliary latitudes,
11  * Technical Report, SRI International, December 2022.
12  * https://arxiv.org/abs/2212.05818
13  * .
14  * Copyright (c) Charles Karney (2022) <charles@karney.com> and licensed under
15  * the MIT/X11 License. For more information, see
16  * https://geographiclib.sourceforge.io/
17  **********************************************************************/
18 
19 #include "AuxLatitude.hpp"
20 
21 /// \cond SKIP
22 
23 #if !defined(AUXLATITUDE_UNOPT)
24 // If 1 use the unoptimized formulas in the paper. Otherwise use the
25 // optimizations to improve roundoff for high eccentricity.
26 #define AUXLATITUDE_UNOPT 0
27 #endif
28 
29 #if !defined(AUXLATITUDE_NATIVE_ELLIPTIC)
30 // If 1 use the templated implementation of RD and RF in this file. Otherwise
31 // use the non-templated versions in the EllipticFunction class.
32 #define AUXLATITUDE_NATIVE_ELLIPTIC 1
33 #endif
34 
35 /// \endcond
36 
37 #if !AUXLATITUDE_NATIVE_ELLIPTIC
38 // Use GeographicLib::EllipticFunction class
40 #endif
41 
42 #if defined(_MSC_VER)
43 // Squelch warnings about constant conditional expressions
44 # pragma warning (disable: 4127)
45 #endif
46 
47 namespace GeographicLib {
48 
49  using namespace std;
50 
51  template<typename T>
53  return AuxAngle(std::numeric_limits<real>::quiet_NaN(),
54  std::numeric_limits<real>::quiet_NaN());
55  }
56 
57  template<typename T>
59  if ( isnan( tan() ) ||
60  (fabs(_y) > std::numeric_limits<real>::max()/2 &&
61  fabs(_x) > std::numeric_limits<real>::max()/2) )
62  // deal with
63  // (0,0), (inf,inf), (nan,nan), (nan,x), (y,nan), (toobig,toobig)
64  return NaN();
65  real r = hypot(_y, _x),
66  y = _y/r, x = _x/r;
67  // deal with r = inf, then one of y,x
68  if (isnan(y)) y = copysign(real(1), _y);
69  if (isnan(x)) x = copysign(real(1), _x);
70  return AuxAngle(y, x);
71  }
72 
73  template<typename T>
75  return AuxAngle(copysign(y(), p.y()), copysign(x(), p.x()));
76  }
77 
78  template<typename T>
80  // Do nothing if p.tan() == 0 to preserve signs of y() and x()
81  if (p.tan() != 0) {
82  real x = _x * p._x - _y * p._y;
83  _y = _y * p._x + _x * p._y;
84  _x = x;
85  }
86  return *this;
87  }
88 
89  template<typename T>
91  : tol_( sqrt(numeric_limits<real>::epsilon()) )
92 #if AUXLATITUDE_UNOPT
93  , bmin_( numeric_limits<real>::min() )
94  , bmax_( numeric_limits<real>::max() )
95 #else
96  , bmin_( log2(numeric_limits<real>::min()) )
97  , bmax_( log2(numeric_limits<real>::max()) )
98 #endif
99  , _f( f )
100  , _fm1( 1 - _f )
101  , _e2( _f * (2 - _f) )
102  , _e2m1( _fm1 * _fm1 )
103  , _e12( _e2/(1 - _e2) )
104  , _e12p1( 1 / _e2m1 )
105  , _n( _f/(2 - _f) )
106  , _e( sqrt(fabs(_e2)) )
107  , _e1( sqrt(fabs(_e12)) )
108  , _n2( _n * _n )
109 #if AUXLATITUDE_UNOPT
110  , _q( _e12p1 + (_f == 0 ? 1 : (_f > 0 ? asinh(_e1) : atan(_e)) / _e) )
111 #else
112  , _q( _e12p1 + (_f == 0 ? 1 : (_f > 0 ? atanh(_e) : atan(_e)) / _e) )
113 #endif
114  {
115  if (!(isfinite(_f) && _f < 1 &&
116  isfinite(_n) && fabs(_n) < 1))
117  throw GeographicErr("Bad ellipsoid parameters");
118  fill(_c, _c + Lmax * AUXNUMBER * AUXNUMBER,
119  numeric_limits<real>::quiet_NaN());
120  }
121 
122  template<typename T>
124  : tol_( sqrt(numeric_limits<real>::epsilon()) )
125 #if AUXLATITUDE_UNOPT
126  , bmin_( numeric_limits<real>::min() )
127  , bmax_( numeric_limits<real>::max() )
128 #else
129  , bmin_( log2(numeric_limits<real>::min()) )
130  , bmax_( log2(numeric_limits<real>::max()) )
131 #endif
132  , _f( (a - b) / a )
133  , _fm1( b / a )
134  , _e2( ((a - b) * (a + b)) / (a * a) )
135  , _e2m1( (b * b) / (a * a) )
136  , _e12( ((a - b) * (a + b)) / (b * b) )
137  , _e12p1( (a * a) / (b * b) )
138  , _n( (a - b) / (a + b) )
139  , _e( sqrt(fabs(a - b) * (a + b)) / a )
140  , _e1( sqrt(fabs(a - b) * (a + b)) / b )
141  , _n2( _n * _n )
142 #if AUXLATITUDE_UNOPT
143  , _q( _e12p1 + (_f == 0 ? 1 : (_f > 0 ? asinh(_e1) : atan(_e)) / _e) )
144 #else
145  , _q( _e12p1 + (_f == 0 ? 1 : (_f > 0 ? atanh(_e) : atan(_e)) / _e) )
146 #endif
147  {
148  if (!(isfinite(_f) && _f < 1 &&
149  isfinite(_n) && fabs(_n) < 1))
150  throw GeographicErr("Bad ellipsoid parameters");
151  fill(_c, _c + Lmax * AUXNUMBER * AUXNUMBER,
152  numeric_limits<real>::quiet_NaN());
153  }
154 
155  template<typename T>
156  AuxAngle<T> AuxLatitude<T>::Parametric(const angle& phi, real* diff) const {
157  if (diff) *diff = _fm1;
158  return angle(phi.y() * _fm1, phi.x());
159  }
160 
161  template<typename T>
162  AuxAngle<T> AuxLatitude<T>::Geocentric(const angle& phi, real* diff) const {
163  if (diff) *diff = _e2m1;
164  return angle(phi.y() * _e2m1, phi.x());
165  }
166 
167  template<typename T>
168  AuxAngle<T> AuxLatitude<T>::Rectifying(const angle& phi, real* diff) const {
169  angle beta(Parametric(phi).normalized());
170  real sbeta = fabs(beta.y()), cbeta = fabs(beta.x());
171  real a = 1, b = _fm1, ka = _e2, kb = -_e12, ka1 = _e2m1, kb1 = _e12p1,
172  smu, cmu, mr;
173  if (_f < 0) {
174  swap(a, b); swap(ka, kb); swap(ka1, kb1); swap(sbeta, cbeta);
175  }
176  // now a,b = larger/smaller semiaxis
177  // beta now measured from larger semiaxis
178  // kb,ka = modulus-squared for distance from beta = 0,pi/2
179  // NB kb <= 0; 0 <= ka <= 1
180  // sa = b*E(beta,sqrt(kb)), sb = a*E(beta',sqrt(ka))
181  // 1 - ka * (1 - sb2) = 1 -ka + ka*sb2
182  real
183  sb2 = sbeta * sbeta,
184  cb2 = cbeta * cbeta,
185  db2 = 1 - kb * sb2,
186  da2 = ka1 + ka * sb2,
187  // DLMF Eq. 19.25.9
188  sa = b * sbeta * ( RF(cb2, db2, 1) - kb * sb2 * RD(cb2, db2, 1) / 3 ),
189  // DLMF Eq. 19.25.10 with complementary angles
190  sb = a * cbeta * ( ka1 * RF(sb2, da2, 1)
191  + ka * ka1 * cb2 * RD(sb2, 1, da2) / 3
192  + ka * sbeta / sqrt(da2) );
193  // sa + sb = 2*EllipticFunction::RG(a*a, b*b) = a*E(e) = b*E(i*e')
194  // mr = a*E(e)*(2/pi) = b*E(i*e')*(2/pi)
195  mr = (2 * (sa + sb)) / Math::pi<real>();
196  smu = sin(sa / mr);
197  cmu = sin(sb / mr);
198  if (_f < 0) { swap(smu, cmu); swap(a, b); }
199  // mu is normalized
200  angle mu(angle(smu, cmu).copyquadrant(phi));
201  if (diff) {
202  real cphi = phi.normalized().x(), tphi = phi.tan();
203  if (isfinite(tphi) || isnan(tphi)) {
204  cmu = mu.x(); cbeta = beta.x();
205  *diff = _fm1 * b/mr * Math::sq(cbeta / cmu) * (cbeta / cphi);
206  } else
207  *diff = _fm1 * mr/a;
208  }
209  return mu;
210  }
211 
212  template<typename T>
213  AuxAngle<T> AuxLatitude<T>::Conformal(const angle& phi, real* diff) const {
214  real tphi = fabs(phi.tan()), tchi = tphi;
215  if ( !( !isfinite(tphi) || tphi == 0 || _f == 0 ) ) {
216  real scphi = sc(tphi),
217  sig = sinh(_e2 * atanhee(tphi) ),
218  scsig = sc(sig);
219 #if AUXLATITUDE_UNOPT
220  tchi = tphi * scsig - sig * scphi;
221 #else
222  if (_f <= 0) {
223  tchi = tphi * scsig - sig * scphi;
224  } else {
225  // The general expression for tchi is
226  // tphi * scsig - sig * scphi
227  // This involves cancellation if f > 0, so change to
228  // (tphi - sig) * (tphi + sig) / (tphi * scsig + sig * scphi)
229  // To control overflow, write as (sigtphi = sig / tphi)
230  // (tphi - sig) * (1 + sigtphi) / (scsig + sigtphi * scphi)
231  real sigtphi = sig / tphi, tphimsig;
232  if (sig < tphi / 2)
233  tphimsig = tphi - sig;
234  else {
235  // Still have possibly dangerous cancellation in tphi - sig.
236  //
237  // Write tphi - sig = (1 - e) * Dg(1, e)
238  // Dg(x, y) = (g(x) - g(y)) / (x - y)
239  // g(x) = sinh(x * atanh(sphi * x))
240  // Note sinh(atanh(sphi)) = tphi
241  // Turn the crank on divided differences, substitute
242  // sphi = tphi/sc(tphi)
243  // atanh(x) = asinh(x/sqrt(1-x^2))
244  real
245  em1 = _e2m1 / (1 + _e), // 1 - e
246  atanhs = asinh(tphi), // atanh(sphi)
247  scbeta = sc(_fm1 * tphi), // sec(beta)
248  scphibeta = sc(tphi) / scbeta, // sec(phi)/sec(beta)
249  atanhes = asinh(_e * tphi / scbeta), // atanh(e * sphi)
250  t1 = (atanhs - _e * atanhes)/2,
251  t2 = asinh(em1 * (tphi * scphibeta)) / em1,
252  Dg = cosh((atanhs + _e * atanhes)/2) * (sinh(t1) / t1)
253  * ((atanhs + atanhes)/2 + (1 + _e)/2 * t2);
254  tphimsig = em1 * Dg; // tphi - sig
255  }
256  tchi = tphimsig * (1 + sigtphi) / (scsig + sigtphi * scphi);
257  }
258 #endif
259  }
260  angle chi(angle(tchi).copyquadrant(phi));
261  if (diff) {
262  if (isfinite(tphi) || isnan(tphi)) {
263  real cchi = chi.normalized().x(),
264  cphi = phi.normalized().x(),
265  cbeta = Parametric(phi).normalized().x();
266  *diff = _e2m1 * (cbeta / cchi) * (cbeta / cphi);
267  } else {
268 #if AUXLATITUDE_UNOPT
269  real ss = _f > 0 ? sinh(_e * asinh(_e1)) : sinh(-_e * atan(_e));
270  *diff = _f > 0 ? 1/( sc(ss) + ss ) : sc(ss) - ss;
271 #else
272  real ss = _f > 0 ? sinh(_e * atanh(_e)) : sinh(-_e * atan(_e));
273  *diff = sc(ss) - ss;
274 #endif
275  }
276  }
277  return chi;
278  }
279 
280  template<typename T>
281  AuxAngle<T> AuxLatitude<T>::Authalic(const angle& phi, real* diff) const {
282  real tphi = fabs(phi.tan());
283  angle xi(phi), phin(phi.normalized());
284  if ( !( !isfinite(tphi) || tphi == 0 || _f == 0 ) ) {
285  real qv = q(tphi),
286  Dqp = Dq(tphi),
287  Dqm = (_q + qv) / (1 + fabs(phin.y())); // Dq(-tphi)
288  xi = angle( copysign(qv, phi.y()), phin.x() * sqrt(Dqp * Dqm) );
289  }
290  if (diff) {
291  if (isfinite(tphi) || isnan(tphi)) {
292  real
293  cbeta = Parametric(phi).normalized().x(),
294  cxi = xi.normalized().x();
295  *diff =
296  (2/_q) * Math::sq(cbeta / cxi) * (cbeta / cxi) * (cbeta / phin.x());
297  } else
298  *diff = _e2m1 * sqrt(_q/2);
299  }
300  return xi;
301  }
302 
303  template<typename T>
305  real* diff) const {
306  switch (auxout) {
307  case GEOGRAPHIC: if (diff) *diff = 1; return phi; break;
308  case PARAMETRIC: return Parametric(phi, diff); break;
309  case GEOCENTRIC: return Geocentric(phi, diff); break;
310  case RECTIFYING: return Rectifying(phi, diff); break;
311  case CONFORMAL : return Conformal (phi, diff); break;
312  case AUTHALIC : return Authalic (phi, diff); break;
313  default:
314  if (diff) *diff = numeric_limits<real>::quiet_NaN();
315  return angle::NaN();
316  break;
317  }
318  }
319 
320  template<typename T>
322  int* niter) const {
323  int n = 0; if (niter) *niter = n;
324  real tphi = _fm1;
325  switch (auxin) {
326  case GEOGRAPHIC: return zeta; break;
327  // case PARAMETRIC: break;
328  case PARAMETRIC: return angle(zeta.y() / _fm1, zeta.x()); break;
329  // case GEOCENTRIC: tphi *= _fm1 ; break;
330  case GEOCENTRIC: return angle(zeta.y() / _e2m1, zeta.x()); break;
331  case RECTIFYING: tphi *= sqrt(_fm1); break;
332  case CONFORMAL : tphi *= _fm1 ; break;
333  case AUTHALIC : tphi *= cbrt(_fm1); break;
334  default: return angle::NaN(); break;
335  }
336 
337  // Drop through to solution by Newton's method
338 #if AUXLATITUDE_UNOPT
339  real tzeta = fabs(zeta.tan());
340  if (!isfinite(tzeta)) return zeta;
341  tphi = tzeta / tphi;
342  real bmin = fmin(tphi, bmin_), bmax = fmax(tphi, bmax_);
343  for (int sign = 0, osign = 0, ntrip = 0; n < numit_;) {
344  ++n;
345  real diff;
346  angle zeta1(ToAuxiliary(auxin, angle(tphi), &diff));
347  real tzeta1 = zeta1.tan();
348  osign = sign;
349  if (tzeta1 == tzeta)
350  break;
351  else if (tzeta1 > tzeta) {
352  sign = 1;
353  bmax = tphi;
354  } else {
355  sign = -1;
356  bmin = tphi;
357  }
358  real dtphi = -(tzeta1 - tzeta) / diff;
359  tphi += dtphi;
360  if (!(fabs(dtphi) >= tol_ * tphi))
361  break;
362  // Skip the Newton guards if UNOPT
363  if (false && ((sign * osign < 0 && n - ntrip > 2) ||
364  tphi >= bmax || tphi <= bmin)) {
365  sign = 0; ntrip = n;
366  tphi = sqrt(bmin * bmax);
367  }
368  }
369 #else
370  real tzeta = fabs(zeta.tan()), ltzeta = log2(tzeta);
371  if (!isfinite(ltzeta)) return zeta;
372  tphi = tzeta / tphi;
373  real ltphi = log2(tphi),
374  bmin = fmin(ltphi, bmin_), bmax = fmax(ltphi, bmax_);
375  for (int sign = 0, osign = 0, ntrip = 0; n < numit_;) {
376  ++n;
377  real diff;
378  angle zeta1(ToAuxiliary(auxin, angle(tphi), &diff));
379  real tzeta1 = zeta1.tan(), ltzeta1 = log2(tzeta1);
380  // Convert derivative from dtan(zeta)/dtan(phi) to
381  // dlog(tan(zeta))/dlog(tan(phi))
382  diff *= tphi/tzeta1;
383  osign = sign;
384  if (tzeta1 == tzeta)
385  break;
386  else if (tzeta1 > tzeta) {
387  sign = 1;
388  bmax = ltphi;
389  } else {
390  sign = -1;
391  bmin = ltphi;
392  }
393  real dltphi = -(ltzeta1 - ltzeta) / diff;
394  ltphi += dltphi;
395  tphi = exp2(ltphi);
396  if (!(fabs(dltphi) >= tol_)) {
397  ++n;
398  // Final Newton iteration without the logs
399  zeta1 = ToAuxiliary(auxin, angle(tphi), &diff);
400  tphi -= (zeta1.tan() - tzeta) / diff;
401  break;
402  }
403  if ((sign * osign < 0 && n - ntrip > 2) ||
404  ltphi >= bmax || ltphi <= bmin) {
405  sign = 0; ntrip = n;
406  ltphi = (bmin + bmax) / 2;
407  tphi = exp2(ltphi);
408  }
409  }
410 #endif
411  if (niter) *niter = n;
412  return angle(tphi).copyquadrant(zeta);
413  }
414 
415  template<typename T>
416  AuxAngle<T> AuxLatitude<T>::Convert(int auxin, int auxout, const angle& zeta,
417  bool series) const {
418  int k = ind(auxout, auxin);
419  if (k < 0) return angle::NaN();
420  if (auxin == auxout) return zeta;
421  if (series) {
422  if ( isnan(_c[Lmax * k]) ) fillcoeff(auxin, auxout, k);
423  angle zetan(zeta.normalized());
424  real d = Clenshaw(zetan.y(), zetan.x(), _c + Lmax * k, Lmax);
425  zetan += angle::radians(d);
426  return zetan;
427  } else {
428  if (auxin < 3 && auxout < 3)
429  return angle(zeta.y() * pow(_fm1, auxout - auxin), zeta.x());
430  else
431  return ToAuxiliary(auxout, FromAuxiliary(auxin, zeta));
432  }
433  }
434 
435  template<typename T>
437 #if AUXLATITUDE_NATIVE_ELLIPTIC
438  // Carlson, eqs 2.28 - 2.34
439  static const real
440  tolRD = pow(real(0.2) * (numeric_limits<real>::epsilon() * real(0.01)),
441  1/real(8));
442  real
443  A0 = (x + y + 3*z)/5,
444  An = A0,
445  Q = fmax(fmax(fabs(A0-x), fabs(A0-y)), fabs(A0-z)) / tolRD,
446  x0 = x,
447  y0 = y,
448  z0 = z,
449  mul = 1,
450  s = 0;
451  while (Q >= mul * fabs(An)) {
452  // Max 7 trips
453  real lam = sqrt(x0)*sqrt(y0) + sqrt(y0)*sqrt(z0) + sqrt(z0)*sqrt(x0);
454  s += 1/(mul * sqrt(z0) * (z0 + lam));
455  An = (An + lam)/4;
456  x0 = (x0 + lam)/4;
457  y0 = (y0 + lam)/4;
458  z0 = (z0 + lam)/4;
459  mul *= 4;
460  }
461  real
462  X = (A0 - x) / (mul * An),
463  Y = (A0 - y) / (mul * An),
464  Z = -(X + Y) / 3,
465  E2 = X*Y - 6*Z*Z,
466  E3 = (3*X*Y - 8*Z*Z)*Z,
467  E4 = 3 * (X*Y - Z*Z) * Z*Z,
468  E5 = X*Y*Z*Z*Z;
469  // https://dlmf.nist.gov/19.36.E2
470  // Polynomial is
471  // (1 - 3*E2/14 + E3/6 + 9*E2^2/88 - 3*E4/22 - 9*E2*E3/52 + 3*E5/26
472  // - E2^3/16 + 3*E3^2/40 + 3*E2*E4/20 + 45*E2^2*E3/272
473  // - 9*(E3*E4+E2*E5)/68)
474  return ((471240 - 540540 * E2) * E5 +
475  (612612 * E2 - 540540 * E3 - 556920) * E4 +
476  E3 * (306306 * E3 + E2 * (675675 * E2 - 706860) + 680680) +
477  E2 * ((417690 - 255255 * E2) * E2 - 875160) + 4084080) /
478  (4084080 * mul * An * sqrt(An)) + 3 * s;
479 #else
481  Math::real(z)));
482 #endif
483  }
484 
485  template<typename T>
486  typename AuxLatitude<T>::real AuxLatitude<T>::RF(real x, real y, real z) {
487 #if AUXLATITUDE_NATIVE_ELLIPTIC
488  // Carlson, eqs 2.2 - 2.7
489  static const real tolRF =
490  pow(3 * numeric_limits<real>::epsilon() * real(0.01), 1/real(8));
491  real
492  A0 = (x + y + z)/3,
493  An = A0,
494  Q = fmax(fmax(fabs(A0-x), fabs(A0-y)), fabs(A0-z)) / tolRF,
495  x0 = x,
496  y0 = y,
497  z0 = z,
498  mul = 1;
499  while (Q >= mul * fabs(An)) {
500  // Max 6 trips
501  real lam = sqrt(x0)*sqrt(y0) + sqrt(y0)*sqrt(z0) + sqrt(z0)*sqrt(x0);
502  An = (An + lam)/4;
503  x0 = (x0 + lam)/4;
504  y0 = (y0 + lam)/4;
505  z0 = (z0 + lam)/4;
506  mul *= 4;
507  }
508  real
509  X = (A0 - x) / (mul * An),
510  Y = (A0 - y) / (mul * An),
511  Z = - (X + Y),
512  E2 = X*Y - Z*Z,
513  E3 = X*Y*Z;
514  // https://dlmf.nist.gov/19.36.E1
515  // Polynomial is
516  // (1 - E2/10 + E3/14 + E2^2/24 - 3*E2*E3/44
517  // - 5*E2^3/208 + 3*E3^2/104 + E2^2*E3/16)
518  // convert to Horner form...
519  return (E3 * (6930 * E3 + E2 * (15015 * E2 - 16380) + 17160) +
520  E2 * ((10010 - 5775 * E2) * E2 - 24024) + 240240) /
521  (240240 * sqrt(An));
522 #else
524  Math::real(z)));
525 #endif
526  }
527 
528  template<typename T>
529  typename AuxLatitude<T>::real AuxLatitude<T>::atanhee(real tphi) const {
530 #if AUXLATITUDE_UNOPT
531  real sphi = sn(tphi);
532  return _f == 0 ? sphi :
533  (_f < 0 ? atan( _e * sphi ) : atanh( _e * sphi )) / _e;
534 #else
535  real s = _f <= 0 ? sn(tphi) : sn(_fm1 * tphi);
536  return _f == 0 ? s :
537  // atanh(e * sphi) = asinh(e' * sbeta)
538  (_f < 0 ? atan( _e * s ) : asinh( _e1 * s )) / _e;
539 #endif
540  }
541 
542  template<typename T>
543  typename AuxLatitude<T>::real AuxLatitude<T>::q(real tphi) const {
544 #if AUXLATITUDE_UNOPT
545  real sphi = sn(tphi);
546  return atanhee(tphi) + sphi / (1 - _e2 * sphi*sphi);
547 #else
548  real scbeta = sc(_fm1 * tphi);
549  return atanhee(tphi) + (tphi / scbeta) * (sc(tphi) / scbeta);
550 #endif
551  }
552 
553  template<typename T>
554  typename AuxLatitude<T>::real AuxLatitude<T>::Dq(real tphi) const {
555 #if AUXLATITUDE_UNOPT
556  real sphi = sn(tphi), d = 1 - sphi;
557  if (tphi <= 0)
558  // This branch is not reached; this case is open-coded in Authalic.
559  return (_q - q(tphi)) / d;
560  else if (d == 0)
561  return 2 / Math::sq(_e2m1);
562  else
563  return (_f == 0 ? 1 / (1 - _e2 * sphi) :
564  (_f < 0 ? atan(_e * d / (1 - _e2 * sphi)) / (_e * d) :
565  atanh(_e * d / (1 - _e2 * sphi)) / (_e * d) )) +
566  (1 + _e2 * sphi) / ((1 - _e2 * sphi*sphi) * _e2m1);
567 #else
568  real scphi = sc(tphi), sphi = sn(tphi),
569  // d = (1 - sphi) can underflow to zero for large tphi
570  d = tphi > 0 ? 1 / (scphi * scphi * (1 + sphi)) : 1 - sphi;
571  if (tphi <= 0)
572  // This branch is not reached; this case is open-coded in Authalic.
573  return (_q - q(tphi)) / d;
574  else if (d == 0)
575  return 2 / Math::sq(_e2m1);
576  else {
577  // General expression for Dq(1, sphi) is
578  // atanh(e * d / (1 - e2 * sphi)) / (e * d) +
579  // (1 + e2 * sphi) / ((1 - e2 * sphi * sphi) * e2m1);
580  // atanh( e * d / (1 - e2 * sphi))
581  // = atanh( e * d * scphi/(scphi - e2 * tphi))
582  // =
583  real scbeta = sc(_fm1 * tphi);
584  return (_f == 0 ? 1 :
585  (_f > 0 ? asinh(_e1 * d * scphi / scbeta) :
586  atan(_e * d / (1 - _e2 * sphi))) / (_e * d) ) +
587  (_f > 0 ?
588  ((scphi + _e2 * tphi) / (_e2m1 * scbeta)) * (scphi / scbeta) :
589  (1 + _e2 * sphi) / ((1 - _e2 * sphi*sphi) * _e2m1) );
590  }
591 #endif
592  }
593 
594  template<typename T>
595  void AuxLatitude<T>::fillcoeff(int auxin, int auxout, int k) const {
596 #if GEOGRAPHICLIB_AUXLATITUDE_ORDER == 4
597  static const real coeffs[] = {
598  // C[phi,phi] skipped
599  // C[phi,beta]; even coeffs only
600  0, 1,
601  0, 1/real(2),
602  1/real(3),
603  1/real(4),
604  // C[phi,theta]; even coeffs only
605  -2, 2,
606  -4, 2,
607  8/real(3),
608  4,
609  // C[phi,mu]; even coeffs only
610  -27/real(32), 3/real(2),
611  -55/real(32), 21/real(16),
612  151/real(96),
613  1097/real(512),
614  // C[phi,chi]
615  116/real(45), -2, -2/real(3), 2,
616  -227/real(45), -8/real(5), 7/real(3),
617  -136/real(35), 56/real(15),
618  4279/real(630),
619  // C[phi,xi]
620  -2582/real(14175), -16/real(35), 4/real(45), 4/real(3),
621  -11966/real(14175), 152/real(945), 46/real(45),
622  3802/real(14175), 3044/real(2835),
623  6059/real(4725),
624  // C[beta,phi]; even coeffs only
625  0, -1,
626  0, 1/real(2),
627  -1/real(3),
628  1/real(4),
629  // C[beta,beta] skipped
630  // C[beta,theta]; even coeffs only
631  0, 1,
632  0, 1/real(2),
633  1/real(3),
634  1/real(4),
635  // C[beta,mu]; even coeffs only
636  -9/real(32), 1/real(2),
637  -37/real(96), 5/real(16),
638  29/real(96),
639  539/real(1536),
640  // C[beta,chi]
641  38/real(45), -1/real(3), -2/real(3), 1,
642  -7/real(9), -14/real(15), 5/real(6),
643  -34/real(21), 16/real(15),
644  2069/real(1260),
645  // C[beta,xi]
646  -1082/real(14175), -46/real(315), 4/real(45), 1/real(3),
647  -338/real(2025), 68/real(945), 17/real(90),
648  1102/real(14175), 461/real(2835),
649  3161/real(18900),
650  // C[theta,phi]; even coeffs only
651  2, -2,
652  -4, 2,
653  -8/real(3),
654  4,
655  // C[theta,beta]; even coeffs only
656  0, -1,
657  0, 1/real(2),
658  -1/real(3),
659  1/real(4),
660  // C[theta,theta] skipped
661  // C[theta,mu]; even coeffs only
662  -23/real(32), -1/real(2),
663  -5/real(96), 5/real(16),
664  1/real(32),
665  283/real(1536),
666  // C[theta,chi]
667  4/real(9), -2/real(3), -2/real(3), 0,
668  -23/real(45), -4/real(15), 1/real(3),
669  -24/real(35), 2/real(5),
670  83/real(126),
671  // C[theta,xi]
672  -2102/real(14175), -158/real(315), 4/real(45), -2/real(3),
673  934/real(14175), -16/real(945), 16/real(45),
674  922/real(14175), -232/real(2835),
675  719/real(4725),
676  // C[mu,phi]; even coeffs only
677  9/real(16), -3/real(2),
678  -15/real(32), 15/real(16),
679  -35/real(48),
680  315/real(512),
681  // C[mu,beta]; even coeffs only
682  3/real(16), -1/real(2),
683  1/real(32), -1/real(16),
684  -1/real(48),
685  -5/real(512),
686  // C[mu,theta]; even coeffs only
687  13/real(16), 1/real(2),
688  33/real(32), -1/real(16),
689  -5/real(16),
690  -261/real(512),
691  // C[mu,mu] skipped
692  // C[mu,chi]
693  41/real(180), 5/real(16), -2/real(3), 1/real(2),
694  557/real(1440), -3/real(5), 13/real(48),
695  -103/real(140), 61/real(240),
696  49561/real(161280),
697  // C[mu,xi]
698  -1609/real(28350), 121/real(1680), 4/real(45), -1/real(6),
699  16463/real(453600), 26/real(945), -29/real(720),
700  449/real(28350), -1003/real(45360),
701  -40457/real(2419200),
702  // C[chi,phi]
703  -82/real(45), 4/real(3), 2/real(3), -2,
704  -13/real(9), -16/real(15), 5/real(3),
705  34/real(21), -26/real(15),
706  1237/real(630),
707  // C[chi,beta]
708  -16/real(45), 0, 2/real(3), -1,
709  19/real(45), -2/real(5), 1/real(6),
710  16/real(105), -1/real(15),
711  17/real(1260),
712  // C[chi,theta]
713  -2/real(9), 2/real(3), 2/real(3), 0,
714  43/real(45), 4/real(15), -1/real(3),
715  2/real(105), -2/real(5),
716  -55/real(126),
717  // C[chi,mu]
718  1/real(360), -37/real(96), 2/real(3), -1/real(2),
719  437/real(1440), -1/real(15), -1/real(48),
720  37/real(840), -17/real(480),
721  -4397/real(161280),
722  // C[chi,chi] skipped
723  // C[chi,xi]
724  -2312/real(14175), -88/real(315), 34/real(45), -2/real(3),
725  6079/real(14175), -184/real(945), 1/real(45),
726  772/real(14175), -106/real(2835),
727  -167/real(9450),
728  // C[xi,phi]
729  538/real(4725), 88/real(315), -4/real(45), -4/real(3),
730  -2482/real(14175), 8/real(105), 34/real(45),
731  -898/real(14175), -1532/real(2835),
732  6007/real(14175),
733  // C[xi,beta]
734  34/real(675), 32/real(315), -4/real(45), -1/real(3),
735  74/real(2025), -4/real(315), -7/real(90),
736  2/real(14175), -83/real(2835),
737  -797/real(56700),
738  // C[xi,theta]
739  778/real(4725), 62/real(105), -4/real(45), 2/real(3),
740  12338/real(14175), -32/real(315), 4/real(45),
741  -1618/real(14175), -524/real(2835),
742  -5933/real(14175),
743  // C[xi,mu]
744  1297/real(18900), -817/real(10080), -4/real(45), 1/real(6),
745  -29609/real(453600), -2/real(35), 49/real(720),
746  -2917/real(56700), 4463/real(90720),
747  331799/real(7257600),
748  // C[xi,chi]
749  2458/real(4725), 46/real(315), -34/real(45), 2/real(3),
750  3413/real(14175), -256/real(315), 19/real(45),
751  -15958/real(14175), 248/real(567),
752  16049/real(28350),
753  // C[xi,xi] skipped
754  };
755  static const int ptrs[] = {
756  0, 0, 6, 12, 18, 28, 38, 44, 44, 50, 56, 66, 76, 82, 88, 88, 94, 104,
757  114, 120, 126, 132, 132, 142, 152, 162, 172, 182, 192, 192, 202, 212,
758  222, 232, 242, 252, 252,
759  };
760 #elif GEOGRAPHICLIB_AUXLATITUDE_ORDER == 6
761  static const real coeffs[] = {
762  // C[phi,phi] skipped
763  // C[phi,beta]; even coeffs only
764  0, 0, 1,
765  0, 0, 1/real(2),
766  0, 1/real(3),
767  0, 1/real(4),
768  1/real(5),
769  1/real(6),
770  // C[phi,theta]; even coeffs only
771  2, -2, 2,
772  6, -4, 2,
773  -8, 8/real(3),
774  -16, 4,
775  32/real(5),
776  32/real(3),
777  // C[phi,mu]; even coeffs only
778  269/real(512), -27/real(32), 3/real(2),
779  6759/real(4096), -55/real(32), 21/real(16),
780  -417/real(128), 151/real(96),
781  -15543/real(2560), 1097/real(512),
782  8011/real(2560),
783  293393/real(61440),
784  // C[phi,chi]
785  -2854/real(675), 26/real(45), 116/real(45), -2, -2/real(3), 2,
786  2323/real(945), 2704/real(315), -227/real(45), -8/real(5), 7/real(3),
787  73814/real(2835), -1262/real(105), -136/real(35), 56/real(15),
788  -399572/real(14175), -332/real(35), 4279/real(630),
789  -144838/real(6237), 4174/real(315),
790  601676/real(22275),
791  // C[phi,xi]
792  28112932/real(212837625), 60136/real(467775), -2582/real(14175),
793  -16/real(35), 4/real(45), 4/real(3),
794  251310128/real(638512875), -21016/real(51975), -11966/real(14175),
795  152/real(945), 46/real(45),
796  -8797648/real(10945935), -94388/real(66825), 3802/real(14175),
797  3044/real(2835),
798  -1472637812/real(638512875), 41072/real(93555), 6059/real(4725),
799  455935736/real(638512875), 768272/real(467775),
800  4210684958LL/real(1915538625),
801  // C[beta,phi]; even coeffs only
802  0, 0, -1,
803  0, 0, 1/real(2),
804  0, -1/real(3),
805  0, 1/real(4),
806  -1/real(5),
807  1/real(6),
808  // C[beta,beta] skipped
809  // C[beta,theta]; even coeffs only
810  0, 0, 1,
811  0, 0, 1/real(2),
812  0, 1/real(3),
813  0, 1/real(4),
814  1/real(5),
815  1/real(6),
816  // C[beta,mu]; even coeffs only
817  205/real(1536), -9/real(32), 1/real(2),
818  1335/real(4096), -37/real(96), 5/real(16),
819  -75/real(128), 29/real(96),
820  -2391/real(2560), 539/real(1536),
821  3467/real(7680),
822  38081/real(61440),
823  // C[beta,chi]
824  -3118/real(4725), -1/real(3), 38/real(45), -1/real(3), -2/real(3), 1,
825  -247/real(270), 50/real(21), -7/real(9), -14/real(15), 5/real(6),
826  17564/real(2835), -5/real(3), -34/real(21), 16/real(15),
827  -49877/real(14175), -28/real(9), 2069/real(1260),
828  -28244/real(4455), 883/real(315),
829  797222/real(155925),
830  // C[beta,xi]
831  7947332/real(212837625), 11824/real(467775), -1082/real(14175),
832  -46/real(315), 4/real(45), 1/real(3),
833  39946703/real(638512875), -16672/real(155925), -338/real(2025),
834  68/real(945), 17/real(90),
835  -255454/real(1563705), -101069/real(467775), 1102/real(14175),
836  461/real(2835),
837  -189032762/real(638512875), 1786/real(18711), 3161/real(18900),
838  80274086/real(638512875), 88868/real(467775),
839  880980241/real(3831077250LL),
840  // C[theta,phi]; even coeffs only
841  -2, 2, -2,
842  6, -4, 2,
843  8, -8/real(3),
844  -16, 4,
845  -32/real(5),
846  32/real(3),
847  // C[theta,beta]; even coeffs only
848  0, 0, -1,
849  0, 0, 1/real(2),
850  0, -1/real(3),
851  0, 1/real(4),
852  -1/real(5),
853  1/real(6),
854  // C[theta,theta] skipped
855  // C[theta,mu]; even coeffs only
856  499/real(1536), -23/real(32), -1/real(2),
857  6565/real(12288), -5/real(96), 5/real(16),
858  -77/real(128), 1/real(32),
859  -4037/real(7680), 283/real(1536),
860  1301/real(7680),
861  17089/real(61440),
862  // C[theta,chi]
863  -3658/real(4725), 2/real(9), 4/real(9), -2/real(3), -2/real(3), 0,
864  61/real(135), 68/real(45), -23/real(45), -4/real(15), 1/real(3),
865  9446/real(2835), -46/real(35), -24/real(35), 2/real(5),
866  -34712/real(14175), -80/real(63), 83/real(126),
867  -2362/real(891), 52/real(45),
868  335882/real(155925),
869  // C[theta,xi]
870  216932/real(2627625), 109042/real(467775), -2102/real(14175),
871  -158/real(315), 4/real(45), -2/real(3),
872  117952358/real(638512875), -7256/real(155925), 934/real(14175),
873  -16/real(945), 16/real(45),
874  -7391576/real(54729675), -25286/real(66825), 922/real(14175),
875  -232/real(2835),
876  -67048172/real(638512875), 268/real(18711), 719/real(4725),
877  46774256/real(638512875), 14354/real(467775),
878  253129538/real(1915538625),
879  // C[mu,phi]; even coeffs only
880  -3/real(32), 9/real(16), -3/real(2),
881  135/real(2048), -15/real(32), 15/real(16),
882  105/real(256), -35/real(48),
883  -189/real(512), 315/real(512),
884  -693/real(1280),
885  1001/real(2048),
886  // C[mu,beta]; even coeffs only
887  -1/real(32), 3/real(16), -1/real(2),
888  -9/real(2048), 1/real(32), -1/real(16),
889  3/real(256), -1/real(48),
890  3/real(512), -5/real(512),
891  -7/real(1280),
892  -7/real(2048),
893  // C[mu,theta]; even coeffs only
894  -15/real(32), 13/real(16), 1/real(2),
895  -1673/real(2048), 33/real(32), -1/real(16),
896  349/real(256), -5/real(16),
897  963/real(512), -261/real(512),
898  -921/real(1280),
899  -6037/real(6144),
900  // C[mu,mu] skipped
901  // C[mu,chi]
902  7891/real(37800), -127/real(288), 41/real(180), 5/real(16), -2/real(3),
903  1/real(2),
904  -1983433/real(1935360), 281/real(630), 557/real(1440), -3/real(5),
905  13/real(48),
906  167603/real(181440), 15061/real(26880), -103/real(140), 61/real(240),
907  6601661/real(7257600), -179/real(168), 49561/real(161280),
908  -3418889/real(1995840), 34729/real(80640),
909  212378941/real(319334400),
910  // C[mu,xi]
911  12674323/real(851350500), -384229/real(14968800), -1609/real(28350),
912  121/real(1680), 4/real(45), -1/real(6),
913  -31621753811LL/real(1307674368000LL), -431/real(17325),
914  16463/real(453600), 26/real(945), -29/real(720),
915  -32844781/real(1751349600), 3746047/real(119750400), 449/real(28350),
916  -1003/real(45360),
917  10650637121LL/real(326918592000LL), 629/real(53460),
918  -40457/real(2419200),
919  205072597/real(20432412000LL), -1800439/real(119750400),
920  -59109051671LL/real(3923023104000LL),
921  // C[chi,phi]
922  4642/real(4725), 32/real(45), -82/real(45), 4/real(3), 2/real(3), -2,
923  -1522/real(945), 904/real(315), -13/real(9), -16/real(15), 5/real(3),
924  -12686/real(2835), 8/real(5), 34/real(21), -26/real(15),
925  -24832/real(14175), -12/real(5), 1237/real(630),
926  109598/real(31185), -734/real(315),
927  444337/real(155925),
928  // C[chi,beta]
929  -998/real(4725), 2/real(5), -16/real(45), 0, 2/real(3), -1,
930  -2/real(27), -22/real(105), 19/real(45), -2/real(5), 1/real(6),
931  116/real(567), -22/real(105), 16/real(105), -1/real(15),
932  2123/real(14175), -8/real(105), 17/real(1260),
933  128/real(4455), -1/real(105),
934  149/real(311850),
935  // C[chi,theta]
936  1042/real(4725), -14/real(45), -2/real(9), 2/real(3), 2/real(3), 0,
937  -712/real(945), -4/real(45), 43/real(45), 4/real(15), -1/real(3),
938  274/real(2835), 124/real(105), 2/real(105), -2/real(5),
939  21068/real(14175), -16/real(105), -55/real(126),
940  -9202/real(31185), -22/real(45),
941  -90263/real(155925),
942  // C[chi,mu]
943  -96199/real(604800), 81/real(512), 1/real(360), -37/real(96), 2/real(3),
944  -1/real(2),
945  1118711/real(3870720), -46/real(105), 437/real(1440), -1/real(15),
946  -1/real(48),
947  -5569/real(90720), 209/real(4480), 37/real(840), -17/real(480),
948  830251/real(7257600), 11/real(504), -4397/real(161280),
949  108847/real(3991680), -4583/real(161280),
950  -20648693/real(638668800),
951  // C[chi,chi] skipped
952  // C[chi,xi]
953  -55271278/real(212837625), 27128/real(93555), -2312/real(14175),
954  -88/real(315), 34/real(45), -2/real(3),
955  106691108/real(638512875), -65864/real(155925), 6079/real(14175),
956  -184/real(945), 1/real(45),
957  5921152/real(54729675), -14246/real(467775), 772/real(14175),
958  -106/real(2835),
959  75594328/real(638512875), -5312/real(467775), -167/real(9450),
960  2837636/real(638512875), -248/real(13365),
961  -34761247/real(1915538625),
962  // C[xi,phi]
963  -44732/real(2837835), 20824/real(467775), 538/real(4725), 88/real(315),
964  -4/real(45), -4/real(3),
965  -12467764/real(212837625), -37192/real(467775), -2482/real(14175),
966  8/real(105), 34/real(45),
967  100320856/real(1915538625), 54968/real(467775), -898/real(14175),
968  -1532/real(2835),
969  -5884124/real(70945875), 24496/real(467775), 6007/real(14175),
970  -839792/real(19348875), -23356/real(66825),
971  570284222/real(1915538625),
972  // C[xi,beta]
973  -70496/real(8513505), 2476/real(467775), 34/real(675), 32/real(315),
974  -4/real(45), -1/real(3),
975  53836/real(212837625), 3992/real(467775), 74/real(2025), -4/real(315),
976  -7/real(90),
977  -661844/real(1915538625), 7052/real(467775), 2/real(14175),
978  -83/real(2835),
979  1425778/real(212837625), 934/real(467775), -797/real(56700),
980  390088/real(212837625), -3673/real(467775),
981  -18623681/real(3831077250LL),
982  // C[xi,theta]
983  -4286228/real(42567525), -193082/real(467775), 778/real(4725),
984  62/real(105), -4/real(45), 2/real(3),
985  -61623938/real(70945875), 92696/real(467775), 12338/real(14175),
986  -32/real(315), 4/real(45),
987  427003576/real(1915538625), 612536/real(467775), -1618/real(14175),
988  -524/real(2835),
989  427770788/real(212837625), -8324/real(66825), -5933/real(14175),
990  -9153184/real(70945875), -320044/real(467775),
991  -1978771378/real(1915538625),
992  // C[xi,mu]
993  -9292991/real(302702400), 7764059/real(239500800), 1297/real(18900),
994  -817/real(10080), -4/real(45), 1/real(6),
995  36019108271LL/real(871782912000LL), 35474/real(467775),
996  -29609/real(453600), -2/real(35), 49/real(720),
997  3026004511LL/real(30648618000LL), -4306823/real(59875200),
998  -2917/real(56700), 4463/real(90720),
999  -368661577/real(4036032000LL), -102293/real(1871100),
1000  331799/real(7257600),
1001  -875457073/real(13621608000LL), 11744233/real(239500800),
1002  453002260127LL/real(7846046208000LL),
1003  // C[xi,chi]
1004  2706758/real(42567525), -55222/real(93555), 2458/real(4725),
1005  46/real(315), -34/real(45), 2/real(3),
1006  -340492279/real(212837625), 516944/real(467775), 3413/real(14175),
1007  -256/real(315), 19/real(45),
1008  4430783356LL/real(1915538625), 206834/real(467775), -15958/real(14175),
1009  248/real(567),
1010  62016436/real(70945875), -832976/real(467775), 16049/real(28350),
1011  -651151712/real(212837625), 15602/real(18711),
1012  2561772812LL/real(1915538625),
1013  // C[xi,xi] skipped
1014  };
1015  static const int ptrs[] = {
1016  0, 0, 12, 24, 36, 57, 78, 90, 90, 102, 114, 135, 156, 168, 180, 180, 192,
1017  213, 234, 246, 258, 270, 270, 291, 312, 333, 354, 375, 396, 396, 417,
1018  438, 459, 480, 501, 522, 522,
1019  };
1020 #elif GEOGRAPHICLIB_AUXLATITUDE_ORDER == 8
1021  static const real coeffs[] = {
1022  // C[phi,phi] skipped
1023  // C[phi,beta]; even coeffs only
1024  0, 0, 0, 1,
1025  0, 0, 0, 1/real(2),
1026  0, 0, 1/real(3),
1027  0, 0, 1/real(4),
1028  0, 1/real(5),
1029  0, 1/real(6),
1030  1/real(7),
1031  1/real(8),
1032  // C[phi,theta]; even coeffs only
1033  -2, 2, -2, 2,
1034  -8, 6, -4, 2,
1035  16, -8, 8/real(3),
1036  40, -16, 4,
1037  -32, 32/real(5),
1038  -64, 32/real(3),
1039  128/real(7),
1040  32,
1041  // C[phi,mu]; even coeffs only
1042  -6607/real(24576), 269/real(512), -27/real(32), 3/real(2),
1043  -155113/real(122880), 6759/real(4096), -55/real(32), 21/real(16),
1044  87963/real(20480), -417/real(128), 151/real(96),
1045  2514467/real(245760), -15543/real(2560), 1097/real(512),
1046  -69119/real(6144), 8011/real(2560),
1047  -5962461/real(286720), 293393/real(61440),
1048  6459601/real(860160),
1049  332287993/real(27525120),
1050  // C[phi,chi]
1051  189416/real(99225), 16822/real(4725), -2854/real(675), 26/real(45),
1052  116/real(45), -2, -2/real(3), 2,
1053  141514/real(8505), -31256/real(1575), 2323/real(945), 2704/real(315),
1054  -227/real(45), -8/real(5), 7/real(3),
1055  -2363828/real(31185), 98738/real(14175), 73814/real(2835),
1056  -1262/real(105), -136/real(35), 56/real(15),
1057  14416399/real(935550), 11763988/real(155925), -399572/real(14175),
1058  -332/real(35), 4279/real(630),
1059  258316372/real(1216215), -2046082/real(31185), -144838/real(6237),
1060  4174/real(315),
1061  -2155215124LL/real(14189175), -115444544/real(2027025),
1062  601676/real(22275),
1063  -170079376/real(1216215), 38341552/real(675675),
1064  1383243703/real(11351340),
1065  // C[phi,xi]
1066  -1683291094/real(37574026875LL), 22947844/real(1915538625),
1067  28112932/real(212837625), 60136/real(467775), -2582/real(14175),
1068  -16/real(35), 4/real(45), 4/real(3),
1069  -14351220203LL/real(488462349375LL), 1228352/real(3007125),
1070  251310128/real(638512875), -21016/real(51975), -11966/real(14175),
1071  152/real(945), 46/real(45),
1072  505559334506LL/real(488462349375LL), 138128272/real(147349125),
1073  -8797648/real(10945935), -94388/real(66825), 3802/real(14175),
1074  3044/real(2835),
1075  973080708361LL/real(488462349375LL), -45079184/real(29469825),
1076  -1472637812/real(638512875), 41072/real(93555), 6059/real(4725),
1077  -1385645336626LL/real(488462349375LL), -550000184/real(147349125),
1078  455935736/real(638512875), 768272/real(467775),
1079  -2939205114427LL/real(488462349375LL), 443810768/real(383107725),
1080  4210684958LL/real(1915538625),
1081  101885255158LL/real(54273594375LL), 387227992/real(127702575),
1082  1392441148867LL/real(325641566250LL),
1083  // C[beta,phi]; even coeffs only
1084  0, 0, 0, -1,
1085  0, 0, 0, 1/real(2),
1086  0, 0, -1/real(3),
1087  0, 0, 1/real(4),
1088  0, -1/real(5),
1089  0, 1/real(6),
1090  -1/real(7),
1091  1/real(8),
1092  // C[beta,beta] skipped
1093  // C[beta,theta]; even coeffs only
1094  0, 0, 0, 1,
1095  0, 0, 0, 1/real(2),
1096  0, 0, 1/real(3),
1097  0, 0, 1/real(4),
1098  0, 1/real(5),
1099  0, 1/real(6),
1100  1/real(7),
1101  1/real(8),
1102  // C[beta,mu]; even coeffs only
1103  -4879/real(73728), 205/real(1536), -9/real(32), 1/real(2),
1104  -86171/real(368640), 1335/real(4096), -37/real(96), 5/real(16),
1105  2901/real(4096), -75/real(128), 29/real(96),
1106  1082857/real(737280), -2391/real(2560), 539/real(1536),
1107  -28223/real(18432), 3467/real(7680),
1108  -733437/real(286720), 38081/real(61440),
1109  459485/real(516096),
1110  109167851/real(82575360),
1111  // C[beta,chi]
1112  -25666/real(99225), 4769/real(4725), -3118/real(4725), -1/real(3),
1113  38/real(45), -1/real(3), -2/real(3), 1,
1114  193931/real(42525), -14404/real(4725), -247/real(270), 50/real(21),
1115  -7/real(9), -14/real(15), 5/real(6),
1116  -1709614/real(155925), -36521/real(14175), 17564/real(2835), -5/real(3),
1117  -34/real(21), 16/real(15),
1118  -637699/real(85050), 2454416/real(155925), -49877/real(14175),
1119  -28/real(9), 2069/real(1260),
1120  48124558/real(1216215), -20989/real(2835), -28244/real(4455),
1121  883/real(315),
1122  -16969807/real(1091475), -2471888/real(184275), 797222/real(155925),
1123  -1238578/real(42525), 2199332/real(225225),
1124  87600385/real(4540536),
1125  // C[beta,xi]
1126  -5946082372LL/real(488462349375LL), 9708931/real(1915538625),
1127  7947332/real(212837625), 11824/real(467775), -1082/real(14175),
1128  -46/real(315), 4/real(45), 1/real(3),
1129  190673521/real(69780335625LL), 164328266/real(1915538625),
1130  39946703/real(638512875), -16672/real(155925), -338/real(2025),
1131  68/real(945), 17/real(90),
1132  86402898356LL/real(488462349375LL), 236067184/real(1915538625),
1133  -255454/real(1563705), -101069/real(467775), 1102/real(14175),
1134  461/real(2835),
1135  110123070361LL/real(488462349375LL), -98401826/real(383107725),
1136  -189032762/real(638512875), 1786/real(18711), 3161/real(18900),
1137  -200020620676LL/real(488462349375LL), -802887278/real(1915538625),
1138  80274086/real(638512875), 88868/real(467775),
1139  -296107325077LL/real(488462349375LL), 66263486/real(383107725),
1140  880980241/real(3831077250LL),
1141  4433064236LL/real(18091198125LL), 37151038/real(127702575),
1142  495248998393LL/real(1302566265000LL),
1143  // C[theta,phi]; even coeffs only
1144  2, -2, 2, -2,
1145  -8, 6, -4, 2,
1146  -16, 8, -8/real(3),
1147  40, -16, 4,
1148  32, -32/real(5),
1149  -64, 32/real(3),
1150  -128/real(7),
1151  32,
1152  // C[theta,beta]; even coeffs only
1153  0, 0, 0, -1,
1154  0, 0, 0, 1/real(2),
1155  0, 0, -1/real(3),
1156  0, 0, 1/real(4),
1157  0, -1/real(5),
1158  0, 1/real(6),
1159  -1/real(7),
1160  1/real(8),
1161  // C[theta,theta] skipped
1162  // C[theta,mu]; even coeffs only
1163  -14321/real(73728), 499/real(1536), -23/real(32), -1/real(2),
1164  -201467/real(368640), 6565/real(12288), -5/real(96), 5/real(16),
1165  2939/real(4096), -77/real(128), 1/real(32),
1166  1155049/real(737280), -4037/real(7680), 283/real(1536),
1167  -19465/real(18432), 1301/real(7680),
1168  -442269/real(286720), 17089/real(61440),
1169  198115/real(516096),
1170  48689387/real(82575360),
1171  // C[theta,chi]
1172  64424/real(99225), 76/real(225), -3658/real(4725), 2/real(9), 4/real(9),
1173  -2/real(3), -2/real(3), 0,
1174  2146/real(1215), -2728/real(945), 61/real(135), 68/real(45),
1175  -23/real(45), -4/real(15), 1/real(3),
1176  -95948/real(10395), 428/real(945), 9446/real(2835), -46/real(35),
1177  -24/real(35), 2/real(5),
1178  29741/real(85050), 4472/real(525), -34712/real(14175), -80/real(63),
1179  83/real(126),
1180  280108/real(13365), -17432/real(3465), -2362/real(891), 52/real(45),
1181  -48965632/real(4729725), -548752/real(96525), 335882/real(155925),
1182  -197456/real(15795), 51368/real(12285),
1183  1461335/real(174636),
1184  // C[theta,xi]
1185  -230886326/real(6343666875LL), -189115382/real(1915538625),
1186  216932/real(2627625), 109042/real(467775), -2102/real(14175),
1187  -158/real(315), 4/real(45), -2/real(3),
1188  -11696145869LL/real(69780335625LL), 288456008/real(1915538625),
1189  117952358/real(638512875), -7256/real(155925), 934/real(14175),
1190  -16/real(945), 16/real(45),
1191  91546732346LL/real(488462349375LL), 478700902/real(1915538625),
1192  -7391576/real(54729675), -25286/real(66825), 922/real(14175),
1193  -232/real(2835),
1194  218929662961LL/real(488462349375LL), -67330724/real(383107725),
1195  -67048172/real(638512875), 268/real(18711), 719/real(4725),
1196  -129039188386LL/real(488462349375LL), -117954842/real(273648375),
1197  46774256/real(638512875), 14354/real(467775),
1198  -178084928947LL/real(488462349375LL), 2114368/real(34827975),
1199  253129538/real(1915538625),
1200  6489189398LL/real(54273594375LL), 13805944/real(127702575),
1201  59983985827LL/real(325641566250LL),
1202  // C[mu,phi]; even coeffs only
1203  57/real(2048), -3/real(32), 9/real(16), -3/real(2),
1204  -105/real(4096), 135/real(2048), -15/real(32), 15/real(16),
1205  -105/real(2048), 105/real(256), -35/real(48),
1206  693/real(16384), -189/real(512), 315/real(512),
1207  693/real(2048), -693/real(1280),
1208  -1287/real(4096), 1001/real(2048),
1209  -6435/real(14336),
1210  109395/real(262144),
1211  // C[mu,beta]; even coeffs only
1212  19/real(2048), -1/real(32), 3/real(16), -1/real(2),
1213  7/real(4096), -9/real(2048), 1/real(32), -1/real(16),
1214  -3/real(2048), 3/real(256), -1/real(48),
1215  -11/real(16384), 3/real(512), -5/real(512),
1216  7/real(2048), -7/real(1280),
1217  9/real(4096), -7/real(2048),
1218  -33/real(14336),
1219  -429/real(262144),
1220  // C[mu,theta]; even coeffs only
1221  509/real(2048), -15/real(32), 13/real(16), 1/real(2),
1222  2599/real(4096), -1673/real(2048), 33/real(32), -1/real(16),
1223  -2989/real(2048), 349/real(256), -5/real(16),
1224  -43531/real(16384), 963/real(512), -261/real(512),
1225  5545/real(2048), -921/real(1280),
1226  16617/real(4096), -6037/real(6144),
1227  -19279/real(14336),
1228  -490925/real(262144),
1229  // C[mu,mu] skipped
1230  // C[mu,chi]
1231  -18975107/real(50803200), 72161/real(387072), 7891/real(37800),
1232  -127/real(288), 41/real(180), 5/real(16), -2/real(3), 1/real(2),
1233  148003883/real(174182400), 13769/real(28800), -1983433/real(1935360),
1234  281/real(630), 557/real(1440), -3/real(5), 13/real(48),
1235  79682431/real(79833600), -67102379/real(29030400), 167603/real(181440),
1236  15061/real(26880), -103/real(140), 61/real(240),
1237  -40176129013LL/real(7664025600LL), 97445/real(49896),
1238  6601661/real(7257600), -179/real(168), 49561/real(161280),
1239  2605413599LL/real(622702080), 14644087/real(9123840),
1240  -3418889/real(1995840), 34729/real(80640),
1241  175214326799LL/real(58118860800LL), -30705481/real(10378368),
1242  212378941/real(319334400),
1243  -16759934899LL/real(3113510400LL), 1522256789/real(1383782400),
1244  1424729850961LL/real(743921418240LL),
1245  // C[mu,xi]
1246  -375027460897LL/real(125046361440000LL),
1247  7183403063LL/real(560431872000LL), 12674323/real(851350500),
1248  -384229/real(14968800), -1609/real(28350), 121/real(1680), 4/real(45),
1249  -1/real(6),
1250  30410873385097LL/real(2000741783040000LL),
1251  1117820213/real(122594472000LL), -31621753811LL/real(1307674368000LL),
1252  -431/real(17325), 16463/real(453600), 26/real(945), -29/real(720),
1253  151567502183LL/real(17863765920000LL),
1254  -116359346641LL/real(3923023104000LL), -32844781/real(1751349600),
1255  3746047/real(119750400), 449/real(28350), -1003/real(45360),
1256  -317251099510901LL/real(8002967132160000LL), -13060303/real(766215450),
1257  10650637121LL/real(326918592000LL), 629/real(53460),
1258  -40457/real(2419200),
1259  -2105440822861LL/real(125046361440000LL),
1260  146875240637LL/real(3923023104000LL), 205072597/real(20432412000LL),
1261  -1800439/real(119750400),
1262  91496147778023LL/real(2000741783040000LL), 228253559/real(24518894400LL),
1263  -59109051671LL/real(3923023104000LL),
1264  126430355893LL/real(13894040160000LL),
1265  -4255034947LL/real(261534873600LL),
1266  -791820407649841LL/real(42682491371520000LL),
1267  // C[chi,phi]
1268  1514/real(1323), -8384/real(4725), 4642/real(4725), 32/real(45),
1269  -82/real(45), 4/real(3), 2/real(3), -2,
1270  142607/real(42525), -2288/real(1575), -1522/real(945), 904/real(315),
1271  -13/real(9), -16/real(15), 5/real(3),
1272  120202/real(51975), 44644/real(14175), -12686/real(2835), 8/real(5),
1273  34/real(21), -26/real(15),
1274  -1097407/real(187110), 1077964/real(155925), -24832/real(14175),
1275  -12/real(5), 1237/real(630),
1276  -12870194/real(1216215), 1040/real(567), 109598/real(31185),
1277  -734/real(315),
1278  -126463/real(72765), -941912/real(184275), 444337/real(155925),
1279  3463678/real(467775), -2405834/real(675675),
1280  256663081/real(56756700),
1281  // C[chi,beta]
1282  1384/real(11025), -34/real(4725), -998/real(4725), 2/real(5),
1283  -16/real(45), 0, 2/real(3), -1,
1284  -12616/real(42525), 1268/real(4725), -2/real(27), -22/real(105),
1285  19/real(45), -2/real(5), 1/real(6),
1286  1724/real(51975), -1858/real(14175), 116/real(567), -22/real(105),
1287  16/real(105), -1/real(15),
1288  115249/real(935550), -26836/real(155925), 2123/real(14175), -8/real(105),
1289  17/real(1260),
1290  140836/real(1216215), -424/real(6237), 128/real(4455), -1/real(105),
1291  210152/real(4729725), -31232/real(2027025), 149/real(311850),
1292  30208/real(6081075), -499/real(225225),
1293  -68251/real(113513400),
1294  // C[chi,theta]
1295  -1738/real(11025), 18/real(175), 1042/real(4725), -14/real(45),
1296  -2/real(9), 2/real(3), 2/real(3), 0,
1297  23159/real(42525), 332/real(945), -712/real(945), -4/real(45),
1298  43/real(45), 4/real(15), -1/real(3),
1299  13102/real(31185), -1352/real(945), 274/real(2835), 124/real(105),
1300  2/real(105), -2/real(5),
1301  -2414843/real(935550), 1528/real(4725), 21068/real(14175), -16/real(105),
1302  -55/real(126),
1303  60334/real(93555), 20704/real(10395), -9202/real(31185), -22/real(45),
1304  40458083/real(14189175), -299444/real(675675), -90263/real(155925),
1305  -3818498/real(6081075), -8962/real(12285),
1306  -4259027/real(4365900),
1307  // C[chi,mu]
1308  -7944359/real(67737600), 5406467/real(38707200), -96199/real(604800),
1309  81/real(512), 1/real(360), -37/real(96), 2/real(3), -1/real(2),
1310  -24749483/real(348364800), -51841/real(1209600), 1118711/real(3870720),
1311  -46/real(105), 437/real(1440), -1/real(15), -1/real(48),
1312  6457463/real(17740800), -9261899/real(58060800), -5569/real(90720),
1313  209/real(4480), 37/real(840), -17/real(480),
1314  -324154477/real(7664025600LL), -466511/real(2494800),
1315  830251/real(7257600), 11/real(504), -4397/real(161280),
1316  -22894433/real(124540416), 8005831/real(63866880), 108847/real(3991680),
1317  -4583/real(161280),
1318  2204645983LL/real(12915302400LL), 16363163/real(518918400),
1319  -20648693/real(638668800),
1320  497323811/real(12454041600LL), -219941297/real(5535129600LL),
1321  -191773887257LL/real(3719607091200LL),
1322  // C[chi,chi] skipped
1323  // C[chi,xi]
1324  -17451293242LL/real(488462349375LL), 308365186/real(1915538625),
1325  -55271278/real(212837625), 27128/real(93555), -2312/real(14175),
1326  -88/real(315), 34/real(45), -2/real(3),
1327  -101520127208LL/real(488462349375LL), 149984636/real(1915538625),
1328  106691108/real(638512875), -65864/real(155925), 6079/real(14175),
1329  -184/real(945), 1/real(45),
1330  10010741462LL/real(37574026875LL), -99534832/real(383107725),
1331  5921152/real(54729675), -14246/real(467775), 772/real(14175),
1332  -106/real(2835),
1333  1615002539/real(75148053750LL), -35573728/real(273648375),
1334  75594328/real(638512875), -5312/real(467775), -167/real(9450),
1335  -3358119706LL/real(488462349375LL), 130601488/real(1915538625),
1336  2837636/real(638512875), -248/real(13365),
1337  46771947158LL/real(488462349375LL), -3196/real(3553875),
1338  -34761247/real(1915538625),
1339  -18696014/real(18091198125LL), -2530364/real(127702575),
1340  -14744861191LL/real(651283132500LL),
1341  // C[xi,phi]
1342  -88002076/real(13956067125LL), -86728/real(16372125),
1343  -44732/real(2837835), 20824/real(467775), 538/real(4725), 88/real(315),
1344  -4/real(45), -4/real(3),
1345  -2641983469LL/real(488462349375LL), -895712/real(147349125),
1346  -12467764/real(212837625), -37192/real(467775), -2482/real(14175),
1347  8/real(105), 34/real(45),
1348  8457703444LL/real(488462349375LL), 240616/real(4209975),
1349  100320856/real(1915538625), 54968/real(467775), -898/real(14175),
1350  -1532/real(2835),
1351  -4910552477LL/real(97692469875LL), -4832848/real(147349125),
1352  -5884124/real(70945875), 24496/real(467775), 6007/real(14175),
1353  9393713176LL/real(488462349375LL), 816824/real(13395375),
1354  -839792/real(19348875), -23356/real(66825),
1355  -4532926649LL/real(97692469875LL), 1980656/real(54729675),
1356  570284222/real(1915538625),
1357  -14848113968LL/real(488462349375LL), -496894276/real(1915538625),
1358  224557742191LL/real(976924698750LL),
1359  // C[xi,beta]
1360  29232878/real(97692469875LL), -18484/real(4343625), -70496/real(8513505),
1361  2476/real(467775), 34/real(675), 32/real(315), -4/real(45), -1/real(3),
1362  -324943819/real(488462349375LL), -4160804/real(1915538625),
1363  53836/real(212837625), 3992/real(467775), 74/real(2025), -4/real(315),
1364  -7/real(90),
1365  -168643106/real(488462349375LL), 237052/real(383107725),
1366  -661844/real(1915538625), 7052/real(467775), 2/real(14175),
1367  -83/real(2835),
1368  113042383/real(97692469875LL), -2915326/real(1915538625),
1369  1425778/real(212837625), 934/real(467775), -797/real(56700),
1370  -558526274/real(488462349375LL), 6064888/real(1915538625),
1371  390088/real(212837625), -3673/real(467775),
1372  155665021/real(97692469875LL), 41288/real(29469825),
1373  -18623681/real(3831077250LL),
1374  504234982/real(488462349375LL), -6205669/real(1915538625),
1375  -8913001661LL/real(3907698795000LL),
1376  // C[xi,theta]
1377  182466964/real(8881133625LL), 53702182/real(212837625),
1378  -4286228/real(42567525), -193082/real(467775), 778/real(4725),
1379  62/real(105), -4/real(45), 2/real(3),
1380  367082779691LL/real(488462349375LL), -32500616/real(273648375),
1381  -61623938/real(70945875), 92696/real(467775), 12338/real(14175),
1382  -32/real(315), 4/real(45),
1383  -42668482796LL/real(488462349375LL), -663111728/real(383107725),
1384  427003576/real(1915538625), 612536/real(467775), -1618/real(14175),
1385  -524/real(2835),
1386  -327791986997LL/real(97692469875LL), 421877252/real(1915538625),
1387  427770788/real(212837625), -8324/real(66825), -5933/real(14175),
1388  74612072536LL/real(488462349375LL), 6024982024LL/real(1915538625),
1389  -9153184/real(70945875), -320044/real(467775),
1390  489898512247LL/real(97692469875LL), -46140784/real(383107725),
1391  -1978771378/real(1915538625),
1392  -42056042768LL/real(488462349375LL), -2926201612LL/real(1915538625),
1393  -2209250801969LL/real(976924698750LL),
1394  // C[xi,mu]
1395  39534358147LL/real(2858202547200LL),
1396  -25359310709LL/real(1743565824000LL), -9292991/real(302702400),
1397  7764059/real(239500800), 1297/real(18900), -817/real(10080), -4/real(45),
1398  1/real(6),
1399  -13216941177599LL/real(571640509440000LL),
1400  -14814966289LL/real(245188944000LL), 36019108271LL/real(871782912000LL),
1401  35474/real(467775), -29609/real(453600), -2/real(35), 49/real(720),
1402  -27782109847927LL/real(250092722880000LL),
1403  99871724539LL/real(1569209241600LL), 3026004511LL/real(30648618000LL),
1404  -4306823/real(59875200), -2917/real(56700), 4463/real(90720),
1405  168979300892599LL/real(1600593426432000LL),
1406  2123926699/real(15324309000LL), -368661577/real(4036032000LL),
1407  -102293/real(1871100), 331799/real(7257600),
1408  1959350112697LL/real(9618950880000LL),
1409  -493031379277LL/real(3923023104000LL), -875457073/real(13621608000LL),
1410  11744233/real(239500800),
1411  -145659994071373LL/real(800296713216000LL),
1412  -793693009/real(9807557760LL), 453002260127LL/real(7846046208000LL),
1413  -53583096419057LL/real(500185445760000LL),
1414  103558761539LL/real(1426553856000LL),
1415  12272105438887727LL/real(128047474114560000LL),
1416  // C[xi,chi]
1417  -64724382148LL/real(97692469875LL), 16676974/real(30405375),
1418  2706758/real(42567525), -55222/real(93555), 2458/real(4725),
1419  46/real(315), -34/real(45), 2/real(3),
1420  85904355287LL/real(37574026875LL), 158999572/real(1915538625),
1421  -340492279/real(212837625), 516944/real(467775), 3413/real(14175),
1422  -256/real(315), 19/real(45),
1423  2986003168LL/real(37574026875LL), -7597644214LL/real(1915538625),
1424  4430783356LL/real(1915538625), 206834/real(467775), -15958/real(14175),
1425  248/real(567),
1426  -375566203/real(39037950), 851209552/real(174139875),
1427  62016436/real(70945875), -832976/real(467775), 16049/real(28350),
1428  5106181018156LL/real(488462349375LL), 3475643362LL/real(1915538625),
1429  -651151712/real(212837625), 15602/real(18711),
1430  34581190223LL/real(8881133625LL), -10656173804LL/real(1915538625),
1431  2561772812LL/real(1915538625),
1432  -5150169424688LL/real(488462349375LL), 873037408/real(383107725),
1433  7939103697617LL/real(1953849397500LL),
1434  // C[xi,xi] skipped
1435  };
1436  static const int ptrs[] = {
1437  0, 0, 20, 40, 60, 96, 132, 152, 152, 172, 192, 228, 264, 284, 304, 304,
1438  324, 360, 396, 416, 436, 456, 456, 492, 528, 564, 600, 636, 672, 672,
1439  708, 744, 780, 816, 852, 888, 888,
1440  };
1441 #else
1442 #error "Unsupported value for GEOGRAPHICLIB_AUXLATITUDE_ORDER"
1443 #endif
1444 
1445  static_assert(sizeof(ptrs) / sizeof(int) == AUXNUMBER*AUXNUMBER+1,
1446  "Mismatch in size of ptrs array");
1447  static_assert(sizeof(coeffs) / sizeof(real) ==
1448  (RECTIFYING+1)*RECTIFYING *
1449  (Lmax * (Lmax + 3) - 2*(Lmax/2))/4 // Even only arrays
1450  + (AUXNUMBER*(AUXNUMBER-1) - (RECTIFYING+1)*RECTIFYING) *
1451  (Lmax * (Lmax + 1))/2,
1452  "Mismatch in size of coeffs array");
1453 
1454  if (k < 0) return; // auxout or auxin out of range
1455  if (auxout == auxin)
1456  fill(_c + Lmax * k, _c + Lmax * (k + 1), 0);
1457  else {
1458  int o = ptrs[k];
1459  real d = _n;
1460  if (auxin <= RECTIFYING && auxout <= RECTIFYING) {
1461  for (int l = 0; l < Lmax; ++l) {
1462  int m = (Lmax - l - 1) / 2; // order of polynomial in n^2
1463  _c[Lmax * k + l] = d * Math::polyval(m, coeffs + o, _n2);
1464  o += m + 1;
1465  d *= _n;
1466  }
1467  } else {
1468  for (int l = 0; l < Lmax; ++l) {
1469  int m = (Lmax - l - 1); // order of polynomial in n
1470  _c[Lmax * k + l] = d * Math::polyval(m, coeffs + o, _n);
1471  o += m + 1;
1472  d *= _n;
1473  }
1474  }
1475  // assert (o == ptrs[AUXNUMBER * auxout + auxin + 1])
1476  }
1477  }
1478 
1479  template<typename T>
1480  typename AuxLatitude<T>::real
1481  AuxLatitude<T>::Clenshaw(real szeta, real czeta,
1482  const real c[], int K, bool alt) {
1483  // Evaluate
1484  // y = sum(c[k] * sin( (2*k+2) * zeta), i, 0, K-1)
1485  // Approx operation count = (K + 5) mult and (2 * K + 2) add
1486  // If alt, use the Reinsch optimizations for small szeta and small czeta.
1487  int k = K;
1488  real u0 = 0; // accumulator for sum
1489  if (alt && 2*fabs(szeta) < 1) {
1490  real X = -4 * szeta * szeta,
1491  d0 = 0; // other accumulator for sum
1492  for (; k > 0;) {
1493  d0 = X * u0 + d0 + c[--k];
1494  u0 = d0 + u0;
1495  }
1496  } else if (alt && 2*fabs(czeta) < 1) {
1497  real X = 4 * czeta * czeta,
1498  d0 = 0; // other accumulator for sum
1499  for (; k > 0;) {
1500  d0 = X * u0 - d0 + c[--k];
1501  u0 = d0 - u0;
1502  }
1503  } else {
1504  real X = 2 * (czeta - szeta) * (czeta + szeta), // 2 * cos(2*zeta)
1505  u1 = 0; // other accumulator for sum
1506  for (; k > 0;) {
1507  real t = X * u0 - u1 + c[--k];
1508  u1 = u0; u0 = t;
1509  }
1510  }
1511  return 2 * szeta * czeta * u0; // sin(2*zeta) * u0
1512  }
1513 
1514  template class AuxAngle<Math::real>;
1515  template class AuxLatitude<Math::real>;
1516 #if GEOGRAPHICLIB_PRECISION != 2
1517  template class AuxAngle<double>;
1518  template class AuxLatitude<double>;
1519 #endif
1520 
1521 } // namespace GeographicLib
Header for the GeographicLib::AuxLatitude and GeographicLib::AuxAngle classes.
Header for GeographicLib::EllipticFunction class.
GeographicLib::Math::real real
Definition: GeodSolve.cpp:31
An accurate representation of angles.
Definition: AuxLatitude.hpp:55
AuxAngle & operator+=(const AuxAngle &p)
Definition: AuxLatitude.cpp:79
AuxAngle copyquadrant(const AuxAngle &p) const
Definition: AuxLatitude.cpp:74
static AuxAngle NaN()
Definition: AuxLatitude.cpp:52
AuxAngle normalized() const
Definition: AuxLatitude.cpp:58
Conversions between auxiliary latitudes.
angle FromAuxiliary(int auxin, const angle &zeta, int *niter=nullptr) const
angle ToAuxiliary(int auxout, const angle &phi, real *diff=nullptr) const
angle Convert(int auxin, int auxout, const angle &zeta, bool series=false) const
static real RD(real x, real y, real z)
static real RF(real x, real y, real z)
Geocentric coordinates
Definition: Geocentric.hpp:67
Exception handling for GeographicLib.
Definition: Constants.hpp:316
static T sq(T x)
Definition: Math.hpp:212
static T polyval(int N, const T p[], T x)
Definition: Math.hpp:271
Namespace for GeographicLib.
Definition: Accumulator.cpp:12
void swap(GeographicLib::NearestNeighbor< dist_t, pos_t, distfun_t > &a, GeographicLib::NearestNeighbor< dist_t, pos_t, distfun_t > &b)