src/Geoid.cpp
Go to the documentation of this file.
1 
10 #include <GeographicLib/Geoid.hpp>
11 // For getenv
12 #include <cstdlib>
14 
15 #if !defined(GEOGRAPHICLIB_DATA)
16 # if defined(_WIN32)
17 # define GEOGRAPHICLIB_DATA "C:/ProgramData/GeographicLib"
18 # else
19 # define GEOGRAPHICLIB_DATA "/usr/local/share/GeographicLib"
20 # endif
21 #endif
22 
23 #if !defined(GEOGRAPHICLIB_GEOID_DEFAULT_NAME)
24 # define GEOGRAPHICLIB_GEOID_DEFAULT_NAME "egm96-5"
25 #endif
26 
27 #if defined(_MSC_VER)
28 // Squelch warnings about unsafe use of getenv
29 # pragma warning (disable: 4996)
30 #endif
31 
32 namespace GeographicLib {
33 
34  using namespace std;
35 
36  // This is the transfer matrix for a 3rd order fit with a 12-point stencil
37  // with weights
38  //
39  // \x -1 0 1 2
40  // y
41  // -1 . 1 1 .
42  // 0 1 2 2 1
43  // 1 1 2 2 1
44  // 2 . 1 1 .
45  //
46  // A algorithm for n-dimensional polynomial fits is described in
47  // F. H. Lesh,
48  // Multi-dimensional least-squares polynomial curve fitting,
49  // CACM 2, 29-30 (1959).
50  //
51  // Here's the Maxima code to generate this matrix:
52  //
53  // /* The stencil and the weights */
54  // xarr:[
55  // 0, 1,
56  // -1, 0, 1, 2,
57  // -1, 0, 1, 2,
58  // 0, 1]$
59  // yarr:[
60  // -1,-1,
61  // 0, 0, 0, 0,
62  // 1, 1, 1, 1,
63  // 2, 2]$
64  // warr:[
65  // 1, 1,
66  // 1, 2, 2, 1,
67  // 1, 2, 2, 1,
68  // 1, 1]$
69  //
70  // /* [x exponent, y exponent] for cubic fit */
71  // pows:[
72  // [0,0],
73  // [1,0],[0,1],
74  // [2,0],[1,1],[0,2],
75  // [3,0],[2,1],[1,2],[0,3]]$
76  //
77  // basisvec(x,y,pows):=map(lambda([ex],(if ex[1]=0 then 1 else x^ex[1])*
78  // (if ex[2]=0 then 1 else y^ex[2])),pows)$
79  // addterm(x,y,f,w,pows):=block([a,b,bb:basisvec(x,y,pows)],
80  // a:w*(transpose(bb).bb),
81  // b:(w*f) * bb,
82  // [a,b])$
83  //
84  // c3row(k):=block([a,b,c,pows:pows,n],
85  // n:length(pows),
86  // a:zeromatrix(n,n),
87  // b:copylist(part(a,1)),
88  // c:[a,b],
89  // for i:1 thru length(xarr) do
90  // c:c+addterm(xarr[i],yarr[i],if i=k then 1 else 0,warr[i],pows),
91  // a:c[1],b:c[2],
92  // part(transpose( a^^-1 . transpose(b)),1))$
93  // c3:[]$
94  // for k:1 thru length(warr) do c3:endcons(c3row(k),c3)$
95  // c3:apply(matrix,c3)$
96  // c0:part(ratsimp(
97  // genmatrix(yc,1,length(warr)).abs(c3).genmatrix(yd,length(pows),1)),2)$
98  // c3:c0*c3$
99 
100  const int Geoid::c0_ = 240; // Common denominator
101  const int Geoid::c3_[stencilsize_ * nterms_] = {
102  9, -18, -88, 0, 96, 90, 0, 0, -60, -20,
103  -9, 18, 8, 0, -96, 30, 0, 0, 60, -20,
104  9, -88, -18, 90, 96, 0, -20, -60, 0, 0,
105  186, -42, -42, -150, -96, -150, 60, 60, 60, 60,
106  54, 162, -78, 30, -24, -90, -60, 60, -60, 60,
107  -9, -32, 18, 30, 24, 0, 20, -60, 0, 0,
108  -9, 8, 18, 30, -96, 0, -20, 60, 0, 0,
109  54, -78, 162, -90, -24, 30, 60, -60, 60, -60,
110  -54, 78, 78, 90, 144, 90, -60, -60, -60, -60,
111  9, -8, -18, -30, -24, 0, 20, 60, 0, 0,
112  -9, 18, -32, 0, 24, 30, 0, 0, -60, 20,
113  9, -18, -8, 0, -24, -30, 0, 0, 60, 20,
114  };
115 
116  // Like c3, but with the coeffs of x, x^2, and x^3 constrained to be zero.
117  // Use this at the N pole so that the height in independent of the longitude
118  // there.
119  //
120  // Here's the Maxima code to generate this matrix (continued from above).
121  //
122  // /* figure which terms to exclude so that fit is indep of x at y=0 */
123  // mask:part(zeromatrix(1,length(pows)),1)+1$
124  // for i:1 thru length(pows) do
125  // if pows[i][1]>0 and pows[i][2]=0 then mask[i]:0$
126  //
127  // /* Same as c3row but with masked pows. */
128  // c3nrow(k):=block([a,b,c,powsa:[],n,d,e],
129  // for i:1 thru length(mask) do if mask[i]>0 then
130  // powsa:endcons(pows[i],powsa),
131  // n:length(powsa),
132  // a:zeromatrix(n,n),
133  // b:copylist(part(a,1)),
134  // c:[a,b],
135  // for i:1 thru length(xarr) do
136  // c:c+addterm(xarr[i],yarr[i],if i=k then 1 else 0,warr[i],powsa),
137  // a:c[1],b:c[2],
138  // d:part(transpose( a^^-1 . transpose(b)),1),
139  // e:[],
140  // for i:1 thru length(mask) do
141  // if mask[i]>0 then (e:endcons(first(d),e),d:rest(d)) else e:endcons(0,e),
142  // e)$
143  // c3n:[]$
144  // for k:1 thru length(warr) do c3n:endcons(c3nrow(k),c3n)$
145  // c3n:apply(matrix,c3n)$
146  // c0n:part(ratsimp(
147  // genmatrix(yc,1,length(warr)).abs(c3n).genmatrix(yd,length(pows),1)),2)$
148  // c3n:c0n*c3n$
149 
150  const int Geoid::c0n_ = 372; // Common denominator
151  const int Geoid::c3n_[stencilsize_ * nterms_] = {
152  0, 0, -131, 0, 138, 144, 0, 0, -102, -31,
153  0, 0, 7, 0, -138, 42, 0, 0, 102, -31,
154  62, 0, -31, 0, 0, -62, 0, 0, 0, 31,
155  124, 0, -62, 0, 0, -124, 0, 0, 0, 62,
156  124, 0, -62, 0, 0, -124, 0, 0, 0, 62,
157  62, 0, -31, 0, 0, -62, 0, 0, 0, 31,
158  0, 0, 45, 0, -183, -9, 0, 93, 18, 0,
159  0, 0, 216, 0, 33, 87, 0, -93, 12, -93,
160  0, 0, 156, 0, 153, 99, 0, -93, -12, -93,
161  0, 0, -45, 0, -3, 9, 0, 93, -18, 0,
162  0, 0, -55, 0, 48, 42, 0, 0, -84, 31,
163  0, 0, -7, 0, -48, -42, 0, 0, 84, 31,
164  };
165 
166  // Like c3n, but y -> 1-y so that h is independent of x at y = 1. Use this
167  // at the S pole so that the height in independent of the longitude there.
168  //
169  // Here's the Maxima code to generate this matrix (continued from above).
170  //
171  // /* Transform c3n to c3s by transforming y -> 1-y */
172  // vv:[
173  // v[11],v[12],
174  // v[7],v[8],v[9],v[10],
175  // v[3],v[4],v[5],v[6],
176  // v[1],v[2]]$
177  // poly:expand(vv.(c3n/c0n).transpose(basisvec(x,1-y,pows)))$
178  // c3sf[i,j]:=coeff(coeff(coeff(poly,v[i]),x,pows[j][1]),y,pows[j][2])$
179  // c3s:genmatrix(c3sf,length(vv),length(pows))$
180  // c0s:part(ratsimp(
181  // genmatrix(yc,1,length(warr)).abs(c3s).genmatrix(yd,length(pows),1)),2)$
182  // c3s:c0s*c3s$
183 
184  const int Geoid::c0s_ = 372; // Common denominator
185  const int Geoid::c3s_[stencilsize_ * nterms_] = {
186  18, -36, -122, 0, 120, 135, 0, 0, -84, -31,
187  -18, 36, -2, 0, -120, 51, 0, 0, 84, -31,
188  36, -165, -27, 93, 147, -9, 0, -93, 18, 0,
189  210, 45, -111, -93, -57, -192, 0, 93, 12, 93,
190  162, 141, -75, -93, -129, -180, 0, 93, -12, 93,
191  -36, -21, 27, 93, 39, 9, 0, -93, -18, 0,
192  0, 0, 62, 0, 0, 31, 0, 0, 0, -31,
193  0, 0, 124, 0, 0, 62, 0, 0, 0, -62,
194  0, 0, 124, 0, 0, 62, 0, 0, 0, -62,
195  0, 0, 62, 0, 0, 31, 0, 0, 0, -31,
196  -18, 36, -64, 0, 66, 51, 0, 0, -102, 31,
197  18, -36, 2, 0, -66, -51, 0, 0, 102, 31,
198  };
199 
200  Geoid::Geoid(const std::string& name, const std::string& path, bool cubic,
201  bool threadsafe)
202  : _name(name)
203  , _dir(path)
204  , _cubic(cubic)
205  , _a( Constants::WGS84_a() )
206  , _e2( (2 - Constants::WGS84_f()) * Constants::WGS84_f() )
207  , _degree( Math::degree() )
208  , _eps( sqrt(numeric_limits<real>::epsilon()) )
209  , _threadsafe(false) // Set after cache is read
210  {
211  GEOGRAPHICLIB_STATIC_ASSERT(sizeof(pixel_t) == pixel_size_,
212  "pixel_t has the wrong size");
213  if (_dir.empty())
215  _filename = _dir + "/" + _name + (pixel_size_ != 4 ? ".pgm" : ".pgm4");
216  _file.open(_filename.c_str(), ios::binary);
217  if (!(_file.good()))
218  throw GeographicErr("File not readable " + _filename);
219  string s;
220  if (!(getline(_file, s) && s == "P5"))
221  throw GeographicErr("File not in PGM format " + _filename);
223  _scale = 0;
224  _maxerror = _rmserror = -1;
225  _description = "NONE";
226  _datetime = "UNKNOWN";
227  while (getline(_file, s)) {
228  if (s.empty())
229  continue;
230  if (s[0] == '#') {
231  istringstream is(s);
232  string commentid, key;
233  if (!(is >> commentid >> key) || commentid != "#")
234  continue;
235  if (key == "Description" || key =="DateTime") {
236  string::size_type p =
237  s.find_first_not_of(" \t", unsigned(is.tellg()));
238  if (p != string::npos)
239  (key == "Description" ? _description : _datetime) = s.substr(p);
240  } else if (key == "Offset") {
241  if (!(is >> _offset))
242  throw GeographicErr("Error reading offset " + _filename);
243  } else if (key == "Scale") {
244  if (!(is >> _scale))
245  throw GeographicErr("Error reading scale " + _filename);
246  } else if (key == (_cubic ? "MaxCubicError" : "MaxBilinearError")) {
247  // It's not an error if the error can't be read
248  is >> _maxerror;
249  } else if (key == (_cubic ? "RMSCubicError" : "RMSBilinearError")) {
250  // It's not an error if the error can't be read
251  is >> _rmserror;
252  }
253  } else {
254  istringstream is(s);
255  if (!(is >> _width >> _height))
256  throw GeographicErr("Error reading raster size " + _filename);
257  break;
258  }
259  }
260  {
261  unsigned maxval;
262  if (!(_file >> maxval))
263  throw GeographicErr("Error reading maxval " + _filename);
264  if (maxval != pixel_max_)
265  throw GeographicErr("Incorrect value of maxval " + _filename);
266  // Add 1 for whitespace after maxval
267  _datastart = (unsigned long long)(_file.tellg()) + 1ULL;
268  _swidth = (unsigned long long)(_width);
269  }
271  throw GeographicErr("Offset not set " + _filename);
272  if (_scale == 0)
273  throw GeographicErr("Scale not set " + _filename);
274  if (_scale < 0)
275  throw GeographicErr("Scale must be positive " + _filename);
276  if (_height < 2 || _width < 2)
277  // Coarsest grid spacing is 180deg.
278  throw GeographicErr("Raster size too small " + _filename);
279  if (_width & 1)
280  // This is so that longitude grids can be extended thru the poles.
281  throw GeographicErr("Raster width is odd " + _filename);
282  if (!(_height & 1))
283  // This is so that latitude grid includes the equator.
284  throw GeographicErr("Raster height is even " + _filename);
285  _file.seekg(0, ios::end);
286  if (!_file.good() ||
287  _datastart + pixel_size_ * _swidth * (unsigned long long)(_height) !=
288  (unsigned long long)(_file.tellg()))
289  // Possibly this test should be "<" because the file contains, e.g., a
290  // second image. However, for now we are more strict.
291  throw GeographicErr("File has the wrong length " + _filename);
292  _rlonres = _width / real(360);
293  _rlatres = (_height - 1) / real(180);
294  _cache = false;
295  _ix = _width;
296  _iy = _height;
297  // Ensure that file errors throw exceptions
298  _file.exceptions(ifstream::eofbit | ifstream::failbit | ifstream::badbit);
299  if (threadsafe) {
300  CacheAll();
301  _file.close();
302  _threadsafe = true;
303  }
304  }
305 
307  lat = Math::LatFix(lat);
308  if (Math::isnan(lat) || Math::isnan(lon)) {
309  return Math::NaN();
310  }
311  lon = Math::AngNormalize(lon);
312  real
313  fx = lon * _rlonres,
314  fy = -lat * _rlatres;
315  int
316  ix = int(floor(fx)),
317  iy = min((_height - 1)/2 - 1, int(floor(fy)));
318  fx -= ix;
319  fy -= iy;
320  iy += (_height - 1)/2;
321  ix += ix < 0 ? _width : (ix >= _width ? -_width : 0);
322  real v00 = 0, v01 = 0, v10 = 0, v11 = 0;
323  real t[nterms_];
324 
325  if (_threadsafe || !(ix == _ix && iy == _iy)) {
326  if (!_cubic) {
327  v00 = rawval(ix , iy );
328  v01 = rawval(ix + 1, iy );
329  v10 = rawval(ix , iy + 1);
330  v11 = rawval(ix + 1, iy + 1);
331  } else {
333  int k = 0;
334  v[k++] = rawval(ix , iy - 1);
335  v[k++] = rawval(ix + 1, iy - 1);
336  v[k++] = rawval(ix - 1, iy );
337  v[k++] = rawval(ix , iy );
338  v[k++] = rawval(ix + 1, iy );
339  v[k++] = rawval(ix + 2, iy );
340  v[k++] = rawval(ix - 1, iy + 1);
341  v[k++] = rawval(ix , iy + 1);
342  v[k++] = rawval(ix + 1, iy + 1);
343  v[k++] = rawval(ix + 2, iy + 1);
344  v[k++] = rawval(ix , iy + 2);
345  v[k++] = rawval(ix + 1, iy + 2);
346 
347  const int* c3x = iy == 0 ? c3n_ : (iy == _height - 2 ? c3s_ : c3_);
348  int c0x = iy == 0 ? c0n_ : (iy == _height - 2 ? c0s_ : c0_);
349  for (unsigned i = 0; i < nterms_; ++i) {
350  t[i] = 0;
351  for (unsigned j = 0; j < stencilsize_; ++j)
352  t[i] += v[j] * c3x[nterms_ * j + i];
353  t[i] /= c0x;
354  }
355  }
356  } else { // same cell; used cached coefficients
357  if (!_cubic) {
358  v00 = _v00;
359  v01 = _v01;
360  v10 = _v10;
361  v11 = _v11;
362  } else
363  copy(_t, _t + nterms_, t);
364  }
365  if (!_cubic) {
366  real
367  a = (1 - fx) * v00 + fx * v01,
368  b = (1 - fx) * v10 + fx * v11,
369  c = (1 - fy) * a + fy * b,
370  h = _offset + _scale * c;
371  if (!_threadsafe) {
372  _ix = ix;
373  _iy = iy;
374  _v00 = v00;
375  _v01 = v01;
376  _v10 = v10;
377  _v11 = v11;
378  }
379  return h;
380  } else {
381  real h = t[0] + fx * (t[1] + fx * (t[3] + fx * t[6])) +
382  fy * (t[2] + fx * (t[4] + fx * t[7]) +
383  fy * (t[5] + fx * t[8] + fy * t[9]));
384  h = _offset + _scale * h;
385  if (!_threadsafe) {
386  _ix = ix;
387  _iy = iy;
388  copy(t, t + nterms_, _t);
389  }
390  return h;
391  }
392  }
393 
394  void Geoid::CacheClear() const {
395  if (!_threadsafe) {
396  _cache = false;
397  try {
398  _data.clear();
399  // Use swap to release memory back to system
400  vector< vector<pixel_t> >().swap(_data);
401  }
402  catch (const exception&) {
403  }
404  }
405  }
406 
407  void Geoid::CacheArea(real south, real west, real north, real east) const {
408  if (_threadsafe)
409  throw GeographicErr("Attempt to change cache of threadsafe Geoid");
410  if (south > north) {
411  CacheClear();
412  return;
413  }
414  south = Math::LatFix(south);
415  north = Math::LatFix(north);
416  west = Math::AngNormalize(west); // west in [-180, 180)
417  east = Math::AngNormalize(east);
418  if (east <= west)
419  east += 360; // east - west in (0, 360]
420  int
421  iw = int(floor(west * _rlonres)),
422  ie = int(floor(east * _rlonres)),
423  in = int(floor(-north * _rlatres)) + (_height - 1)/2,
424  is = int(floor(-south * _rlatres)) + (_height - 1)/2;
425  in = max(0, min(_height - 2, in));
426  is = max(0, min(_height - 2, is));
427  is += 1;
428  ie += 1;
429  if (_cubic) {
430  in -= 1;
431  is += 1;
432  iw -= 1;
433  ie += 1;
434  }
435  if (ie - iw >= _width - 1) {
436  // Include entire longitude range
437  iw = 0;
438  ie = _width - 1;
439  } else {
440  ie += iw < 0 ? _width : (iw >= _width ? -_width : 0);
441  iw += iw < 0 ? _width : (iw >= _width ? -_width : 0);
442  }
443  int oysize = int(_data.size());
444  _xsize = ie - iw + 1;
445  _ysize = is - in + 1;
446  _xoffset = iw;
447  _yoffset = in;
448 
449  try {
450  _data.resize(_ysize, vector<pixel_t>(_xsize));
451  for (int iy = min(oysize, _ysize); iy--;)
452  _data[iy].resize(_xsize);
453  }
454  catch (const bad_alloc&) {
455  CacheClear();
456  throw GeographicErr("Insufficient memory for caching " + _filename);
457  }
458 
459  try {
460  for (int iy = in; iy <= is; ++iy) {
461  int iy1 = iy, iw1 = iw;
462  if (iy < 0 || iy >= _height) {
463  // Allow points "beyond" the poles to support interpolation
464  iy1 = iy1 < 0 ? -iy1 : 2 * (_height - 1) - iy1;
465  iw1 += _width/2;
466  if (iw1 >= _width)
467  iw1 -= _width;
468  }
469  int xs1 = min(_width - iw1, _xsize);
470  filepos(iw1, iy1);
471  Utility::readarray<pixel_t, pixel_t, true>
472  (_file, &(_data[iy - in][0]), xs1);
473  if (xs1 < _xsize) {
474  // Wrap around longitude = 0
475  filepos(0, iy1);
476  Utility::readarray<pixel_t, pixel_t, true>
477  (_file, &(_data[iy - in][xs1]), _xsize - xs1);
478  }
479  }
480  _cache = true;
481  }
482  catch (const exception& e) {
483  CacheClear();
484  throw GeographicErr(string("Error filling cache ") + e.what());
485  }
486  }
487 
488  std::string Geoid::DefaultGeoidPath() {
489  string path;
490  char* geoidpath = getenv("GEOGRAPHICLIB_GEOID_PATH");
491  if (geoidpath)
492  path = string(geoidpath);
493  if (!path.empty())
494  return path;
495  char* datapath = getenv("GEOGRAPHICLIB_DATA");
496  if (datapath)
497  path = string(datapath);
498  return (!path.empty() ? path : string(GEOGRAPHICLIB_DATA)) + "/geoids";
499  }
500 
501  std::string Geoid::DefaultGeoidName() {
502  string name;
503  char* geoidname = getenv("GEOGRAPHICLIB_GEOID_NAME");
504  if (geoidname)
505  name = string(geoidname);
506  return !name.empty() ? name : string(GEOGRAPHICLIB_GEOID_DEFAULT_NAME);
507  }
508 
509 } // namespace GeographicLib
static T AngNormalize(T x)
Definition: Math.hpp:440
const gtsam::Symbol key('X', 0)
static T NaN()
Definition: Math.hpp:830
#define max(a, b)
Definition: datatypes.h:20
const bool _cubic
Definition: Geoid.hpp:104
#define GEOGRAPHICLIB_GEOID_DEFAULT_NAME
Definition: src/Geoid.cpp:24
std::string _description
Definition: Geoid.hpp:108
Scalar * b
Definition: benchVecAdd.cpp:17
real _t[nterms_]
Definition: Geoid.hpp:121
static const double lat
#define GEOGRAPHICLIB_DATA
Definition: src/Geoid.cpp:19
Header for GeographicLib::Utility class.
static const int c3n_[stencilsize_ *nterms_]
Definition: Geoid.hpp:100
#define min(a, b)
Definition: datatypes.h:19
static bool isnan(T x)
Definition: Math.hpp:853
static const int c0s_
Definition: Geoid.hpp:98
static const unsigned nterms_
Definition: Geoid.hpp:95
static T LatFix(T x)
Definition: Math.hpp:467
void CacheArea(real south, real west, real north, real east) const
Definition: src/Geoid.cpp:407
Scalar Scalar * c
Definition: benchVecAdd.cpp:17
static const int c3_[stencilsize_ *nterms_]
Definition: Geoid.hpp:99
Mathematical functions needed by GeographicLib.
Definition: Math.hpp:102
Definition: BFloat16.h:88
static std::string DefaultGeoidPath()
Definition: src/Geoid.cpp:488
Geoid(const Geoid &)
std::string _dir
Definition: Geoid.hpp:103
static double epsilon
Definition: testRot3.cpp:37
void filepos(int ix, int iy) const
Definition: Geoid.hpp:122
unsigned short pixel_t
Definition: Geoid.hpp:86
EIGEN_DEVICE_FUNC const FloorReturnType floor() const
const double degree
int EIGEN_BLAS_FUNC() swap(int *n, RealScalar *px, int *incx, RealScalar *py, int *incy)
Definition: level1_impl.h:130
const double fy
real height(real lat, real lon) const
Definition: src/Geoid.cpp:306
std::string _datetime
Definition: Geoid.hpp:108
real rawval(int ix, int iy) const
Definition: Geoid.hpp:131
static const int c0n_
Definition: Geoid.hpp:97
Array< int, Dynamic, 1 > v
static const int c0_
Definition: Geoid.hpp:96
Definition: main.h:100
Namespace for GeographicLib.
Array< double, 1, 3 > e(1./3., 0.5, 2.)
RealScalar s
static const unsigned pixel_max_
Definition: Geoid.hpp:88
std::vector< std::vector< pixel_t > > _data
Definition: Geoid.hpp:114
std::string _filename
Definition: Geoid.hpp:103
static const unsigned stencilsize_
Definition: Geoid.hpp:94
Constants needed by GeographicLib
Definition: Constants.hpp:131
const double h
unsigned long long _swidth
Definition: Geoid.hpp:111
Exception handling for GeographicLib.
Definition: Constants.hpp:389
static const int c3s_[stencilsize_ *nterms_]
Definition: Geoid.hpp:101
Math::real real
Definition: Geoid.hpp:84
static EIGEN_DEPRECATED const end_t end
static std::string DefaultGeoidName()
Definition: src/Geoid.cpp:501
float * p
std::ifstream _file
Definition: Geoid.hpp:106
unsigned long long _datastart
Definition: Geoid.hpp:111
static const double lon
static const unsigned pixel_size_
Definition: Geoid.hpp:87
Jet< T, N > sqrt(const Jet< T, N > &f)
Definition: jet.h:418
Annotation for function names.
Definition: attr.h:48
int EIGEN_BLAS_FUNC() copy(int *n, RealScalar *px, int *incx, RealScalar *py, int *incy)
Definition: level1_impl.h:29
const double fx
std::string _name
Definition: Geoid.hpp:103
void CacheClear() const
Definition: src/Geoid.cpp:394
Header for GeographicLib::Geoid class.
v resize(3)
void CacheAll() const
Definition: Geoid.hpp:262
std::ptrdiff_t j
Point2 t(10, 10)


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:34:18