svn-gvsig-desktop / tags / v1_10_0_Build_1258 / libraries / libjni-proj4 / src / geod_inv.c
History | View | Annotate | Download (1.51 KB)
1 |
#ifndef lint
|
---|---|
2 |
static const char SCCSID[]="@(#)geod_inv.c 4.5 95/09/23 GIE REL"; |
3 |
#endif
|
4 |
# include "projects.h" |
5 |
# include "geodesic.h" |
6 |
# define DTOL 1e-12 |
7 |
void
|
8 |
geod_inv(void) {
|
9 |
double th1,th2,thm,dthm,dlamm,dlam,sindlamm,costhm,sinthm,cosdthm,
|
10 |
sindthm,L,E,cosd,d,X,Y,T,sind,tandlammp,u,v,D,A,B; |
11 |
|
12 |
if (ellipse) {
|
13 |
th1 = atan(onef * tan(phi1)); |
14 |
th2 = atan(onef * tan(phi2)); |
15 |
} else {
|
16 |
th1 = phi1; |
17 |
th2 = phi2; |
18 |
} |
19 |
thm = .5 * (th1 + th2);
|
20 |
dthm = .5 * (th2 - th1);
|
21 |
dlamm = .5 * ( dlam = adjlon(lam2 - lam1) );
|
22 |
if (fabs(dlam) < DTOL && fabs(dthm) < DTOL) {
|
23 |
al12 = al21 = geod_S = 0.;
|
24 |
return;
|
25 |
} |
26 |
sindlamm = sin(dlamm); |
27 |
costhm = cos(thm); sinthm = sin(thm); |
28 |
cosdthm = cos(dthm); sindthm = sin(dthm); |
29 |
L = sindthm * sindthm + (cosdthm * cosdthm - sinthm * sinthm) |
30 |
* sindlamm * sindlamm; |
31 |
d = acos(cosd = 1 - L - L);
|
32 |
if (ellipse) {
|
33 |
E = cosd + cosd; |
34 |
sind = sin( d ); |
35 |
Y = sinthm * cosdthm; |
36 |
Y *= (Y + Y) / (1. - L);
|
37 |
T = sindthm * costhm; |
38 |
T *= (T + T) / L; |
39 |
X = Y + T; |
40 |
Y -= T; |
41 |
T = d / sind; |
42 |
D = 4. * T * T;
|
43 |
A = D * E; |
44 |
B = D + D; |
45 |
geod_S = geod_a * sind * (T - f4 * (T * X - Y) + |
46 |
f64 * (X * (A + (T - .5 * (A - E)) * X) -
|
47 |
Y * (B + E * Y) + D * X * Y)); |
48 |
tandlammp = tan(.5 * (dlam - .25 * (Y + Y - E * (4. - X)) * |
49 |
(f2 * T + f64 * (32. * T - (20. * T - A) |
50 |
* X - (B + 4.) * Y)) * tan(dlam)));
|
51 |
} else {
|
52 |
geod_S = geod_a * d; |
53 |
tandlammp = tan(dlamm); |
54 |
} |
55 |
u = atan2(sindthm , (tandlammp * costhm)); |
56 |
v = atan2(cosdthm , (tandlammp * sinthm)); |
57 |
al12 = adjlon(TWOPI + v - u); |
58 |
al21 = adjlon(TWOPI - v - u); |
59 |
} |