diff options
| -rw-r--r-- | examples/pj_obs_api_test.c | 132 | ||||
| -rw-r--r-- | src/Makefile.am | 2 | ||||
| -rw-r--r-- | src/PJ_cart.c | 336 | ||||
| -rw-r--r-- | src/lib_proj.cmake | 3 | ||||
| -rw-r--r-- | src/makefile.vc | 2 | ||||
| -rw-r--r-- | src/pj_ctx.c | 2 | ||||
| -rw-r--r-- | src/pj_list.h | 1 | ||||
| -rw-r--r-- | src/pj_malloc.c | 128 | ||||
| -rw-r--r-- | src/pj_obs_api.c | 248 | ||||
| -rw-r--r-- | src/pj_observation.c | 114 | ||||
| -rw-r--r-- | src/proj.def | 18 | ||||
| -rw-r--r-- | src/proj.h | 199 | ||||
| -rw-r--r-- | src/proj_api.h | 102 | ||||
| -rw-r--r-- | src/projects.h | 16 |
14 files changed, 1072 insertions, 231 deletions
diff --git a/examples/pj_obs_api_test.c b/examples/pj_obs_api_test.c new file mode 100644 index 00000000..49d1224e --- /dev/null +++ b/examples/pj_obs_api_test.c @@ -0,0 +1,132 @@ +/******************************************************************************* + Tiny test of an evolving new API, demonstrating simple examples of + 2D and 3D transformations. + + The main transformation setup object is PJ, well known from the two + former proj APIs (projects.h and proj_api.h) + + The main data object PJ_OBS is new, but encapsulates the older + LP, XY etc. objects in a framework for storing a 3D observation taken at + a 4D point in space-time, and including some additional metadata (e.g. + a serial number or an epsg code). + + PJ_OBS uses unions to enforce explicit statement of what kind of + coordinates are expected at a given spot in the code. + + The primary elements of the API are: + + pj_create (char *desc): + Create a new PJ object from a text description. + pj_apply (PJ *P, int direction, PJ_OBS obs): + Forward or inverse transformation of obs. + pj_error (PJ *P): + Check error status of P. + pj_free (PJ *P): + Clean up when done with the transformation object. + + These 4 functions should cover most of the common use cases. Additional + functionality for handling thread contexts are provided by the functions + pj_debug_set, pj_error_set, pj_log_set, pj_context_renew, + pj_context_inherit, pj_context_free, pj_fileapi_set. + + The proj thread contexts have not seen widespread use, so one of the + intentions with this new API is to make them less visible on the API + surface: Contexts do not have a life by themselves, they are only + visible through their associated PJs, and the number of functions + supporting them is limited. + + Compilation example: + gcc -pedantic -W -Wall -Wextra -I../src -o pj_proj_test pj_proj_test.c -L../../../build/proj.4/lib -lproj_4_9 + + Thomas Knudsen, 2016-10-30/11-03 +*******************************************************************************/ +#include <proj.h> + + +int main (void) { + PJ *p; + PJ_OBS a, b; + int err; + double dist; + + /* Log everything libproj offers to log for you */ + pj_debug_set (0, PJ_LOG_DEBUG_MINOR); + + /* An utm projection on the GRS80 ellipsoid */ + p = pj_create ("+proj=utm +zone=32 +ellps=GRS80"); + if (0==p) + return puts ("Oops"), 0; + + /* zero initialize everything, then set (longitude, latitude) to (12, 55) */ + a = pj_obs_null; + /* a.coo.lp: The coordinate part of a, interpreted as a classic LP pair */ + a.coo.lp.lam = TORAD(12); + a.coo.lp.phi = TORAD(55); + + /* Forward projection */ + b = pj_apply (p, PJ_FWD, a); + printf ("FWD: %15.9f %15.9f\n", b.coo.enh.e, b.coo.enh.n); + + /* Inverse projection */ + a = pj_apply (p, PJ_INV, b); + printf ("INV: %15.9f %15.9f\n", TODEG(a.coo.lpz.lam), TODEG(a.coo.lpz.phi)); + + /* Null projection */ + a = pj_apply (p, PJ_IDENT, a); + printf ("IDENT: %15.9f %15.9f\n", TODEG(a.coo.lpz.lam), TODEG(a.coo.lpz.phi)); + + /* Forward again, to get two linear items for comparison */ + a = pj_apply (p, PJ_FWD, a); + printf ("FWD: %15.9f %15.9f\n", b.coo.enh.e, b.coo.enh.n); + + dist = pj_obs_dist_2d (a, b); + printf ("Roundtrip deviation, (nm): %15.9f\n", dist*1e9); + + /* Invalid projection */ + a = pj_apply (p, 42, a); + if (a.coo.lpz.lam!=DBL_MAX) + printf ("%15.9f %15.9f\n", a.coo.lpz.lam, a.coo.lpz.phi); + err = pj_error (p); + printf ("pj_error: %d\n", err); + + /* Clean up */ + pj_free (p); + + /* Now do some 3D transformations */ + p = pj_create ("+proj=cart +ellps=GRS80"); + if (0==p) + return puts ("Oops"), 0; + + /* zero initialize everything, then set (longitude, latitude, height) to (12, 55, 100) */ + a = b = pj_obs_null; + a.coo.lpz.lam = TORAD(12); + a.coo.lpz.phi = TORAD(55); + a.coo.lpz.z = 100; + + /* Forward projection: 3D-Cartesian-to-Ellipsoidal */ + b = pj_apply (p, PJ_FWD, a); + printf ("FWD: %15.9f %15.9f %15.9f\n", b.coo.xyz.x, b.coo.xyz.y, b.coo.xyz.z); + + /* Check roundtrip precision for 10000 iterations each way */ + dist = pj_roundtrip (p, PJ_FWD, 10000, a); + printf ("Roundtrip deviation, fwd (nm): %15.9f\n", dist*1e9*1e5); + dist = pj_roundtrip (p, PJ_INV, 10000, b); + printf ("Roundtrip deviation, inv (nm): %15.9f\n", dist*1e9); + + /* Inverse projection: Ellipsoidal-to-3D-Cartesian */ + b = pj_apply (p, PJ_INV, b); + printf ("INV: %15.9f %15.9f %15.9f\n", TODEG(b.coo.lpz.lam), TODEG(b.coo.lpz.phi), b.coo.lpz.z); + + /* Move p to another context */ + pj_context_renew (p); + b = pj_apply (p, PJ_FWD, b); + printf ("CTX1: %15.9f %15.9f %15.9f\n", b.coo.xyz.x, b.coo.xyz.y, b.coo.xyz.z); + + /* Move it back to the default context */ + pj_context_free (p); + b = pj_apply (p, PJ_INV, b); + printf ("CTX0: %15.9f %15.9f %15.9f\n", TODEG(b.coo.lpz.lam), TODEG(b.coo.lpz.phi), b.coo.lpz.z); + + pj_free (p); + return 0; +} diff --git a/src/Makefile.am b/src/Makefile.am index 65944111..1d53375a 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -74,7 +74,7 @@ libproj_la_SOURCES = \ jniproj.c pj_mutex.c pj_initcache.c pj_apply_vgridshift.c geodesic.c \ pj_strtod.c \ \ - pj_observation.c + pj_obs_api.c PJ_cart.c install-exec-local: rm -f $(DESTDIR)$(bindir)/invproj$(EXEEXT) diff --git a/src/PJ_cart.c b/src/PJ_cart.c new file mode 100644 index 00000000..5df95c31 --- /dev/null +++ b/src/PJ_cart.c @@ -0,0 +1,336 @@ +/****************************************************************************** + * Project: PROJ.4 + * Purpose: Convert between ellipsoidal, geodetic coordinates and + * cartesian, geocentric coordinates. + * + * Formally, this functionality is also found in the PJ_geocent.c + * code. + * + * Actually, however, the PJ_geocent transformations are carried + * out in concert between 2D stubs in PJ_geocent.c and 3D code + * placed in pj_transform.c. + * + * For pipeline-style datum shifts, we do need direct access + * to the full 3D interface for this functionality. + * + * Hence this code, which may look like "just another PJ_geocent" + * but really is something substantially different. + * + * Author: Thomas Knudsen, thokn@sdfe.dk + * + ****************************************************************************** + * Copyright (c) 2016, Thomas Knudsen / SDFE + * + * 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. + *****************************************************************************/ + +#define PJ_LIB__ +#include <proj.h> +#include <projects.h> +#include <assert.h> +#include <stddef.h> +#include <math.h> +#include <errno.h> +PROJ_HEAD(cart, "Geodetic/cartesian conversions"); + + +/************************************************************** + CARTESIAN / GEODETIC CONVERSIONS +*************************************************************** + This material follows: + + Bernhard Hofmann-Wellenhof & Helmut Moritz: + Physical Geodesy, 2nd edition. + Springer, 2005. + + chapter 5.6: Coordinate transformations + (HM, below), + + and + + Wikipedia: Geographic Coordinate Conversion, + https://en.m.wikipedia.org/wiki/Geographic_coordinate_conversion + + (WP, below). + + The cartesian-to-geodetic conversion is based on Bowring's + celebrated method: + + B. R. Bowring: + Transformation from spatial to geographical coordinates + Survey Review 23(181), pp. 323-327, 1976 + + (BB, below), + + but could probably use some TLC from a newer and faster + algorithm: + + Toshio Fukushima: + Transformation from Cartesian to Geodetic Coordinates + Accelerated by Halley’s Method + Journal of Geodesy, February 2006 + + (TF, below). + + + These routines are probably not as robust at those in + geocent.c, at least thay haven't been through as heavy + use as their geocent sisters. Some care has been taken + to avoid singularities, but extreme cases (e.g. setting + es, the squared eccentricity, to 1), will cause havoc. + +**************************************************************/ + + +static void freeup (PJ *P) { + pj_freeup_plain (P); + return; +} + + +/*********************************************************************/ +static double semiminor_axis (double a, double es) { +/*********************************************************************/ + return a * sqrt (1 - es); +} + +/*********************************************************************/ +static double normal_radius_of_curvature (double a, double es, double phi) { +/*********************************************************************/ + double s = sin(phi); + if (es==0) + return a; + /* This is from WP. HM formula 2-149 gives an a,b version */ + return a / sqrt (1 - es*s*s); +} + + +/*********************************************************************/ +static double second_eccentricity_squared (double a, double es) { +/********************************************************************** + Follows the definition, e'2 = (a2-b2)/b2, but uses the identity + (a2-b2) = (a-b)(a+b), for improved numerical precision. +***********************************************************************/ + double b = semiminor_axis (a, es); + return (a - b) * (a + b) / (b*b); +} + + +/*********************************************************************/ +static XYZ cartesian (LPZ geodetic, PJ *P) { +/*********************************************************************/ + double N, h, lam, phi, cosphi = cos(geodetic.phi); + XYZ xyz; + + N = normal_radius_of_curvature(P->a, P->es, geodetic.phi); + phi = geodetic.phi; + lam = geodetic.lam; + h = geodetic.z; + + /* HM formula 5-27 (z formula follows WP) */ + xyz.x = (N + h) * cosphi * cos(lam); + xyz.y = (N + h) * cosphi * sin(lam); + xyz.z = (N * (1 - P->es) + h) * sin(phi); + + /*********************************************************************/ + /* */ + /* For historical reasons, in proj, plane coordinates are measured */ + /* in units of the semimajor axis. Since 3D handling is grafted on */ + /* later, this is not the case for heights. And even though this */ + /* coordinate really is 3D cartesian, the z-part looks like a height */ + /* to proj. Hence, we have the somewhat unusual situation of having */ + /* a point coordinate with differing units between dimensions. */ + /* */ + /* The scaling and descaling is handled by the pj_fwd/inv functions. */ + /* */ + /*********************************************************************/ + xyz.x /= P->a; + xyz.y /= P->a; + + return xyz; +} + + +/*********************************************************************/ +static LPZ geodetic (XYZ cartesian, PJ *P) { +/*********************************************************************/ + double N, b, p, theta, c, s, e2s; + LPZ lpz; + + cartesian.x *= P->a; + cartesian.y *= P->a; + + /* Perpendicular distance from point to Z-axis (HM eq. 5-28) */ + p = hypot (cartesian.x, cartesian.y); + + /* Ancillary ellipsoidal parameters */ + b = semiminor_axis (P->a, P->es); + e2s = second_eccentricity_squared (P->a, P->es); + + /* HM eq. (5-37) */ + theta = atan2 (cartesian.z * P->a, p * b); + + /* HM eq. (5-36) (from BB, 1976) */ + c = cos(theta); + s = sin(theta); + lpz.phi = atan2 (cartesian.z + e2s*b*s*s*s, p - P->es*P->a*c*c*c); + lpz.lam = atan2 (cartesian.y, cartesian.x); + N = normal_radius_of_curvature (P->a, P->es, lpz.phi); + + + c = cos(lpz.phi); + if (fabs(c) < 1e-6) { + /* poleward of 89.99994 deg, we avoid division by zero */ + /* by computing the height as the cartesian z value */ + /* minus the semiminor axis length */ + lpz.z = cartesian.z - (cartesian.z > 0? b: -b); + } + else + lpz.z = p / c - N; + + return lpz; +} + + + +/* In effect, 2 cartesian coordinates of a point on the ellipsoid. Rather pointless, but... */ +static XY cart_forward (LP lp, PJ *P) { + PJ_TRIPLET point; + point.lp = lp; + point.lpz.z = 0; + + point.xyz = cartesian (point.lpz, P); + return point.xy; +} + +/* And the other way round. Still rather pointless, but... */ +static LP cart_reverse (XY xy, PJ *P) { + PJ_TRIPLET point; + point.xy = xy; + point.xyz.z = 0; + + point.lpz = geodetic (point.xyz, P); + return point.lp; +} + + + +/*********************************************************************/ +PJ *PROJECTION(cart) { +/*********************************************************************/ + P->fwd3d = cartesian; + P->inv3d = geodetic; + P->fwd = cart_forward; + P->inv = cart_reverse; + return P; +} + +#ifndef PJ_SELFTEST +/* selftest stub */ +int pj_cart_selftest (void) {return 0;} +#else +/* Testing quite a bit of the pj_obs_api as a side effect (inspired by pj_obs_api_test.c) */ +int pj_cart_selftest (void) { + PJ *p; + PJ_OBS a, b; + int err; + double dist; + + /* Log everything libproj offers to log for you */ + pj_debug_set (0, PJ_LOG_DEBUG_MINOR); + + /* An utm projection on the GRS80 ellipsoid */ + p = pj_create ("+proj=utm +zone=32 +ellps=GRS80"); + if (0==p) + return 1; + + /* Turn off logging */ + pj_debug_set (0, PJ_LOG_NONE); + + /* zero initialize everything, then set (longitude, latitude) to (12, 55) */ + a = pj_obs_null; + /* a.coo.lp: The coordinate part of a, interpreted as a classic LP pair */ + a.coo.lp.lam = TORAD(12); + a.coo.lp.phi = TORAD(55); + + /* Forward projection */ + b = pj_apply (p, PJ_FWD, a); + + /* Inverse projection */ + a = pj_apply (p, PJ_INV, b); + + /* Null projection */ + a = pj_apply (p, PJ_IDENT, a); + + /* Forward again, to get two linear items for comparison */ + a = pj_apply (p, PJ_FWD, a); + + dist = pj_obs_dist_2d (a, b); + if (dist > 2e-9) + return 2; + + /* Invalid projection */ + a = pj_apply (p, 42, a); + if (a.coo.lpz.lam!=DBL_MAX) + return 3; + err = pj_error (p); + if (0==err) + return 4; + + /* Clear error */ + pj_error_set (p, 0); + + /* Clean up */ + pj_free (p); + + /* Now do some 3D transformations */ + p = pj_create ("+proj=cart +ellps=GRS80"); + if (0==p) + return 5; + + /* zero initialize everything, then set (longitude, latitude, height) to (12, 55, 100) */ + a = b = pj_obs_null; + a.coo.lpz.lam = TORAD(12); + a.coo.lpz.phi = TORAD(55); + a.coo.lpz.z = 100; + + /* Forward projection: 3D-Cartesian-to-Ellipsoidal */ + b = pj_apply (p, PJ_FWD, a); + + /* Check roundtrip precision for 10000 iterations each way */ + dist = pj_roundtrip (p, PJ_FWD, 10000, a); + dist = pj_roundtrip (p, PJ_INV, 10000, b); + if (dist > 2e-9) + return 6; + + /* Inverse projection: Ellipsoidal-to-3D-Cartesian */ + b = pj_apply (p, PJ_INV, b); + + /* Move p to another context */ + pj_context_renew (p); + b = pj_apply (p, PJ_FWD, b); + + /* Move it back to the default context */ + pj_context_free (p); + b = pj_apply (p, PJ_INV, b); + + pj_free (p); + return 0; +} +#endif diff --git a/src/lib_proj.cmake b/src/lib_proj.cmake index 54e59ae7..2386a102 100644 --- a/src/lib_proj.cmake +++ b/src/lib_proj.cmake @@ -50,6 +50,7 @@ SET(SRC_LIBPROJ_PJ PJ_boggs.c PJ_bonne.c PJ_calcofi.c + PJ_cart.c PJ_cass.c PJ_cc.c PJ_cea.c @@ -193,7 +194,7 @@ SET(SRC_LIBPROJ_CORE pj_mlfn.c pj_msfn.c pj_mutex.c - pj_observation.c + pj_obs_api.c pj_open_lib.c pj_param.c pj_phi2.c diff --git a/src/makefile.vc b/src/makefile.vc index f06a769a..51b4320b 100644 --- a/src/makefile.vc +++ b/src/makefile.vc @@ -59,7 +59,7 @@ support = \ pj_strtod.obj pj_run_selftests.obj pj_generic_selftest.obj pipeline = \ - pj_observation.obj + pj_obs_api.obj PJ_cart.obj geodesic = geodesic.obj diff --git a/src/pj_ctx.c b/src/pj_ctx.c index 43d12255..9ae3d4e0 100644 --- a/src/pj_ctx.c +++ b/src/pj_ctx.c @@ -192,5 +192,3 @@ projFileAPI *pj_ctx_get_fileapi( projCtx ctx ) { return ctx->fileapi; } - - diff --git a/src/pj_list.h b/src/pj_list.h index a6806556..3569e335 100644 --- a/src/pj_list.h +++ b/src/pj_list.h @@ -17,6 +17,7 @@ PROJ_HEAD(bipc, "Bipolar conic of western hemisphere") PROJ_HEAD(boggs, "Boggs Eumorphic") PROJ_HEAD(bonne, "Bonne (Werner lat_1=90)") PROJ_HEAD(calcofi, "Cal Coop Ocean Fish Invest Lines/Stations") +PROJ_HEAD(cart, "Geodetic/cartesian conversions") PROJ_HEAD(cass, "Cassini") PROJ_HEAD(cc, "Central Cylindrical") PROJ_HEAD(cea, "Equal Area Cylindrical") diff --git a/src/pj_malloc.c b/src/pj_malloc.c index aab69e99..b9d8f825 100644 --- a/src/pj_malloc.c +++ b/src/pj_malloc.c @@ -1,28 +1,68 @@ +/****************************************************************************** + * Project: PROJ.4 + * Purpose: Memory management for proj.4. + * This version includes an implementation of generic destructors, + * for memory deallocation for the large majority of PJ-objects + * that do not allocate anything alse than the PJ-object itself, + * and its associated opaque object - i.e. no additional malloc'ed + * memory inside the opaque object. + * + * Author: Gerald I. Evenden (Original proj.4 author), + * Frank Warmerdam (2000) pj_malloc? + * Thomas Knudsen (2016) - freeup/dealloc parts + * + ****************************************************************************** + * Copyright (c) 2000, Frank Warmerdam + * Copyright (c) 2016, Thomas Knudsen / SDFE + * + * 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. + *****************************************************************************/ + + /* allocate and deallocate memory */ /* These routines are used so that applications can readily replace ** projection system memory allocation/deallocation call with custom ** application procedures. */ + #include <projects.h> #include <errno.h> - void * -pj_malloc(size_t size) { -/* -/ Currently, pj_malloc is a hack to solve an errno problem. -/ The problem is described in more details at -/ https://bugzilla.redhat.com/bugzilla/show_bug.cgi?id=86420. -/ It seems, that pj_init and similar functions incorrectly -/ (under debian/glibs-2.3.2) assume that pj_malloc resets -/ errno after success. pj_malloc tries to mimic this. -*/ +/**********************************************************************/ +void *pj_malloc(size_t size) { +/*********************************************************************** +Currently, pj_malloc is a hack to solve an errno problem. +The problem is described in more details at +https://bugzilla.redhat.com/bugzilla/show_bug.cgi?id=86420. +It seems, that pj_init and similar functions incorrectly +(under debian/glibs-2.3.2) assume that pj_malloc resets +errno after success. pj_malloc tries to mimic this. +***********************************************************************/ int old_errno = errno; void *res = malloc(size); if ( res && !old_errno ) errno = 0; return res; } - void -pj_dalloc(void *ptr) { + +/**********************************************************************/ +void pj_dalloc(void *ptr) { +/**********************************************************************/ free(ptr); } @@ -70,3 +110,67 @@ pointer" to signal an error in a multi level allocation: pj_dalloc (ptr); return 0; } + + + + +/*****************************************************************************/ +void *pj_freeup_msg_plain (PJ *P, int errlev) { /* Destructor */ +/***************************************************************************** + Does memory deallocation for "plain" PJ objects, i.e. that vast majority + of PJs where the opaque object does not contain any additionally + allocated memory. +******************************************************************************/ + if (0==P) + return 0; + + if (0!=errlev) + pj_ctx_set_errno (P->ctx, errlev); + + if (0==P->opaque) + return pj_dealloc (P); + + pj_dealloc (P->opaque); + return pj_dealloc(P); +} + + +/*****************************************************************************/ +void pj_freeup_plain (PJ *P) { +/***************************************************************************** + Adapts pj_freeup_msg_plain to the format defined for the callback in + the PJ object. + + i.e. reduces most instances of projection deallocation code to: + + static void freeup (PJ *P) { + pj_freeup_plain (P); + return; + } + + rather than: + + static void *freeup_msg_add (PJ *P, int errlev) { + if (0==P) + return 0; + pj_ctx_set_errno (P->ctx, errlev); + + if (0==P->opaque) + return pj_dealloc (P); + + (* projection specific deallocation goes here *) + + pj_dealloc (P->opaque); + return pj_dealloc(P); + } + + (* Adapts pipeline_freeup to the format defined for the PJ object *) + static void freeup_msg_add (PJ *P) { + freeup_new_add (P, 0); + return; + } + + ******************************************************************************/ + pj_freeup_msg_plain (P, 0); + return; + } diff --git a/src/pj_obs_api.c b/src/pj_obs_api.c new file mode 100644 index 00000000..ccec8a64 --- /dev/null +++ b/src/pj_obs_api.c @@ -0,0 +1,248 @@ +/****************************************************************************** + * Project: PROJ.4 + * Purpose: Implement a (currently minimalistic) proj API based primarily + * on the PJ_OBS generic geodetic data type. + * + * proj thread contexts have not seen widespread use, so one of the + * intentions with this new API is to make them less visible on the + * API surface: Contexts do not have a life by themselves, they are + * visible only through their associated PJs, and the number of + * functions supporting them is limited. + * + * Author: Thomas Knudsen, thokn@sdfe.dk, 2016-06-09/2016-11-06 + * + ****************************************************************************** + * Copyright (c) 2016, Thomas Knudsen/SDFE + * + * 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. + *****************************************************************************/ +#define PJ_OBS_C +#include <proj.h> +#include <projects.h> +#include <float.h> +#include <math.h> + + +/* Used as return value in case of errors */ +const PJ_OBS pj_obs_error = { + /* Cannot use HUGE_VAL here: MSVC misimplements HUGE_VAL as something that is not compile time constant */ + {{DBL_MAX,DBL_MAX,DBL_MAX,DBL_MAX}}, + {{DBL_MAX,DBL_MAX,DBL_MAX}}, + 0, 0 +}; + +/* Used for zero-initializing new objects */ +const PJ_OBS pj_obs_null = { + {{0, 0, 0, 0}}, + {{0, 0, 0}}, + 0, 0 +}; + +/* Magic object signaling proj system shutdown mode to routines taking a PJ * arg */ +const PJ *pj_shutdown = (PJ *) &pj_shutdown; + +/* Euclidean distance between two 2D coordinates stored in PJ_OBSs */ +double pj_obs_dist_2d (PJ_OBS a, PJ_OBS b) { + double *A = a.coo.v, *B = b.coo.v; + return hypot (A[0] - B[0], A[1] - B[1]); +} + +/* Euclidean distance between two 3D coordinates stored in PJ_OBSs */ +double pj_obs_dist_3d (PJ_OBS a, PJ_OBS b) { + double *A = a.coo.v, *B = b.coo.v; + return hypot (hypot (A[0] - B[0], A[1] - B[1]), A[2] - B[2]); +} + + +PJ_OBS pj_fwdobs (PJ_OBS obs, PJ *P) { + if (0!=P->fwd3d) { + obs.coo.xyz = pj_fwd3d (obs.coo.lpz, P); + return obs; + } + if (0!=P->fwd) { + obs.coo.xy = pj_fwd (obs.coo.lp, P); + return obs; + } + pj_error_set (P, EINVAL); + return pj_obs_error; +} + + +PJ_OBS pj_invobs (PJ_OBS obs, PJ *P) { + if (0!=P->inv3d) { + obs.coo.lpz = pj_inv3d (obs.coo.xyz, P); + return obs; + } + if (0!=P->inv) { + obs.coo.lp = pj_inv (obs.coo.xy, P); + return obs; + } + pj_error_set (P, EINVAL); + return pj_obs_error; +} + + +/* Apply the transformation P to the observation obs */ +PJ_OBS pj_apply (PJ *P, enum pj_direction direction, PJ_OBS obs) { + if (0==P) + return obs; + + switch (direction) { + case PJ_FWD: + return pj_fwdobs (obs, P); + case PJ_INV: + return pj_invobs (obs, P); + case PJ_IDENT: + return obs; + default: + break; + } + + pj_error_set (P, EINVAL); + return pj_obs_error; +} + + +/* Measure numerical deviation after n roundtrips fwd-inv (or inv-fwd) */ +double pj_roundtrip (PJ *P, enum pj_direction direction, int n, PJ_OBS obs) { + int i; + PJ_OBS o, u; + + if (0==P) + return HUGE_VAL; + + if (n < 1) { + pj_error_set (P, EINVAL); + return HUGE_VAL; + } + + o.coo = obs.coo; + + for (i = 0; i < n; i++) { + switch (direction) { + case PJ_FWD: + u.coo.xyz = pj_fwd3d (o.coo.lpz, P); + o.coo.lpz = pj_inv3d (u.coo.xyz, P); + break; + case PJ_INV: + u.coo.lpz = pj_inv3d (o.coo.xyz, P); + o.coo.xyz = pj_fwd3d (u.coo.lpz, P); + break; + default: + pj_error_set (P, EINVAL); + return HUGE_VAL; + } + } + + return pj_obs_dist_3d (o, obs); +} + + +PJ *pj_create (const char *definition) { + return pj_init_plus (definition); +} + + +PJ *pj_create_argv (int argc, char **argv) { + return pj_init (argc, argv); +} + + +/* From here: Minimum viable support for contexts. The first four functions */ +/* relate to error reporting, debugging, and logging, hence being generically */ +/* useful. The remaining is a compact implementation of the more low level */ +/* proj_api.h thread contexts, which may or may not be useful */ + +int pj_error (PJ *P) { + return pj_ctx_get_errno (pj_get_ctx(P)); +} + + +void pj_error_set (PJ *P, int err) { + pj_ctx_set_errno (pj_get_ctx(P), err); +} + + +/* Set debug level 0-3. Higher number means more debug info. 0 turns it off */ +void pj_debug_set (PJ *P, enum pj_debug_level debuglevel) { + PJ_CONTEXT *ctx; + if (0==P) + ctx = pj_get_default_ctx(); + else + ctx = pj_get_ctx (P); + ctx->debug_level = debuglevel; +} + + +/* Put a new logging function into P's context. The opaque object app_data is passed as first arg at each call to the logger */ +void pj_log_set (PJ *P, void *app_data, void (*log)(void *, int, const char *)) { + PJ_CONTEXT *ctx = pj_get_ctx (P); + ctx->app_data = app_data; + if (0!=log) + ctx->logger = log; +} + + +/* Move P to a new context - initially a copy of the default context */ +int pj_context_renew (PJ *P) { + PJ_CONTEXT *ctx = pj_ctx_alloc (); + if (0==ctx) { + pj_error_set (P, ENOMEM); + return 1; + } + + pj_set_ctx (P, ctx); + return 0; +} + +/* Move daughter to mother's context */ +void pj_context_inherit (PJ *mother, PJ *daughter) { + pj_set_ctx (daughter, pj_get_ctx (mother)); +} + + +void pj_context_free (const PJ *P) { + PJ_CONTEXT *ctx; + + if (0==P) + return; + + /* During shutdown it should be OK to free the default context - but this needs more work */ + if (pj_shutdown==P) { + /* The obvious solution "pj_ctx_free (pj_get_default_ctx ())" fails */ + return; + } + + ctx = pj_get_ctx ((PJ *) P); + + /* Otherwise, trying to free the default context is a no-op */ + if (pj_get_default_ctx ()==ctx) + return; + + /* The common (?) case: free the context and move the PJ to the default context */ + pj_ctx_free (ctx); + pj_set_ctx ((PJ *) P, pj_get_default_ctx ()); +} + + +/* Minimum support for fileapi - which may never have been used... */ +void pj_fileapi_set (PJ *P, void *fileapi) { + PJ_CONTEXT *ctx = pj_get_ctx (P); + ctx->fileapi = (projFileAPI *) fileapi; +} diff --git a/src/pj_observation.c b/src/pj_observation.c deleted file mode 100644 index 0d1f7acf..00000000 --- a/src/pj_observation.c +++ /dev/null @@ -1,114 +0,0 @@ -/****************************************************************************** - * Project: PROJ.4 - * Purpose: Implement some service routines for the PJ_OBSERVATION generic - * geodetic data type - * Author: Thomas Knudsen, thokn@sdfe.dk, 2016-06-09/2016-10-11 - * - ****************************************************************************** - * Copyright (c) 2016, Thomas Knudsen/SDFE - * - * 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 <proj.h> -#include <projects.h> - - -/* Constructor for the OBSERVATION object */ -PJ_OBSERVATION pj_observation ( - double x, double y, double z, double t, - double o, double p, double k, - int id, unsigned int flags -) { - PJ_OBSERVATION g; - g.coo.xyzt.x = x; - g.coo.xyzt.y = y; - g.coo.xyzt.z = z; - g.coo.xyzt.t = t; - g.anc.opk.o = o; - g.anc.opk.p = p; - g.anc.opk.k = k; - g.id = id; - g.flags = flags; - return g; -} - - -PJ_OBSERVATION pj_apply (PJ *P, enum pj_direction direction, PJ_OBSERVATION obs) { - - switch (direction) { - case 1: - obs.coo.xyz = pj_fwd3d (obs.coo.lpz, P); - return obs; - case -1: - obs.coo.lpz = pj_inv3d (obs.coo.xyz, P); - return obs; - case 0: - return obs; - default: - pj_ctx_set_errno (P->ctx, EINVAL); - } - - return pj_observation (HUGE_VAL,HUGE_VAL,HUGE_VAL,HUGE_VAL,HUGE_VAL,HUGE_VAL,HUGE_VAL,0,0); -} - - -/* initial attempt, following a suggestion by Kristian Evers */ -double pj_roundtrip(PJ *P, enum pj_direction direction, int n, PJ_OBSERVATION obs) { - double d; - PJ_OBSERVATION o; - - /* multiple roundtrips not supported yet */ - if (n > 1) { - pj_ctx_set_errno (P->ctx, EINVAL); - return HUGE_VAL; - } - - switch (direction) { - case 1: - o.coo.xyz = pj_fwd3d (obs.coo.lpz, P); - break; - case -1: - o.coo.lpz = pj_inv3d (obs.coo.xyz, P); - break; - default: - pj_ctx_set_errno (P->ctx, EINVAL); - return HUGE_VAL; - } - - o.coo.xyz = pj_fwd3d (obs.coo.lpz, P); - o.coo.lpz = pj_inv3d (o.coo.xyz, P); - - /* Should add "both ways" here */ - d = hypot (hypot (o.coo.v[0] - obs.coo.v[0], o.coo.v[1] - obs.coo.v[1]), o.coo.v[2] - obs.coo.v[2]); - return d; -} - - -int pj_show_triplet (FILE *stream, const char *banner, PJ_TRIPLET point) { - int i = 0; - if (banner) - i += fprintf (stream, "%s", banner); - - i += fprintf(stream, "%16.10f ", point.xyz.x); - i += fprintf(stream, "%16.10f ", point.xyz.y); - i += fprintf(stream, "%16.4f", point.xyz.z); - if (banner) - i += fprintf(stream, "\n"); - return i; -} diff --git a/src/proj.def b/src/proj.def index dbf27f96..09b32c32 100644 --- a/src/proj.def +++ b/src/proj.def @@ -89,7 +89,17 @@ EXPORTS geod_polygon_testpoint @87 geod_polygon_clear @88 pj_run_selftests @89 - pj_observation @90 - pj_show_triplet @91 - pj_apply @92 - + pj_error @90 + pj_create @91 + pj_create_argv @92 + pj_apply @93 + pj_fwdobs @94 + pj_invobs @95 + pj_roundtrip @96 + pj_debug_set @97 + pj_error_set @98 + pj_log_set @99 + pj_context_renew @100 + pj_context_inherit @101 + pj_context_free @102 + pj_fileapi_set @103 @@ -1,23 +1,88 @@ /****************************************************************************** * Project: PROJ.4 - * Purpose: Revised, experimental, "bare essentials" API for PROJ.4. - * Intended as the foundation for added geodetic functionality. + * Purpose: Revised, experimental API for PROJ.4, intended as the foundation + * for added geodetic functionality. * - * Introduces the OBSERVATION data type, for generic coordinate - * and ancillary data handling. + * The original proj API (defined in projects.h) has grown organically + * over the years, but it has also grown somewhat messy. * - * Also introduces the PJ_SPATIOTEMPORAL and PJ_TRIPLET unions - * making it possible to make explicit the previously used - * "implicit type punning", where a XY is turned into a LP by - * re#defining both as UV, behind the back of the user. + * The same has happened with the newer high level API (defined in + * proj_api.h): To support various historical objectives, proj_api.h + * contains a rather complex combination of conditional defines and + * typedefs. Probably for good (historical) reasons, which are not + * always evident from today's perspective. + * + * This is an evolving attempt at creating a re-rationalized API + * with primary design goals focused on sanitizing the namespaces. + * Hence, all symbols exposed are being moved to the pj_ namespace, + * while all data types are being moved to the PJ_ namespace. + * + * Please note that this API is *orthogonal* to the previous APIs: + * Apart from some inclusion guards, projects.h and proj_api.h are not + * touched - if you do not include proj.h, the projects and proj_api + * APIs should work as they always have. + * + * A few implementation details: + * + * Apart from the namespacing efforts, I'm trying to eliminate three + * proj_api elements, which I have found especially confusing. + * + * FIRST and foremost, I try to avoid typedef'ing away pointer + * semantics. I agree that it can be occasionally useful, but I + * prefer having the pointer nature of function arguments being + * explicitly visible. + * + * Hence, projCtx has been replaced by PJ_CONTEXT *. + * and projPJ has been replaced by PJ * + * + * SECOND, I try to eliminate cases of information hiding implemented + * by redefining data types to void pointers. + * + * I prefer using a combination of forward declarations and typedefs. + * Hence: + * typedef void *projCtx; + * Has been replaced by: + * struct projCtx_t; + * typedef struct projCtx_t PJ_CONTEXT; + * This makes it possible for the calling program to know that the + * PJ_CONTEXT data type exists, and handle pointers to that data type + * without having any idea about its internals. + * + * (obviously, in this example, struct projCtx_t should also be + * renamed struct pj_ctx some day...) + * + * THIRD, I try to eliminate implicit type punning. Hence this API + * introduces the PJ_OBS ("observation") data type, for generic + * coordinate and handling of ancillary data. + * + * It includes the PJ_COORD and PJ_TRIPLET unions making it possible + * to make explicit the previously used "implicit type punning", where + * a XY is turned into a LP by re#defining both as UV, behind the back + * of the user. + * + * The PJ_COORD union is used for storing 1D, 2D, 3D and 4D coordinates. + * The PJ_TRIPLET union is used for storing any set of up to 3 related + * observations. At the application code level, the names of these + * unions will usually not be used - they will only be accessed via + * their tag names in the PJ_OBS data type. * * The bare essentials API presented here follows the PROJ.4 * convention of sailing the coordinate to be reprojected, up on * the stack ("call by value"), and symmetrically returning the - * result on the stack. Although the OBSERVATION object is 4 times + * result on the stack. Although the PJ_OBS object is 4 times * as large as the traditional XY and LP objects, timing results * have shown the overhead to be very reasonable. * + * In its current incarnation, the API is focused on compactness: + * Currently it consists of only nine functions. + * + * Hence, due to the proj_ctx subsystem being little used, it has + * been replaced by a tiny set of 3 functions, making it possible, + * but not very convenient, to do almost everything possible with + * the original ctx API. + * + * See pj_proj_test.c for an example of how to use the API. + * * Author: Thomas Knudsen, <thokn@sdfe.dk> * ****************************************************************************** @@ -35,7 +100,7 @@ * * 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 + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO COORD 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 @@ -46,11 +111,21 @@ #endif #include <math.h> +#include <float.h> +#include <limits.h> + #include <stdio.h> #include <stdlib.h> #include <string.h> #include <errno.h> +#ifdef PROJECTS_H +#error proj.h must be included before projects.h +#endif +#ifdef PROJ_API_H +#error proj.h must be included before proj_api.h +#endif + #ifndef PROJ_H #define PROJ_H #ifdef __cplusplus @@ -58,33 +133,26 @@ extern "C" { #endif -/****************************************************************************** - proj.h is included by projects.h in order to determine the size of the - PJ_OBSERVATION object. +/* We need to access the PJ_VERSION symbol from proj_api.h */ +#define PROJ_API_INCLUDED_FOR_PJ_VERSION_ONLY +#include <proj_api.h> +#undef PROJ_API_INCLUDED_FOR_PJ_VERSION_ONLY - In order to stomp as little as possible on the traditional proj.4 name - space, proj.h is littered with inclusion guards, leaving only the minimum - possible implementation when included from projects.h - The PJ_OBSERVATION object is fully defined if proj.h is included alone or - in connection with *but before* projects.h (the latter may be needed in - some cases, where it is necessary to access this "bare essentials" API, - while still having direct access to PJ object internals) -******************************************************************************/ +/* first forward declare everything needed */ /* Data type for generic geodetic observations */ -struct PJ_OBSERVATION; -typedef struct PJ_OBSERVATION PJ_OBSERVATION; +struct PJ_OBS; +typedef struct PJ_OBS PJ_OBS; /* Data type for generic geodetic 3D data */ union PJ_TRIPLET; typedef union PJ_TRIPLET PJ_TRIPLET; /* Data type for generic geodetic 3D data plus epoch information */ -union PJ_SPATIOTEMPORAL; -typedef union PJ_SPATIOTEMPORAL PJ_SPATIOTEMPORAL; +union PJ_COORD; +typedef union PJ_COORD PJ_COORD; -#ifndef PROJECTS_H /* Data type for projection/transformation information */ struct PJconsts; typedef struct PJconsts PJ; /* the PJ object herself */ @@ -115,21 +183,18 @@ typedef struct { double lam, phi, z; } LPZ; /* Degrees, minutes, and seconds */ typedef struct { double d, m, s; } PJ_DMS; -/* Geoid undulation (N) and deflections of the vertical (z, e) */ -typedef struct { double N, z, e; } PJ_NZE; +/* Geoid undulation (N) and deflections of the vertical (eta, zeta) */ +typedef struct { double e, z, N; } PJ_EZN; /* Ellipsoidal parameters */ typedef struct { double a, f; } PJ_AF; -#endif /* ndef PROJECTS_H */ -union PJ_SPATIOTEMPORAL { -#ifndef PROJECTS_H +union PJ_COORD { PJ_XYZT xyzt; PJ_UVWT uvwt; PJ_ENHT enht; PJ_LPZT lpzt; PJ_ENH enh; -#endif double v[4]; /* Who cares - it's just a vector! */ XYZ xyz; UVW uvw; @@ -143,12 +208,10 @@ union PJ_SPATIOTEMPORAL { /* Avoid preprocessor renaming and implicit type-punning: Use a union to make it explicit */ union PJ_TRIPLET { -#ifndef PROJECTS_H PJ_OPK opk; PJ_ENH enh; - PJ_NZE nze; + PJ_EZN ezn; PJ_AF af; -#endif double v[3]; /* Who cares - it's just a vector! */ XYZ xyz; LPZ lpz; @@ -158,39 +221,68 @@ union PJ_TRIPLET { UV uv; }; -struct PJ_OBSERVATION { - PJ_SPATIOTEMPORAL coo; /* coordinate data */ +struct PJ_OBS { + PJ_COORD coo; /* coordinate data */ PJ_TRIPLET anc; /* ancillary data */ int id; /* integer ancillary data - e.g. observation number, EPSG code... */ unsigned int flags; /* additional data, intended for flags */ }; +/* The context type - properly namespaced synonym for projCtx */ +struct projCtx_t; +typedef struct projCtx_t PJ_CONTEXT; +typedef int *PJ_FILE; -#ifndef PROJECTS_H +/* Manage the transformation definition object PJ */ +PJ *pj_create (const char *definition); +PJ *pj_create_argv (int argc, char **argv); +void pj_free (PJ *P); +int pj_error (PJ *P); -/* Direction: "+" forward, "-" reverse, 0: do nothing */ -enum pj_direction { - PJ_FWD = -1, - PJ_IDENT = 0, - PJ_INV = 1 +/* High level functionality for handling thread contexts */ +enum pj_debug_level { + PJ_LOG_NONE = 0, + PJ_LOG_ERROR = 1, + PJ_LOG_DEBUG_MAJOR = 2, + PJ_LOG_DEBUG_MINOR = 3 }; +void pj_debug_set (PJ *P, enum pj_debug_level debuglevel); +void pj_error_set (PJ *P, int err); +void pj_log_set (PJ *P, void *app_data, void (*log)(void *, int, const char *)); + +/* Lower level functionality for handling thread contexts */ +int pj_context_renew (PJ *P); +void pj_context_inherit (PJ *mother, PJ *daughter); +void pj_context_free (const PJ *P); + +/* Lowest level: Minimum support for fileapi */ +void pj_fileapi_set (PJ *P, void *fileapi); /* Apply transformation to observation - in forward or inverse direction */ -PJ_OBSERVATION pj_apply (PJ *P, enum pj_direction direction, PJ_OBSERVATION obs); +enum pj_direction { + PJ_FWD = 1, /* Forward */ + PJ_IDENT = 0, /* Do nothing */ + PJ_INV = -1 /* Inverse */ +}; +PJ_OBS pj_apply (PJ *P, enum pj_direction direction, PJ_OBS obs); /* Measure internal consistency - in forward or inverse direction */ -double pj_roundtrip(PJ *P, enum pj_direction direction, int n, PJ_OBSERVATION obs); +double pj_roundtrip (PJ *P, enum pj_direction direction, int n, PJ_OBS obs); + +/* Euclidean distance between two 2D coordinates stored in PJ_OBSs */ +double pj_obs_dist_2d (PJ_OBS a, PJ_OBS b); +/* Euclidean distance between two 3D coordinates stored in PJ_OBSs */ +double pj_obs_dist_3d (PJ_OBS a, PJ_OBS b); -int pj_show_triplet (FILE *stream, const char *banner, PJ_TRIPLET point); -/* Constructor for the OBSERVATION object */ -PJ_OBSERVATION pj_observation ( - double x, double y, double z, double t, - double o, double p, double k, - int id, unsigned int flags -); + +#ifndef PJ_OBS_C +extern const PJ_OBS pj_obs_error; +extern const PJ_OBS pj_obs_null; +extern const PJ *pj_shutdown; +#endif #ifndef TODEG #define TODEG(rad) ((rad)*180.0/M_PI) @@ -199,10 +291,9 @@ PJ_OBSERVATION pj_observation ( #define TORAD(deg) ((deg)*M_PI/180.0) #endif -#endif /* ndef PROJECTS_H */ - #ifdef __cplusplus } #endif + #endif /* ndef PROJ_H */ diff --git a/src/proj_api.h b/src/proj_api.h index 24a6f053..706e35ab 100644 --- a/src/proj_api.h +++ b/src/proj_api.h @@ -25,45 +25,77 @@ * DEALINGS IN THE SOFTWARE. *****************************************************************************/ + + + /* + * This version number should be updated with every release! The format of + * PJ_VERSION is + * + * * Before version 4.10.0: PJ_VERSION=MNP where M, N, and P are the major, + * minor, and patch numbers; e.g., PJ_VERSION=493 for version 4.9.3. + * + * * Version 4.10.0 and later: PJ_VERSION=MMMNNNPP later where MMM, NNN, PP + * are the major, minor, and patch numbers (the minor and patch numbers + * are padded with leading zeros if necessary); e.g., PJ_VERSION=401000 + * for version 4.10.0. + */ + #ifndef PJ_VERSION + #define PJ_VERSION 493 + #endif + + +/* If we're not asked for PJ_VERSION only, give them everything */ +#ifndef PROJ_API_INCLUDED_FOR_PJ_VERSION_ONLY /* General projections header file */ #ifndef PROJ_API_H #define PROJ_API_H -/* standard inclusions */ -#include <math.h> -#include <stdlib.h> - #ifdef __cplusplus extern "C" { #endif -/* - * This version number should be updated with every release! The format of - * PJ_VERSION is - * - * * Before version 4.10.0: PJ_VERSION=MNP where M, N, and P are the major, - * minor, and patch numbers; e.g., PJ_VERSION=493 for version 4.9.3. - * - * * Version 4.10.0 and later: PJ_VERSION=MMMNNNPP later where MMM, NNN, PP - * are the major, minor, and patch numbers (the minor and patch numbers - * are padded with leading zeros if necessary); e.g., PJ_VERSION=401000 - * for version 4.10.0. - */ -#define PJ_VERSION 493 + +/* standard inclusions */ +#include <math.h> +#include <stdlib.h> + /* pj_init() and similar functions can be used with a non-C locale */ /* Can be detected too at runtime if the symbol pj_atof exists */ #define PJ_LOCALE_SAFE 1 extern char const pj_release[]; /* global release id string */ +extern int pj_errno; /* global error return code */ #define RAD_TO_DEG 57.295779513082321 #define DEG_TO_RAD .017453292519943296 -extern int pj_errno; /* global error return code */ +#if defined(PROJECTS_H) || defined(PROJ_H) +#define PROJ_API_H_NOT_INVOKED_AS_PRIMARY_API +#endif + + + +#ifndef PROJ_H +/* In proj.h these macros are replaced by the enumeration pj_debug_level */ +#define PJ_LOG_NONE 0 +#define PJ_LOG_ERROR 1 +#define PJ_LOG_DEBUG_MAJOR 2 +#define PJ_LOG_DEBUG_MINOR 3 +#endif + -#if !defined(PROJECTS_H) +#ifdef PROJ_API_H_NOT_INVOKED_AS_PRIMARY_API + /* These make the function declarations below conform with classic proj */ + typedef PJ *projPJ; /* projPJ is a pointer to PJ */ + typedef projCtx_t *projCtx; /* projCtx is a pointer to projCtx_t */ +# define projXY XY +# define projLP LP +# define projXYZ XYZ +# define projLPZ LPZ +#else + /* i.e. proj_api invoked as primary API */ typedef struct { double u, v; } projUV; typedef struct { double u, v, w; } projUVW; typedef void *projPJ; @@ -72,15 +104,12 @@ extern int pj_errno; /* global error return code */ #define projXYZ projUVW #define projLPZ projUVW typedef void *projCtx; -#else - typedef PJ *projPJ; - typedef projCtx_t *projCtx; -# define projXY XY -# define projLP LP -# define projXYZ XYZ -# define projLPZ LPZ #endif + + + +/* If included *after* proj.h finishes, we have alternative names */ /* file reading api, like stdio */ typedef int *PAFile; typedef struct projFileAPI_t { @@ -91,14 +120,20 @@ typedef struct projFileAPI_t { void (*FClose)(PAFile); } projFileAPI; + + /* procedure prototypes */ +projCtx pj_get_default_ctx(void); +projCtx pj_get_ctx( projPJ ); + projXY pj_fwd(projLP, projPJ); projLP pj_inv(projXY, projPJ); projXYZ pj_fwd3d(projLPZ, projPJ); projLPZ pj_inv3d(projXYZ, projPJ); + int pj_transform( projPJ src, projPJ dst, long point_count, int point_offset, double *x, double *y, double *z ); int pj_datum_transform( projPJ src, projPJ dst, long point_count, int point_offset, @@ -128,6 +163,8 @@ projPJ pj_init_ctx( projCtx, int, char ** ); projPJ pj_init_plus_ctx( projCtx, const char * ); char *pj_get_def(projPJ, int); projPJ pj_latlong_from_proj( projPJ ); + + void *pj_malloc(size_t); void pj_dalloc(void *); void *pj_calloc (size_t n, size_t size); @@ -139,8 +176,8 @@ void pj_acquire_lock(void); void pj_release_lock(void); void pj_cleanup_lock(void); -projCtx pj_get_default_ctx(void); -projCtx pj_get_ctx( projPJ ); +int pj_run_selftests (int verbosity); + void pj_set_ctx( projPJ, projCtx ); projCtx pj_ctx_alloc(void); void pj_ctx_free( projCtx ); @@ -168,17 +205,10 @@ 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 -#define PJ_LOG_DEBUG_MINOR 3 #ifdef __cplusplus } #endif #endif /* ndef PROJ_API_H */ - +#endif /* ndef PROJ_API_INCLUDED_FOR_PJ_VERSION_ONLY */ diff --git a/src/projects.h b/src/projects.h index f6ac6acf..d07c43c7 100644 --- a/src/projects.h +++ b/src/projects.h @@ -168,6 +168,7 @@ typedef struct { double u, v, w; } UVW; /* Forward declarations and typedefs for stuff needed inside the PJ object */ struct PJconsts; +struct PJ_OBS; struct pj_opaque; struct ARG_list; struct FACTORS; @@ -181,6 +182,7 @@ enum pj_io_units { }; #ifndef PROJ_H typedef struct PJconsts PJ; /* the PJ object herself */ +typedef struct PJ_OBS PJ_OBS; #endif struct PJ_REGION_S { @@ -191,11 +193,9 @@ struct PJ_REGION_S { }; -#include "proj.h" /* Need this for sizeof(PJ_OBSERVATION) */ struct projCtx_t; typedef struct projCtx_t projCtx_t; - /* base projection data structure */ struct PJconsts { @@ -230,12 +230,13 @@ struct PJconsts { **************************************************************************************/ + XY (*fwd)(LP, PJ *); LP (*inv)(XY, PJ *); XYZ (*fwd3d)(LPZ, PJ *); LPZ (*inv3d)(XYZ, PJ *); - PJ_OBSERVATION (*fwdobs)(PJ_OBSERVATION, PJ *); - PJ_OBSERVATION (*invobs)(PJ_OBSERVATION, PJ *); + PJ_OBS (*fwdobs)(PJ_OBS, PJ *); + PJ_OBS (*invobs)(PJ_OBS, PJ *); void (*spc)(LP, PJ *, struct FACTORS *); @@ -280,8 +281,8 @@ struct PJconsts { double one_es; /* 1 - e^2 */ double rone_es; /* 1/one_es */ - /* The flattenings */ + /* The flattenings */ double f; /* first flattening */ double f2; /* second flattening */ double n; /* third flattening */ @@ -589,7 +590,6 @@ typedef struct _PJ_GridCatalog { struct _PJ_GridCatalog *next; } PJ_GridCatalog; - /* procedure prototypes */ double dmstor(const char *, char **); double dmstor_ctx(projCtx ctx, const char *, char **); @@ -722,9 +722,13 @@ struct PJ_PRIME_MERIDIANS *pj_get_prime_meridians_ref( void ); double pj_atof( const char* nptr ); double pj_strtod( const char *nptr, char **endptr ); +void pj_freeup_plain (PJ *P); #ifdef __cplusplus } #endif +#ifndef PROJECTS_H_ATEND +#define PROJECTS_H_ATEND +#endif #endif /* end of basic projections header */ |