From e2ce502a6eebd1250329f602e79bab1171f57403 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Tue, 3 Apr 2018 17:02:25 +0200 Subject: Remove proj.spec file The file seem to be unused by the packagers and we don't use ourselves. Instead of trying to keep it up to date with each release we might as well get rid of the thing. --- HOWTO-RELEASE | 18 ++++++------------ proj.spec | 51 --------------------------------------------------- 2 files changed, 6 insertions(+), 63 deletions(-) delete mode 100644 proj.spec diff --git a/HOWTO-RELEASE b/HOWTO-RELEASE index 3053523b..d020d486 100644 --- a/HOWTO-RELEASE +++ b/HOWTO-RELEASE @@ -46,33 +46,27 @@ 7) Update NEWS file with summary of release changes. -8) Update the proj.spec file for the new version: - - PACKAGE_VERSION - - Source0 - - %setup line - - rm -rf line - -9) If this is a major release, prepare a branch. +8) If this is a major release, prepare a branch. git checkout -b 5.0 -10) Tag the release +9) Tag the release git tag 5.0.0 -11) Do a "make dist-all" in the proj root directory. After some grinding +10) Do a "make dist-all" in the proj root directory. After some grinding this should result in files like proj-4.4.3.tar.gz and proj-4.4.3.zip being created. These are full source distributions. -12) Put these in the proj ftp area on /osgeo/download/proj on +11) Put these in the proj ftp area on /osgeo/download/proj on download.osgeo.org. This can be done via scp using a command like the following. scp proj-4.4.3.{tar.gz,zip} warmerdam@download.osgeo.org:/osgeo/download/proj -13) Announce the new release on the PROJ.4 and MetaCRS mailing lists. +12) Announce the new release on the PROJ.4 and MetaCRS mailing lists. -14) Issue a new release report on GitHub: +13) Issue a new release report on GitHub: https://github.com/OSGeo/proj.4/releases/new diff --git a/proj.spec b/proj.spec deleted file mode 100644 index b0854077..00000000 --- a/proj.spec +++ /dev/null @@ -1,51 +0,0 @@ -%define PACKAGE_NAME proj -%define PACKAGE_VERSION 5.0.0 -%define PACKAGE_URL http://proj4.org -%define _prefix /usr - -Summary: Cartographic projection and geodetic transformation software -Name: %PACKAGE_NAME -Version: %PACKAGE_VERSION -Release: 1 -Source0: proj-5.0.0.tar.gz -License: MIT, Copyright (c) 2000, Frank Warmerdam -Group: Applications/GIS -Vendor: Intevation GmbH -Distribution: FreeGIS CD - -BuildRoot: %{_builddir}/%{name}-root -Prefix: %{_prefix} -BuildRequires: libtool -Conflicts: PROJ.4 - -%description -PROJ is a generic coordinate transformation software, that transforms -coordinates from one coordinate reference system (CRS) to another. This -includes cartographic projections as well as geodetic transformations. - -%prep -%setup -D -n proj-5.0.0 -%{_builddir}/%{name}-%{version}/autogen.sh -%configure - -%build -make - -%install -rm -rf $RPM_BUILD_ROOT -%makeinstall - -%clean -rm -rf %{_builddir}/proj-5.0.0 -rm -rf $RPM_BUILD_ROOT - -%files -%defattr(-,root,root) -%{_libdir}/* -%{_bindir}/* -%{_includedir}/* -%{_datadir}/proj/* -%{_mandir}/man1/* -%{_mandir}/man3/* - -%doc AUTHORS COPYING ChangeLog NEWS README.md -- cgit v1.2.3 From 51967aa03df0bbe134c38475a2633430dddefc82 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Tue, 3 Apr 2018 17:45:23 +0200 Subject: Add --skip-lines option to cct With this it is possible to skip the header when transforming coordinates from a file. --- docs/source/apps/cct.rst | 7 +++++++ src/cct.c | 24 ++++++++++++++++++++---- 2 files changed, 27 insertions(+), 4 deletions(-) diff --git a/docs/source/apps/cct.rst b/docs/source/apps/cct.rst index dc2944f2..20aec06c 100644 --- a/docs/source/apps/cct.rst +++ b/docs/source/apps/cct.rst @@ -44,6 +44,13 @@ The following control parameters can appear in any order: Specify a fixed observation height to be used for all input data. +.. option:: s , --skip-lines= + + .. versionadded:: 5.1.0 + + Skip the first *n* lines of input. This applies to any kind of input, wether + it comes from ``STDIN``, a file or interactive user input. + .. option:: -v, --verbose Write non-essential, but potentially useful, information to stderr. diff --git a/src/cct.c b/src/cct.c index 6b8f35f3..2952f2b7 100644 --- a/src/cct.c +++ b/src/cct.c @@ -105,6 +105,7 @@ static const char usage[] = { " -z value Provide a fixed z value for all input data (e.g. -z 0)\n" " -t value Provide a fixed t value for all input data (e.g. -t 0)\n" " -I Do the inverse transformation\n" + " -s n Skip n first lines of a infile\n" " -v Verbose: Provide non-essential informational output.\n" " Repeat -v for more verbosity (e.g. -vv)\n" "--------------------------------------------------------------------------------\n" @@ -116,6 +117,7 @@ static const char usage[] = { " --time Alias for -t\n" " --verbose Alias for -v\n" " --inverse Alias for -I\n" + " --skip-lines Alias for -s\n" " --help Alias for -h\n" " --version Print version number\n" "--------------------------------------------------------------------------------\n" @@ -156,13 +158,19 @@ int main(int argc, char **argv) { OPTARGS *o; FILE *fout = stdout; char *buf; - int nfields = 4, direction = 1, verbose; + int nfields = 4, direction = 1, skip_lines = 0, verbose; double fixed_z = HUGE_VAL, fixed_time = HUGE_VAL; int columns_xyzt[] = {1, 2, 3, 4}; const char *longflags[] = {"v=verbose", "h=help", "I=inverse", "version", 0}; - const char *longkeys[] = {"o=output", "c=columns", "z=height", "t=time", 0}; - - o = opt_parse (argc, argv, "hvI", "cozt", longflags, longkeys); + const char *longkeys[] = { + "o=output", + "c=columns", + "z=height", + "t=time", + "s=skip-lines", + 0}; + + o = opt_parse (argc, argv, "hvI", "cozts", longflags, longkeys); if (0==o) return 0; @@ -202,6 +210,10 @@ int main(int argc, char **argv) { nfields--; } + if (opt_given (o, "s")) { + skip_lines = atoi (opt_arg(o, "s")); + } + if (opt_given (o, "c")) { /* cppcheck-suppress invalidscanf */ int ncols = sscanf (opt_arg (o, "c"), "%d,%d,%d,%d", columns_xyzt, columns_xyzt+1, columns_xyzt+2, columns_xyzt+3); @@ -265,6 +277,10 @@ int main(int argc, char **argv) { continue; } point = parse_input_line (buf, columns_xyzt, fixed_z, fixed_time); + if (skip_lines > 0) { + skip_lines--; + continue; + } if (HUGE_VAL==point.xyzt.x) { char *c = column (buf, 1); -- cgit v1.2.3 From 64922f1694fba917c92116b0f595e76d2304046b Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Sun, 1 Apr 2018 19:50:34 +0200 Subject: Add a CITATION file (fixes #309) --- CITATION | 14 ++++++++++++++ HOWTO-RELEASE | 14 ++++++++------ Makefile.am | 2 +- README.md | 4 ++++ 4 files changed, 27 insertions(+), 7 deletions(-) create mode 100644 CITATION diff --git a/CITATION b/CITATION new file mode 100644 index 00000000..82643d84 --- /dev/null +++ b/CITATION @@ -0,0 +1,14 @@ +To cite PROJ in publications use: + + PROJ contributors (2018). PROJ coordinate transformation software + library. Open Source Geospatial Foundation. URL http://proj4.org/. + +A BibTeX entry for LaTeX users is + + @Manual{, + title = {{PROJ} coordinate transformation software library}, + author = {{PROJ contributors}}, + organization = {Open Source Geospatial Foundation}, + year = {2018}, + url = {http://proj4.org/}, + } diff --git a/HOWTO-RELEASE b/HOWTO-RELEASE index d020d486..fbf99670 100644 --- a/HOWTO-RELEASE +++ b/HOWTO-RELEASE @@ -46,27 +46,29 @@ 7) Update NEWS file with summary of release changes. -8) If this is a major release, prepare a branch. +8) Update the year in CITATION + +9) If this is a major release, prepare a branch. git checkout -b 5.0 -9) Tag the release +10) Tag the release git tag 5.0.0 -10) Do a "make dist-all" in the proj root directory. After some grinding +11) Do a "make dist-all" in the proj root directory. After some grinding this should result in files like proj-4.4.3.tar.gz and proj-4.4.3.zip being created. These are full source distributions. -11) Put these in the proj ftp area on /osgeo/download/proj on +12) Put these in the proj ftp area on /osgeo/download/proj on download.osgeo.org. This can be done via scp using a command like the following. scp proj-4.4.3.{tar.gz,zip} warmerdam@download.osgeo.org:/osgeo/download/proj -12) Announce the new release on the PROJ.4 and MetaCRS mailing lists. +13) Announce the new release on the PROJ.4 and MetaCRS mailing lists. -13) Issue a new release report on GitHub: +14) Issue a new release report on GitHub: https://github.com/OSGeo/proj.4/releases/new diff --git a/Makefile.am b/Makefile.am index f4a82f9c..39a1201b 100644 --- a/Makefile.am +++ b/Makefile.am @@ -1,6 +1,6 @@ SUBDIRS = src man nad jniwrap cmake test -EXTRA_DIST = makefile.vc nmake.opt CMakeLists.txt +EXTRA_DIST = makefile.vc nmake.opt CMakeLists.txt CITATION pkgconfigdir = $(libdir)/pkgconfig pkgconfig_DATA = proj.pc diff --git a/README.md b/README.md index 32e3df66..981be967 100644 --- a/README.md +++ b/README.md @@ -144,3 +144,7 @@ where PROJ will find them after installation. The default location of `PROJ_LIB` on UNIX-based systems is `/usr/local/share/proj` but it may be changed to a different directory. On Windows you have to define `PROJ_LIB` yourself. + +### Citing PROJ in publications + +See [CITATION](CITATION) -- cgit v1.2.3 From ac78787bf3bfa15dabbbad0042ca762020542365 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Sat, 7 Apr 2018 15:35:08 +0200 Subject: Pipeline: make sure geocentric/cartesian space transform is done with original ellipsoid parameters (before any projection mess with them) --- src/PJ_pipeline.c | 2 ++ src/proj_4D_api.c | 4 ++-- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/src/PJ_pipeline.c b/src/PJ_pipeline.c index b60e05a8..aa7d76f8 100644 --- a/src/PJ_pipeline.c +++ b/src/PJ_pipeline.c @@ -306,6 +306,8 @@ static void set_ellipsoid(PJ *P) { /* the PJ you provided". */ proj_errno_reset (P); } + P->a_orig = P->a; + P->es_orig = P->es; pj_calc_ellipsoid_params (P, P->a, P->es); diff --git a/src/proj_4D_api.c b/src/proj_4D_api.c index 11a56ac0..aed4685b 100644 --- a/src/proj_4D_api.c +++ b/src/proj_4D_api.c @@ -485,7 +485,7 @@ Returns 1 on success, 0 on failure if (0==d[0] && 0==d[1] && 0==d[2] && 0==d[3] && 0==d[4] && 0==d[5] && 0==d[6]) { /* If the current ellipsoid is not WGS84, then make sure the */ /* change in ellipsoid is still done. */ - if (!(fabs(P->a - 6378137.0) < 1e-8 && fabs(P->f - 1./ 298.257223563) < 1e-15)) { + if (!(fabs(P->a_orig - 6378137.0) < 1e-8 && fabs(P->es_orig - 0.0066943799901413) < 1e-15)) { do_cart = 1; } break; @@ -512,7 +512,7 @@ Returns 1 on success, 0 on failure /* geocentric/cartesian space or we need to do a Helmert transform. */ if (P->is_geocent || P->helmert || do_cart) { char def[150]; - sprintf (def, "break_cs2cs_recursion proj=cart a=%40.20g f=%40.20g", P->a, P->f); + sprintf (def, "break_cs2cs_recursion proj=cart a=%40.20g es=%40.20g", P->a_orig, P->es_orig); Q = proj_create (P->ctx, def); if (0==Q) return 0; -- cgit v1.2.3 From 3bcd7bb02392d6ecda5b2a49854d855c92bfcbd4 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Sun, 8 Apr 2018 17:26:55 +0200 Subject: Override style of Sphinx versionmodified directives so they stand out more clearly in the text --- docs/source/_static/theme_overrides.css | 13 +++++++++++++ 1 file changed, 13 insertions(+) diff --git a/docs/source/_static/theme_overrides.css b/docs/source/_static/theme_overrides.css index b82d2881..4e612939 100644 --- a/docs/source/_static/theme_overrides.css +++ b/docs/source/_static/theme_overrides.css @@ -55,3 +55,16 @@ div.math { overflow-x: auto; overflow-y: hidden; } + +/* Override style of versionmodified directives */ +.versionadded { + font-style: italic; +} + +.versionchanged { + font-style: italic; +} + +.deprecated { + font-style: italic; +} -- cgit v1.2.3 From 89af2753651f8c21eae9df631a71c0a3b03bcc3b Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Mon, 9 Apr 2018 02:02:26 +0200 Subject: Fix bad example in eqc docs --- docs/source/operations/projections/eqc.rst | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/docs/source/operations/projections/eqc.rst b/docs/source/operations/projections/eqc.rst index be431977..5833bd61 100644 --- a/docs/source/operations/projections/eqc.rst +++ b/docs/source/operations/projections/eqc.rst @@ -64,9 +64,8 @@ Example using EPSG 32662 (WGS84 Plate Carrée):: Example using Plate Carrée projection with true scale at latitude 30° and central meridian 90°W:: - $ echo -88 30 | proj +proj=eqc +lat_ts=30 +lat_0=90w - -8483684.61 13358338.90 - + $ echo -88 30 | proj +proj=eqc +lat_ts=30 +lon_0=90w + 192811.01 3339584.72 Mathematical definition ####################### -- cgit v1.2.3 From 3a2ddb6c6efbccf388ea89e177ca51fd25946ecf Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Sat, 7 Apr 2018 13:41:12 +0200 Subject: Add webmerc projection --- docs/source/operations/projections/index.rst | 1 + docs/source/operations/projections/webmerc.rst | 78 ++++++++++++++++++++++++++ src/PJ_merc.c | 17 ++++++ src/pj_list.h | 1 + test/gie/4D-API_cs2cs-style.gie | 13 +++++ 5 files changed, 110 insertions(+) create mode 100644 docs/source/operations/projections/webmerc.rst diff --git a/docs/source/operations/projections/index.rst b/docs/source/operations/projections/index.rst index 9255d9fd..77decd6c 100644 --- a/docs/source/operations/projections/index.rst +++ b/docs/source/operations/projections/index.rst @@ -141,6 +141,7 @@ Projections map the spherical 3D space to a flat 2D space. wag5 wag6 wag7 + webmerc weren wink1 wink2 diff --git a/docs/source/operations/projections/webmerc.rst b/docs/source/operations/projections/webmerc.rst new file mode 100644 index 00000000..0c465afc --- /dev/null +++ b/docs/source/operations/projections/webmerc.rst @@ -0,0 +1,78 @@ +.. _webmerc: + +******************************************************************************** +Web Mercator / Pseudo Mercator +******************************************************************************** + +The Web Mercator / Pseudo Mercator projection is a cylindrical map projection. +This is a variant of the regular :ref:`merc` projection, except that the computation +is done on a sphere, using the semi-major axis of the ellipsoid. + +From `Wikipedia `_: + + This projection is widely used by the Web Mercator, Google Web Mercator, + Spherical Mercator, WGS 84 Web Mercator[1] or WGS 84/Pseudo-Mercator is a + variant of the Mercator projection and is the de facto standard for Web + mapping applications. [...] + It is used by virtually all major online map providers [...] + Its official EPSG identifier is EPSG:3857, although others have been used + historically. + + ++---------------------+----------------------------------------------------------+ +| **Classification** | Cylindrical (non conformant if used with ellipsoid) | ++---------------------+----------------------------------------------------------+ +| **Version** | Added in proj 5.1.0 | ++---------------------+----------------------------------------------------------+ +| **Available forms** | Forward and inverse, spherical projection | ++---------------------+----------------------------------------------------------+ +| **Defined area** | Global, but best used near the equator | ++---------------------+----------------------------------------------------------+ +| **Options** | Neither lat_0, lon_0, lat_ts or k_0 should be used. | ++---------------------+----------------------------------------------------------+ + +Usage +######## + +Example:: + + $ echo 2 49 | proj +proj=webmerc +datum=WGS84 + 222638.98 6274861.39 + +Mathematical definition +####################### + +The formulas describing the Mercator projection are all taken from G. Evenden's libproj manuals [Evenden2005]_. + +Forward projection +================== + +.. math:: + + x = \lambda + +.. math:: + + y = \ln \left[ \tan \left(\frac{\pi}{4} + \frac{\phi}{2} \right) \right] + + +Inverse projection +================== + +.. math:: + + \lambda = {x} + +.. math:: + + \phi = \frac{\pi}{2} - 2 \arctan \left[ e^{-y} \right] + + + +Further reading +############### + +#. `Wikipedia `_ + + + diff --git a/src/PJ_merc.c b/src/PJ_merc.c index 18303c45..a3e5e846 100644 --- a/src/PJ_merc.c +++ b/src/PJ_merc.c @@ -1,8 +1,10 @@ #define PJ_LIB__ +#include "proj_internal.h" #include "proj.h" #include "projects.h" PROJ_HEAD(merc, "Mercator") "\n\tCyl, Sph&Ell\n\tlat_ts="; +PROJ_HEAD(webmerc, "Web Mercator / Pseudo Mercator") "\n\tCyl, Sph\n\t"; #define EPS10 1.e-10 @@ -76,3 +78,18 @@ PJ *PROJECTION(merc) { return P; } +PJ *PROJECTION(webmerc) { + + /* Overriding k_0, lat_0 and lon_0 with fixed parameters */ + P->k0 = 1.0; + P->phi0 = 0.0; + P->lam0 = 0.0; + + P->b = P->a; + /* Clean up the ellipsoidal parameters to reflect the sphere */ + P->es = P->e = P->f = 0; + pj_calc_ellipsoid_params (P, P->a, 0); + P->inv = s_inverse; + P->fwd = s_forward; + return P; +} diff --git a/src/pj_list.h b/src/pj_list.h index 09a66034..440f7972 100644 --- a/src/pj_list.h +++ b/src/pj_list.h @@ -156,6 +156,7 @@ PROJ_HEAD(wag4, "Wagner IV") PROJ_HEAD(wag5, "Wagner V") PROJ_HEAD(wag6, "Wagner VI") PROJ_HEAD(wag7, "Wagner VII") +PROJ_HEAD(webmerc, "Web Mercator / Pseudo Mercator") PROJ_HEAD(weren, "Werenskiold I") PROJ_HEAD(wink1, "Winkel I") PROJ_HEAD(wink2, "Winkel II") diff --git a/test/gie/4D-API_cs2cs-style.gie b/test/gie/4D-API_cs2cs-style.gie index bbd9ee39..6eed4faa 100644 --- a/test/gie/4D-API_cs2cs-style.gie +++ b/test/gie/4D-API_cs2cs-style.gie @@ -215,6 +215,19 @@ accept 487147.594520173 4934316.46263998 0 expect -10370728.80 5552839.74 0 ------------------------------------------------------------------------------- +------------------------------------------------------------------------------- +Test that Google's Web Mercator with +proj=webmerc +------------------------------------------------------------------------------- +operation proj=pipeline step init=epsg:26915 inv step proj=webmerc datum=WGS84 +------------------------------------------------------------------------------- +tolerance 20 cm +accept 487147.594520173 4934316.46263998 +expect -10370728.80 5552839.74 + +accept 487147.594520173 4934316.46263998 0 +expect -10370728.80 5552839.74 0 +------------------------------------------------------------------------------- + ------------------------------------------------------------------------------- Test that +datum parameters are handled correctly in pipelines. See #872 for details. -- cgit v1.2.3 From bc8026a46606abf48679d81d80991623d5b91ca7 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Mon, 9 Apr 2018 22:27:59 +0200 Subject: Add description of projection and +lat_b parameter --- docs/source/operations/projections/airy.rst | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/docs/source/operations/projections/airy.rst b/docs/source/operations/projections/airy.rst index d2730592..eefda4e9 100644 --- a/docs/source/operations/projections/airy.rst +++ b/docs/source/operations/projections/airy.rst @@ -3,6 +3,11 @@ ******************************************************************************** Airy ******************************************************************************** + +The Airy projection is an azimuthal minimum error projection for the region +within the small or great circle defined by an angular distance, +:math:`\phi_b`, from the tangency point of the plane :math:`( \lambda_0, \phi_0 )`. + +---------------------+----------------------------------------------------------+ | **Classification** | | +---------------------+----------------------------------------------------------+ @@ -14,7 +19,10 @@ Airy +---------------------+----------------------------------------------------------+ | **Options** | +---------------------+----------------------------------------------------------+ -| `+lat_b` | | +| `+lat_b` | Angular distance from tangency point of the plane | +| | :math:`( \lambda_0, \phi_0 )` where the error is kept at | +| | minimum. Defaults to :math:`90^{\circ}` which is | +| | suitable for hemispherical maps. | +---------------------+----------------------------------------------------------+ | `+no_cut` | Do not cut at hemisphere limit | +---------------------+----------------------------------------------------------+ -- cgit v1.2.3 From b7c59722cb0058391e47d2899d510961a0985c59 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Mon, 9 Apr 2018 22:28:20 +0200 Subject: Add more test cases for airy projection --- test/gie/builtins.gie | 68 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 68 insertions(+) diff --git a/test/gie/builtins.gie b/test/gie/builtins.gie index 4d744135..84d5c832 100644 --- a/test/gie/builtins.gie +++ b/test/gie/builtins.gie @@ -318,6 +318,74 @@ expect -189109.886908621 94583.752387504 accept -2 -1 expect -189109.886908621 -94583.752387504 +------------------------------------------------------------------------------- +Test north polar aspect +------------------------------------------------------------------------------- +operation +proj=airy +R=1 +lat_0=90 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 0 +expect 0 -1.3863 +accept 0 90 +expect 0 0 +accept 0 -90 +expect failure errno tolerance_condition + + +------------------------------------------------------------------------------- +Test south polar aspect +------------------------------------------------------------------------------- +operation +proj=airy +R=1 +lat_0=-90 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 0 +expect 0 1.3863 +accept 0 -90 +expect 0 0 +accept 0 90 +expect failure errno tolerance_condition + +------------------------------------------------------------------------------- +Test oblique aspect +------------------------------------------------------------------------------- +operation +proj=airy +R=1 +lon_0=45 +lat_0=45 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 45 45 +expect 0 0 +accept 0 0 +expext -0.7336 -0.5187 +accept -45 -45 +expect failure errno tolerance_condition + +------------------------------------------------------------------------------- +Test that coordinates on the opposing hemisphere are projected when using ++no_cut. +------------------------------------------------------------------------------- +operation +proj=airy +R=1 +lat_0=-90 +no_cut +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 10 +expect 0 1.5677 + + +------------------------------------------------------------------------------- +Test the +lat_b parameter +------------------------------------------------------------------------------- +operation +proj=airy +R=1 +lat_b=89.99999999 # check tolerance +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 0 +expect 0 0 +------------------------------------------------------------------------------- +operation +proj=airy +R=1 +lat_b=30 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 0 +expect 0 0 +accept 25 25 +expect 0.3821 0.4216 + =============================================================================== Aitoff -- cgit v1.2.3 From 8bd2fc699e61eb3665b7d16d450886a914258933 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Wed, 11 Apr 2018 00:03:03 +0200 Subject: Add more tests for laea --- test/gie/builtins.gie | 151 +++++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 149 insertions(+), 2 deletions(-) diff --git a/test/gie/builtins.gie b/test/gie/builtins.gie index 4d744135..47d63bf1 100644 --- a/test/gie/builtins.gie +++ b/test/gie/builtins.gie @@ -2143,7 +2143,7 @@ Lambert Azimuthal Equal Area =============================================================================== ------------------------------------------------------------------------------- -operation +proj=laea +ellps=GRS80 +lat_1=0.5 +lat_2=2 +operation +proj=laea +ellps=GRS80 ------------------------------------------------------------------------------- tolerance 0.1 mm accept 2 1 @@ -2155,6 +2155,9 @@ expect -222602.471450095 110589.827224410 accept -2 -1 expect -222602.471450095 -110589.827224409 +accept 180 0 +expect failure errno tolerance_condition + direction inverse accept 200 100 expect 0.001796631 0.000904369 @@ -2166,7 +2169,7 @@ accept -200 -100 expect -0.001796631 -0.000904369 ------------------------------------------------------------------------------- -operation +proj=laea +R=6400000 +lat_1=0.5 +lat_2=2 +operation +proj=laea +R=6400000 ------------------------------------------------------------------------------- tolerance 0.1 mm accept 2 1 @@ -2178,6 +2181,9 @@ expect -223365.281370125 111716.668072916 accept -2 -1 expect -223365.281370125 -111716.668072916 +accept 180 0 +expect failure errno tolerance_condition + direction inverse accept 200 100 expect 0.001790493 0.000895247 @@ -2188,6 +2194,147 @@ expect -0.001790493 0.000895247 accept -200 -100 expect -0.001790493 -0.000895247 +------------------------------------------------------------------------------- +Test oblique aspect of the spherical form +------------------------------------------------------------------------------- +operation +proj=laea +R=1 +lat_0=45 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 45 +expect 0 0 +accept 0 0 +expect 0 -0.7654 +accept 0 90 +expect 0 0.7654 +accept 0 -45 +expect 0 -1.4142 +accept 45 45 +expect 0.5194 0.1521 + +tolerance 0.1 mm +accept 45 45 +roundtrip 100 + +# error when waaay outside the sphere +direction inverse +accept 0 10 +expect failure errno tolerance_condition + +------------------------------------------------------------------------------- +Test oblique aspect of the ellipsoidal form +------------------------------------------------------------------------------- +operation +proj=laea +ellps=GRS80 +lat_0=45 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 45 +expect 0 0 +accept 0 0 +expect 0 -4860248.8602 +accept 0 90 +expect 0 4886594.2207 +accept 0 -45 +expect 0 -8984728.0442 +accept 45 45 +expect 3318800.8682 968788.2336 + +tolerance 10 cm +accept 45 45 +roundtrip 100 + +# test rho Date: Wed, 11 Apr 2018 17:56:08 +0200 Subject: Enhance the precision of Spherical Mercator projection near the equator This uses a linear approximation of tan(x+pi/4) for better precision at small latitudes. As a result points of latitude 0 maintain a 0 Y coordinate and 0,0 is transformed to 0,0 --- src/PJ_merc.c | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/src/PJ_merc.c b/src/PJ_merc.c index a3e5e846..3801b828 100644 --- a/src/PJ_merc.c +++ b/src/PJ_merc.c @@ -8,6 +8,13 @@ PROJ_HEAD(webmerc, "Web Mercator / Pseudo Mercator") "\n\tCyl, Sph\n\t"; #define EPS10 1.e-10 +static double _tan_near_fort_pi(double x) { + if (fabs(x) <= __DBL_EPSILON__) { + return 2*x + 1.0; + } + return tan(M_FORTPI + x); +} + static XY e_forward (LP lp, PJ *P) { /* Ellipsoidal, forward */ XY xy = {0.0,0.0}; if (fabs(fabs(lp.phi) - M_HALFPI) <= EPS10) { @@ -27,7 +34,7 @@ static XY s_forward (LP lp, PJ *P) { /* Spheroidal, forward */ return xy; } xy.x = P->k0 * lp.lam; - xy.y = P->k0 * log(tan(M_FORTPI + .5 * lp.phi)); + xy.y = P->k0 * log(_tan_near_fort_pi(.5 * lp.phi)); return xy; } -- cgit v1.2.3 From 60d86d4914b816143b4145ccae0def0420f2a70d Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Wed, 11 Apr 2018 23:11:23 +0200 Subject: Update docs for nsper --- docs/source/operations/projections/nsper.rst | 23 ++++++++++++++++++++++- 1 file changed, 22 insertions(+), 1 deletion(-) diff --git a/docs/source/operations/projections/nsper.rst b/docs/source/operations/projections/nsper.rst index 4a8651fc..3f46f2b4 100644 --- a/docs/source/operations/projections/nsper.rst +++ b/docs/source/operations/projections/nsper.rst @@ -4,7 +4,28 @@ Near-sided perspective ******************************************************************************** +The near-sided perspective projection simulates a view from a height +:math:`h` similar to how a satellite in orbit would see it. + ++---------------------+-------------------------------------------------------------+ +| **Classification** | Azimuthal. Neither conformal nor equal area. | ++---------------------+-------------------------------------------------------------+ +| **Available forms** | Forward and inverse spherical projection | ++---------------------+-------------------------------------------------------------+ +| **Defined area** | Global, although for one hemisphere at a time. | ++---------------------+-------------------------------------------------------------+ +| **Options** | ++---------------------+-------------------------------------------------------------+ +| `+h` | Height of the view point above the Earth and must be in | +| | the same units as the radius of the sphere. Required. | ++---------------------+-------------------------------------------------------------| +| `+lat_0` | Latitude (in degrees) of the view position. Defaults to 0. | ++---------------------+-------------------------------------------------------------+ +| `+lon_0` | Longitude (in degrees) of the view position. Defaults to 0. | ++---------------------+-------------------------------------------------------------+ + + .. image:: ./images/nsper.png :scale: 50% - :alt: Near-sided perspective + :alt: Near-sided perspective -- cgit v1.2.3 From 38443f7545dad1d849a889921976fa94e9b76264 Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Wed, 11 Apr 2018 23:11:41 +0200 Subject: Add more tests for nsper --- test/gie/builtins.gie | 54 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 54 insertions(+) diff --git a/test/gie/builtins.gie b/test/gie/builtins.gie index 4f6d6a44..fb8a53c3 100644 --- a/test/gie/builtins.gie +++ b/test/gie/builtins.gie @@ -3347,6 +3347,60 @@ expect -0.001790493 0.000895247 accept -200 -100 expect -0.001790493 -0.000895247 +------------------------------------------------------------------------------- +Test north polar aspect +------------------------------------------------------------------------------- +operation +proj=nsper +R=1 +h=3 +lat_0=90 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 90 +expect 0 0 +accept 45 45 +expect 0.4555 -0.4555 +roundtrip 100 +accept 0 0 +expect failure errno tolerance_condition + +direction inverse +accept 0 0 +expect 0 90 +accept 0 2 # projected coordinate is outside the sphere +expect failure errno tolerance_condition + +------------------------------------------------------------------------------- +Test south polar aspect +------------------------------------------------------------------------------- +operation +proj=nsper +R=1 +h=3 +lat_0=-90 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 -90 +expect 0 0 +accept -45 -45 +expect -0.4555 0.4555 +roundtrip 100 +accept 0 0 +expect failure errno tolerance_condition + +------------------------------------------------------------------------------- +operation +proj=nsper +R=1 +h=3 +lat_0=45 +------------------------------------------------------------------------------- +tolerance 0.1 mm +accept 0 45 +expect 0 0 +accept 0 90 +expect 0 0.6442 +accept 0 0 +expect 0 -0.6442 +accept 45 45 +expect 0.4767 0.1396 +roundtrip 100 + + +------------------------------------------------------------------------------- +operation +proj=nsper +R=1 +h=0 +------------------------------------------------------------------------------- +expect failure errno h_less_than_zero + =============================================================================== New Zealand Map Grid -- cgit v1.2.3 From bd2e1363b021e1a1c6a4ebbc0aa55e8dc98e2587 Mon Sep 17 00:00:00 2001 From: Javier Goizueta Date: Thu, 12 Apr 2018 08:50:59 +0200 Subject: Fix: use proper double epsilon declaration --- src/PJ_merc.c | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/PJ_merc.c b/src/PJ_merc.c index 3801b828..2c89fbc6 100644 --- a/src/PJ_merc.c +++ b/src/PJ_merc.c @@ -2,6 +2,7 @@ #include "proj_internal.h" #include "proj.h" #include "projects.h" +#include PROJ_HEAD(merc, "Mercator") "\n\tCyl, Sph&Ell\n\tlat_ts="; PROJ_HEAD(webmerc, "Web Mercator / Pseudo Mercator") "\n\tCyl, Sph\n\t"; @@ -9,7 +10,7 @@ PROJ_HEAD(webmerc, "Web Mercator / Pseudo Mercator") "\n\tCyl, Sph\n\t"; #define EPS10 1.e-10 static double _tan_near_fort_pi(double x) { - if (fabs(x) <= __DBL_EPSILON__) { + if (fabs(x) <= DBL_EPSILON) { return 2*x + 1.0; } return tan(M_FORTPI + x); -- cgit v1.2.3 From f2b3604a94349697828baa7de63a30a97ac15b2b Mon Sep 17 00:00:00 2001 From: Javier Goizueta Date: Thu, 12 Apr 2018 09:47:43 +0200 Subject: Use log1p in forward spherical mercator --- src/PJ_merc.c | 28 ++++++++++++++++++++++++---- 1 file changed, 24 insertions(+), 4 deletions(-) diff --git a/src/PJ_merc.c b/src/PJ_merc.c index 2c89fbc6..0bf98625 100644 --- a/src/PJ_merc.c +++ b/src/PJ_merc.c @@ -9,11 +9,31 @@ PROJ_HEAD(webmerc, "Web Mercator / Pseudo Mercator") "\n\tCyl, Sph\n\t"; #define EPS10 1.e-10 -static double _tan_near_fort_pi(double x) { +#if !defined(HAVE_C99_MATH) +#define HAVE_C99_MATH 0 +#endif + +#if HAVE_C99_MATH +#define log1px log1p +#else +static double log1px(double x) { + volatile double + y = 1 + x, + z = y - 1; + /* Here's the explanation for this magic: y = 1 + z, exactly, and z + * approx x, thus log(y)/z (which is nearly constant near z = 0) returns + * a good approximation to the true log(1 + x)/x. The multiplication x * + * (log(y)/z) introduces little additional error. */ + return z == 0 ? x : x * log(y) / z; +} +#endif + +static double logtanpfpim1(double x) { /* log(tan(x/2 + M_FORTPI)) */ if (fabs(x) <= DBL_EPSILON) { - return 2*x + 1.0; + /* tan(M_FORTPI + .5 * x) can be approximated by 1.0 + x */ + return log1px(x); } - return tan(M_FORTPI + x); + return log(tan(M_FORTPI + .5 * x)); } static XY e_forward (LP lp, PJ *P) { /* Ellipsoidal, forward */ @@ -35,7 +55,7 @@ static XY s_forward (LP lp, PJ *P) { /* Spheroidal, forward */ return xy; } xy.x = P->k0 * lp.lam; - xy.y = P->k0 * log(_tan_near_fort_pi(.5 * lp.phi)); + xy.y = P->k0 * logtanpfpim1(lp.phi); return xy; } -- cgit v1.2.3 From be27d12a028d07ed656080004c949c0da7a6012e Mon Sep 17 00:00:00 2001 From: Kristian Evers Date: Fri, 13 Apr 2018 18:11:29 +0200 Subject: Always ignore out-commented lines in cct Previous to this commit cct would return the following $ cct -c 2,3,4 -t 0 -I +proj=cart +ellps=GRS80 BLAH 3579685.56545 508396.50343 5236837.50646 8.0832413787 55.5578176654 99.9833 0.0000 8.0832413787 55.5578176654 99.9833 0.0000 where the second input should not be parsed as a valid coordinate. With this commit that no longer happens and the input is returned verbatim back to the user. Closes #932 --- src/cct.c | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/src/cct.c b/src/cct.c index 2952f2b7..00fc2127 100644 --- a/src/cct.c +++ b/src/cct.c @@ -271,6 +271,7 @@ int main(int argc, char **argv) { while (opt_input_loop (o, optargs_file_format_text)) { int err; void *ret = fgets (buf, 10000, o->input); + char *c = column (buf, 1); opt_eof_handler (o); if (0==ret) { fprintf (stderr, "Read error in record %d\n", (int) o->record_index); @@ -281,15 +282,14 @@ int main(int argc, char **argv) { skip_lines--; continue; } - if (HUGE_VAL==point.xyzt.x) { - char *c = column (buf, 1); - /* if it's a comment or blank line, we reflect it */ - if (c && ((*c=='\0') || (*c=='#'))) { - fprintf (fout, "%s\n", buf); - continue; - } + /* if it's a comment or blank line, we reflect it */ + if (c && ((*c=='\0') || (*c=='#'))) { + fprintf (fout, "%s", buf); + continue; + } + if (HUGE_VAL==point.xyzt.x) { /* otherwise, it must be a syntax error */ fprintf (fout, "# Record %d UNREADABLE: %s", (int) o->record_index, buf); if (verbose) -- cgit v1.2.3