svn-gvsig-desktop / branches / v2_0_0_prep / libraries / libjni-proj4 / src / PJ_eqdc.c @ 29018
History | View | Annotate | Download (2.28 KB)
1 |
#ifndef lint
|
---|---|
2 |
static const char SCCSID[]="@(#)PJ_eqdc.c 4.2 94/03/16 GIE REL"; |
3 |
#endif
|
4 |
#define PROJ_PARMS__ \
|
5 |
double phi1; \
|
6 |
double phi2; \
|
7 |
double n; \
|
8 |
double rho; \
|
9 |
double rho0; \
|
10 |
double c; \
|
11 |
double *en; \
|
12 |
int ellips;
|
13 |
#define PJ_LIB__
|
14 |
#include <projects.h> |
15 |
PROJ_HEAD(eqdc, "Equidistant Conic")
|
16 |
"\n\tConic, Sph&Ell\n\tlat_1= lat_2=";
|
17 |
# define EPS10 1.e-10 |
18 |
FORWARD(e_forward); /* sphere & ellipsoid */
|
19 |
P->rho = P->c - (P->ellips ? pj_mlfn(lp.phi, sin(lp.phi), |
20 |
cos(lp.phi), P->en) : lp.phi); |
21 |
xy.x = P->rho * sin( lp.lam *= P->n ); |
22 |
xy.y = P->rho0 - P->rho * cos(lp.lam); |
23 |
return (xy);
|
24 |
} |
25 |
INVERSE(e_inverse); /* sphere & ellipsoid */
|
26 |
if ((P->rho = hypot(xy.x, xy.y = P->rho0 - xy.y)) != 0.0 ) { |
27 |
if (P->n < 0.) { |
28 |
P->rho = -P->rho; |
29 |
xy.x = -xy.x; |
30 |
xy.y = -xy.y; |
31 |
} |
32 |
lp.phi = P->c - P->rho; |
33 |
if (P->ellips)
|
34 |
lp.phi = pj_inv_mlfn(lp.phi, P->es, P->en); |
35 |
lp.lam = atan2(xy.x, xy.y) / P->n; |
36 |
} else {
|
37 |
lp.lam = 0.;
|
38 |
lp.phi = P->n > 0. ? HALFPI : - HALFPI;
|
39 |
} |
40 |
return (lp);
|
41 |
} |
42 |
SPECIAL(fac) { |
43 |
double sinphi, cosphi;
|
44 |
|
45 |
sinphi = sin(lp.phi); |
46 |
cosphi = cos(lp.phi); |
47 |
fac->code |= IS_ANAL_HK; |
48 |
fac->h = 1.;
|
49 |
fac->k = P->n * (P->c - (P->ellips ? pj_mlfn(lp.phi, sinphi, |
50 |
cosphi, P->en) : lp.phi)) / pj_msfn(sinphi, cosphi, P->es); |
51 |
} |
52 |
FREEUP; if (P) { if (P->en) pj_dalloc(P->en); pj_dalloc(P); } } |
53 |
ENTRY1(eqdc, en) |
54 |
double cosphi, sinphi;
|
55 |
int secant;
|
56 |
|
57 |
P->phi1 = pj_param(P->params, "rlat_1").f;
|
58 |
P->phi2 = pj_param(P->params, "rlat_2").f;
|
59 |
if (fabs(P->phi1 + P->phi2) < EPS10) E_ERROR(-21); |
60 |
if (!(P->en = pj_enfn(P->es)))
|
61 |
E_ERROR_0; |
62 |
P->n = sinphi = sin(P->phi1); |
63 |
cosphi = cos(P->phi1); |
64 |
secant = fabs(P->phi1 - P->phi2) >= EPS10; |
65 |
if( (P->ellips = (P->es > 0.)) ) { |
66 |
double ml1, m1;
|
67 |
|
68 |
m1 = pj_msfn(sinphi, cosphi, P->es); |
69 |
P->en = pj_enfn(P->es); |
70 |
ml1 = pj_mlfn(P->phi1, sinphi, cosphi, P->en); |
71 |
if (secant) { /* secant cone */ |
72 |
sinphi = sin(P->phi2); |
73 |
cosphi = cos(P->phi2); |
74 |
P->n = (m1 - pj_msfn(sinphi, cosphi, P->es)) / |
75 |
(pj_mlfn(P->phi2, sinphi, cosphi, P->en) - ml1); |
76 |
} |
77 |
P->c = ml1 + m1 / P->n; |
78 |
P->rho0 = P->c - pj_mlfn(P->phi0, sin(P->phi0), |
79 |
cos(P->phi0), P->en); |
80 |
} else {
|
81 |
if (secant)
|
82 |
P->n = (cosphi - cos(P->phi2)) / (P->phi2 - P->phi1); |
83 |
P->c = P->phi1 + cos(P->phi1) / P->n; |
84 |
P->rho0 = P->c - P->phi0; |
85 |
} |
86 |
P->inv = e_inverse; |
87 |
P->fwd = e_forward; |
88 |
P->spc = fac; |
89 |
ENDENTRY(P) |