OSSIM - Open Source Software Image Map  Version 1.9.0 (20180803)
Public Member Functions | Protected Attributes | List of all members
ossimGeodeticEvaluator Class Reference

#include <ossimGeodeticEvaluator.h>

Public Member Functions

 ossimGeodeticEvaluator (const ossimEllipsoid &ell=ossimEllipsoid())
 constructor. More...
 
 ~ossimGeodeticEvaluator ()
 virtual destructor. More...
 
bool direct (const ossimGpt &p1, const double &az1, const double &d, ossimGpt &p2, double &az2)
 Evaluate direct problem. More...
 
bool inverse (const ossimGpt &p1, const ossimGpt &p2, double &d, double &az12, double &az21)
 Evaluate Vincenty inverse problem. More...
 

Protected Attributes

double m_A
 
double m_B
 
double m_F
 
double m_A2
 
double m_B2
 
double m_Ecc2
 
double m_2ndEcc2
 

Detailed Description

Definition at line 17 of file ossimGeodeticEvaluator.h.

Constructor & Destructor Documentation

◆ ossimGeodeticEvaluator()

ossimGeodeticEvaluator::ossimGeodeticEvaluator ( const ossimEllipsoid ell = ossimEllipsoid())

constructor.

Defaults to WGS-84

Definition at line 29 of file ossimGeodeticEvaluator.cpp.

References m_2ndEcc2, m_A, m_A2, m_B, and m_B2.

30  :
31  m_A(ell.getA()),
32  m_B(ell.getB()),
33  m_F(ell.getFlattening()),
35 {
36  m_A2 = m_A*m_A;
37  m_B2 = m_B*m_B;
38  m_2ndEcc2 = (m_A2 - m_B2) / m_B2;
39 }
const double & getA() const
double eccentricitySquared() const
const double & getB() const
const double & getFlattening() const

◆ ~ossimGeodeticEvaluator()

ossimGeodeticEvaluator::~ossimGeodeticEvaluator ( )

virtual destructor.

Definition at line 48 of file ossimGeodeticEvaluator.cpp.

49 {
50  if (traceExec()) ossimNotify(ossimNotifyLevel_DEBUG)
51  << "DEBUG: ~ossimGeodeticEvaluator(): returning..." << std::endl;
52 }
OSSIMDLLEXPORT std::ostream & ossimNotify(ossimNotifyLevel level=ossimNotifyLevel_WARN)

Member Function Documentation

◆ direct()

bool ossimGeodeticEvaluator::direct ( const ossimGpt p1,
const double &  az1,
const double &  d,
ossimGpt p2,
double &  az2 
)

Evaluate direct problem.

Given: Point 1 position, azimuth & distance to point 2. Find: Point 2 position, azimuth from point 2 to point 1.

Parameters
p1Point 1.
az1Azimuth from point 1 to point 2.
dDistance between points 1 & 2.
p2Point 2.
az2Azimuth from point 2 to point 1.

Definition at line 65 of file ossimGeodeticEvaluator.cpp.

70 {
71  if (traceDebug())
72  {
74  << "\nossimGeodeticEvaluator::direct DEBUG:" << std::endl;
75  }
76 
77  bool computationOK = false;
78 
79  double phi1 = p1.latr();
80  // double lam1 = p1.lonr();
81 
82  double cosalpha1 = cos(az12);
83  double sinalpha1 = sin(az12);
84  double tanU1 = (1.0 - m_F) * tan(phi1);
85  double cosU1 = 1.0 / sqrt( 1.0 + tanU1 * tanU1 );
86  double sinU1 = tanU1 * cosU1;
87 
88  // eq. 1
89  double sigma1 = atan2(tanU1, cosalpha1);
90 
91  // eq. 2
92  double sinalpha = cosU1 * sinalpha1;
93 
94  double sin2alpha = sinalpha * sinalpha;
95  double cos2alpha = 1.0 - sin2alpha;
96 
97  double u2 = cos2alpha * m_2ndEcc2;
98 
99  // Eq. 3
100  double A = 1 + u2 / 16384 * (4096 + u2 * (-768 + u2 * (320 - 175*u2)));
101 
102  // Eq. 4
103  double B = u2 / 1024 * (256 + u2 * (-128 + u2 * (74 - 47*u2)));
104 
105  // iterate until there is a negligible change in sigma
106  double dSig;
107  double d_bA = d / (m_B * A);
108  double sigma = d_bA;
109  double sinsig;
110  double cossig;
111  double prevSigma = M_PI;
112  double sigmaM2;
113  double cos2sigm;
114  double cos2sigm2;
115 
116  while (abs(sigma - prevSigma) > 1.e-12)
117  {
118  // eq. 5
119  sigmaM2 = 2.0*sigma1 + sigma;
120  cos2sigm = cos(sigmaM2);
121  cos2sigm2 = cos2sigm * cos2sigm;
122  sinsig = sin(sigma);
123  cossig = cos(sigma);
124  double sin2sig = sinsig*sinsig;
125 
126  // Eq. 6
127  dSig = B * sinsig * (cos2sigm + B / 4 * (cossig * (-1 + 2*cos2sigm2) - B / 6 * cos2sigm * (-3 + 4*sin2sig) * (-3 + 4*cos2sigm2)));
128 
129  // eq. 7
130  prevSigma = sigma;
131  sigma = d_bA + dSig;
132  }
133 
134  sigmaM2 = 2.0*sigma1 + sigma;
135  cos2sigm = cos(sigmaM2);
136  cos2sigm2 = cos2sigm * cos2sigm;
137  cossig = cos(sigma);
138  sinsig = sin(sigma);
139 
140  // eq. 8
141  double tmp = sinU1*sinsig - cosU1*cossig*cosalpha1;
142  double phi2 = atan2( sinU1*cossig + cosU1*sinsig*cosalpha1,
143  (1.0-m_F) * sqrt( sin2alpha + tmp*tmp));
144 
145  // eq. 9
146  double lam = atan2(sinsig*sinalpha1, cosU1*cossig - sinU1*sinsig*cosalpha1);
147 
148  // eq. 10
149  double C = m_F/16 * cos2alpha * (4 + m_F * (4 - 3*cos2alpha));
150 
151  // eq. 11
152  double L = lam - (1-C) * m_F * sinalpha * (sigma + C * sinsig * (cos2sigm + C * cossig * (-1 + 2*cos2sigm2)));
153 
154  // eq. 12
155  double alpha2 = atan2(sinalpha, -tmp);
156 
157 
158  p2.latr(phi2);
159  p2.lonr(p1.lonr() + L);
160  az21 = alpha2 + M_PI;
161  if (az21 < 0.0)
162  az21 += 2.0*M_PI;
163  if (az21 > 2.0*M_PI)
164  az21 -= 2.0*M_PI;
165 
166  return computationOK;
167 }
#define abs(a)
Definition: auxiliary.h:74
#define A(r, c)
#define M_PI
double lonr() const
Returns the longitude in radian measure.
Definition: ossimGpt.h:76
double latr() const
latr().
Definition: ossimGpt.h:66
OSSIMDLLEXPORT std::ostream & ossimNotify(ossimNotifyLevel level=ossimNotifyLevel_WARN)

◆ inverse()

bool ossimGeodeticEvaluator::inverse ( const ossimGpt p1,
const ossimGpt p2,
double &  d,
double &  az12,
double &  az21 
)

Evaluate Vincenty inverse problem.

Given: Point 1 & point 2 positions. Find: Distance and azimuths between points.

Parameters
p1Point 1.
p2Point 2.
dDistance between points 1 & 2.
az12Azimuth from point 1 to point 2.
az21Azimuth from point 2 to point 1.

Definition at line 181 of file ossimGeodeticEvaluator.cpp.

186 {
187  if (traceDebug())
188  {
190  << "\nossimGeodeticEvaluator::inverse DEBUG:" << std::endl;
191  }
192 
193  bool computationOK = false;
194 
195  double phi1 = p1.latr();
196  double lam1 = p1.lonr();
197  double phi2 = p2.latr();
198  double lam2 = p2.lonr();
199 
200  double U1 = atan((1.0-m_F) * tan(phi1));
201  double sinU1 = sin(U1);
202  double cosU1 = cos(U1);
203 
204  double U2 = atan((1.0-m_F) * tan(phi2));
205  double sinU2 = sin(U2);
206  double cosU2 = cos(U2);
207 
208  double sinU1sinU2 = sinU1 * sinU2;
209  double cosU1sinU2 = cosU1 * sinU2;
210  double sinU1cosU2 = sinU1 * cosU2;
211  double cosU1cosU2 = cosU1 * cosU2;
212 
213  // Eq. 13
214  // 1st approximation
215  double L = lam2 - lam1;
216  double lam = L;
217 
218 
219  // Auxiliaries
220  double sig = 0.0;
221  double dSig = 0.0;
222  double lam0;
223  double cos2alpha;
224  double cos2sigm;
225  double cos2sigm2;
226  double sinsig;
227  double sin2sig;
228  double cossig;
229  bool converged = false;
230 
231 
232  int iter = 0;
233  int iterMax = 100;
234  int iterMin = 4;
235 
236  // Iterative loop
237  do
238  {
239  iter++;
240 
241  double sinlam = sin(lam);
242  double coslam = cos(lam);
243 
244  // Eq. 14
245  double tmp = cosU1sinU2 - sinU1cosU2 * coslam;
246  sin2sig = cosU2*sinlam * cosU2*sinlam + tmp * tmp;
247  sinsig = sqrt(sin2sig);
248 
249  // Eq. 15
250  cossig = sinU1sinU2 + (cosU1cosU2 * coslam);
251 
252  // Eq. 16
253  sig = atan2(sinsig, cossig);
254 
255  // Eq. 17
256  double sinalpha;
257  if (sin2sig == 0.0)
258  sinalpha = 0.0;
259  else
260  sinalpha = cosU1cosU2 * sinlam / sinsig;
261  cos2alpha = 1.0 - sinalpha*sinalpha;
262 
263  // Eq. 18
264  if (cos2alpha == 0.0)
265  cos2sigm = 0.0;
266  else
267  cos2sigm = cossig - 2 * sinU1sinU2 / cos2alpha;
268  cos2sigm2 = cos2sigm * cos2sigm;
269 
270  // Eq. 10
271  double C = m_F/16 * cos2alpha * (4 + m_F * (4 - 3*cos2alpha));
272 
273  // Eq. 11
274  lam0 = lam;
275  lam = L + (1-C) * m_F * sinalpha * (sig + C * sinsig * (cos2sigm + C * cossig * (-1 + 2*cos2sigm2)));
276 
277 
278  // Check delta lambda for convergence
279  double delta = abs(lam - lam0);
280  if (delta < 1e-13 && iter>=iterMin)
281  {
282  converged = true;
283  }
284 
285  } while (!converged && iter<=iterMax);
286 
287 
288  double u2 = cos2alpha * m_2ndEcc2;
289 
290  // Eq. 3
291  double A = 1 + u2 / 16384 * (4096 + u2 * (-768 + u2 * (320 - 175*u2)));
292 
293  // Eq. 4
294  double B = u2 / 1024 * (256 + u2 * (-128 + u2 * (74 - 47*u2)));
295 
296  // Eq. 6
297  dSig = B * sinsig * (cos2sigm + B / 4 * (cossig * (-1 + 2*cos2sigm2) - B / 6 * cos2sigm * (-3 + 4*sin2sig) * (-3 + 4*cos2sigm2)));
298 
299 
300  // Compute distance
301  // Eq. 19
302  d = m_B * A * (sig - dSig);
303 
304 
305  // Compute azimuths
306  // Diverged, same meridian or antipodal
307  if (!converged)
308  {
309  if (phi1 > phi2)
310  {
311  az12 = M_PI;
312  az21 = 0.0;
313  }
314  else if (phi1 < phi2)
315  {
316  az12 = 0.0;
317  az21 = M_PI;
318  }
319  else
320  {
321  az12 = ossim::nan();
322  az21 = ossim::nan();
323  }
324  }
325 
326  // Converged
327  else
328  {
329  // Eq. 20
330  az12 = atan2(cosU2 * sin(lam), (cosU1sinU2 - sinU1cosU2 * cos(lam)));
331  if (az12 < 0.0) az12 += 2.0*M_PI;
332 
333  // Eq. 21
334  az21 = atan2(cosU1 * sin(lam), (-sinU1cosU2 + cosU1sinU2 * cos(lam))) + M_PI;
335  if (az21 < 0.0) az21 += 2.0*M_PI;
336 
337  computationOK = true;
338  }
339 
340  return computationOK;
341  }
double nan()
Method to return ieee floating point double precision NAN.
Definition: ossimCommon.h:135
#define abs(a)
Definition: auxiliary.h:74
#define A(r, c)
#define M_PI
double lonr() const
Returns the longitude in radian measure.
Definition: ossimGpt.h:76
double latr() const
latr().
Definition: ossimGpt.h:66
OSSIMDLLEXPORT std::ostream & ossimNotify(ossimNotifyLevel level=ossimNotifyLevel_WARN)

Member Data Documentation

◆ m_2ndEcc2

double ossimGeodeticEvaluator::m_2ndEcc2
protected

Definition at line 69 of file ossimGeodeticEvaluator.h.

Referenced by ossimGeodeticEvaluator().

◆ m_A

double ossimGeodeticEvaluator::m_A
protected

Definition at line 63 of file ossimGeodeticEvaluator.h.

Referenced by ossimGeodeticEvaluator().

◆ m_A2

double ossimGeodeticEvaluator::m_A2
protected

Definition at line 66 of file ossimGeodeticEvaluator.h.

Referenced by ossimGeodeticEvaluator().

◆ m_B

double ossimGeodeticEvaluator::m_B
protected

Definition at line 64 of file ossimGeodeticEvaluator.h.

Referenced by ossimGeodeticEvaluator().

◆ m_B2

double ossimGeodeticEvaluator::m_B2
protected

Definition at line 67 of file ossimGeodeticEvaluator.h.

Referenced by ossimGeodeticEvaluator().

◆ m_Ecc2

double ossimGeodeticEvaluator::m_Ecc2
protected

Definition at line 68 of file ossimGeodeticEvaluator.h.

◆ m_F

double ossimGeodeticEvaluator::m_F
protected

Definition at line 65 of file ossimGeodeticEvaluator.h.


The documentation for this class was generated from the following files: