summaryrefslogtreecommitdiff
path: root/navit/tools/gpx2navit_txt/src/geod_inv.c
diff options
context:
space:
mode:
Diffstat (limited to 'navit/tools/gpx2navit_txt/src/geod_inv.c')
-rw-r--r--navit/tools/gpx2navit_txt/src/geod_inv.c73
1 files changed, 73 insertions, 0 deletions
diff --git a/navit/tools/gpx2navit_txt/src/geod_inv.c b/navit/tools/gpx2navit_txt/src/geod_inv.c
new file mode 100644
index 000000000..d938da8a3
--- /dev/null
+++ b/navit/tools/gpx2navit_txt/src/geod_inv.c
@@ -0,0 +1,73 @@
+#ifndef lint
+static const char SCCSID[] =
+ "@(#)geod_inv.c 4.5 95/09/23 GIE REL";
+#endif
+# include "projects.h"
+# include "geodesic.h"
+# define DTOL 1e-12
+void geod_inv(void)
+{
+ double th1,
+ th2,
+ thm,
+ dthm,
+ dlamm,
+ dlam,
+ sindlamm,
+ costhm,
+ sinthm,
+ cosdthm,
+ sindthm, L, E, cosd, d, X, Y, T, sind, tandlammp, u, v, D, A, B;
+
+ if (ellipse) {
+ th1 = atan(onef * tan(phi1));
+ th2 = atan(onef * tan(phi2));
+ } else {
+ th1 = phi1;
+ th2 = phi2;
+ }
+ thm = .5 * (th1 + th2);
+ dthm = .5 * (th2 - th1);
+ dlamm = .5 * (dlam = adjlon(lam2 - lam1));
+ if (fabs(dlam) < DTOL && fabs(dthm) < DTOL) {
+ al12 = al21 = geod_S = 0.;
+ return;
+ }
+ sindlamm = sin(dlamm);
+ costhm = cos(thm);
+ sinthm = sin(thm);
+ cosdthm = cos(dthm);
+ sindthm = sin(dthm);
+ L = sindthm * sindthm + (cosdthm * cosdthm - sinthm * sinthm)
+ * sindlamm * sindlamm;
+ d = acos(cosd = 1 - L - L);
+ if (ellipse) {
+ E = cosd + cosd;
+ sind = sin(d);
+ Y = sinthm * cosdthm;
+ Y *= (Y + Y) / (1. - L);
+ T = sindthm * costhm;
+ T *= (T + T) / L;
+ X = Y + T;
+ Y -= T;
+ T = d / sind;
+ D = 4. * T * T;
+ A = D * E;
+ B = D + D;
+ geod_S = geod_a * sind * (T - f4 * (T * X - Y) +
+ f64 * (X * (A + (T - .5 * (A - E)) * X) -
+ Y * (B + E * Y) + D * X * Y));
+ tandlammp = tan(.5 * (dlam - .25 * (Y + Y - E * (4. - X)) *
+ (f2 * T + f64 * (32. * T - (20. * T - A)
+ * X - (B +
+ 4.) * Y)) *
+ tan(dlam)));
+ } else {
+ geod_S = geod_a * d;
+ tandlammp = tan(dlamm);
+ }
+ u = atan2(sindthm, (tandlammp * costhm));
+ v = atan2(cosdthm, (tandlammp * sinthm));
+ al12 = adjlon(TWOPI + v - u);
+ al21 = adjlon(TWOPI - v - u);
+}