aboutsummaryrefslogtreecommitdiff
path: root/src/projections/adams.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/projections/adams.cpp')
-rw-r--r--src/projections/adams.cpp150
1 files changed, 148 insertions, 2 deletions
diff --git a/src/projections/adams.cpp b/src/projections/adams.cpp
index 389ca054..4150fff2 100644
--- a/src/projections/adams.cpp
+++ b/src/projections/adams.cpp
@@ -6,9 +6,21 @@
* PROJ by Kristian Evers. Original code found in file src/proj_guyou.c, see
* https://github.com/rouault/libproj4/blob/master/libproject-1.01/src/proj_guyou.c
* for reference.
+* Fix for Peirce Quincuncial projection to diamond or square by Toby C. Wilkinson to
+* correctly flip out southern hemisphere into the four triangles of Peirce's
+* quincunx. The fix inspired by a similar rotate and translate solution applied
+* by Jonathan Feinberg for cartopy, see
+* https://github.com/jonathf/cartopy/blob/8172cac7fc45cafc86573d408ce85b74258a9c28/lib/cartopy/peircequincuncial.py
+* Added original code for horizontal and vertical arrangement of hemispheres by Toby
+* C. Wilkinson to allow creation of lateral quincuncial projections, such as Grieger's
+* Triptychial, see, e.g.:
+* - Grieger, B. (2020). “Optimized global map projections for specific applications:
+* the triptychial projection and the Spilhaus projection”. EGU2020-9885.
+* https://doi.org/10.5194/egusphere-egu2020-9885
*
* Copyright (c) 2005, 2006, 2009 Gerald I. Evenden
* Copyright (c) 2020 Kristian Evers
+* Copyright (c) 2021 Toby C Wilkinson
*
* Related material
* ----------------
@@ -48,8 +60,20 @@ enum projection_type {
ADAMS_WS2,
};
+enum peirce_type {
+ PEIRCE_Q_SQUARE,
+ PEIRCE_Q_DIAMOND,
+ PEIRCE_Q_NHEMISPHERE,
+ PEIRCE_Q_SHEMISPHERE,
+ PEIRCE_Q_HORIZONTAL,
+ PEIRCE_Q_VERTICAL,
+};
+
struct pj_opaque {
projection_type mode;
+ peirce_type pqtype;
+ double scrollx = 0.0;
+ double scrolly = 0.0;
};
} // anonymous namespace
@@ -118,9 +142,18 @@ static PJ_XY adams_forward(PJ_LP lp, PJ *P) {
}
break;
case PEIRCE_Q: {
- if( lp.phi < -TOL ) {
+ /* lam0 - note that the original Peirce model used a central meridian of around -70deg, but the default within proj is +lon0=0 */
+ if (Q->pqtype == PEIRCE_Q_NHEMISPHERE) {
+ if( lp.phi < -TOL ) {
+ proj_errno_set(P, PROJ_ERR_COORD_TRANSFM_OUTSIDE_PROJECTION_DOMAIN);
+ return proj_coord_error().xy;
+ }
+ }
+ if (Q->pqtype == PEIRCE_Q_SHEMISPHERE) {
+ if( lp.phi > -TOL ) {
proj_errno_set(P, PROJ_ERR_COORD_TRANSFM_OUTSIDE_PROJECTION_DOMAIN);
return proj_coord_error().xy;
+ }
}
const double sl = sin(lp.lam);
const double cl = cos(lp.lam);
@@ -173,7 +206,71 @@ static PJ_XY adams_forward(PJ_LP lp, PJ *P) {
xy.x = ell_int_5(m);
xy.y = ell_int_5(n);
- if (Q->mode == ADAMS_HEMI || Q->mode == ADAMS_WS2) { /* rotate by 45deg. */
+ if (Q->mode == PEIRCE_Q) {
+ /* Constant complete elliptic integral of the first kind with m=0.5, calculated using https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.ellipk.html . Used as basic as scaled shift distance */
+ constexpr double shd = 1.8540746773013719 * 2;
+
+ /* For square and diamond Quincuncial projections, spin out southern hemisphere to triangular segments of quincunx (before rotation for square)*/
+ if( Q->pqtype == PEIRCE_Q_SQUARE || ( Q->pqtype == PEIRCE_Q_DIAMOND )) {
+ if (lp.phi < 0.) { /* fold out segments */
+ if (lp.lam < ( -0.75 * M_PI )) xy.y = shd - xy.y; /* top left segment, shift up and reflect y */
+ if ( (lp.lam < (-0.25 * M_PI)) && (lp.lam >= ( -0.75 * M_PI ))) xy.x = - shd - xy.x; /* left segment, shift left and reflect x */
+ if ( (lp.lam < (0.25 * M_PI)) && (lp.lam >= ( -0.25 * M_PI ))) xy.y = - shd - xy.y; /* bottom segment, shift down and reflect y */
+ if ( (lp.lam < (0.75 * M_PI)) && (lp.lam >= ( 0.25 * M_PI ))) xy.x = shd - xy.x; /* right segment, shift right and reflect x */
+ if (lp.lam >= (0.75 * M_PI)) xy.y = shd - xy.y; /* top right segment, shift up and reflect y */
+ }
+ }
+
+ /* For square types rotate xy by 45 deg */
+ if( Q->pqtype == PEIRCE_Q_SQUARE ) {
+ const double temp = xy.x;
+ xy.x = RSQRT2 * (xy.x - xy.y);
+ xy.y = RSQRT2 * (temp + xy.y);
+ }
+
+ /* For rectangle Quincuncial projs, spin out southern hemisphere to east (horizontal) or north (vertical) after rotation */
+ if( Q->pqtype == PEIRCE_Q_HORIZONTAL ) {
+ if (lp.phi < 0.) {
+ xy.x = shd - xy.x; /* reflect x to east */
+ }
+ xy.x = xy.x - (shd / 2); /* shift everything so origin is in middle of two hemispheres */
+ }
+ if( Q->pqtype == PEIRCE_Q_VERTICAL ) {
+ if (lp.phi < 0.) {
+ xy.y = shd - xy.y; /* reflect y to north */
+ }
+ xy.y = xy.y - (shd / 2); /* shift everything so origin is in middle of two hemispheres */
+ }
+
+ //if o_scrollx param present, scroll x
+ if (!(Q->scrollx == 0.0) && (Q->pqtype == PEIRCE_Q_HORIZONTAL) ) {
+ double xscale = 2.0;
+ double xthresh = shd / 2;
+ xy.x = xy.x + (Q->scrollx * (xthresh * 2 * xscale)); /*shift relative to proj width*/
+ if(xy.x >= (xthresh * xscale)) {
+ xy.x = xy.x - (shd * xscale);
+ }
+ else if (xy.x < -(xthresh * xscale)) {
+ xy.x = xy.x + (shd * xscale);
+ }
+ }
+
+ //if o_scrolly param present, scroll y
+ if (!(Q->scrolly == 0.0) && (Q->pqtype == PEIRCE_Q_VERTICAL)) {
+ double yscale = 2.0;
+ double ythresh = shd / 2;
+ xy.y = xy.y + (Q->scrolly * (ythresh * 2 * yscale)); /*shift relative to proj height*/
+ if(xy.y >= (ythresh * yscale)) {
+ xy.y = xy.y - (shd * yscale);
+ }
+ else if (xy.y < -(ythresh * yscale)) {
+ xy.y = xy.y + (shd * yscale);
+ }
+ }
+
+ }
+
+ if (Q->mode == ADAMS_HEMI || Q->mode == ADAMS_WS2 ) { /* rotate by 45deg. */
const double temp = xy.x;
xy.x = RSQRT2 * (xy.x - xy.y);
xy.y = RSQRT2 * (temp + xy.y);
@@ -218,6 +315,55 @@ static PJ *setup(PJ *P, projection_type mode) {
if( mode == ADAMS_WS2 )
P->inv = adams_inverse;
+ if( mode == PEIRCE_Q) {
+ // Quincuncial projections type options: square, diamond, hemisphere, horizontal (rectangle) or vertical (rectangle)
+ const char* pqtype = pj_param (P->ctx, P->params, "stype").s;
+
+ if (!pqtype) pqtype = "diamond"; /* default if type value not supplied */
+
+ if (strcmp(pqtype, "square") == 0) {
+ Q->pqtype = PEIRCE_Q_SQUARE;
+ }
+ else if (strcmp(pqtype, "diamond") == 0) {
+ Q->pqtype = PEIRCE_Q_DIAMOND;
+ }
+ else if (strcmp(pqtype, "nhemisphere") == 0) {
+ Q->pqtype = PEIRCE_Q_NHEMISPHERE;
+ }
+ else if (strcmp(pqtype, "shemisphere") == 0) {
+ Q->pqtype = PEIRCE_Q_SHEMISPHERE;
+ }
+ else if (strcmp(pqtype, "horizontal") == 0) {
+ Q->pqtype = PEIRCE_Q_HORIZONTAL;
+ if (pj_param(P->ctx, P->params, "tscrollx").i) {
+ double scrollx;
+ scrollx = pj_param(P->ctx, P->params, "dscrollx").f;
+ if (scrollx > 1 || scrollx < -1) {
+ proj_log_error(P, _("Invalid value for scrollx: |scrollx| should between -1 and 1"));
+ return pj_default_destructor (P, PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
+ }
+ Q->scrollx = scrollx;
+ }
+ }
+ else if (strcmp(pqtype, "vertical") == 0) {
+ Q->pqtype = PEIRCE_Q_VERTICAL;
+ if (pj_param(P->ctx, P->params, "tscrolly").i) {
+ double scrolly;
+ scrolly = pj_param(P->ctx, P->params, "dscrolly").f;
+ if (scrolly > 1 || scrolly < -1) {
+ proj_log_error(P, _("Invalid value for scrolly: |scrolly| should between -1 and 1"));
+ return pj_default_destructor (P, PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
+ }
+ Q->scrolly = scrolly;
+ }
+ }
+ else {
+ proj_log_error (P, _("peirce_q: invalid value for 'type' parameter"));
+ return pj_default_destructor (P, PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
+ }
+
+ }
+
return P;
}
generated by cgit v1.2.3 (git 2.25.1) at 2025-12-15 20:40:28 +0000