From db5fe58b875664cf852fe45eab961cdae8332970 Mon Sep 17 00:00:00 2001 From: Thomas Knudsen Date: Thu, 22 Feb 2018 08:45:18 +0100 Subject: horner: support westings/southings in complex case --- src/PJ_horner.c | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) (limited to 'src') diff --git a/src/PJ_horner.c b/src/PJ_horner.c index d1146aa1..a132eff6 100644 --- a/src/PJ_horner.c +++ b/src/PJ_horner.c @@ -96,6 +96,8 @@ static HORNER *horner_alloc (size_t order, int complex_polynomia); static void horner_free (HORNER *h); struct horner { + int uneg; /* u axis negated? */ + int vneg; /* v axis negated? */ int order; /* maximum degree of polynomium */ int coefs; /* number of coefficients for each polynomium */ double range; /* radius of the region of validity */ @@ -340,11 +342,19 @@ polynomial evaluation engine. c = cb + sz; e = position.u - transformation->fwd_origin->u; n = position.v - transformation->fwd_origin->v; + if (transformation->uneg) + e = -e; + if (transformation->vneg) + n = -n; } else { /* inverse */ cb = transformation->inv_c; c = cb + sz; e = position.u - transformation->inv_origin->u; n = position.v - transformation->inv_origin->v; + if (transformation->uneg) + e = -e; + if (transformation->vneg) + n = -n; } if ((fabs(n) > range) || (fabs(e) > range)) { @@ -454,6 +464,10 @@ PJ *PROJECTION(horner) { P->opaque = (void *) Q; if (complex_horner) { + /* Westings and/or southings? */ + Q->uneg = pj_param_exists (P->params, "uneg") ? 1 : 0; + Q->vneg = pj_param_exists (P->params, "vneg") ? 1 : 0; + n = 2*degree + 2; if (0==parse_coefs (P, Q->fwd_c, "fwd_c", n)) return horner_freeup (P, PJD_ERR_MISSING_ARGS); @@ -461,8 +475,8 @@ PJ *PROJECTION(horner) { return horner_freeup (P, PJD_ERR_MISSING_ARGS); P->fwd4d = complex_horner_forward_4d; P->inv4d = complex_horner_reverse_4d; - } + else { n = horner_number_of_coefficients (degree); if (0==parse_coefs (P, Q->fwd_u, "fwd_u", n)) -- cgit v1.2.3 From a4537c65ea6acea0b01af78c4c01e952d5d2e23c Mon Sep 17 00:00:00 2001 From: Thomas Knudsen Date: Thu, 22 Feb 2018 09:02:41 +0100 Subject: gie expect: ignore unspecd dims, improve reporting/builtins In pipelines including a Helmert shift, we need to run the test through the 4D transformation interface, even though the input coordinate system may be 2D. This can be enforced by appending "0 0" to the 2D coordinate pair in the accept instruction, which is sufficiently recognizable to be considered an idiom for selecting 4D. On return, however, (i.e. in the expect instruction), the last dimensions will contain garbage, and this garbage will be compared with "0 0" when computing the deviation. This obviously leads to nonsensical results, which this commit repairs by zeroing all dimensions *not given* in expect, before computing the deviation. Additionally, the test tolerance for geo/cartesian roundtrip precision has been relaxed from picometer to nanometer level. These tests have shown to intermittently bomb, and as the pm level tolerance is probably a leftover from when deviation was computed in degrees, not meter (and hence a factor of 111000 more tight than intended at its introduction) relaxing it by a factor of 1000 makes ample sense. Also, two new features, introduced while debugging this case has been left in the code: - improved reporting, for verbosity levels higher than 2 - a "skip" instruction, forcing all remaining work to be skipped (i.e. run until something strange happens - then stop to handle debugging, while avoiding additional garbage) --- src/gie.c | 133 +++++++++++++++++++++++++++++++++++++++++++++++++------------- 1 file changed, 105 insertions(+), 28 deletions(-) (limited to 'src') diff --git a/src/gie.c b/src/gie.c index 577ad7cf..d9323bac 100644 --- a/src/gie.c +++ b/src/gie.c @@ -146,7 +146,7 @@ static ffio *ffio_create (const char **tags, size_t n_tags, size_t max_record_si static const char *gie_tags[] = { "", "operation", "accept", "expect", "roundtrip", "banner", "verbose", - "direction", "tolerance", "ignore", "builtins", "echo", "" + "direction", "tolerance", "ignore", "builtins", "echo", "skip", "" }; static const size_t n_gie_tags = sizeof gie_tags / sizeof gie_tags[0]; @@ -177,6 +177,7 @@ typedef struct { PJ_COORD a, b, c, e; PJ_DIRECTION dir; int verbosity; + int skip; int op_id; int op_ok, op_ko, op_skip; int total_ok, total_ko, total_skip; @@ -192,6 +193,8 @@ typedef struct { ffio *F = 0; static gie_ctx T; +int tests=0, succs=0, succ_fails=0, fail_fails=0, succ_rtps=0, fail_rtps=0; +int succ_builtins=0, fail_builtins=0; static const char delim[] = {"-------------------------------------------------------------------------------\n"}; @@ -299,6 +302,13 @@ int main (int argc, char **argv) { fprintf (T.fout, "%sGrand total: %d. Success: %d, Skipped: %d, Failure: %d\n", delim, T.grand_ok+T.grand_ko+T.grand_skip, T.grand_ok, T.grand_skip, T.grand_ko); fprintf (T.fout, "%s", delim); + if (T.verbosity > 1) { + fprintf (T.fout, "Failing roundtrips: %4d, Succeeding roundtrips: %4d\n", fail_rtps, succ_rtps); + fprintf (T.fout, "Failing failures: %4d, Succeeding failures: %4d\n", fail_fails, succ_fails); + fprintf (T.fout, "Failing builtins: %4d, Succeeding builtins: %4d\n", fail_builtins, succ_builtins); + fprintf (T.fout, "Internal counters: %4.4d(%4.4d)\n", tests, succs); + fprintf (T.fout, "%s", delim); + } } else if (T.grand_ko) @@ -315,21 +325,54 @@ int main (int argc, char **argv) { static int another_failure (void) { T.op_ko++; T.total_ko++; + proj_errno_reset (T.P); return 0; } static int another_skip (void) { T.op_skip++; T.total_skip++; + proj_errno_reset (T.P); return 0; } static int another_success (void) { T.op_ok++; T.total_ok++; + proj_errno_reset (T.P); return 0; } +static int another_succeeding_failure (void) { + succ_fails++; + return another_success (); +} + +static int another_failing_failure (void) { + fail_fails++; + return another_failure (); +} + +static int another_succeeding_roundtrip (void) { + succ_rtps++; + return another_success (); +} + +static int another_failing_roundtrip (void) { + fail_rtps++; + return another_failure (); +} + +static int another_succeeding_builtin (void) { + succ_builtins++; + return another_success (); +} + +static int another_failing_builtin (void) { + fail_builtins++; + return another_failure (); +} + static int process_file (const char *fname) { FILE *f; @@ -339,6 +382,9 @@ static int process_file (const char *fname) { T.op_ko = T.total_ko = 0; T.op_skip = T.total_skip = 0; + if (T.skip) + return proj_destroy (T.P), T.P = 0, 0; + f = fopen (fname, "rt"); if (0==f) { if (T.verbosity > 0) { @@ -556,27 +602,27 @@ using the "builtins" command verb. i = pj_unitconvert_selftest (); if (i!=0) { fprintf (T.fout, "pj_unitconvert_selftest fails with %d\n", i); - another_failure(); + another_failing_builtin(); } else - another_success (); + another_succeeding_builtin (); i = pj_cart_selftest (); if (i!=0) { fprintf (T.fout, "pj_cart_selftest fails with %d\n", i); - another_failure(); + another_failing_builtin(); } else - another_success (); + another_succeeding_builtin (); i = pj_horner_selftest (); if (i!=0) { fprintf (T.fout, "pj_horner_selftest fails with %d\n", i); - another_failure(); + another_failing_builtin(); } else - another_success (); + another_succeeding_builtin (); return 0; } @@ -674,7 +720,7 @@ back/forward transformation pairs. r = proj_roundtrip (T.P, T.dir, ntrips, &coo); if (r <= d) - return another_success (); + return another_succeeding_roundtrip (); if (T.verbosity > -1) { if (0==T.op_ko && T.verbosity < 2) @@ -683,7 +729,7 @@ back/forward transformation pairs. fprintf (T.fout, " FAILURE in %s(%d):\n", opt_strip_path (T.curr_file), (int) F->lineno); fprintf (T.fout, " roundtrip deviation: %.6f mm, expected: %.6f mm\n", 1000*r, 1000*d); } - return another_failure (); + return another_failing_roundtrip (); } @@ -723,7 +769,7 @@ static int expect_message_cannot_parse (const char *args) { } static int expect_failure_with_errno_message (int expected, int got) { - another_failure (); + another_failing_failure (); if (T.verbosity < 0) return 1; @@ -779,14 +825,16 @@ Tell GIE what to expect, when transforming the ACCEPTed input if (expect_failure_with_errno && proj_errno (T.P)!=expect_failure_with_errno) return expect_failure_with_errno_message (expect_failure_with_errno, proj_errno(T.P)); - return another_success (); + return another_succeeding_failure (); } /* Otherwise, it's a true failure */ banner (T.operation); - errmsg(3, "%sInvalid operation definition in line no. %d: %s\n", - delim, (int) T.operation_lineno, pj_strerrno(proj_errno(T.P))); - return another_failure (); + errmsg (3, "%sInvalid operation definition in line no. %d:\n %s (errno=%s/%d)\n", + delim, (int) T.operation_lineno, pj_strerrno(proj_errno(T.P)), + err_const_from_errno (proj_errno(T.P)), proj_errno(T.P) + ); + return another_failing_failure (); } /* We may still successfully fail even if the proj_create succeeded */ @@ -797,20 +845,24 @@ Tell GIE what to expect, when transforming the ACCEPTed input ci = proj_angular_input (T.P, T.dir)? torad_coord (T.P, T.dir, T.a): T.a; co = expect_trans_n_dim (ci); - /* Failed to fail? - that's a failure */ - if (co.xyz.x!=HUGE_VAL) - return another_failure (); - if (expect_failure_with_errno) { - printf ("errno=%d, expected=%d\n", proj_errno (T.P), expect_failure_with_errno); if (proj_errno (T.P)==expect_failure_with_errno) - return another_success (); - - return another_failure (); + return another_succeeding_failure (); + printf ("errno=%d, expected=%d\n", proj_errno (T.P), expect_failure_with_errno); + return another_failing_failure (); } - /* Yes, we failed successfully */ - return another_success (); + + /* Succeeded in failing? - that's a success */ + if (co.xyz.x==HUGE_VAL) + return another_succeeding_failure (); + + /* Failed to fail? - that's a failure */ + banner (T.operation); + errmsg (3, "%sFailed to fail. Operation definition in line no. %d\n", + delim, (int) T.operation_lineno + ); + return another_failing_failure (); } @@ -822,10 +874,12 @@ Tell GIE what to expect, when transforming the ACCEPTed input printf ("left: %d right: %d\n", T.P->left, T.P->right); } + tests++; T.e = parse_coord (args); if (HUGE_VAL==T.e.v[0]) return expect_message_cannot_parse (args); + /* expected angular values, probably in degrees */ ce = proj_angular_output (T.P, T.dir)? torad_coord (T.P, T.dir, T.e): T.e; if (T.verbosity > 3) @@ -836,8 +890,14 @@ Tell GIE what to expect, when transforming the ACCEPTed input if (T.verbosity > 3) printf ("ACCEPTS %.12f %.12f %.12f %.12f\n", ci.v[0],ci.v[1],ci.v[2],ci.v[3]); - /* angular output from proj_trans comes in radians */ + /* do the transformation, but mask off dimensions not given in expect-ation */ co = expect_trans_n_dim (ci); + if (T.dimensions_given < 4) + co.v[3] = 0; + if (T.dimensions_given < 3) + co.v[2] = 0; + + /* angular output from proj_trans comes in radians */ T.b = proj_angular_output (T.P, T.dir)? todeg_coord (T.P, T.dir, co): co; if (T.verbosity > 3) printf ("GOT %.12f %.12f %.12f %.12f\n", co.v[0],co.v[1],co.v[2],co.v[3]); @@ -856,6 +916,7 @@ Tell GIE what to expect, when transforming the ACCEPTed input if (d > T.tolerance) return expect_message (d, args); + succs++; another_success (); return 0; @@ -893,7 +954,20 @@ fprintf (T.fout, "%s\n", args); +/*****************************************************************************/ +static int skip (const char *args) { +/***************************************************************************** +Indicate that the remaining material should be skipped. Mostly for debugging. +******************************************************************************/ + T.skip = 1; + (void) args; + return 0; +} + + static int dispatch (const char *cmnd, const char *args) { + if (T.skip) + return SKIP; if (0==strcmp (cmnd, "operation")) return operation ((char *) args); if (0==strcmp (cmnd, "accept")) return accept (args); if (0==strcmp (cmnd, "expect")) return expect (args); @@ -902,9 +976,10 @@ static int dispatch (const char *cmnd, const char *args) { if (0==strcmp (cmnd, "verbose")) return verbose (args); if (0==strcmp (cmnd, "direction")) return direction (args); if (0==strcmp (cmnd, "tolerance")) return tolerance (args); - if (0==strcmp (cmnd, "ignore")) return ignore (args); + if (0==strcmp (cmnd, "ignore")) return ignore (args); if (0==strcmp (cmnd, "builtins")) return builtins (args); if (0==strcmp (cmnd, "echo")) return echo (args); + if (0==strcmp (cmnd, "skip")) return skip (args); return 0; } @@ -1223,6 +1298,8 @@ static int nextline (ffio *G) { Read next line of input file. Returns 1 on success, 0 on failure. ****************************************************************************************/ G->next_args[0] = 0; + if (T.skip) + return 0; if (0==fgets (G->next_args, (int) G->next_args_size - 1, G->f)) return 0; if (feof (G->f)) @@ -1570,7 +1647,7 @@ static int pj_cart_selftest (void) { /* Forward projection: Ellipsoidal-to-3D-Cartesian */ dist = proj_roundtrip (P, PJ_FWD, 1, &a); - if (dist > 1e-12) + if (dist > 1e-9) return 8; /* Test at the South Pole */ @@ -1582,7 +1659,7 @@ static int pj_cart_selftest (void) { /* Forward projection: Ellipsoidal-to-3D-Cartesian */ dist = proj_roundtrip (P, PJ_FWD, 1, &a); - if (dist > 1e-12) + if (dist > 1e-9) return 9; -- cgit v1.2.3 From 7e2b86e5d1598944626152c89fda0d1d9d5c1b0c Mon Sep 17 00:00:00 2001 From: Mike Toews Date: Fri, 23 Feb 2018 14:10:48 +1300 Subject: Tidy a few typos --- src/PJ_calcofi.c | 2 +- src/cct.c | 4 ++-- src/pj_mlfn.c | 2 +- 3 files changed, 4 insertions(+), 4 deletions(-) (limited to 'src') diff --git a/src/PJ_calcofi.c b/src/PJ_calcofi.c index 5c0b64ab..955e82c2 100644 --- a/src/PJ_calcofi.c +++ b/src/PJ_calcofi.c @@ -13,7 +13,7 @@ PROJ_HEAD(calcofi, /* Conversions for the California Cooperative Oceanic Fisheries Investigations Line/Station coordinate system following the algorithm of: -Eber, L.E., and R.P. Hewitt. 1979. Conversion algorithms for the CALCOFI +Eber, L.E., and R.P. Hewitt. 1979. Conversion algorithms for the CalCOFI station grid. California Cooperative Oceanic Fisheries Investigations Reports 20:135-137. (corrected for typographical errors). http://www.calcofi.org/publications/calcofireports/v20/Vol_20_Eber___Hewitt.pdf diff --git a/src/cct.c b/src/cct.c index dc68122d..25e97af9 100644 --- a/src/cct.c +++ b/src/cct.c @@ -24,7 +24,7 @@ professor of Geodesy at the University of Copenhagen, mentor and advisor for a generation of Danish geodesists, colleague and collaborator for two generations of global geodesists, Secretary General for the International Association of Geodesy, IAG (1995--2007), fellow of the -Amercan Geophysical Union (1991), recipient of the IAG Levallois Medal +American Geophysical Union (1991), recipient of the IAG Levallois Medal (2007), the European Geosciences Union Vening Meinesz Medal (2008), and of numerous other honours. @@ -34,7 +34,7 @@ the development of geodesy - both through his scientific contributions, comprising more than 250 publications, and by his mentoring and teaching of the next generations of geodesists. -As Christian was an avid Fortran programmer, and a keen Unix connoiseur, +As Christian was an avid Fortran programmer, and a keen Unix connoisseur, he would have enjoyed to know that his initials would be used to name a modest Unix style transformation filter, hinting at the tireless aspect of his personality, which was certainly one of the reasons he accomplished diff --git a/src/pj_mlfn.c b/src/pj_mlfn.c index e00f2bf1..84282d2b 100644 --- a/src/pj_mlfn.c +++ b/src/pj_mlfn.c @@ -1,5 +1,5 @@ #include -/* meridinal distance for ellipsoid and inverse +/* meridional distance for ellipsoid and inverse ** 8th degree - accurate to < 1e-5 meters when used in conjunction ** with typical major axis values. ** Inverse determines phi to EPS (1e-11) radians, about 1e-6 seconds. -- cgit v1.2.3