diff options
| -rwxr-xr-x | nad/testvarious | 4 | ||||
| -rw-r--r-- | nad/tv_out.dist | 16 | ||||
| -rw-r--r-- | src/PJ_aeqd.c | 17 | ||||
| -rw-r--r-- | src/geod_interface.c | 22 | ||||
| -rw-r--r-- | src/proj_api.h | 4 |
5 files changed, 28 insertions, 35 deletions
diff --git a/nad/testvarious b/nad/testvarious index 694d0d19..988b6d7b 100755 --- a/nad/testvarious +++ b/nad/testvarious @@ -656,7 +656,7 @@ EOF echo "##############################################################" >> ${OUT} echo "Test bug 244 " >> ${OUT} # -$EXE -f '%.12f' \ +$EXE -f '%.9f' \ +init=epsg:4326 \ +to +proj=aeqd +lon_0=130.0 +lat_0=40.0 +a=6378137 +b=6378137 +units=m +no_defs \ -E >>${OUT} <<EOF @@ -669,7 +669,7 @@ $EXE -f '%.12f' \ +proj=aeqd +lon_0=130.0 +lat_0=40.0 +a=6378137 +b=6378137 +units=m +no_defs \ +to +init=epsg:4326 \ -E >>${OUT} <<EOF -987122.4183302754536 -14429896.5395309105515 +987122.418330284 -14429896.539530909 EOF echo "##############################################################" >> ${OUT} echo "Test bug 245 (use +datum=carthage)" >> ${OUT} diff --git a/nad/tv_out.dist b/nad/tv_out.dist index d36ab5b2..8cdd7074 100644 --- a/nad/tv_out.dist +++ b/nad/tv_out.dist @@ -20,8 +20,8 @@ Test simple prime meridian handling. ############################################################## Test support for the lon_wrap switch. 1d00'00.000"W 10d00'00.000"N 0.0 359dE 10dN 0.000 -0d00'00.000"W 10d00'00.000"N 0.0 360dE 10dN 0.000 -0d00'00.000"E 10d00'00.000"N 0.0 360dE 10dN 0.000 +0d00'00.000"W 10d00'00.000"N 0.0 0dE 10dN 0.000 +0d00'00.000"E 10d00'00.000"N 0.0 0dE 10dN 0.000 1d00'00.000"E 45d00'00.000"N 0.0 1dE 45dN 0.000 179d00'00.000"E 45d00'00.000"N 0.0 179dE 45dN 0.000 181d00'00.000"E 45d00'00.000"N 0.0 181dE 45dN 0.000 @@ -101,8 +101,8 @@ Test healpix forward projection on sphere 90.0 0 7.85398 0.00000 0.00000 -90.0 0 -7.85398 0.00000 0.00000 -180 0 -15.70796 0.00000 0.00000 --180 90.0 -19.63495 7.85398 0.00000 --180 -90.0 -19.63495 -7.85398 0.00000 +-180 90.0 -11.78097 7.85398 0.00000 +-180 -90.0 -11.78097 -7.85398 0.00000 0 60.0 1.43738 5.36437 0.00000 0 -60.0 1.43738 -5.36437 0.00000 Test healpix forward projection on ellipsoid @@ -114,8 +114,8 @@ Test healpix forward projection on ellipsoid 90.0 0 6.78898 0.00000 0.00000 -90.0 0 -6.78898 0.00000 0.00000 -180 0 -13.57797 0.00000 0.00000 --180 90.0 -16.97246 6.78898 0.00000 --180 -90.0 -16.97246 -6.78898 0.00000 +-180 90.0 -10.18348 6.78898 0.00000 +-180 -90.0 -10.18348 -6.78898 0.00000 0 60.0 0.00000 3.35128 0.00000 0 -60.0 0.00000 -3.35128 0.00000 Test healpix inverse projection on ellipsoid @@ -307,10 +307,10 @@ Test bug 229 (2) 13 -10 13.0000000000000 -10.0000000000000 0.0000000000000 ############################################################## Test bug 244 --140.100000 -87.000000 987122.418330275454 -14429896.539530910552 0.000000000000 +-140.100000 -87.000000 987122.418330284 -14429896.539530909 0.000000000 ############################################################## Test bug 244 (2) -987122.4183302754536 -14429896.5395309105515 -140.100000000000 -87.000000000000 0.000000000000 +987122.418330284 -14429896.539530909 -140.100000000000 -87.000000000000 0.000000000000 ############################################################## Test bug 245 (use +datum=carthage) 10 34 592302.9819462 3762148.7340610 -30.3110170 diff --git a/src/PJ_aeqd.c b/src/PJ_aeqd.c index 739dfdd6..22a75ac8 100644 --- a/src/PJ_aeqd.c +++ b/src/PJ_aeqd.c @@ -44,7 +44,6 @@ PROJ_HEAD(aeqd, "Azimuthal Equidistant") "\n\tAzi, Sph&Ell\n\tlat_0 guam"; #define EPS10 1.e-10 #define TOL 1.e-14 -#define RHO 57.295779513082320876798154814105 #define N_POLE 0 #define S_POLE 1 @@ -85,11 +84,11 @@ FORWARD(e_forward); /* elliptical */ break; } - phi1 = P->phi0*RHO; lam1 = P->lam0*RHO; - phi2 = lp.phi*RHO; lam2 = (lp.lam+P->lam0)*RHO; + phi1 = P->phi0 / DEG_TO_RAD; lam1 = P->lam0 / DEG_TO_RAD; + phi2 = lp.phi / DEG_TO_RAD; lam2 = (lp.lam+P->lam0) / DEG_TO_RAD; geod_inverse(&P->g, phi1, lam1, phi2, lam2, &s12, &azi1, &azi2); - azi1 /= RHO; + azi1 *= DEG_TO_RAD; xy.x = s12 * sin(azi1) / P->a; xy.y = s12 * cos(azi1) / P->a; break; @@ -160,13 +159,13 @@ INVERSE(e_inverse); /* elliptical */ x2 = xy.x * P->a; y2 = xy.y * P->a; - lat1 = P->phi0 * RHO; - lon1 = P->lam0 * RHO; - azi1 = atan2(x2, y2) * RHO; + lat1 = P->phi0 / DEG_TO_RAD; + lon1 = P->lam0 / DEG_TO_RAD; + azi1 = atan2(x2, y2) / DEG_TO_RAD; s12 = sqrt(x2 * x2 + y2 * y2); geod_direct(&P->g, lat1, lon1, azi1, s12, &lat2, &lon2, &azi2); - lp.phi = lat2 / RHO; - lp.lam = lon2 / RHO; + lp.phi = lat2 * DEG_TO_RAD; + lp.lam = lon2 * DEG_TO_RAD; lp.lam -= P->lam0; } else { /* Polar */ lp.phi = pj_inv_mlfn(P->ctx, P->mode == N_POLE ? P->Mp - c : P->Mp + c, diff --git a/src/geod_interface.c b/src/geod_interface.c index 63b16b1e..a30377ac 100644 --- a/src/geod_interface.c +++ b/src/geod_interface.c @@ -1,20 +1,14 @@ #include "projects.h" #include "geod_interface.h" -/* DEG_IN is a crock to work around the problem that dmstor.c uses the wrong - * value for pi/180 (namely .0174532925199433 which is an inaccurately - * truncated version of DEG_TO_RAD). - */ -#define DEG_IN .0174532925199433 -#define DEG_OUT DEG_TO_RAD; - void geod_ini(void) { geod_init(&GlobalGeodesic, geod_a, geod_f); } void geod_pre(void) { double - lat1 = phi1 / DEG_IN, lon1 = lam1 / DEG_IN, azi1 = al12 / DEG_IN; + lat1 = phi1 / DEG_TO_RAD, lon1 = lam1 / DEG_TO_RAD, + azi1 = al12 / DEG_TO_RAD; geod_lineinit(&GlobalGeodesicLine, &GlobalGeodesic, lat1, lon1, azi1, 0U); } @@ -23,17 +17,17 @@ void geod_for(void) { s12 = geod_S, lat2, lon2, azi2; geod_position(&GlobalGeodesicLine, s12, &lat2, &lon2, &azi2); azi2 += azi2 >= 0 ? -180 : 180; /* Compute back azimuth */ - phi2 = lat2 * DEG_OUT; - lam2 = lon2 * DEG_OUT; - al21 = azi2 * DEG_OUT; + phi2 = lat2 * DEG_TO_RAD; + lam2 = lon2 * DEG_TO_RAD; + al21 = azi2 * DEG_TO_RAD; } void geod_inv(void) { double - lat1 = phi1 / DEG_IN, lon1 = lam1 / DEG_IN, - lat2 = phi2 / DEG_IN, lon2 = lam2 / DEG_IN, + lat1 = phi1 / DEG_TO_RAD, lon1 = lam1 / DEG_TO_RAD, + lat2 = phi2 / DEG_TO_RAD, lon2 = lam2 / DEG_TO_RAD, azi1, azi2, s12; geod_inverse(&GlobalGeodesic, lat1, lon1, lat2, lon2, &s12, &azi1, &azi2); azi2 += azi2 >= 0 ? -180 : 180; /* Compute back azimuth */ - al12 = azi1 * DEG_OUT; al21 = azi2 * DEG_OUT; geod_S = s12; + al12 = azi1 * DEG_TO_RAD; al21 = azi2 * DEG_TO_RAD; geod_S = s12; } diff --git a/src/proj_api.h b/src/proj_api.h index 99faeafb..9e1e49a8 100644 --- a/src/proj_api.h +++ b/src/proj_api.h @@ -46,8 +46,8 @@ extern "C" { extern char const pj_release[]; /* global release id string */ -#define RAD_TO_DEG 57.29577951308232 -#define DEG_TO_RAD .0174532925199432958 +#define RAD_TO_DEG 57.295779513082321 +#define DEG_TO_RAD .017453292519943296 extern int pj_errno; /* global error return code */ |