aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/PJ_aea.c316
-rw-r--r--src/PJ_cass.c85
-rw-r--r--src/lib_proj.cmake2
-rw-r--r--src/pj_generic_selftest.c200
-rw-r--r--src/pj_init.c2
-rw-r--r--src/pj_run_selftests.c63
-rw-r--r--src/proj.c15
-rw-r--r--src/proj_api.h3
-rw-r--r--src/projects.h22
9 files changed, 589 insertions, 119 deletions
diff --git a/src/PJ_aea.c b/src/PJ_aea.c
index 0881ba78..8285e68f 100644
--- a/src/PJ_aea.c
+++ b/src/PJ_aea.c
@@ -1,7 +1,9 @@
/******************************************************************************
* Project: PROJ.4
* Purpose: Implementation of the aea (Albers Equal Area) projection.
- * Author: Gerald Evenden
+ * and the leac (Lambert Equal Area Conic) projection
+ * Author: Gerald Evenden (1995)
+ * Thomas Knudsen (2016) - revise/add regression tests
*
******************************************************************************
* Copyright (c) 1995, Gerald Evenden
@@ -25,35 +27,21 @@
* DEALINGS IN THE SOFTWARE.
*****************************************************************************/
-#define PROJ_PARMS__ \
- double ec; \
- double n; \
- double c; \
- double dd; \
- double n2; \
- double rho0; \
- double rho; \
- double phi1; \
- double phi2; \
- double *en; \
- int ellips;
-
#define PJ_LIB__
#include <projects.h>
# define EPS10 1.e-10
# define TOL7 1.e-7
-PROJ_HEAD(aea, "Albers Equal Area")
- "\n\tConic Sph&Ell\n\tlat_1= lat_2=";
-PROJ_HEAD(leac, "Lambert Equal Area Conic")
- "\n\tConic, Sph&Ell\n\tlat_1= south";
+PROJ_HEAD(aea, "Albers Equal Area") "\n\tConic Sph&Ell\n\tlat_1= lat_2=";
+PROJ_HEAD(leac, "Lambert Equal Area Conic") "\n\tConic, Sph&Ell\n\tlat_1= south";
+
+
/* determine latitude angle phi-1 */
# define N_ITER 15
# define EPSILON 1.0e-7
# define TOL 1.0e-10
- static double
-phi1_(double qs, double Te, double Tone_es) {
+static double phi1_(double qs, double Te, double Tone_es) {
int i;
double Phi, sinpi, cospi, con, com, dphi;
@@ -73,86 +61,276 @@ phi1_(double qs, double Te, double Tone_es) {
} while (fabs(dphi) > TOL && --i);
return( i ? Phi : HUGE_VAL );
}
-FORWARD(e_forward); /* ellipsoid & spheroid */
- if ((P->rho = P->c - (P->ellips ? P->n * pj_qsfn(sin(lp.phi),
- P->e, P->one_es) : P->n2 * sin(lp.phi))) < 0.) F_ERROR
- P->rho = P->dd * sqrt(P->rho);
- xy.x = P->rho * sin( lp.lam *= P->n );
- xy.y = P->rho0 - P->rho * cos(lp.lam);
- return (xy);
+
+
+struct pj_opaque {
+ double ec;
+ double n;
+ double c;
+ double dd;
+ double n2;
+ double rho0;
+ double rho;
+ double phi1;
+ double phi2;
+ double *en;
+ int ellips;
+};
+
+
+
+static XY e_forward (LP lp, PJ *P) { /* Ellipsoid/spheroid, forward */
+ XY xy = {0.0,0.0};
+ struct pj_opaque *Q = P->opaque;
+ if ((Q->rho = Q->c - (Q->ellips ? Q->n * pj_qsfn(sin(lp.phi),
+ P->e, P->one_es) : Q->n2 * sin(lp.phi))) < 0.) F_ERROR
+ Q->rho = Q->dd * sqrt(Q->rho);
+ xy.x = Q->rho * sin( lp.lam *= Q->n );
+ xy.y = Q->rho0 - Q->rho * cos(lp.lam);
+ return xy;
}
-INVERSE(e_inverse) /* ellipsoid & spheroid */;
- if( (P->rho = hypot(xy.x, xy.y = P->rho0 - xy.y)) != 0.0 ) {
- if (P->n < 0.) {
- P->rho = -P->rho;
+
+
+static LP e_inverse (XY xy, PJ *P) { /* Ellipsoid/spheroid, inverse */
+ LP lp = {0.0,0.0};
+ struct pj_opaque *Q = P->opaque;
+ if( (Q->rho = hypot(xy.x, xy.y = Q->rho0 - xy.y)) != 0.0 ) {
+ if (Q->n < 0.) {
+ Q->rho = -Q->rho;
xy.x = -xy.x;
xy.y = -xy.y;
}
- lp.phi = P->rho / P->dd;
- if (P->ellips) {
- lp.phi = (P->c - lp.phi * lp.phi) / P->n;
- if (fabs(P->ec - fabs(lp.phi)) > TOL7) {
+ lp.phi = Q->rho / Q->dd;
+ if (Q->ellips) {
+ lp.phi = (Q->c - lp.phi * lp.phi) / Q->n;
+ if (fabs(Q->ec - fabs(lp.phi)) > TOL7) {
if ((lp.phi = phi1_(lp.phi, P->e, P->one_es)) == HUGE_VAL)
I_ERROR
} else
lp.phi = lp.phi < 0. ? -HALFPI : HALFPI;
- } else if (fabs(lp.phi = (P->c - lp.phi * lp.phi) / P->n2) <= 1.)
+ } else if (fabs(lp.phi = (Q->c - lp.phi * lp.phi) / Q->n2) <= 1.)
lp.phi = asin(lp.phi);
else
lp.phi = lp.phi < 0. ? -HALFPI : HALFPI;
- lp.lam = atan2(xy.x, xy.y) / P->n;
+ lp.lam = atan2(xy.x, xy.y) / Q->n;
} else {
lp.lam = 0.;
- lp.phi = P->n > 0. ? HALFPI : - HALFPI;
+ lp.phi = Q->n > 0. ? HALFPI : - HALFPI;
}
- return (lp);
+ return lp;
}
-FREEUP; if (P) { if (P->en) pj_dalloc(P->en); pj_dalloc(P); } }
- static PJ *
-setup(PJ *P) {
+
+
+static void *freeup_new (PJ *P) { /* Destructor */
+ if (0==P)
+ return 0;
+
+ if (0==P->opaque)
+ return pj_dealloc (P);
+
+ pj_dealloc (P->opaque->en);
+ pj_dealloc (P->opaque);
+ return pj_dealloc(P);
+}
+
+
+static void freeup (PJ *P) {
+ freeup_new (P);
+ return;
+}
+
+
+static PJ *setup(PJ *P) {
double cosphi, sinphi;
int secant;
+ struct pj_opaque *Q = P->opaque;
- if (fabs(P->phi1 + P->phi2) < EPS10) E_ERROR(-21);
- P->n = sinphi = sin(P->phi1);
- cosphi = cos(P->phi1);
- secant = fabs(P->phi1 - P->phi2) >= EPS10;
- if( (P->ellips = (P->es > 0.))) {
+ P->inv = e_inverse;
+ P->fwd = e_forward;
+
+ if (fabs(Q->phi1 + Q->phi2) < EPS10) E_ERROR(-21);
+ Q->n = sinphi = sin(Q->phi1);
+ cosphi = cos(Q->phi1);
+ secant = fabs(Q->phi1 - Q->phi2) >= EPS10;
+ if( (Q->ellips = (P->es > 0.))) {
double ml1, m1;
- if (!(P->en = pj_enfn(P->es))) E_ERROR_0;
+ if (!(Q->en = pj_enfn(P->es))) E_ERROR_0;
m1 = pj_msfn(sinphi, cosphi, P->es);
ml1 = pj_qsfn(sinphi, P->e, P->one_es);
if (secant) { /* secant cone */
double ml2, m2;
- sinphi = sin(P->phi2);
- cosphi = cos(P->phi2);
+ sinphi = sin(Q->phi2);
+ cosphi = cos(Q->phi2);
m2 = pj_msfn(sinphi, cosphi, P->es);
ml2 = pj_qsfn(sinphi, P->e, P->one_es);
- P->n = (m1 * m1 - m2 * m2) / (ml2 - ml1);
+ Q->n = (m1 * m1 - m2 * m2) / (ml2 - ml1);
}
- P->ec = 1. - .5 * P->one_es * log((1. - P->e) /
+ Q->ec = 1. - .5 * P->one_es * log((1. - P->e) /
(1. + P->e)) / P->e;
- P->c = m1 * m1 + P->n * ml1;
- P->dd = 1. / P->n;
- P->rho0 = P->dd * sqrt(P->c - P->n * pj_qsfn(sin(P->phi0),
+ Q->c = m1 * m1 + Q->n * ml1;
+ Q->dd = 1. / Q->n;
+ Q->rho0 = Q->dd * sqrt(Q->c - Q->n * pj_qsfn(sin(P->phi0),
P->e, P->one_es));
} else {
- if (secant) P->n = .5 * (P->n + sin(P->phi2));
- P->n2 = P->n + P->n;
- P->c = cosphi * cosphi + P->n2 * sinphi;
- P->dd = 1. / P->n;
- P->rho0 = P->dd * sqrt(P->c - P->n2 * sin(P->phi0));
+ if (secant) Q->n = .5 * (Q->n + sin(Q->phi2));
+ Q->n2 = Q->n + Q->n;
+ Q->c = cosphi * cosphi + Q->n2 * sinphi;
+ Q->dd = 1. / Q->n;
+ Q->rho0 = Q->dd * sqrt(Q->c - Q->n2 * sin(P->phi0));
}
- P->inv = e_inverse; P->fwd = e_forward;
+
return P;
}
-ENTRY1(aea,en)
- P->phi1 = pj_param(P->ctx, P->params, "rlat_1").f;
- P->phi2 = pj_param(P->ctx, P->params, "rlat_2").f;
-ENDENTRY(setup(P))
-ENTRY1(leac,en)
- P->phi2 = pj_param(P->ctx, P->params, "rlat_1").f;
- P->phi1 = pj_param(P->ctx, P->params, "bsouth").i ? - HALFPI: HALFPI;
-ENDENTRY(setup(P))
+
+
+PJ *PROJECTION(aea) {
+ struct pj_opaque *Q = pj_calloc (1, sizeof (struct pj_opaque));
+ if (0==Q)
+ return freeup_new (P);
+ P->opaque = Q;
+
+ P->pfree = freeup;
+ P->descr = des_aea;
+
+ Q->phi1 = pj_param(P->ctx, P->params, "rlat_1").f;
+ Q->phi2 = pj_param(P->ctx, P->params, "rlat_2").f;
+ setup(P);
+ return P;
+}
+
+
+PJ *PROJECTION(leac) {
+ struct pj_opaque *Q = pj_calloc (1, sizeof (struct pj_opaque));
+ if (0==Q)
+ return freeup_new (P);
+ P->opaque = Q;
+
+ P->pfree = freeup;
+ P->descr = des_leac;
+ Q->phi2 = pj_param(P->ctx, P->params, "rlat_1").f;
+ Q->phi1 = pj_param(P->ctx, P->params, "bsouth").i ? - HALFPI: HALFPI;
+ setup(P);
+ return P;
+}
+
+
+#ifdef PJ_OMIT_SELFTEST
+int pj_aea_selftest (void) {return 0;}
+#else
+
+int pj_aea_selftest (void) {
+ double tolerance_lp = 1e-10;
+ double tolerance_xy = 1e-7;
+
+ char e_args[] = {"+proj=aea +ellps=GRS80 +lat_1=0 +lat_2=2"};
+ char s_args[] = {"+proj=aea +ellps=GRS80 +lat_1=0 +lat_2=2"};
+
+ LP fwd_in[] = {
+ { 2, 1},
+ { 2,-1},
+ {-2, 1},
+ {-2,-1}
+ };
+
+ XY e_fwd_expect[] = {
+ {222571.60875710563, 110653.32674302977},
+ {222706.30650839131, -110484.26714439997},
+ {-222571.60875710563, 110653.32674302977},
+ {-222706.30650839131, -110484.26714439997},
+ };
+
+ XY s_fwd_expect[] = {
+ {223334.08517088494, 111780.43188447191},
+ {223470.15499168713, -111610.33943099028},
+ {-223334.08517088494, 111780.43188447191},
+ {-223470.15499168713, -111610.33943099028},
+ };
+
+ XY inv_in[] = {
+ { 200, 100},
+ { 200,-100},
+ {-200, 100},
+ {-200,-100}
+ };
+
+ LP e_inv_expect[] = {
+ {0.0017966310597749514, 0.00090436885862202158},
+ {0.0017966300767030448, -0.00090437009538581453},
+ {-0.0017966310597749514, 0.00090436885862202158},
+ {-0.0017966300767030448, -0.00090437009538581453},
+ };
+
+ LP s_inv_expect[] = {
+ {0.0017904935979658752, 0.00089524594491375306},
+ {0.0017904926216016812, -0.00089524716502493225},
+ {-0.0017904935979658752, 0.00089524594491375306},
+ {-0.0017904926216016812, -0.00089524716502493225},
+ };
+ return pj_generic_selftest (e_args, s_args, tolerance_xy, tolerance_lp, 4, 4, fwd_in, e_fwd_expect, s_fwd_expect, inv_in, e_inv_expect, s_inv_expect);
+}
+
+
+#endif
+
+
+
+#ifdef PJ_OMIT_SELFTEST
+int pj_leac_selftest (void) {return 0;}
+#else
+
+int pj_leac_selftest (void) {
+ double tolerance_lp = 1e-10;
+ double tolerance_xy = 1e-7;
+
+ char e_args[] = {"+proj=leac +ellps=GRS80 +lat_1=0 +lat_2=2"};
+ char s_args[] = {"+proj=leac +ellps=GRS80 +lat_1=0 +lat_2=2"};
+
+ LP fwd_in[] = {
+ { 2, 1},
+ { 2,-1},
+ {-2, 1},
+ {-2,-1}
+ };
+
+ XY e_fwd_expect[] = {
+ {220685.14054297868, 112983.50088939646},
+ {224553.31227982609, -108128.63674487274},
+ {-220685.14054297868, 112983.50088939646},
+ {-224553.31227982609, -108128.63674487274},
+ };
+
+ XY s_fwd_expect[] = {
+ {221432.86859285168, 114119.45452653214},
+ {225331.72412711097, -109245.82943505641},
+ {-221432.86859285168, 114119.45452653214},
+ {-225331.72412711097, -109245.82943505641},
+ };
+
+ XY inv_in[] = {
+ { 200, 100},
+ { 200,-100},
+ {-200, 100},
+ {-200,-100}
+ };
+
+ LP e_inv_expect[] = {
+ {0.0017966446840328458, 0.00090435171340223211},
+ {0.0017966164523713021, -0.00090438724081843625},
+ {-0.0017966446840328458, 0.00090435171340223211},
+ {-0.0017966164523713021, -0.00090438724081843625},
+ };
+
+ LP s_inv_expect[] = {
+ {0.0017905070979748127, 0.00089522906964877795},
+ {0.001790479121519977, -0.00089526404022281043},
+ {-0.0017905070979748127, 0.00089522906964877795},
+ {-0.001790479121519977, -0.00089526404022281043},
+ };
+
+ return pj_generic_selftest (e_args, s_args, tolerance_xy, tolerance_lp, 4, 4, fwd_in, e_fwd_expect, s_fwd_expect, inv_in, e_inv_expect, s_inv_expect);
+}
+
+
+#endif
diff --git a/src/PJ_cass.c b/src/PJ_cass.c
index 9e26ad8f..5f885376 100644
--- a/src/PJ_cass.c
+++ b/src/PJ_cass.c
@@ -11,7 +11,7 @@ PROJ_HEAD(cass, "Cassini") "\n\tCyl, Sph&Ell";
# define C5 .06666666666666666666
-struct opaque {
+struct pj_opaque {
/* These are the only opaque members actually initialized */
double *en;
double m0;
@@ -32,17 +32,17 @@ struct opaque {
static XY e_forward (LP lp, PJ *P) { /* Ellipsoidal, forward */
XY xy = {0.0,0.0};
- struct opaque *O = P->opaq;
- xy.y = pj_mlfn(lp.phi, O->n = sin(lp.phi), O->c = cos(lp.phi), O->en);
- O->n = 1./sqrt(1. - P->es * O->n * O->n);
- O->tn = tan(lp.phi); O->t = O->tn * O->tn;
- O->a1 = lp.lam * O->c;
- O->c *= P->es * O->c / (1 - P->es);
- O->a2 = O->a1 * O->a1;
- xy.x = O->n * O->a1 * (1. - O->a2 * O->t *
- (C1 - (8. - O->t + 8. * O->c) * O->a2 * C2));
- xy.y -= O->m0 - O->n * O->tn * O->a2 *
- (.5 + (5. - O->t + 6. * O->c) * O->a2 * C3);
+ struct pj_opaque *Q = P->opaque;
+ xy.y = pj_mlfn(lp.phi, Q->n = sin(lp.phi), Q->c = cos(lp.phi), Q->en);
+ Q->n = 1./sqrt(1. - P->es * Q->n * Q->n);
+ Q->tn = tan(lp.phi); Q->t = Q->tn * Q->tn;
+ Q->a1 = lp.lam * Q->c;
+ Q->c *= P->es * Q->c / (1 - P->es);
+ Q->a2 = Q->a1 * Q->a1;
+ xy.x = Q->n * Q->a1 * (1. - Q->a2 * Q->t *
+ (C1 - (8. - Q->t + 8. * Q->c) * Q->a2 * C2));
+ xy.y -= Q->m0 - Q->n * Q->tn * Q->a2 *
+ (.5 + (5. - Q->t + 6. * Q->c) * Q->a2 * C3);
return xy;
}
@@ -57,45 +57,48 @@ static XY s_forward (LP lp, PJ *P) { /* Spheroidal, forward */
static LP e_inverse (XY xy, PJ *P) { /* Ellipsoidal, inverse */
LP lp = {0.0,0.0};
- struct opaque *O = P->opaq;
+ struct pj_opaque *Q = P->opaque;
double ph1;
- ph1 = pj_inv_mlfn(P->ctx, O->m0 + xy.y, P->es, O->en);
- O->tn = tan(ph1); O->t = O->tn * O->tn;
- O->n = sin(ph1);
- O->r = 1. / (1. - P->es * O->n * O->n);
- O->n = sqrt(O->r);
- O->r *= (1. - P->es) * O->n;
- O->dd = xy.x / O->n;
- O->d2 = O->dd * O->dd;
- lp.phi = ph1 - (O->n * O->tn / O->r) * O->d2 *
- (.5 - (1. + 3. * O->t) * O->d2 * C3);
- lp.lam = O->dd * (1. + O->t * O->d2 *
- (-C4 + (1. + 3. * O->t) * O->d2 * C5)) / cos(ph1);
+ ph1 = pj_inv_mlfn(P->ctx, Q->m0 + xy.y, P->es, Q->en);
+ Q->tn = tan(ph1); Q->t = Q->tn * Q->tn;
+ Q->n = sin(ph1);
+ Q->r = 1. / (1. - P->es * Q->n * Q->n);
+ Q->n = sqrt(Q->r);
+ Q->r *= (1. - P->es) * Q->n;
+ Q->dd = xy.x / Q->n;
+ Q->d2 = Q->dd * Q->dd;
+ lp.phi = ph1 - (Q->n * Q->tn / Q->r) * Q->d2 *
+ (.5 - (1. + 3. * Q->t) * Q->d2 * C3);
+ lp.lam = Q->dd * (1. + Q->t * Q->d2 *
+ (-C4 + (1. + 3. * Q->t) * Q->d2 * C5)) / cos(ph1);
return lp;
}
static LP s_inverse (XY xy, PJ *P) { /* Spheroidal, inverse */
LP lp = {0.0,0.0};
- lp.phi = asin(sin(P->opaq->dd = xy.y + P->phi0) * cos(xy.x));
- lp.lam = atan2(tan(xy.x), cos(P->opaq->dd));
+ lp.phi = asin(sin(P->opaque->dd = xy.y + P->phi0) * cos(xy.x));
+ lp.lam = atan2(tan(xy.x), cos(P->opaque->dd));
return lp;
}
-static void freeup(PJ *P) { /* Destructor */
+static void *freeup_new(PJ *P) { /* Destructor */
if (0==P)
- return;
- if (0==P->opaq) {
- pj_dealloc (P);
- return;
- }
- pj_dealloc(P->opaq->en);
- pj_dealloc(P->opaq);
- pj_dealloc(P);
+ return 0;
+ if (0==P->opaque)
+ return pj_dealloc (P);
+
+ pj_dealloc(P->opaque->en);
+ pj_dealloc(P->opaque);
+ return pj_dealloc(P);
}
+static void freeup(PJ *P) { /* Destructor */
+ freeup_new (P);
+ return;
+}
PJ *PROJECTION(cass) {
@@ -107,19 +110,19 @@ PJ *PROJECTION(cass) {
}
/* Otherwise it's ellipsoidal */
- P->opaq = pj_calloc (1, sizeof (struct opaque));
- if (0==P->opaq) {
+ P->opaque = pj_calloc (1, sizeof (struct pj_opaque));
+ if (0==P->opaque) {
freeup (P);
return 0;
}
- P->opaq->en = pj_enfn(P->es);
- if (0==P->opaq->en) {
+ P->opaque->en = pj_enfn(P->es);
+ if (0==P->opaque->en) {
freeup (P);
return 0;
}
- P->opaq->m0 = pj_mlfn(P->phi0, sin(P->phi0), cos(P->phi0), P->opaq->en);
+ P->opaque->m0 = pj_mlfn(P->phi0, sin(P->phi0), cos(P->phi0), P->opaque->en);
P->inv = e_inverse;
P->fwd = e_forward;
diff --git a/src/lib_proj.cmake b/src/lib_proj.cmake
index a8797a74..eeb0582d 100644
--- a/src/lib_proj.cmake
+++ b/src/lib_proj.cmake
@@ -172,6 +172,7 @@ SET(SRC_LIBPROJ_CORE
pj_fwd3d.c
pj_gauss.c
pj_gc_reader.c
+ pj_generic_selftest.c
pj_geocent.c
pj_gridcatalog.c
pj_gridinfo.c
@@ -195,6 +196,7 @@ SET(SRC_LIBPROJ_CORE
pj_pr_list.c
pj_qsfn.c
pj_release.c
+ pj_run_selftests.c
pj_strerrno.c
pj_transform.c
pj_tsfn.c
diff --git a/src/pj_generic_selftest.c b/src/pj_generic_selftest.c
new file mode 100644
index 00000000..e2d24baf
--- /dev/null
+++ b/src/pj_generic_selftest.c
@@ -0,0 +1,200 @@
+/******************************************************************************
+ * Project: PROJ.4
+ * Purpose: Generic regression test for PROJ.4 projection algorithms.
+ * Author: Thomas Knudsen
+ *
+ ******************************************************************************
+ * Copyright (c) 2016, Thomas Knudsen
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a
+ * copy of this software and associated documentation files (the "Software"),
+ * to deal in the Software without restriction, including without limitation
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,
+ * and/or sell copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included
+ * in all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+ * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+ * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+ * DEALINGS IN THE SOFTWARE.
+ *****************************************************************************/
+
+
+#include <stdio.h>
+#define PJ_LIB__
+#include <projects.h>
+
+
+static int deviates_xy (XY expected, XY got, double tolerance);
+static int deviates_lp (LP expected, LP got, double tolerance);
+static XY pj_fwd_deg (LP in, PJ *P);
+
+
+/**********************************************************************/
+int pj_generic_selftest (
+/**********************************************************************/
+ char *e_args,
+ char *s_args,
+ double tolerance_xy,
+ double tolerance_lp,
+ int n_fwd,
+ int n_inv,
+ LP *fwd_in,
+ XY *e_fwd_expect,
+ XY *s_fwd_expect,
+ XY *inv_in,
+ LP *e_inv_expect,
+ LP *s_inv_expect
+) {
+/***********************************************************************
+
+Generic regression test for PROJ.4 projection algorithms, testing both
+ellipsoidal ("e_") and spheroidal ("s_") versions of the projection
+algorithms in both forward ("_fwd_") and inverse ("_inv_") mode.
+
+Compares the "known good" results in <e_fwd_expect> and <s_fwd_expect>
+with the actual results obtained by transforming the forward input data
+set in <fwd_in> with pj_fwd() using setup arguments <e_args> and
+<s_args>, respectively.
+
+Then
+
+Compares the "known good" results in <e_inv_expect> and <s_inv_expect>
+with the actual results obtained by transforming the inverse input data
+set in <inv_in> with pj_inv() using setup arguments <e_args> and
+<s_args>, respectively.
+
+Any of the pointers passed may be set to 0, indicating "don't test this
+part".
+
+Returns 0 if all data agree to within the accuracy specified in
+<tolerance_xy> and <tolerance_lp>. Non-zero otherwise.
+
+***********************************************************************/
+ int i;
+
+ PJ *P;
+
+puts ("testing");
+ if (e_args) {
+
+ puts ("e_args");
+ puts (e_args);
+
+ P = pj_init_plus(e_args);
+ if (0==P)
+ return 2;
+
+ /* Test forward ellipsoidal */
+ if (e_fwd_expect) { puts ("e_fwd");
+ for (i = 0; i < n_fwd; i++)
+ if (deviates_xy (e_fwd_expect[i], pj_fwd_deg ( fwd_in[i], P ), tolerance_xy))
+ break;
+ if ( i != n_fwd )
+ return 100 + i;
+ }
+
+ /* Test inverse ellipsoidal */
+ if (e_inv_expect) { puts ("e_inv");
+ for (i = 0; i < n_inv; i++)
+ if (deviates_lp (e_inv_expect[i], pj_inv ( inv_in[i], P ), tolerance_lp))
+ break;
+ if ( i != n_inv )
+ return 200 + i;
+ }
+
+ pj_free (P);
+ }
+
+
+ if (s_args) {
+ puts ("s_args");
+ puts (s_args);
+ P = pj_init_plus(s_args);
+ if (0==P)
+ return 3;
+
+ /* Test forward spherical */
+ if (s_fwd_expect) { puts ("s_fwd");
+ for (i = 0; i < n_fwd; i++)
+ if (deviates_xy (s_fwd_expect[i], pj_fwd_deg ( fwd_in[i], P ), tolerance_xy))
+ break;
+ if ( i != n_fwd )
+ return 300 + i;
+ }
+
+ /* Test inverse spherical */
+ if (s_inv_expect) { puts ("s_inv");
+ for (i = 0; i < n_inv; i++)
+ if (deviates_lp (s_inv_expect[i], pj_inv ( inv_in[i], P ), tolerance_lp))
+ break;
+ if ( i != n_inv )
+ return 400 + i;
+ }
+
+ pj_free (P);
+ }
+
+ return 0;
+}
+
+
+
+/**********************************************************************/
+static int deviates_xy (XY expected, XY got, double tolerance) {
+/***********************************************************************
+
+ Determine whether two XYs deviate by more than <tolerance>.
+
+***********************************************************************/
+ if (HUGE_VAL== expected.x)
+ return 0;
+ if (HUGE_VAL== expected.y)
+ return 0;
+ if (hypot ( expected.x - got.x, expected.y - got.y ) > tolerance)
+ return 1;
+ return 0;
+}
+
+
+/**********************************************************************/
+static int deviates_lp (LP expected, LP got, double tolerance) {
+/***********************************************************************
+
+ Determine whether two LPs deviate by more than <tolerance>.
+
+ This one is slightly tricky, since the <expected> LP is
+ supposed to be represented as degrees (since it was at some
+ time written down by a real human), whereas the <got> LP is
+ represented in radians (since it is supposed to be the result
+ output from pj_inv)
+
+***********************************************************************/
+ if (HUGE_VAL== expected.lam)
+ return 0;
+ if (HUGE_VAL== expected.phi)
+ return 0;
+ if (hypot ( DEG_TO_RAD * expected.lam - got.lam, DEG_TO_RAD * expected.phi - got.phi ) > tolerance)
+ return 1;
+ return 0;
+}
+
+
+/**********************************************************************/
+static XY pj_fwd_deg (LP in, PJ *P) {
+/***********************************************************************
+
+ Wrapper for pj_fwd, accepting input in degrees.
+
+***********************************************************************/
+ LP in_rad;
+ in_rad.lam = DEG_TO_RAD * in.lam;
+ in_rad.phi = DEG_TO_RAD * in.phi;
+ return pj_fwd (in_rad, P);
+}
diff --git a/src/pj_init.c b/src/pj_init.c
index ca7b9e4a..79b64bb3 100644
--- a/src/pj_init.c
+++ b/src/pj_init.c
@@ -661,5 +661,5 @@ pj_free(PJ *P) {
void pj_prepare (PJ *P, const char *description, void (*freeup)(struct PJconsts *), size_t sizeof_struct_opaque) {
P->descr = description;
P->pfree = freeup;
- P->opaq = pj_calloc (1, sizeof_struct_opaque);
+ P->opaque = pj_calloc (1, sizeof_struct_opaque);
}
diff --git a/src/pj_run_selftests.c b/src/pj_run_selftests.c
new file mode 100644
index 00000000..44573b2b
--- /dev/null
+++ b/src/pj_run_selftests.c
@@ -0,0 +1,63 @@
+/******************************************************************************
+ * Project: PROJ.4
+ * Purpose: Generic regression test for PROJ.4 projection algorithms.
+ * Author: Thomas Knudsen
+ *
+ ******************************************************************************
+ * Copyright (c) 2016, Thomas Knudsen
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a
+ * copy of this software and associated documentation files (the "Software"),
+ * to deal in the Software without restriction, including without limitation
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,
+ * and/or sell copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included
+ * in all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+ * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+ * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+ * DEALINGS IN THE SOFTWARE.
+ *****************************************************************************/
+
+
+#include <stdio.h>
+#define PJ_LIB__
+#include <projects.h>
+
+
+extern int pj_aea_selftest(void);
+extern int pj_leac_selftest(void);
+
+
+
+static void run_one_test (const char *mnemonic, int (testfunc)(void), int verbosity, int *n_ok, int *n_ko) {
+ int ret = testfunc ();
+ if (ret)
+ (*n_ko)++;
+ else
+ (*n_ok)++;
+
+ if (verbosity)
+ printf ("Testing: %8s - return code: %d\n", mnemonic, ret);
+}
+
+
+int pj_run_selftests (int verbosity) {
+ int n_ok = 0, n_ko = 0;
+
+ if (verbosity)
+ printf ("Running internal regression tests\n");
+ run_one_test ("aea", pj_aea_selftest, verbosity, &n_ok, &n_ko);
+ run_one_test ("leac", pj_leac_selftest, verbosity, &n_ok, &n_ko);
+
+ if (0==verbosity)
+ printf ("Internal regression tests done. ");
+ printf ("Total: %d, Failure: %d, Success: %d\n", n_ok+n_ko, n_ko, n_ok);
+ return n_ko;
+}
diff --git a/src/proj.c b/src/proj.c
index 23f9cca9..40c6df94 100644
--- a/src/proj.c
+++ b/src/proj.c
@@ -42,7 +42,7 @@ postscale = 0;
*oform = (char *)0, /* output format for x-y or decimal degrees */
*oterr = "*\t*", /* output line for unprojectable input */
*usage =
-"%s\nusage: %s [ -beEfiIlormsStTvVwW [args] ] [ +opts[=arg] ] [ files ]\n";
+"%s\nusage: %s [ -bCeEfiIlormsStTvVwW [args] ] [ +opts[=arg] ] [ files ]\n";
static struct FACTORS
facs;
static double
@@ -158,7 +158,7 @@ vprocess(FILE *fid) {
if (!oform)
oform = "%.3f";
if (bin_in || bin_out)
- emess(1,"binary I/O not available in -V option");
+ emess(1,"binary I/O not available in -V option");
for (;;) {
++emess_dat.File_line;
if (!(s = fgets(line, MAX_LINE, fid)))
@@ -267,6 +267,9 @@ int main(int argc, char **argv) {
case 'b': /* binary I/O */
bin_in = bin_out = 1;
continue;
+ case 'C': /* Check - run internal regression tests */
+ pj_run_selftests (very_verby);
+ continue;
case 'v': /* monitor dump of initialization */
mon = 1;
continue;
@@ -300,8 +303,8 @@ int main(int argc, char **argv) {
char *str;
for (lp = pj_get_list_ref() ; lp->id ; ++lp) {
- if( strcmp(lp->id,"latlong") == 0
- || strcmp(lp->id,"longlat") == 0
+ if( strcmp(lp->id,"latlong") == 0
+ || strcmp(lp->id,"longlat") == 0
|| strcmp(lp->id,"geocent") == 0 )
continue;
@@ -353,7 +356,7 @@ int main(int argc, char **argv) {
continue; /* artificial */
case 'e': /* error line alternative */
if (--argc <= 0)
- noargument:
+ noargument:
emess(1,"missing argument for -%c",*arg);
oterr = *++argv;
continue;
@@ -364,7 +367,7 @@ int main(int argc, char **argv) {
case 'm': /* cartesian multiplier */
if (--argc <= 0) goto noargument;
postscale = 1;
- if (!strncmp("1/",*++argv,2) ||
+ if (!strncmp("1/",*++argv,2) ||
!strncmp("1:",*argv,2)) {
if((fscale = atof((*argv)+2)) == 0.)
goto badscale;
diff --git a/src/proj_api.h b/src/proj_api.h
index dd25d61a..24a6f053 100644
--- a/src/proj_api.h
+++ b/src/proj_api.h
@@ -168,6 +168,9 @@ char *pj_ctx_fgets(projCtx ctx, char *line, int size, PAFile file);
PAFile pj_open_lib(projCtx, const char *, const char *);
+int pj_run_selftests (int verbosity);
+
+
#define PJ_LOG_NONE 0
#define PJ_LOG_ERROR 1
#define PJ_LOG_DEBUG_MAJOR 2
diff --git a/src/projects.h b/src/projects.h
index 3106eb76..2b638e94 100644
--- a/src/projects.h
+++ b/src/projects.h
@@ -229,7 +229,7 @@ typedef struct ARG_list {
#ifdef PJ_LIB__
/* we need this forward declaration in order to be able to add a
pointer to struct opaque to the typedef struct PJconsts below */
- struct opaque;
+ struct pj_opaque;
#endif
typedef struct PJconsts {
@@ -290,7 +290,7 @@ typedef struct PJconsts {
double last_after_date;
#ifdef PJ_LIB__
- struct opaque *opaq;
+ struct pj_opaque *opaque;
#endif
#ifdef PROJ_PARMS__
@@ -381,6 +381,24 @@ PJ *pj_projection_specific_setup_##name (PJ *P)
#endif
+int pj_generic_selftest (
+ char *e_args,
+ char *s_args,
+ double tolerance_xy,
+ double tolerance_lp,
+ int n_fwd,
+ int n_inv,
+ LP *fwd_in,
+ XY *e_fwd_expect,
+ XY *s_fwd_expect,
+ XY *inv_in,
+ LP *e_inv_expect,
+ LP *s_inv_expect
+);
+
+
+
+
#define MAX_TAB_ID 80
typedef struct { float lam, phi; } FLP;
typedef struct { int lam, phi; } ILP;
generated by cgit v1.2.3 (git 2.25.1) at 2025-12-31 04:42:28 +0000