svn-gvsig-desktop / tags / v2_0_0_Build_2027 / libraries / libjni-proj4 / src / PJ_oea.c @ 39128
History | View | Annotate | Download (1.69 KB)
1 |
#ifndef lint
|
---|---|
2 |
static const char SCCSID[]="@(#)PJ_oea.c 4.1 94/02/15 GIE REL"; |
3 |
#endif
|
4 |
#define PROJ_PARMS__ \
|
5 |
double theta; \
|
6 |
double m, n; \
|
7 |
double two_r_m, two_r_n, rm, rn, hm, hn; \
|
8 |
double cp0, sp0;
|
9 |
#define PJ_LIB__
|
10 |
#include <projects.h> |
11 |
PROJ_HEAD(oea, "Oblated Equal Area") "\n\tMisc Sph\n\tn= m= theta="; |
12 |
FORWARD(s_forward); /* sphere */
|
13 |
double Az, M, N, cp, sp, cl, shz;
|
14 |
|
15 |
cp = cos(lp.phi); |
16 |
sp = sin(lp.phi); |
17 |
cl = cos(lp.lam); |
18 |
Az = aatan2(cp * sin(lp.lam), P->cp0 * sp - P->sp0 * cp * cl) + P->theta; |
19 |
shz = sin(0.5 * aacos(P->sp0 * sp + P->cp0 * cp * cl)); |
20 |
M = aasin(shz * sin(Az)); |
21 |
N = aasin(shz * cos(Az) * cos(M) / cos(M * P->two_r_m)); |
22 |
xy.y = P->n * sin(N * P->two_r_n); |
23 |
xy.x = P->m * sin(M * P->two_r_m) * cos(N) / cos(N * P->two_r_n); |
24 |
return (xy);
|
25 |
} |
26 |
INVERSE(s_inverse); /* sphere */
|
27 |
double N, M, xp, yp, z, Az, cz, sz, cAz;
|
28 |
|
29 |
N = P->hn * aasin(xy.y * P->rn); |
30 |
M = P->hm * aasin(xy.x * P->rm * cos(N * P->two_r_n) / cos(N)); |
31 |
xp = 2. * sin(M);
|
32 |
yp = 2. * sin(N) * cos(M * P->two_r_m) / cos(M);
|
33 |
cAz = cos(Az = aatan2(xp, yp) - P->theta); |
34 |
z = 2. * aasin(0.5 * hypot(xp, yp)); |
35 |
sz = sin(z); |
36 |
cz = cos(z); |
37 |
lp.phi = aasin(P->sp0 * cz + P->cp0 * sz * cAz); |
38 |
lp.lam = aatan2(sz * sin(Az), |
39 |
P->cp0 * cz - P->sp0 * sz * cAz); |
40 |
return (lp);
|
41 |
} |
42 |
FREEUP; if (P) pj_dalloc(P); }
|
43 |
ENTRY0(oea) |
44 |
if (((P->n = pj_param(P->params, "dn").f) <= 0.) || |
45 |
((P->m = pj_param(P->params, "dm").f) <= 0.)) |
46 |
E_ERROR(-39)
|
47 |
else {
|
48 |
P->theta = pj_param(P->params, "rtheta").f;
|
49 |
P->sp0 = sin(P->phi0); |
50 |
P->cp0 = cos(P->phi0); |
51 |
P->rn = 1./ P->n;
|
52 |
P->rm = 1./ P->m;
|
53 |
P->two_r_n = 2. * P->rn;
|
54 |
P->two_r_m = 2. * P->rm;
|
55 |
P->hm = 0.5 * P->m; |
56 |
P->hn = 0.5 * P->n; |
57 |
P->fwd = s_forward; |
58 |
P->inv = s_inverse; |
59 |
P->es = 0.;
|
60 |
} |
61 |
ENDENTRY(P) |