From 911852371a6a109445cd79e5176a89c539b2e768 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Sun, 24 Feb 2019 15:58:16 +0100 Subject: Robinson: fix wrong values for forward path for latitudes >= 87.5 (fixes #1172), and fix inaccurate inverse method --- src/projections/robin.cpp | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) (limited to 'src/projections/robin.cpp') diff --git a/src/projections/robin.cpp b/src/projections/robin.cpp index 8f142aad..9f7908f6 100644 --- a/src/projections/robin.cpp +++ b/src/projections/robin.cpp @@ -7,7 +7,7 @@ PROJ_HEAD(robin, "Robinson") "\n\tPCyl, Sph"; #define V(C,z) (C.c0 + z * (C.c1 + z * (C.c2 + z * C.c3))) -#define DV(C,z) (C.c1 + z * (C.c2 + C.c2 + z * 3. * C.c3)) +#define DV(C,z) (C.c1 + 2 * z * C.c2 + z * z * 3. * C.c3) /* note: following terms based upon 5 deg. intervals in degrees. @@ -74,7 +74,7 @@ static const struct COEFS Y[] = { #define RC1 0.08726646259971647884 #define NODES 18 #define ONEEPS 1.000001 -#define EPS 1e-8 +#define EPS 1e-10 /* Not sure at all of the appropriate number for MAX_ITER... */ #define MAX_ITER 100 @@ -90,7 +90,7 @@ static PJ_XY s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - if (i >= NODES) i = NODES - 1; + if (i >= NODES) i = NODES; dphi = RAD_TO_DEG * (dphi - RC1 * i); xy.x = V(X[i], dphi) * FXC * lp.lam; xy.y = V(Y[i], dphi) * FYC; @@ -133,10 +133,8 @@ static PJ_LP s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ T = Y[i]; /* first guess, linear interp */ t = 5. * (lp.phi - T.c0)/(Y[i+1].c0 - T.c0); - /* make into root */ - T.c0 = (float)(T.c0 - lp.phi); for (iters = MAX_ITER; iters ; --iters) { /* Newton-Raphson */ - t -= t1 = V(T,t) / DV(T,t); + t -= t1 = (V(T,t) - lp.phi) / DV(T,t); if (fabs(t1) < EPS) break; } -- cgit v1.2.3 From 095d2204f8bb05d172936aebbb1e9e44852c049f Mon Sep 17 00:00:00 2001 From: Chris Mayo Date: Fri, 29 Mar 2019 19:17:37 +0000 Subject: Remove duplicate instances of #include "proj_internal.h" Introduced by "Merge projects.h into proj_internal.h" 8ab6f683. --- src/projections/robin.cpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) (limited to 'src/projections/robin.cpp') diff --git a/src/projections/robin.cpp b/src/projections/robin.cpp index 9f7908f6..dfb750dd 100644 --- a/src/projections/robin.cpp +++ b/src/projections/robin.cpp @@ -1,8 +1,7 @@ #define PJ_LIB__ -#include "proj_math.h" -#include "proj_internal.h" #include "proj.h" #include "proj_internal.h" +#include "proj_math.h" PROJ_HEAD(robin, "Robinson") "\n\tPCyl, Sph"; -- cgit v1.2.3 From 00980bf63fae6d350f425c44a648f33d7c09a931 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Fri, 26 Apr 2019 18:18:30 +0200 Subject: Prefix inverse and forward functions by their projection names This is mostly to have better OSSFuzz report. Currently a lot of bug summaries are like `proj4/standard_fuzzer: Divide-by-zero in s_inverse` By prefixing the projection name, we will get better reports, like `Divide-by-zero in airy_s_inverse` This also makes it slightly easier to set a breakpoint by function name. --- src/projections/robin.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) (limited to 'src/projections/robin.cpp') diff --git a/src/projections/robin.cpp b/src/projections/robin.cpp index dfb750dd..c08ac0e2 100644 --- a/src/projections/robin.cpp +++ b/src/projections/robin.cpp @@ -77,7 +77,7 @@ static const struct COEFS Y[] = { /* Not sure at all of the appropriate number for MAX_ITER... */ #define MAX_ITER 100 -static PJ_XY s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ +static PJ_XY robin_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; long i; double dphi; @@ -99,7 +99,7 @@ static PJ_XY s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ } -static PJ_LP s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ +static PJ_LP robin_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; long i; double t, t1; @@ -149,8 +149,8 @@ static PJ_LP s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ *PROJECTION(robin) { P->es = 0.; - P->inv = s_inverse; - P->fwd = s_forward; + P->inv = robin_s_inverse; + P->fwd = robin_s_forward; return P; } -- cgit v1.2.3