From 6b0567b3c01ada8412a446310af9784124453099 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Philippe=20Rivi=C3=A8re?= Date: Fri, 21 Sep 2018 17:55:08 +0200 Subject: the Bertin 1953 projection --- src/Makefile.am | 2 +- src/PJ_bertin1953.c | 108 ++++++++++++++++++++++++++++++++++++++++++++++++++++ src/lib_proj.cmake | 1 + src/makefile.vc | 2 +- src/pj_list.h | 1 + 5 files changed, 112 insertions(+), 2 deletions(-) create mode 100644 src/PJ_bertin1953.c (limited to 'src') diff --git a/src/Makefile.am b/src/Makefile.am index 84846c1c..e8bea068 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -55,7 +55,7 @@ libproj_la_SOURCES = \ PJ_mill.c PJ_ocea.c PJ_omerc.c PJ_somerc.c \ PJ_tcc.c PJ_tcea.c PJ_times.c PJ_tmerc.c \ PJ_airy.c PJ_aitoff.c PJ_august.c PJ_bacon.c \ - PJ_chamb.c PJ_hammer.c PJ_lagrng.c PJ_larr.c \ + PJ_bertin1953.c PJ_chamb.c PJ_hammer.c PJ_lagrng.c PJ_larr.c \ PJ_lask.c PJ_latlong.c PJ_nocol.c PJ_ob_tran.c PJ_oea.c \ PJ_tpeqd.c PJ_vandg.c PJ_vandg2.c PJ_vandg4.c \ PJ_wag7.c PJ_lcca.c PJ_geos.c proj_etmerc.c \ diff --git a/src/PJ_bertin1953.c b/src/PJ_bertin1953.c new file mode 100644 index 00000000..19aa4d53 --- /dev/null +++ b/src/PJ_bertin1953.c @@ -0,0 +1,108 @@ +/* + Created by Jacques Bertin in 1953, this projection was the go-to choice + of the French cartographic school when they wished to represent phenomena + on a global scale. + + Formula designed by Philippe Rivière, 2017. + https://visionscarto.net/bertin-projection-1953 + + Port to PROJ by Philippe Rivière, 21 September 2018 +*/ + +#define PJ_LIB__ + +#include +#include + +#include "proj.h" +#include "projects.h" + +PROJ_HEAD(bertin1953, "Bertin 1953") + "\n\tMisc Sph no inv."; + +struct pj_opaque { + double w; + double m, rm; + double cosDeltaPhi, sinDeltaPhi, cosDeltaGamma, sinDeltaGamma, deltaLambda; +}; + + +static XY s_forward (LP lp, PJ *P) { + XY xy = {0.0,0.0}; + struct pj_opaque *Q = P->opaque; + + // Projection constants + double fu = 1.4, k = 12., w = 1.68; + + // Aliases + double lambda = lp.lam, phi = lp.phi; + + // Variable + double d; + + // Apply rotation + lambda += Q->deltaLambda; + + double cosphi = cos(phi), + x = cos(lambda) * cosphi, + y = sin(lambda) * cosphi, + z = sin(phi), + z0 = z * Q->cosDeltaPhi + x * Q->sinDeltaPhi; + lambda = atan2(y * Q->cosDeltaGamma - z0 * Q->sinDeltaGamma, + x * Q->cosDeltaPhi - z * Q->sinDeltaPhi); + z0 = z0 * Q->cosDeltaGamma + y * Q->sinDeltaGamma; + phi = asin(z0); + + lambda = adjlon(lambda); + + // Adjust pre-projection + if (lambda + phi < -fu) { + double u = (lambda - phi + 1.6) * (lambda + phi + fu) / 8.; + lambda += u; + phi -= 0.8 * u * sin(phi + M_PI / 2.); + } + + // Project with Hammer (1.68,2) + cosphi = cos(phi); + d = sqrt(2./(1. + cosphi * cos(lambda / 2.))); + xy.x = w * d * cosphi * sin(lambda / 2.); + xy.y = d * sin(phi); + + // Adjust post-projection + d = (1. - cos(lambda * phi)) / k; + if (xy.y < 0.) { + xy.x *= 1. + d; + } + if (xy.y > 0.) { + xy.x *= 1. + d / 1.5 * xy.x * xy.x; + } + + return xy; +} + + +PJ *PROJECTION(bertin1953) { + struct pj_opaque *Q = pj_calloc (1, sizeof (struct pj_opaque)); + if (0==Q) + return pj_default_destructor (P, ENOMEM); + P->opaque = Q; + + // force +lon_0 = -16.5 + double deltaLambda = P->lam0 -16.5 / 180. * M_PI; + + // force +lat_0=-42 + double deltaPhi = -42. / 180. * M_PI, + deltaGamma = 0. / 180. * M_PI; + + Q->deltaLambda = deltaLambda; + Q->cosDeltaPhi = cos(deltaPhi); + Q->sinDeltaPhi = sin(deltaPhi); + Q->cosDeltaGamma = cos(deltaGamma); + Q->sinDeltaGamma = sin(deltaGamma); + + P->es = 0.; + P->fwd = s_forward; + // P->inv = s_inverse; + + return P; +} diff --git a/src/lib_proj.cmake b/src/lib_proj.cmake index 0779568a..29f3e94f 100644 --- a/src/lib_proj.cmake +++ b/src/lib_proj.cmake @@ -45,6 +45,7 @@ SET(SRC_LIBPROJ_PJ PJ_august.c PJ_axisswap.c PJ_bacon.c + PJ_bertin1953.c PJ_bipc.c PJ_boggs.c PJ_bonne.c diff --git a/src/makefile.vc b/src/makefile.vc index cf9d878c..fc6a6861 100644 --- a/src/makefile.vc +++ b/src/makefile.vc @@ -21,7 +21,7 @@ cylinder = \ PJ_gstmerc.obj proj_etmerc.obj PJ_comill.obj misc = \ - PJ_airy.obj PJ_aitoff.obj PJ_august.obj PJ_bacon.obj \ + PJ_airy.obj PJ_aitoff.obj PJ_august.obj PJ_bacon.obj PJ_bertin1953.obj \ PJ_chamb.obj PJ_hammer.obj PJ_lagrng.obj PJ_larr.obj \ PJ_lask.obj PJ_nocol.obj PJ_ob_tran.obj PJ_oea.obj \ PJ_sch.obj PJ_tpeqd.obj PJ_vandg.obj PJ_vandg2.obj \ diff --git a/src/pj_list.h b/src/pj_list.h index 0a8d9d1b..b4fc2306 100644 --- a/src/pj_list.h +++ b/src/pj_list.h @@ -14,6 +14,7 @@ PROJ_HEAD(apian, "Apian Globular I") PROJ_HEAD(august, "August Epicycloidal") PROJ_HEAD(axisswap, "Axis ordering") PROJ_HEAD(bacon, "Bacon Globular") +PROJ_HEAD(bertin1953, "Bertin 1953") PROJ_HEAD(bipc, "Bipolar conic of western hemisphere") PROJ_HEAD(boggs, "Boggs Eumorphic") PROJ_HEAD(bonne, "Bonne (Werner lat_1=90)") -- cgit v1.2.3 From c9962986f4c39219b614872ef5af9c8c7b8c7503 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Philippe=20Rivi=C3=A8re?= Date: Fri, 21 Sep 2018 22:42:58 +0200 Subject: coding style (https://travis-ci.com/OSGeo/proj.4/jobs/147225960 & https://travis-ci.com/OSGeo/proj.4/jobs/147274068) --- src/PJ_bertin1953.c | 45 +++++++++++++++++++++++---------------------- 1 file changed, 23 insertions(+), 22 deletions(-) (limited to 'src') diff --git a/src/PJ_bertin1953.c b/src/PJ_bertin1953.c index 19aa4d53..f273f0e9 100644 --- a/src/PJ_bertin1953.c +++ b/src/PJ_bertin1953.c @@ -31,23 +31,23 @@ static XY s_forward (LP lp, PJ *P) { XY xy = {0.0,0.0}; struct pj_opaque *Q = P->opaque; - // Projection constants + /* Projection constants */ double fu = 1.4, k = 12., w = 1.68; - // Aliases + /* Aliases */ double lambda = lp.lam, phi = lp.phi; - // Variable + /* Variable */ double d; - // Apply rotation + /* Apply rotation */ + double cosphi, x, y, z, z0; lambda += Q->deltaLambda; - - double cosphi = cos(phi), - x = cos(lambda) * cosphi, - y = sin(lambda) * cosphi, - z = sin(phi), - z0 = z * Q->cosDeltaPhi + x * Q->sinDeltaPhi; + cosphi = cos(phi); + x = cos(lambda) * cosphi; + y = sin(lambda) * cosphi; + z = sin(phi); + z0 = z * Q->cosDeltaPhi + x * Q->sinDeltaPhi; lambda = atan2(y * Q->cosDeltaGamma - z0 * Q->sinDeltaGamma, x * Q->cosDeltaPhi - z * Q->sinDeltaPhi); z0 = z0 * Q->cosDeltaGamma + y * Q->sinDeltaGamma; @@ -55,20 +55,20 @@ static XY s_forward (LP lp, PJ *P) { lambda = adjlon(lambda); - // Adjust pre-projection + /* Adjust pre-projection */ if (lambda + phi < -fu) { - double u = (lambda - phi + 1.6) * (lambda + phi + fu) / 8.; - lambda += u; - phi -= 0.8 * u * sin(phi + M_PI / 2.); + d = (lambda - phi + 1.6) * (lambda + phi + fu) / 8.; + lambda += d; + phi -= 0.8 * d * sin(phi + M_PI / 2.); } - // Project with Hammer (1.68,2) + /* Project with Hammer (1.68,2) */ cosphi = cos(phi); d = sqrt(2./(1. + cosphi * cos(lambda / 2.))); xy.x = w * d * cosphi * sin(lambda / 2.); xy.y = d * sin(phi); - // Adjust post-projection + /* Adjust post-projection */ d = (1. - cos(lambda * phi)) / k; if (xy.y < 0.) { xy.x *= 1. + d; @@ -82,17 +82,18 @@ static XY s_forward (LP lp, PJ *P) { PJ *PROJECTION(bertin1953) { + double deltaLambda, deltaPhi, deltaGamma; struct pj_opaque *Q = pj_calloc (1, sizeof (struct pj_opaque)); if (0==Q) return pj_default_destructor (P, ENOMEM); P->opaque = Q; - // force +lon_0 = -16.5 - double deltaLambda = P->lam0 -16.5 / 180. * M_PI; + /* force +lon_0 = -16.5 */ + deltaLambda = P->lam0 -16.5 / 180. * M_PI; - // force +lat_0=-42 - double deltaPhi = -42. / 180. * M_PI, - deltaGamma = 0. / 180. * M_PI; + /* force +lat_0=-42 */ + deltaPhi = -42. / 180. * M_PI; + deltaGamma = 0. / 180. * M_PI; Q->deltaLambda = deltaLambda; Q->cosDeltaPhi = cos(deltaPhi); @@ -102,7 +103,7 @@ PJ *PROJECTION(bertin1953) { P->es = 0.; P->fwd = s_forward; - // P->inv = s_inverse; + /* P->inv = s_inverse; */ return P; } -- cgit v1.2.3 From 5f4fa9ccba63ac89e09035688588ced7649ff43f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Philippe=20Rivi=C3=A8re?= Date: Mon, 24 Sep 2018 13:00:31 +0200 Subject: Req. changes for Bertin1953: - classification - tests - coding style --- src/PJ_bertin1953.c | 67 +++++++++++++++++++++-------------------------------- 1 file changed, 26 insertions(+), 41 deletions(-) (limited to 'src') diff --git a/src/PJ_bertin1953.c b/src/PJ_bertin1953.c index f273f0e9..a7c60903 100644 --- a/src/PJ_bertin1953.c +++ b/src/PJ_bertin1953.c @@ -14,6 +14,7 @@ #include #include +#include "proj_internal.h" #include "proj.h" #include "projects.h" @@ -21,8 +22,6 @@ PROJ_HEAD(bertin1953, "Bertin 1953") "\n\tMisc Sph no inv."; struct pj_opaque { - double w; - double m, rm; double cosDeltaPhi, sinDeltaPhi, cosDeltaGamma, sinDeltaGamma, deltaLambda; }; @@ -31,45 +30,38 @@ static XY s_forward (LP lp, PJ *P) { XY xy = {0.0,0.0}; struct pj_opaque *Q = P->opaque; - /* Projection constants */ - double fu = 1.4, k = 12., w = 1.68; + double fu = 1.4, k = 12., w = 1.68, d; - /* Aliases */ - double lambda = lp.lam, phi = lp.phi; - - /* Variable */ - double d; - - /* Apply rotation */ + /* Rotate */ double cosphi, x, y, z, z0; - lambda += Q->deltaLambda; - cosphi = cos(phi); - x = cos(lambda) * cosphi; - y = sin(lambda) * cosphi; - z = sin(phi); + lp.lam += PJ_TORAD(-16.5); + cosphi = cos(lp.phi); + x = cos(lp.lam) * cosphi; + y = sin(lp.lam) * cosphi; + z = sin(lp.phi); z0 = z * Q->cosDeltaPhi + x * Q->sinDeltaPhi; - lambda = atan2(y * Q->cosDeltaGamma - z0 * Q->sinDeltaGamma, + lp.lam = atan2(y * Q->cosDeltaGamma - z0 * Q->sinDeltaGamma, x * Q->cosDeltaPhi - z * Q->sinDeltaPhi); z0 = z0 * Q->cosDeltaGamma + y * Q->sinDeltaGamma; - phi = asin(z0); + lp.phi = asin(z0); - lambda = adjlon(lambda); + lp.lam = adjlon(lp.lam); /* Adjust pre-projection */ - if (lambda + phi < -fu) { - d = (lambda - phi + 1.6) * (lambda + phi + fu) / 8.; - lambda += d; - phi -= 0.8 * d * sin(phi + M_PI / 2.); + if (lp.lam + lp.phi < -fu) { + d = (lp.lam - lp.phi + 1.6) * (lp.lam + lp.phi + fu) / 8.; + lp.lam += d; + lp.phi -= 0.8 * d * sin(lp.phi + M_PI / 2.); } /* Project with Hammer (1.68,2) */ - cosphi = cos(phi); - d = sqrt(2./(1. + cosphi * cos(lambda / 2.))); - xy.x = w * d * cosphi * sin(lambda / 2.); - xy.y = d * sin(phi); + cosphi = cos(lp.phi); + d = sqrt(2./(1. + cosphi * cos(lp.lam / 2.))); + xy.x = w * d * cosphi * sin(lp.lam / 2.); + xy.y = d * sin(lp.phi); /* Adjust post-projection */ - d = (1. - cos(lambda * phi)) / k; + d = (1. - cos(lp.lam * lp.phi)) / k; if (xy.y < 0.) { xy.x *= 1. + d; } @@ -82,28 +74,21 @@ static XY s_forward (LP lp, PJ *P) { PJ *PROJECTION(bertin1953) { - double deltaLambda, deltaPhi, deltaGamma; struct pj_opaque *Q = pj_calloc (1, sizeof (struct pj_opaque)); if (0==Q) return pj_default_destructor (P, ENOMEM); P->opaque = Q; - /* force +lon_0 = -16.5 */ - deltaLambda = P->lam0 -16.5 / 180. * M_PI; - - /* force +lat_0=-42 */ - deltaPhi = -42. / 180. * M_PI; - deltaGamma = 0. / 180. * M_PI; + P->lam0 = 0; + P->phi0 = PJ_TORAD(-42.); - Q->deltaLambda = deltaLambda; - Q->cosDeltaPhi = cos(deltaPhi); - Q->sinDeltaPhi = sin(deltaPhi); - Q->cosDeltaGamma = cos(deltaGamma); - Q->sinDeltaGamma = sin(deltaGamma); + Q->cosDeltaPhi = cos(P->phi0); + Q->sinDeltaPhi = sin(P->phi0); + Q->cosDeltaGamma = 1.; + Q->sinDeltaGamma = 0.; P->es = 0.; P->fwd = s_forward; - /* P->inv = s_inverse; */ return P; } -- cgit v1.2.3 From e649eb7369caa87f0e5e8076fb7bffe72485b646 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Philippe=20Rivi=C3=A8re?= Date: Mon, 24 Sep 2018 13:05:45 +0200 Subject: snake case --- src/PJ_bertin1953.c | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) (limited to 'src') diff --git a/src/PJ_bertin1953.c b/src/PJ_bertin1953.c index a7c60903..5a027da2 100644 --- a/src/PJ_bertin1953.c +++ b/src/PJ_bertin1953.c @@ -22,7 +22,7 @@ PROJ_HEAD(bertin1953, "Bertin 1953") "\n\tMisc Sph no inv."; struct pj_opaque { - double cosDeltaPhi, sinDeltaPhi, cosDeltaGamma, sinDeltaGamma, deltaLambda; + double cos_delta_phi, sin_delta_phi, cos_delta_gamma, sin_delta_gamma, deltaLambda; }; @@ -39,10 +39,10 @@ static XY s_forward (LP lp, PJ *P) { x = cos(lp.lam) * cosphi; y = sin(lp.lam) * cosphi; z = sin(lp.phi); - z0 = z * Q->cosDeltaPhi + x * Q->sinDeltaPhi; - lp.lam = atan2(y * Q->cosDeltaGamma - z0 * Q->sinDeltaGamma, - x * Q->cosDeltaPhi - z * Q->sinDeltaPhi); - z0 = z0 * Q->cosDeltaGamma + y * Q->sinDeltaGamma; + z0 = z * Q->cos_delta_phi + x * Q->sin_delta_phi; + lp.lam = atan2(y * Q->cos_delta_gamma - z0 * Q->sin_delta_gamma, + x * Q->cos_delta_phi - z * Q->sin_delta_phi); + z0 = z0 * Q->cos_delta_gamma + y * Q->sin_delta_gamma; lp.phi = asin(z0); lp.lam = adjlon(lp.lam); @@ -82,10 +82,10 @@ PJ *PROJECTION(bertin1953) { P->lam0 = 0; P->phi0 = PJ_TORAD(-42.); - Q->cosDeltaPhi = cos(P->phi0); - Q->sinDeltaPhi = sin(P->phi0); - Q->cosDeltaGamma = 1.; - Q->sinDeltaGamma = 0.; + Q->cos_delta_phi = cos(P->phi0); + Q->sin_delta_phi = sin(P->phi0); + Q->cos_delta_gamma = 1.; + Q->sin_delta_gamma = 0.; P->es = 0.; P->fwd = s_forward; -- cgit v1.2.3