diff options
| author | Even Rouault <even.rouault@spatialys.com> | 2020-03-06 20:38:38 +0100 |
|---|---|---|
| committer | Even Rouault <even.rouault@spatialys.com> | 2020-03-08 19:58:35 +0100 |
| commit | 60d3df673ca224107eb63e459073fc11ab5f4f16 (patch) | |
| tree | 52c871c17f0c98c90ce6415cdb79dc1e73c14f4d | |
| parent | 38ec5e662a74d40e02e38dc5dca553c3ecb04356 (diff) | |
| download | PROJ-60d3df673ca224107eb63e459073fc11ab5f4f16.tar.gz PROJ-60d3df673ca224107eb63e459073fc11ab5f4f16.zip | |
src/projections/: remove assignments in expression and multiple statements per line
Should hopefully result in no change in results, and hopefully more
readable code...
64 files changed, 530 insertions, 397 deletions
diff --git a/src/projections/aea.cpp b/src/projections/aea.cpp index c9d24245..6ffb4fd6 100644 --- a/src/projections/aea.cpp +++ b/src/projections/aea.cpp @@ -111,7 +111,8 @@ static PJ_XY aea_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoid/spheroid, forward return xy; } Q->rho = Q->dd * sqrt(Q->rho); - xy.x = Q->rho * sin( lp.lam *= Q->n ); + lp.lam *= Q->n; + xy.x = Q->rho * sin( lp.lam ); xy.y = Q->rho0 - Q->rho * cos(lp.lam); return xy; } @@ -120,7 +121,9 @@ static PJ_XY aea_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoid/spheroid, forward static PJ_LP aea_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoid/spheroid, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - if( (Q->rho = hypot(xy.x, xy.y = Q->rho0 - xy.y)) != 0.0 ) { + xy.y = Q->rho0 - xy.y; + Q->rho = hypot(xy.x, xy.y); + if( Q->rho != 0.0 ) { if (Q->n < 0.) { Q->rho = -Q->rho; xy.x = -xy.x; @@ -134,16 +137,20 @@ static PJ_LP aea_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoid/spheroid, inverse proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } - if ((lp.phi = phi1_(lp.phi, P->e, P->one_es)) == HUGE_VAL) { + lp.phi = phi1_(lp.phi, P->e, P->one_es); + if (lp.phi == HUGE_VAL) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } } else lp.phi = lp.phi < 0. ? -M_HALFPI : M_HALFPI; - } else if (fabs(lp.phi = (Q->c - lp.phi * lp.phi) / Q->n2) <= 1.) - lp.phi = asin(lp.phi); - else - lp.phi = lp.phi < 0. ? -M_HALFPI : M_HALFPI; + } else { + lp.phi = (Q->c - lp.phi * lp.phi) / Q->n2; + if (fabs(lp.phi) <= 1.) + lp.phi = asin(lp.phi); + else + lp.phi = lp.phi < 0. ? -M_HALFPI : M_HALFPI; + } lp.lam = atan2(xy.x, xy.y) / Q->n; } else { lp.lam = 0.; @@ -155,8 +162,6 @@ static PJ_LP aea_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoid/spheroid, inverse static PJ *setup(PJ *P) { - double cosphi, sinphi; - int secant; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); P->inv = aea_e_inverse; @@ -166,13 +171,16 @@ static PJ *setup(PJ *P) { return destructor(P, PJD_ERR_LAT_LARGER_THAN_90); if (fabs(Q->phi1 + Q->phi2) < EPS10) return destructor(P, PJD_ERR_CONIC_LAT_EQUAL); - Q->n = sinphi = sin(Q->phi1); - cosphi = cos(Q->phi1); - secant = fabs(Q->phi1 - Q->phi2) >= EPS10; - if( (Q->ellips = (P->es > 0.))) { + double sinphi = sin(Q->phi1); + Q->n = sinphi; + double cosphi = cos(Q->phi1); + const int secant = fabs(Q->phi1 - Q->phi2) >= EPS10; + Q->ellips = (P->es > 0.); + if( Q->ellips ) { double ml1, m1; - if (!(Q->en = pj_enfn(P->es))) + Q->en = pj_enfn(P->es); + if (Q->en == nullptr) return destructor(P, 0); m1 = pj_msfn(sinphi, cosphi, P->es); ml1 = pj_qsfn(sinphi, P->e, P->one_es); diff --git a/src/projections/aeqd.cpp b/src/projections/aeqd.cpp index e0fa5c37..ad8c289e 100644 --- a/src/projections/aeqd.cpp +++ b/src/projections/aeqd.cpp @@ -106,8 +106,8 @@ static PJ_XY aeqd_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward coslam = - coslam; /*-fallthrough*/ case S_POLE: - xy.x = (rho = fabs(Q->Mp - pj_mlfn(lp.phi, sinphi, cosphi, Q->en))) * - sin(lp.lam); + rho = fabs(Q->Mp - pj_mlfn(lp.phi, sinphi, cosphi, Q->en)); + xy.x = rho * sin(lp.lam); xy.y = rho * coslam; break; case EQUIT: @@ -117,8 +117,10 @@ static PJ_XY aeqd_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward break; } - phi1 = P->phi0 / DEG_TO_RAD; lam1 = P->lam0 / DEG_TO_RAD; - phi2 = lp.phi / DEG_TO_RAD; lam2 = (lp.lam+P->lam0) / DEG_TO_RAD; + phi1 = P->phi0 / DEG_TO_RAD; + lam1 = P->lam0 / DEG_TO_RAD; + phi2 = lp.phi / DEG_TO_RAD; + lam2 = (lp.lam+P->lam0) / DEG_TO_RAD; geod_inverse(&Q->g, phi1, lam1, phi2, lam2, &s12, &azi1, &azi2); azi1 *= DEG_TO_RAD; @@ -169,7 +171,8 @@ oblcon: proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - xy.x = (xy.y = (M_HALFPI + lp.phi)) * sin(lp.lam); + xy.y = (M_HALFPI + lp.phi); + xy.x = xy.y * sin(lp.lam); xy.y *= coslam; break; } @@ -187,8 +190,9 @@ static PJ_LP e_guam_inv(PJ_XY xy, PJ *P) { /* Guam elliptical */ lp.phi = P->phi0; for (i = 0; i < 3; ++i) { t = P->e * sin(lp.phi); + t = sqrt(1. - t * t); lp.phi = pj_inv_mlfn(P->ctx, Q->M1 + xy.y - - x2 * tan(lp.phi) * (t = sqrt(1. - t * t)), P->es, Q->en); + x2 * tan(lp.phi) * t, P->es, Q->en); } lp.lam = xy.x * t / cos(lp.phi); return lp; @@ -232,7 +236,8 @@ static PJ_LP aeqd_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); double cosc, c_rh, sinc; - if ((c_rh = hypot(xy.x, xy.y)) > M_PI) { + c_rh = hypot(xy.x, xy.y); + if (c_rh > M_PI) { if (c_rh - EPS10 > M_PI) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; @@ -247,7 +252,7 @@ static PJ_LP aeqd_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse sinc = sin(c_rh); cosc = cos(c_rh); if (Q->mode == EQUIT) { - lp.phi = aasin(P->ctx, xy.y * sinc / c_rh); + lp.phi = aasin(P->ctx, xy.y * sinc / c_rh); xy.x *= sinc; xy.y = cosc * c_rh; } else { @@ -311,7 +316,8 @@ PJ *PROJECTION(aeqd) { case EQUIT: case OBLIQ: Q->N1 = 1. / sqrt(1. - P->es * Q->sinph0 * Q->sinph0); - Q->G = Q->sinph0 * (Q->He = P->e / sqrt(P->one_es)); + Q->He = P->e / sqrt(P->one_es); + Q->G = Q->sinph0 * Q->He; Q->He *= Q->cosph0; break; } diff --git a/src/projections/airy.cpp b/src/projections/airy.cpp index 91b4b7a2..f7602e53 100644 --- a/src/projections/airy.cpp +++ b/src/projections/airy.cpp @@ -77,7 +77,8 @@ static PJ_XY airy_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - if (fabs(s = 1. - cosz) > EPS) { + s = 1. - cosz; + if (fabs(s) > EPS) { t = 0.5 * (1. + cosz); if(t == 0) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); @@ -88,8 +89,7 @@ static PJ_XY airy_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward Krho = 0.5 - Q->Cb; xy.x = Krho * cosphi * sinlam; if (Q->mode == OBLIQ) - xy.y = Krho * (Q->cosph0 * sinphi - - Q->sinph0 * cosphi * coslam); + xy.y = Krho * (Q->cosph0 * sinphi - Q->sinph0 * cosphi * coslam); else xy.y = Krho * sinphi; break; @@ -100,7 +100,8 @@ static PJ_XY airy_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - if ((lp.phi *= 0.5) > EPS) { + lp.phi *= 0.5; + if (lp.phi > EPS) { t = tan(lp.phi); Krho = -2.*(log(cos(lp.phi)) / t + t * Q->Cb); xy.x = Krho * sinlam; diff --git a/src/projections/aitoff.cpp b/src/projections/aitoff.cpp index dadae12d..7920309c 100644 --- a/src/projections/aitoff.cpp +++ b/src/projections/aitoff.cpp @@ -67,7 +67,9 @@ static PJ_XY aitoff_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwa struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); double c, d; - if((d = acos(cos(lp.phi) * cos(c = 0.5 * lp.lam))) != 0.0) {/* basic Aitoff */ + c = 0.5 * lp.lam; + d = acos(cos(lp.phi) * cos(c)); + if(d != 0.0) {/* basic Aitoff */ xy.x = 2. * d * cos(lp.phi) * sin(c) * (xy.y = 1. / sin(d)); xy.y *= d * sin(lp.phi); } else @@ -113,8 +115,10 @@ static PJ_LP aitoff_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver do { iter = 0; do { - sl = sin(lp.lam * 0.5); cl = cos(lp.lam * 0.5); - sp = sin(lp.phi); cp = cos(lp.phi); + sl = sin(lp.lam * 0.5); + cl = cos(lp.lam * 0.5); + sp = sin(lp.phi); + cp = cos(lp.phi); D = cp * cl; C = 1. - D * D; const double denom = pow(C, 1.5); @@ -137,11 +141,14 @@ static PJ_LP aitoff_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver f2p = 0.5 * (f2p + 1.); f2l *= 0.5; } - f1 -= xy.x; f2 -= xy.y; - dl = (f2 * f1p - f1 * f2p) / (dp = f1p * f2l - f2p * f1l); + f1 -= xy.x; + f2 -= xy.y; + dp = f1p * f2l - f2p * f1l; + dl = (f2 * f1p - f1 * f2p) / dp; dp = (f1 * f2l - f2 * f1l) / dp; dl = fmod(dl, M_PI); /* set to interval [-M_PI, M_PI] */ - lp.phi -= dp; lp.lam -= dl; + lp.phi -= dp; + lp.lam -= dl; } while ((fabs(dp) > EPSILON || fabs(dl) > EPSILON) && (iter++ < MAXITER)); if (lp.phi > M_PI_2) lp.phi -= 2.*(lp.phi-M_PI_2); /* correct if symmetrical solution for Aitoff */ if (lp.phi < -M_PI_2) lp.phi -= 2.*(lp.phi+M_PI_2); /* correct if symmetrical solution for Aitoff */ @@ -149,7 +156,8 @@ static PJ_LP aitoff_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver /* calculate x,y coordinates with solution obtained */ if((D = acos(cos(lp.phi) * cos(C = 0.5 * lp.lam))) != 0.0) {/* Aitoff */ - x = 2. * D * cos(lp.phi) * sin(C) * (y = 1. / sin(D)); + y = 1. / sin(D); + x = 2. * D * cos(lp.phi) * sin(C) * y; y *= D * sin(lp.phi); } else x = y = 0.; diff --git a/src/projections/august.cpp b/src/projections/august.cpp index 2104c3cc..049c6d90 100644 --- a/src/projections/august.cpp +++ b/src/projections/august.cpp @@ -11,17 +11,20 @@ PROJ_HEAD(august, "August Epicycloidal") "\n\tMisc Sph, no inv"; static PJ_XY august_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double t, c1, c, x1, x12, y1, y12; - (void) P; - - t = tan(.5 * lp.phi); - c1 = sqrt(1. - t * t); - c = 1. + c1 * cos(lp.lam *= .5); - x1 = sin(lp.lam) * c1 / c; - y1 = t / c; - xy.x = M * x1 * (3. + (x12 = x1 * x1) - 3. * (y12 = y1 * y1)); - xy.y = M * y1 * (3. + 3. * x12 - y12); - return (xy); + double t, c1, c, x1, x12, y1, y12; + (void) P; + + t = tan(.5 * lp.phi); + c1 = sqrt(1. - t * t); + lp.lam *= .5; + c = 1. + c1 * cos(lp.lam); + x1 = sin(lp.lam) * c1 / c; + y1 = t / c; + x12 = x1 * x1; + y12 = y1 * y1; + xy.x = M * x1 * (3. + x12 - 3. * y12); + xy.y = M * y1 * (3. + 3. * x12 - y12); + return (xy); } diff --git a/src/projections/bacon.cpp b/src/projections/bacon.cpp index 5db2d854..3efd4dbe 100644 --- a/src/projections/bacon.cpp +++ b/src/projections/bacon.cpp @@ -23,20 +23,21 @@ PROJ_HEAD(bacon, "Bacon Globular") "\n\tMisc Sph, no inv"; static PJ_XY bacon_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double ax, f; - - xy.y = Q->bacn ? M_HALFPI * sin(lp.phi) : lp.phi; - if ((ax = fabs(lp.lam)) >= EPS) { - if (Q->ortl && ax >= M_HALFPI) - xy.x = sqrt(HLFPI2 - lp.phi * lp.phi + EPS) + ax - M_HALFPI; - else { - f = 0.5 * (HLFPI2 / ax + ax); - xy.x = ax - f + sqrt(f * f - xy.y * xy.y); - } - if (lp.lam < 0.) xy.x = - xy.x; - } else - xy.x = 0.; - return (xy); + double ax, f; + + xy.y = Q->bacn ? M_HALFPI * sin(lp.phi) : lp.phi; + ax = fabs(lp.lam); + if (ax >= EPS) { + if (Q->ortl && ax >= M_HALFPI) + xy.x = sqrt(HLFPI2 - lp.phi * lp.phi + EPS) + ax - M_HALFPI; + else { + f = 0.5 * (HLFPI2 / ax + ax); + xy.x = ax - f + sqrt(f * f - xy.y * xy.y); + } + if (lp.lam < 0.) xy.x = - xy.x; + } else + xy.x = 0.; + return (xy); } diff --git a/src/projections/bipc.cpp b/src/projections/bipc.cpp index 9e991ecc..bf4ba834 100644 --- a/src/projections/bipc.cpp +++ b/src/projections/bipc.cpp @@ -54,7 +54,8 @@ static PJ_XY bipc_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward Az = atan2(sdlam , C45 * (tphi - cdlam)); } if( (tag = (Az > Azba)) ) { - cdlam = cos(sdlam = lp.lam + R110); + sdlam = lp.lam + R110; + cdlam = cos(sdlam); sdlam = sin(sdlam); z = S20 * sphi + C20 * cphi * cdlam; if (fabs(z) > 1.) { @@ -86,7 +87,8 @@ static PJ_XY bipc_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - r = F * (t = pow(tan(.5 * z), n)); + t = pow(tan(.5 * z), n); + r = F * t; if ((al = .5 * (R104 - z)) < 0.) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; @@ -100,7 +102,8 @@ static PJ_XY bipc_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward else al = al < 0. ? -1. : 1.; } else al = acos(al); - if (fabs(t = n * (Av - Az)) < al) + t = n * (Av - Az); + if (fabs(t) < al) r /= cos(al + (tag ? t : -t)); xy.x = r * sin(t); xy.y += (tag ? -r : r) * cos(t); @@ -135,8 +138,10 @@ static PJ_LP bipc_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse c = C45; Av = Azba; } - rl = rp = r = hypot(xy.x, xy.y); - fAz = fabs(Az = atan2(xy.x, xy.y)); + r = hypot(xy.x, xy.y); + rl = rp = r; + Az = atan2(xy.x, xy.y); + fAz = fabs(Az); for (i = NITER; i ; --i) { z = 2. * atan(pow(r / F,1 / n)); al = acos((pow(tan(.5 * z), n) + diff --git a/src/projections/boggs.cpp b/src/projections/boggs.cpp index e7278904..1900879b 100644 --- a/src/projections/boggs.cpp +++ b/src/projections/boggs.cpp @@ -14,25 +14,25 @@ PROJ_HEAD(boggs, "Boggs Eumorphic") "\n\tPCyl, no inv, Sph"; static PJ_XY boggs_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double theta, th1, c; - int i; - (void) P; - - theta = lp.phi; - if (fabs(fabs(lp.phi) - M_HALFPI) < EPS) - xy.x = 0.; - else { - c = sin(theta) * M_PI; - for (i = NITER; i; --i) { - theta -= th1 = (theta + sin(theta) - c) / - (1. + cos(theta)); - if (fabs(th1) < EPS) break; - } - theta *= 0.5; - xy.x = FXC * lp.lam / (1. / cos(lp.phi) + FXC2 / cos(theta)); - } - xy.y = FYC * (lp.phi + M_SQRT2 * sin(theta)); - return (xy); + double theta, th1, c; + int i; + (void) P; + + theta = lp.phi; + if (fabs(fabs(lp.phi) - M_HALFPI) < EPS) + xy.x = 0.; + else { + c = sin(theta) * M_PI; + for (i = NITER; i; --i) { + th1 = (theta + sin(theta) - c) / (1. + cos(theta)); + theta -= th1; + if (fabs(th1) < EPS) break; + } + theta *= 0.5; + xy.x = FXC * lp.lam / (1. / cos(lp.phi) + FXC2 / cos(theta)); + } + xy.y = FYC * (lp.phi + M_SQRT2 * sin(theta)); + return (xy); } diff --git a/src/projections/bonne.cpp b/src/projections/bonne.cpp index 89f69e6d..c94764cf 100644 --- a/src/projections/bonne.cpp +++ b/src/projections/bonne.cpp @@ -25,7 +25,9 @@ static PJ_XY bonne_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forwar struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); double rh, E, c; - rh = Q->am1 + Q->m1 - pj_mlfn(lp.phi, E = sin(lp.phi), c = cos(lp.phi), Q->en); + E = sin(lp.phi); + c = cos(lp.phi); + rh = Q->am1 + Q->m1 - pj_mlfn(lp.phi, E, c, Q->en); if (fabs(rh) > EPS10) { E = c * lp.lam / (rh * sqrt(1. - P->es * E * E)); xy.x = rh * sin(E); @@ -45,7 +47,8 @@ static PJ_XY bonne_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar rh = Q->cphi1 + Q->phi1 - lp.phi; if (fabs(rh) > EPS10) { - xy.x = rh * sin(E = lp.lam * cos(lp.phi) / rh); + E = lp.lam * cos(lp.phi) / rh; + xy.x = rh * sin(E); xy.y = Q->cphi1 - rh * cos(E); } else xy.x = xy.y = 0.; @@ -58,7 +61,8 @@ static PJ_LP bonne_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); double rh; - rh = hypot(xy.x, xy.y = Q->cphi1 - xy.y); + xy.y = Q->cphi1 - xy.y; + rh = hypot(xy.x, xy.y); lp.phi = Q->cphi1 + Q->phi1 - rh; if (fabs(lp.phi) > M_HALFPI) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); @@ -77,7 +81,8 @@ static PJ_LP bonne_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, invers struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); double s, rh; - rh = hypot(xy.x, xy.y = Q->am1 - xy.y); + xy.y = Q->am1 - xy.y; + rh = hypot(xy.x, xy.y); lp.phi = pj_inv_mlfn(P->ctx, Q->am1 + Q->m1 - rh, P->es, Q->en); if ((s = fabs(lp.phi)) < M_HALFPI) { s = sin(lp.phi); @@ -122,8 +127,9 @@ PJ *PROJECTION(bonne) { Q->en = pj_enfn(P->es); if (nullptr==Q->en) return destructor(P, ENOMEM); - Q->m1 = pj_mlfn(Q->phi1, Q->am1 = sin(Q->phi1), - c = cos(Q->phi1), Q->en); + Q->am1 = sin(Q->phi1); + c = cos(Q->phi1); + Q->m1 = pj_mlfn(Q->phi1, Q->am1, c, Q->en); Q->am1 = c / (sqrt(1. - P->es * Q->am1 * Q->am1) * Q->am1); P->inv = bonne_e_inverse; P->fwd = bonne_e_forward; diff --git a/src/projections/cass.cpp b/src/projections/cass.cpp index 9eea10c5..e253cafc 100644 --- a/src/projections/cass.cpp +++ b/src/projections/cass.cpp @@ -30,10 +30,13 @@ static PJ_XY cass_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward PJ_XY xy = {0.0, 0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - xy.y = pj_mlfn (lp.phi, n = sin (lp.phi), c = cos (lp.phi), Q->en); + n = sin (lp.phi); + c = cos (lp.phi); + xy.y = pj_mlfn (lp.phi, n, c, Q->en); n = 1./sqrt(1. - P->es * n*n); - tn = tan(lp.phi); t = tn * tn; + tn = tan(lp.phi); + t = tn * tn; a1 = lp.lam * c; c *= P->es * c / (1 - P->es); a2 = a1 * a1; @@ -61,7 +64,8 @@ static PJ_LP cass_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); ph1 = pj_inv_mlfn (P->ctx, Q->m0 + xy.y, P->es, Q->en); - tn = tan (ph1); t = tn*tn; + tn = tan (ph1); + t = tn*tn; n = sin (ph1); r = 1. / (1. - P->es * n * n); n = sqrt (r); diff --git a/src/projections/cea.cpp b/src/projections/cea.cpp index a7e5fd04..7e6d3212 100644 --- a/src/projections/cea.cpp +++ b/src/projections/cea.cpp @@ -43,9 +43,10 @@ static PJ_LP cea_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse static PJ_LP cea_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double t; - if ((t = fabs(xy.y *= P->k0)) - EPS <= 1.) { + xy.y *= P->k0; + const double t = fabs(xy.y); + if (t - EPS <= 1.) { if (t >= 1.) lp.phi = xy.y < 0. ? -M_HALFPI : M_HALFPI; else @@ -88,7 +89,8 @@ PJ *PROJECTION(cea) { t = sin(t); P->k0 /= sqrt(1. - P->es * t * t); P->e = sqrt(P->es); - if (!(Q->apa = pj_authset(P->es))) + Q->apa = pj_authset(P->es); + if (!(Q->apa)) return pj_default_destructor(P, ENOMEM); Q->qp = pj_qsfn(1., P->e, P->one_es); diff --git a/src/projections/chamb.cpp b/src/projections/chamb.cpp index e6bac22c..36609e79 100644 --- a/src/projections/chamb.cpp +++ b/src/projections/chamb.cpp @@ -129,10 +129,14 @@ PJ *PROJECTION(chamb) { Q->beta_0 = lc(P->ctx,Q->c[0].v.r, Q->c[2].v.r, Q->c[1].v.r); Q->beta_1 = lc(P->ctx,Q->c[0].v.r, Q->c[1].v.r, Q->c[2].v.r); Q->beta_2 = M_PI - Q->beta_0; - Q->p.y = 2. * (Q->c[0].p.y = Q->c[1].p.y = Q->c[2].v.r * sin(Q->beta_0)); + Q->c[0].p.y = Q->c[2].v.r * sin(Q->beta_0); + Q->c[1].p.y = Q->c[0].p.y; + Q->p.y = 2. * Q->c[0].p.y; Q->c[2].p.y = 0.; - Q->c[0].p.x = - (Q->c[1].p.x = 0.5 * Q->c[0].v.r); - Q->p.x = Q->c[2].p.x = Q->c[0].p.x + Q->c[2].v.r * cos(Q->beta_0); + Q->c[1].p.x = 0.5 * Q->c[0].v.r; + Q->c[0].p.x = -Q->c[1].p.x; + Q->c[2].p.x = Q->c[0].p.x + Q->c[2].v.r * cos(Q->beta_0); + Q->p.x = Q->c[2].p.x; P->es = 0.; P->fwd = chamb_s_forward; diff --git a/src/projections/collg.cpp b/src/projections/collg.cpp index 40319a51..1b9d2da7 100644 --- a/src/projections/collg.cpp +++ b/src/projections/collg.cpp @@ -14,7 +14,8 @@ PROJ_HEAD(collg, "Collignon") "\n\tPCyl, Sph"; static PJ_XY collg_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; (void) P; - if ((xy.y = 1. - sin(lp.phi)) <= 0.) + xy.y = 1. - sin(lp.phi); + if (xy.y <= 0.) xy.y = 0.; else xy.y = sqrt(xy.y); @@ -27,7 +28,8 @@ static PJ_XY collg_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar static PJ_LP collg_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; lp.phi = xy.y / FYC - 1.; - if (fabs(lp.phi = 1. - lp.phi * lp.phi) < 1.) + lp.phi = 1. - lp.phi * lp.phi; + if (fabs(lp.phi) < 1.) lp.phi = asin(lp.phi); else if (fabs(lp.phi) > ONEEPS) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); @@ -36,7 +38,8 @@ static PJ_LP collg_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers lp.phi = lp.phi < 0. ? -M_HALFPI : M_HALFPI; } - if ((lp.lam = 1. - sin(lp.phi)) <= 0.) + lp.lam = 1. - sin(lp.phi); + if (lp.lam <= 0.) lp.lam = 0.; else lp.lam = xy.x / (FXC * sqrt(lp.lam)); diff --git a/src/projections/comill.cpp b/src/projections/comill.cpp index afcfbf7f..189e583e 100644 --- a/src/projections/comill.cpp +++ b/src/projections/comill.cpp @@ -59,7 +59,8 @@ static PJ_LP comill_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver y2 = yc * yc; f = (yc * (K1 + y2 * (K2 + K3 * y2))) - xy.y; fder = C1 + y2 * (C2 + C3 * y2); - yc -= tol = f / fder; + tol = f / fder; + yc -= tol; if (fabs(tol) < EPS) { break; } diff --git a/src/projections/eck2.cpp b/src/projections/eck2.cpp index 27b94aed..0d9fd5fa 100644 --- a/src/projections/eck2.cpp +++ b/src/projections/eck2.cpp @@ -29,7 +29,8 @@ static PJ_LP eck2_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse PJ_LP lp = {0.0,0.0}; (void) P; - lp.lam = xy.x / (FXC * ( lp.phi = 2. - fabs(xy.y) / FYC) ); + lp.phi = 2. - fabs(xy.y) / FYC; + lp.lam = xy.x / (FXC * lp.phi); lp.phi = (4. - lp.phi * lp.phi) * C13; if (fabs(lp.phi) >= 1.) { if (fabs(lp.phi) > ONEEPS) { diff --git a/src/projections/eck3.cpp b/src/projections/eck3.cpp index dd04fb39..0deeb4f3 100644 --- a/src/projections/eck3.cpp +++ b/src/projections/eck3.cpp @@ -90,7 +90,8 @@ PJ *PROJECTION(wag6) { return pj_default_destructor (P, ENOMEM); P->opaque = Q; - Q->C_x = Q->C_y = 0.94745; + Q->C_x = 0.94745; + Q->C_y = 0.94745; Q->A = 0.0; Q->B = 0.30396355092701331433; diff --git a/src/projections/eck4.cpp b/src/projections/eck4.cpp index df2caf42..76474401 100644 --- a/src/projections/eck4.cpp +++ b/src/projections/eck4.cpp @@ -28,8 +28,8 @@ static PJ_XY eck4_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward for (i = NITER; i ; --i) { c = cos(lp.phi); s = sin(lp.phi); - lp.phi -= V = (lp.phi + s * (c + 2.) - p) / - (1. + c * (c + 2.) - s * s); + V = (lp.phi + s * (c + 2.) - p) / (1. + c * (c + 2.) - s * s); + lp.phi -= V; if (fabs(V) < EPS) break; } @@ -46,10 +46,10 @@ static PJ_XY eck4_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward static PJ_LP eck4_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double c; lp.phi = aasin(P->ctx,xy.y * RC_y); - lp.lam = xy.x / (C_x * (1. + (c = cos(lp.phi)))); + const double c = cos(lp.phi); + lp.lam = xy.x / (C_x * (1. + c)); lp.phi = aasin(P->ctx,(lp.phi + sin(lp.phi) * (c + 2.)) * RC_p); return lp; } diff --git a/src/projections/eck5.cpp b/src/projections/eck5.cpp index bf0e6a42..b39c0ea3 100644 --- a/src/projections/eck5.cpp +++ b/src/projections/eck5.cpp @@ -26,7 +26,8 @@ static PJ_XY eck5_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward static PJ_LP eck5_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; (void) P; - lp.lam = RXF * xy.x / (1. + cos( lp.phi = RYF * xy.y)); + lp.phi = RYF * xy.y; + lp.lam = RXF * xy.x / (1. + cos(lp.phi)); return lp; } diff --git a/src/projections/eqdc.cpp b/src/projections/eqdc.cpp index 84e37910..659488b1 100644 --- a/src/projections/eqdc.cpp +++ b/src/projections/eqdc.cpp @@ -31,8 +31,9 @@ static PJ_XY eqdc_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward Q->rho = Q->c - (Q->ellips ? pj_mlfn(lp.phi, sin(lp.phi), cos(lp.phi), Q->en) : lp.phi); - xy.x = Q->rho * sin( lp.lam *= Q->n ); - xy.y = Q->rho0 - Q->rho * cos(lp.lam); + const double lam_mul_n = lp.lam * Q->n; + xy.x = Q->rho * sin(lam_mul_n); + xy.y = Q->rho0 - Q->rho * cos(lam_mul_n); return xy; } @@ -93,10 +94,12 @@ PJ *PROJECTION(eqdc) { if (!(Q->en = pj_enfn(P->es))) return destructor(P, ENOMEM); - Q->n = sinphi = sin(Q->phi1); + sinphi = sin(Q->phi1); + Q->n = sinphi; cosphi = cos(Q->phi1); secant = fabs(Q->phi1 - Q->phi2) >= EPS10; - if( (Q->ellips = (P->es > 0.)) ) { + Q->ellips = (P->es > 0.); + if( Q->ellips ) { double ml1, m1; m1 = pj_msfn(sinphi, cosphi, P->es); diff --git a/src/projections/fouc_s.cpp b/src/projections/fouc_s.cpp index 29d6437d..c5989514 100644 --- a/src/projections/fouc_s.cpp +++ b/src/projections/fouc_s.cpp @@ -33,14 +33,14 @@ static PJ_XY fouc_s_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwa static PJ_LP fouc_s_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double V; int i; if (Q->n != 0.0) { lp.phi = xy.y; for (i = MAX_ITER; i ; --i) { - lp.phi -= V = (Q->n * lp.phi + Q->n1 * sin(lp.phi) - xy.y ) / - (Q->n + Q->n1 * cos(lp.phi)); + const double V = (Q->n * lp.phi + Q->n1 * sin(lp.phi) - xy.y ) / + (Q->n + Q->n1 * cos(lp.phi)); + lp.phi -= V; if (fabs(V) < LOOP_TOL) break; } @@ -48,7 +48,7 @@ static PJ_LP fouc_s_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver lp.phi = xy.y < 0. ? -M_HALFPI : M_HALFPI; } else lp.phi = aasin(P->ctx,xy.y); - V = cos(lp.phi); + const double V = cos(lp.phi); lp.lam = xy.x * (Q->n + Q->n1 * V) / V; return lp; } diff --git a/src/projections/geos.cpp b/src/projections/geos.cpp index fcd7d4ee..338f07c2 100644 --- a/src/projections/geos.cpp +++ b/src/projections/geos.cpp @@ -121,7 +121,7 @@ static PJ_XY geos_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward static PJ_LP geos_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double Vx, Vy, Vz, a, b, det, k; + double Vx, Vy, Vz, a, b, k; /* Setting three components of vector from satellite to position.*/ Vx = -1.0; @@ -136,7 +136,8 @@ static PJ_LP geos_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse /* Calculation of terms in cubic equation and determinant.*/ a = Vy * Vy + Vz * Vz + Vx * Vx; b = 2 * Q->radius_g * Vx; - if ((det = (b * b) - 4 * a * Q->C) < 0.) { + const double det = (b * b) - 4 * a * Q->C; + if (det < 0.) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } @@ -158,7 +159,7 @@ static PJ_LP geos_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse static PJ_LP geos_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double Vx, Vy, Vz, a, b, det, k; + double Vx, Vy, Vz, a, b, k; /* Setting three components of vector from satellite to position.*/ Vx = -1.0; @@ -175,7 +176,8 @@ static PJ_LP geos_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse a = Vz / Q->radius_p; a = Vy * Vy + a * a + Vx * Vx; b = 2 * Q->radius_g * Vx; - if ((det = (b * b) - 4 * a * Q->C) < 0.) { + const double det = (b * b) - 4 * a * Q->C; + if (det < 0.) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } diff --git a/src/projections/gn_sinu.cpp b/src/projections/gn_sinu.cpp index 84883cbc..815de8be 100644 --- a/src/projections/gn_sinu.cpp +++ b/src/projections/gn_sinu.cpp @@ -25,9 +25,10 @@ struct pj_opaque { static PJ_XY gn_sinu_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward */ PJ_XY xy = {0.0,0.0}; - double s, c; - xy.y = pj_mlfn(lp.phi, s = sin(lp.phi), c = cos(lp.phi), static_cast<struct pj_opaque*>(P->opaque)->en); + const double s = sin(lp.phi); + const double c = cos(lp.phi); + xy.y = pj_mlfn(lp.phi, s, c, static_cast<struct pj_opaque*>(P->opaque)->en); xy.x = lp.lam * c / sqrt(1. - P->es * s * s); return xy; } @@ -37,7 +38,9 @@ static PJ_LP gn_sinu_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inve PJ_LP lp = {0.0,0.0}; double s; - if ((s = fabs(lp.phi = pj_inv_mlfn(P->ctx, xy.y, P->es, static_cast<struct pj_opaque*>(P->opaque)->en))) < M_HALFPI) { + lp.phi = pj_inv_mlfn(P->ctx, xy.y, P->es, static_cast<struct pj_opaque*>(P->opaque)->en); + s = fabs(lp.phi); + if (s < M_HALFPI) { s = sin(lp.phi); lp.lam = xy.x * sqrt(1. - P->es * s * s) / cos(lp.phi); } else if ((s - EPS10) < M_HALFPI) { @@ -57,13 +60,13 @@ static PJ_XY gn_sinu_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forw if (Q->m == 0.0) lp.phi = Q->n != 1. ? aasin(P->ctx,Q->n * sin(lp.phi)): lp.phi; else { - double k, V; int i; - k = Q->n * sin(lp.phi); + const double k = Q->n * sin(lp.phi); for (i = MAX_ITER; i ; --i) { - lp.phi -= V = (Q->m * lp.phi + sin(lp.phi) - k) / - (Q->m + cos(lp.phi)); + const double V = (Q->m * lp.phi + sin(lp.phi) - k) / + (Q->m + cos(lp.phi)); + lp.phi -= V; if (fabs(V) < LOOP_TOL) break; } @@ -112,7 +115,8 @@ static void setup(PJ *P) { P->inv = gn_sinu_s_inverse; P->fwd = gn_sinu_s_forward; - Q->C_x = (Q->C_y = sqrt((Q->m + 1.) / Q->n))/(Q->m + 1.); + Q->C_y = sqrt((Q->m + 1.) / Q->n); + Q->C_x = Q->C_y/(Q->m + 1.); } diff --git a/src/projections/goode.cpp b/src/projections/goode.cpp index c716649d..fdace387 100644 --- a/src/projections/goode.cpp +++ b/src/projections/goode.cpp @@ -69,12 +69,16 @@ PJ *PROJECTION(goode) { P->destructor = destructor; P->es = 0.; - if (!(Q->sinu = pj_sinu(nullptr)) || !(Q->moll = pj_moll(nullptr))) + Q->sinu = pj_sinu(nullptr); + Q->moll = pj_moll(nullptr); + if (Q->sinu == nullptr || Q->moll == nullptr) return destructor (P, ENOMEM); Q->sinu->es = 0.; Q->sinu->ctx = P->ctx; Q->moll->ctx = P->ctx; - if (!(Q->sinu = pj_sinu(Q->sinu)) || !(Q->moll = pj_moll(Q->moll))) + Q->sinu = pj_sinu(Q->sinu); + Q->moll = pj_moll(Q->moll); + if (Q->sinu == nullptr || Q->moll == nullptr) return destructor (P, ENOMEM); P->fwd = goode_s_forward; diff --git a/src/projections/hammer.cpp b/src/projections/hammer.cpp index 56bdf74e..8d6d9408 100644 --- a/src/projections/hammer.cpp +++ b/src/projections/hammer.cpp @@ -63,12 +63,14 @@ PJ *PROJECTION(hammer) { P->opaque = Q; if (pj_param(P->ctx, P->params, "tW").i) { - if ((Q->w = fabs(pj_param(P->ctx, P->params, "dW").f)) <= 0.) + Q->w = fabs(pj_param(P->ctx, P->params, "dW").f); + if (Q->w <= 0.) return pj_default_destructor (P, PJD_ERR_W_OR_M_ZERO_OR_LESS); } else Q->w = .5; if (pj_param(P->ctx, P->params, "tM").i) { - if ((Q->m = fabs(pj_param(P->ctx, P->params, "dM").f)) <= 0.) + Q->m = fabs(pj_param(P->ctx, P->params, "dM").f); + if (Q->m <= 0.) return pj_default_destructor (P, PJD_ERR_W_OR_M_ZERO_OR_LESS); } else Q->m = 1.; diff --git a/src/projections/hatano.cpp b/src/projections/hatano.cpp index 29a01f80..c10c4e35 100644 --- a/src/projections/hatano.cpp +++ b/src/projections/hatano.cpp @@ -24,13 +24,13 @@ PROJ_HEAD(hatano, "Hatano Asymmetrical Equal Area") "\n\tPCyl, Sph"; static PJ_XY hatano_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double th1, c; int i; (void) P; - c = sin(lp.phi) * (lp.phi < 0. ? CSz : CN); + const double c = sin(lp.phi) * (lp.phi < 0. ? CSz : CN); for (i = NITER; i; --i) { - lp.phi -= th1 = (lp.phi + sin(lp.phi) - c) / (1. + cos(lp.phi)); + const double th1 = (lp.phi + sin(lp.phi) - c) / (1. + cos(lp.phi));; + lp.phi -= th1; if (fabs(th1) < EPS) break; } xy.x = FXC * lp.lam * cos(lp.phi *= .5); diff --git a/src/projections/isea.cpp b/src/projections/isea.cpp index 801f2b82..1ebbeebb 100644 --- a/src/projections/isea.cpp +++ b/src/projections/isea.cpp @@ -339,7 +339,7 @@ static int isea_snyder_forward(struct isea_geo * ll, struct isea_pt * out) double theta; /* additional variables from snyder */ - double q, Rprime, H, Ag, Azprime, Az, dprime, f, rho, + double q, H, Ag, Azprime, Az, dprime, f, rho, x, y; /* variables used to store intermediate results */ @@ -427,7 +427,7 @@ static int isea_snyder_forward(struct isea_geo * ll, struct isea_pt * out) /* eq 5 */ /* Rprime = 0.9449322893 * R; */ /* R' in the paper is for the truncated */ - Rprime = 0.91038328153090290025; + const double Rprime = 0.91038328153090290025; /* eq 6 */ H = acos(sin(Az) * sin(G) * cos(g) - cos(Az) * cos(G)); diff --git a/src/projections/krovak.cpp b/src/projections/krovak.cpp index 98f09199..aef44d42 100644 --- a/src/projections/krovak.cpp +++ b/src/projections/krovak.cpp @@ -198,7 +198,8 @@ PJ *PROJECTION(krovak) { /* we want Bessel as fixed ellipsoid */ P->a = 6377397.155; - P->e = sqrt(P->es = 0.006674372230614); + P->es = 0.006674372230614; + P->e = sqrt(P->es); /* if latitude of projection center is not set, use 49d30'N */ if (!pj_param(P->ctx, P->params, "tlat_0").i) diff --git a/src/projections/labrd.cpp b/src/projections/labrd.cpp index 21d9099f..c9dfdfc6 100644 --- a/src/projections/labrd.cpp +++ b/src/projections/labrd.cpp @@ -27,8 +27,10 @@ static PJ_XY labrd_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forwar V2 = .5 * P->e * Q->A * log ((1. + t)/(1. - t)); ps = 2. * (atan(exp(V1 - V2 + Q->C)) - M_FORTPI); I1 = ps - Q->p0s; - cosps = cos(ps); cosps2 = cosps * cosps; - sinps = sin(ps); sinps2 = sinps * sinps; + cosps = cos(ps); + cosps2 = cosps * cosps; + sinps = sin(ps); + sinps2 = sinps * sinps; I4 = Q->A * cosps; I2 = .5 * Q->A * I4 * sinps; I3 = I2 * Q->A * Q->A * (5. * cosps2 - sinps2) / 12.; @@ -125,7 +127,8 @@ PJ *PROJECTION(labrd) { - Q->A * log( tan(M_FORTPI + .5 * P->phi0)) + log( tan(M_FORTPI + .5 * Q->p0s)); t = Az + Az; - Q->Ca = (1. - cos(t)) * ( Q->Cb = 1. / (12. * Q->kRg * Q->kRg) ); + Q->Cb = 1. / (12. * Q->kRg * Q->kRg); + Q->Ca = (1. - cos(t)) * Q->Cb; Q->Cb *= sin(t); Q->Cc = 3. * (Q->Ca * Q->Ca - Q->Cb * Q->Cb); Q->Cd = 6. * Q->Ca * Q->Cb; diff --git a/src/projections/lcc.cpp b/src/projections/lcc.cpp index 02530d5b..b78b439c 100644 --- a/src/projections/lcc.cpp +++ b/src/projections/lcc.cpp @@ -128,7 +128,8 @@ PJ *PROJECTION(lcc) { } Q->n /= denom; } - Q->c = (Q->rho0 = m1 * pow(ml1, -Q->n) / Q->n); + Q->rho0 = m1 * pow(ml1, -Q->n) / Q->n; + Q->c = Q->rho0; Q->rho0 *= (fabs(fabs(P->phi0) - M_HALFPI) < EPS10) ? 0. : pow(pj_tsfn(P->phi0, sin(P->phi0), P->e), Q->n); } else { diff --git a/src/projections/lcca.cpp b/src/projections/lcca.cpp index 11ecb29c..51dd28aa 100644 --- a/src/projections/lcca.cpp +++ b/src/projections/lcca.cpp @@ -88,8 +88,9 @@ static PJ_XY lcca_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward S = pj_mlfn(lp.phi, sin(lp.phi), cos(lp.phi), Q->en) - Q->M0; dr = fS(S, Q->C); r = Q->r0 - dr; - xy.x = P->k0 * (r * sin( lp.lam *= Q->l ) ); - xy.y = P->k0 * (Q->r0 - r * cos(lp.lam) ); + const double lam_mul_l = lp.lam * Q->l; + xy.x = P->k0 * (r * sin(lam_mul_l)); + xy.y = P->k0 * (Q->r0 - r * cos(lam_mul_l) ); return xy; } diff --git a/src/projections/mbt_fps.cpp b/src/projections/mbt_fps.cpp index 9ce2aa36..85051dc7 100644 --- a/src/projections/mbt_fps.cpp +++ b/src/projections/mbt_fps.cpp @@ -18,19 +18,18 @@ PROJ_HEAD(mbt_fps, "McBryde-Thomas Flat-Pole Sine (No. 2)") "\n\tCyl, Sph"; static PJ_XY mbt_fps_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double k, V, t; - int i; (void) P; - k = C3 * sin(lp.phi); - for (i = MAX_ITER; i ; --i) { - t = lp.phi / C2; - lp.phi -= V = (C1 * sin(t) + sin(lp.phi) - k) / - (C1_2 * cos(t) + cos(lp.phi)); + const double k = C3 * sin(lp.phi); + for (int i = MAX_ITER; i ; --i) { + const double t = lp.phi / C2; + const double V = (C1 * sin(t) + sin(lp.phi) - k) / + (C1_2 * cos(t) + cos(lp.phi)); + lp.phi -= V; if (fabs(V) < LOOP_TOL) break; } - t = lp.phi / C2; + const double t = lp.phi / C2; xy.x = C_x * lp.lam * (1. + 3. * cos(lp.phi)/cos(t) ); xy.y = C_y * sin(t); return xy; @@ -39,9 +38,9 @@ static PJ_XY mbt_fps_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forw static PJ_LP mbt_fps_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double t; - lp.phi = C2 * (t = aasin(P->ctx,xy.y / C_y)); + const double t = aasin(P->ctx,xy.y / C_y); + lp.phi = C2 * t; lp.lam = xy.x / (C_x * (1. + 3. * cos(lp.phi)/cos(t))); lp.phi = aasin(P->ctx,(C1 * sin(t) + sin(lp.phi)) / C3); return (lp); diff --git a/src/projections/mbtfpp.cpp b/src/projections/mbtfpp.cpp index 6cc7c466..cc01cb40 100644 --- a/src/projections/mbtfpp.cpp +++ b/src/projections/mbtfpp.cpp @@ -40,8 +40,10 @@ static PJ_LP mbtfpp_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver } else lp.phi = asin(lp.phi); - lp.lam = xy.x / ( FXC * (2. * cos(C23 * (lp.phi *= 3.)) - 1.) ); - if (fabs(lp.phi = sin(lp.phi) / CSy) >= 1.) { + lp.phi *= 3.; + lp.lam = xy.x / ( FXC * (2. * cos(C23 * lp.phi) - 1.) ); + lp.phi = sin(lp.phi) / CSy; + if (fabs(lp.phi) >= 1.) { if (fabs(lp.phi) > ONEEPS) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; diff --git a/src/projections/mbtfpq.cpp b/src/projections/mbtfpq.cpp index 9a419790..5c7f8ca6 100644 --- a/src/projections/mbtfpq.cpp +++ b/src/projections/mbtfpq.cpp @@ -20,14 +20,13 @@ PROJ_HEAD(mbtfpq, "McBryde-Thomas Flat-Polar Quartic") "\n\tCyl, Sph"; static PJ_XY mbtfpq_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double th1, c; - int i; (void) P; - c = C * sin(lp.phi); - for (i = NITER; i; --i) { - lp.phi -= th1 = (sin(.5*lp.phi) + sin(lp.phi) - c) / - (.5*cos(.5*lp.phi) + cos(lp.phi)); + const double c = C * sin(lp.phi); + for (int i = NITER; i; --i) { + const double th1 = (sin(.5*lp.phi) + sin(lp.phi) - c) / + (.5*cos(.5*lp.phi) + cos(lp.phi)); + lp.phi -= th1; if (fabs(th1) < EPS) break; } xy.x = FXC * lp.lam * (1.0 + 2. * cos(lp.phi)/cos(0.5 * lp.phi)); @@ -46,10 +45,18 @@ static PJ_LP mbtfpq_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inver proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } - else if (lp.phi < 0.) { t = -1.; lp.phi = -M_PI; } - else { t = 1.; lp.phi = M_PI; } - } else - lp.phi = 2. * asin(t = lp.phi); + else if (lp.phi < 0.) { + t = -1.; + lp.phi = -M_PI; + } + else { + t = 1.; + lp.phi = M_PI; + } + } else { + t = lp.phi; + lp.phi = 2. * asin(lp.phi); + } lp.lam = RXC * xy.x / (1. + 2. * cos(lp.phi)/cos(0.5 * lp.phi)); lp.phi = RC * (t + sin(lp.phi)); if (fabs(lp.phi) > 1.) diff --git a/src/projections/misrsom.cpp b/src/projections/misrsom.cpp index 09e2d8f3..71116e1e 100644 --- a/src/projections/misrsom.cpp +++ b/src/projections/misrsom.cpp @@ -53,7 +53,8 @@ static void seraz0(double lam, double mult, PJ *P) { h = sqrt((1. + Q->q * sdsq) / (1. + Q->w * sdsq)) * ((1. + Q->w * sdsq) / (d__1 * d__1) - Q->p22 * Q->ca); sq = sqrt(Q->xj * Q->xj + s * s); - Q->b += fc = mult * (h * Q->xj - s * s) / sq; + fc = mult * (h * Q->xj - s * s) / sq; + Q->b += fc; Q->a2 += fc * cos(lam + lam); Q->a4 += fc * cos(lam * 4.); fc = mult * s * (h + Q->xj) / sq; @@ -89,7 +90,8 @@ static PJ_XY misrsom_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forw fac = lampp - sin(lampp) * M_HALFPI; for (l = 50; l; --l) { lamt = lp.lam + Q->p22 * sav; - if (fabs(c = cos(lamt)) < TOL) + c = cos(lamt); + if (fabs(c) < TOL) lamt -= TOL; xlam = (P->one_es * tanphi * Q->sa + sin(lamt) * Q->ca) / c; lamdp = atan(xlam) + fac; diff --git a/src/projections/moll.cpp b/src/projections/moll.cpp index d511c6e0..5d4f7825 100644 --- a/src/projections/moll.cpp +++ b/src/projections/moll.cpp @@ -23,13 +23,13 @@ struct pj_opaque { static PJ_XY moll_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double k, V; int i; - k = Q->C_p * sin(lp.phi); + const double k = Q->C_p * sin(lp.phi); for (i = MAX_ITER; i ; --i) { - lp.phi -= V = (lp.phi + sin(lp.phi) - k) / - (1. + cos(lp.phi)); + const double V = (lp.phi + sin(lp.phi) - k) / + (1. + cos(lp.phi)); + lp.phi -= V; if (fabs(V) < LOOP_TOL) break; } diff --git a/src/projections/natearth.cpp b/src/projections/natearth.cpp index ff8aef6a..5c096605 100644 --- a/src/projections/natearth.cpp +++ b/src/projections/natearth.cpp @@ -57,7 +57,7 @@ static PJ_XY natearth_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, for static PJ_LP natearth_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double yc, tol, y2, y4, f, fder; + double yc, y2, y4, f, fder; int i; (void) P; @@ -75,7 +75,8 @@ static PJ_LP natearth_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inv y4 = y2 * y2; f = (yc * (B0 + y2 * (B1 + y4 * (B2 + B3 * y2 + B4 * y4)))) - xy.y; fder = C0 + y2 * (C1 + y4 * (C2 + C3 * y2 + C4 * y4)); - yc -= tol = f / fder; + const double tol = f / fder; + yc -= tol; if (fabs(tol) < EPS) { break; } diff --git a/src/projections/natearth2.cpp b/src/projections/natearth2.cpp index 95d73c36..d149ca85 100644 --- a/src/projections/natearth2.cpp +++ b/src/projections/natearth2.cpp @@ -51,7 +51,7 @@ static PJ_XY natearth2_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, fo static PJ_LP natearth2_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double yc, tol, y2, y4, y6, f, fder; + double yc, y2, y4, y6, f, fder; int i; (void) P; @@ -69,7 +69,8 @@ static PJ_LP natearth2_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, in y4 = y2 * y2; f = (yc * (B0 + y4 * y4 * (B1 + B2 * y2 + B3 * y4))) - xy.y; fder = C0 + y4 * y4 * (C1 + C2 * y2 + C3 * y4); - yc -= tol = f / fder; + const double tol = f / fder; + yc -= tol; if (fabs(tol) < EPS) { break; } diff --git a/src/projections/nell.cpp b/src/projections/nell.cpp index 63a0eec1..faa63e59 100644 --- a/src/projections/nell.cpp +++ b/src/projections/nell.cpp @@ -13,16 +13,16 @@ PROJ_HEAD(nell, "Nell") "\n\tPCyl, Sph"; static PJ_XY nell_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double k, V; int i; (void) P; - k = 2. * sin(lp.phi); - V = lp.phi * lp.phi; - lp.phi *= 1.00371 + V * (-0.0935382 + V * -0.011412); + const double k = 2. * sin(lp.phi); + const double phi_pow_2 = lp.phi * lp.phi; + lp.phi *= 1.00371 + phi_pow_2 * (-0.0935382 + phi_pow_2 * -0.011412); for (i = MAX_ITER; i ; --i) { - lp.phi -= V = (lp.phi + sin(lp.phi) - k) / - (1. + cos(lp.phi)); + const double V = (lp.phi + sin(lp.phi) - k) / + (1. + cos(lp.phi)); + lp.phi -= V; if (fabs(V) < LOOP_TOL) break; } diff --git a/src/projections/nell_h.cpp b/src/projections/nell_h.cpp index 63d12391..a0057cda 100644 --- a/src/projections/nell_h.cpp +++ b/src/projections/nell_h.cpp @@ -24,14 +24,14 @@ static PJ_XY nell_h_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwa static PJ_LP nell_h_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double V, c, p; int i; (void) P; - p = 0.5 * xy.y; + const double p = 0.5 * xy.y; for (i = NITER; i ; --i) { - c = cos(0.5 * lp.phi); - lp.phi -= V = (lp.phi - tan(lp.phi/2) - p)/(1. - 0.5/(c*c)); + const double c = cos(0.5 * lp.phi); + const double V = (lp.phi - tan(lp.phi/2) - p)/(1. - 0.5/(c*c)); + lp.phi -= V; if (fabs(V) < EPS) break; } diff --git a/src/projections/nsper.cpp b/src/projections/nsper.cpp index 852d0c98..903946b9 100644 --- a/src/projections/nsper.cpp +++ b/src/projections/nsper.cpp @@ -202,8 +202,10 @@ PJ *PROJECTION(tpers) { omega = pj_param(P->ctx, P->params, "rtilt").f; gamma = pj_param(P->ctx, P->params, "razi").f; Q->tilt = 1; - Q->cg = cos(gamma); Q->sg = sin(gamma); - Q->cw = cos(omega); Q->sw = sin(omega); + Q->cg = cos(gamma); + Q->sg = sin(gamma); + Q->cw = cos(omega); + Q->sw = sin(omega); return setup(P); } diff --git a/src/projections/nzmg.cpp b/src/projections/nzmg.cpp index 2f1a897e..099bd178 100644 --- a/src/projections/nzmg.cpp +++ b/src/projections/nzmg.cpp @@ -65,8 +65,12 @@ static PJ_XY nzmg_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward int i; lp.phi = (lp.phi - P->phi0) * RAD_TO_SEC5; - for (p.r = *(C = tpsi + (i = Ntpsi)); i ; --i) - p.r = *--C + lp.phi * p.r; + i = Ntpsi; + C = tpsi + i; + for (p.r = *C; i ; --i) { + --C; + p.r = *C + lp.phi * p.r; + } p.r *= lp.phi; p.i = lp.lam; p = pj_zpoly1(p, bf, Nbf); @@ -91,15 +95,21 @@ static PJ_LP nzmg_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse f.r -= xy.y; f.i -= xy.x; den = fp.r * fp.r + fp.i * fp.i; - p.r += dp.r = -(f.r * fp.r + f.i * fp.i) / den; - p.i += dp.i = -(f.i * fp.r - f.r * fp.i) / den; + dp.r = -(f.r * fp.r + f.i * fp.i) / den; + dp.i = -(f.i * fp.r - f.r * fp.i) / den; + p.r += dp.r; + p.i += dp.i; if ((fabs(dp.r) + fabs(dp.i)) <= EPSLN) break; } if (nn) { lp.lam = p.i; - for (lp.phi = *(C = tphi + (i = Ntphi)); i ; --i) - lp.phi = *--C + p.r * lp.phi; + i = Ntphi; + C = tphi + i; + for (lp.phi = *C; i ; --i) { + --C; + lp.phi = *C + p.r * lp.phi; + } lp.phi = P->phi0 + p.r * lp.phi * SEC5_TO_RAD; } else lp.lam = lp.phi = HUGE_VAL; diff --git a/src/projections/ob_tran.cpp b/src/projections/ob_tran.cpp index badc6dd8..7cf1cf98 100644 --- a/src/projections/ob_tran.cpp +++ b/src/projections/ob_tran.cpp @@ -61,7 +61,8 @@ static PJ_LP o_inverse(PJ_XY xy, PJ *P) { /* spheroid */ PJ_LP lp = Q->link->inv(xy, Q->link); if (lp.lam != HUGE_VAL) { - coslam = cos(lp.lam -= Q->lamp); + lp.lam -= Q->lamp; + coslam = cos(lp.lam); sinphi = sin(lp.phi); cosphi = cos(lp.phi); /* Formula (5-9) */ @@ -222,7 +223,8 @@ PJ *PROJECTION(ob_tran) { phi1 = pj_param(P->ctx, P->params, "ro_lat_1").f; lam2 = pj_param(P->ctx, P->params, "ro_lon_2").f; phi2 = pj_param(P->ctx, P->params, "ro_lat_2").f; - if (fabs(phi1 - phi2) <= TOL || (con = fabs(phi1)) <= TOL || + con = fabs(phi1); + if (fabs(phi1 - phi2) <= TOL || con <= TOL || fabs(con - M_HALFPI) <= TOL || fabs(fabs(phi2) - M_HALFPI) <= TOL) return destructor(P, PJD_ERR_LAT_1_OR_2_ZERO_OR_90); diff --git a/src/projections/ocea.cpp b/src/projections/ocea.cpp index 646b8638..4fc5391d 100644 --- a/src/projections/ocea.cpp +++ b/src/projections/ocea.cpp @@ -37,12 +37,12 @@ static PJ_XY ocea_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward static PJ_LP ocea_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double t, s; xy.y /= Q->rok; xy.x /= Q->rtk; - t = sqrt(1. - xy.y * xy.y); - lp.phi = asin(xy.y * Q->sinphi + t * Q->cosphi * (s = sin(xy.x))); + const double t = sqrt(1. - xy.y * xy.y); + const double s = sin(xy.x); + lp.phi = asin(xy.y * Q->sinphi + t * Q->cosphi * s); lp.lam = atan2(t * Q->sinphi * s - xy.y * Q->cosphi, t * cos(xy.x)); return lp; diff --git a/src/projections/oea.cpp b/src/projections/oea.cpp index cbedf0c0..61fb0647 100644 --- a/src/projections/oea.cpp +++ b/src/projections/oea.cpp @@ -19,15 +19,14 @@ struct pj_opaque { static PJ_XY oea_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double Az, M, N, cp, sp, cl, shz; - - cp = cos(lp.phi); - sp = sin(lp.phi); - cl = cos(lp.lam); - Az = aatan2(cp * sin(lp.lam), Q->cp0 * sp - Q->sp0 * cp * cl) + Q->theta; - shz = sin(0.5 * aacos(P->ctx, Q->sp0 * sp + Q->cp0 * cp * cl)); - M = aasin(P->ctx, shz * sin(Az)); - N = aasin(P->ctx, shz * cos(Az) * cos(M) / cos(M * Q->two_r_m)); + + const double cp = cos(lp.phi); + const double sp = sin(lp.phi); + const double cl = cos(lp.lam); + const double Az = aatan2(cp * sin(lp.lam), Q->cp0 * sp - Q->sp0 * cp * cl) + Q->theta; + const double shz = sin(0.5 * aacos(P->ctx, Q->sp0 * sp + Q->cp0 * cp * cl)); + const double M = aasin(P->ctx, shz * sin(Az)); + const double N = aasin(P->ctx, shz * cos(Az) * cos(M) / cos(M * Q->two_r_m)); xy.y = Q->n * sin(N * Q->two_r_n); xy.x = Q->m * sin(M * Q->two_r_m) * cos(N) / cos(N * Q->two_r_n); @@ -38,16 +37,16 @@ static PJ_XY oea_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward static PJ_LP oea_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double N, M, xp, yp, z, Az, cz, sz, cAz; - - N = Q->hn * aasin(P->ctx,xy.y * Q->rn); - M = Q->hm * aasin(P->ctx,xy.x * Q->rm * cos(N * Q->two_r_n) / cos(N)); - xp = 2. * sin(M); - yp = 2. * sin(N) * cos(M * Q->two_r_m) / cos(M); - cAz = cos(Az = aatan2(xp, yp) - Q->theta); - z = 2. * aasin(P->ctx, 0.5 * hypot(xp, yp)); - sz = sin(z); - cz = cos(z); + + const double N = Q->hn * aasin(P->ctx,xy.y * Q->rn); + const double M = Q->hm * aasin(P->ctx,xy.x * Q->rm * cos(N * Q->two_r_n) / cos(N)); + const double xp = 2. * sin(M); + const double yp = 2. * sin(N) * cos(M * Q->two_r_m) / cos(M); + const double Az = aatan2(xp, yp) - Q->theta; + const double cAz = cos(Az); + const double z = 2. * aasin(P->ctx, 0.5 * hypot(xp, yp)); + const double sz = sin(z); + const double cz = cos(z); lp.phi = aasin(P->ctx, Q->sp0 * cz + Q->cp0 * sz * cAz); lp.lam = aatan2(sz * sin(Az), Q->cp0 * cz - Q->sp0 * sz * cAz); diff --git a/src/projections/omerc.cpp b/src/projections/omerc.cpp index 954023df..e9f3b833 100644 --- a/src/projections/omerc.cpp +++ b/src/projections/omerc.cpp @@ -48,26 +48,26 @@ struct pj_opaque { static PJ_XY omerc_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double S, T, U, V, W, temp, u, v; + double u, v; if (fabs(fabs(lp.phi) - M_HALFPI) > EPS) { - W = Q->E / pow(pj_tsfn(lp.phi, sin(lp.phi), P->e), Q->B); - temp = 1. / W; - S = .5 * (W - temp); - T = .5 * (W + temp); - V = sin(Q->B * lp.lam); - U = (S * Q->singam - V * Q->cosgam) / T; + const double W = Q->E / pow(pj_tsfn(lp.phi, sin(lp.phi), P->e), Q->B); + const double one_div_W = 1. / W; + const double S = .5 * (W - one_div_W); + const double T = .5 * (W + one_div_W); + const double V = sin(Q->B * lp.lam); + const double U = (S * Q->singam - V * Q->cosgam) / T; if (fabs(fabs(U) - 1.0) < EPS) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } v = 0.5 * Q->ArB * log((1. - U)/(1. + U)); - temp = cos(Q->B * lp.lam); - if(fabs(temp) < TOL) { - u = Q->A * lp.lam; - } else { - u = Q->ArB * atan2((S * Q->cosgam + V * Q->singam), temp); - } + const double temp = cos(Q->B * lp.lam); + if(fabs(temp) < TOL) { + u = Q->A * lp.lam; + } else { + u = Q->ArB * atan2((S * Q->cosgam + V * Q->singam), temp); + } } else { v = lp.phi > 0 ? Q->v_pole_n : Q->v_pole_s; u = Q->ArB * lp.phi; diff --git a/src/projections/ortho.cpp b/src/projections/ortho.cpp index 34b6b689..5f9366c3 100644 --- a/src/projections/ortho.cpp +++ b/src/projections/ortho.cpp @@ -48,7 +48,8 @@ static PJ_XY ortho_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar xy.y = sin(lp.phi); break; case OBLIQ: - if (Q->sinph0 * (sinphi = sin(lp.phi)) + Q->cosph0 * cosphi * coslam < - EPS10) + sinphi = sin(lp.phi); + if (Q->sinph0 * sinphi + Q->cosph0 * cosphi * coslam < - EPS10) return forward_error(P, lp, xy); xy.y = Q->cosph0 * sinphi - Q->sinph0 * cosphi * coslam; break; @@ -69,11 +70,13 @@ static PJ_XY ortho_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar static PJ_LP ortho_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double rh, cosc, sinc; + double sinc; lp.lam = HUGE_VAL; lp.phi = HUGE_VAL; - if ((sinc = (rh = hypot(xy.x, xy.y))) > 1.) { + const double rh = hypot(xy.x, xy.y); + sinc = rh; + if (sinc > 1.) { if ((sinc - 1.) > EPS10) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); proj_log_trace(P, "Point (%.3f, %.3f) is outside the projection boundary"); @@ -81,7 +84,7 @@ static PJ_LP ortho_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers } sinc = 1.; } - cosc = sqrt(1. - sinc * sinc); /* in this range OK */ + const double cosc = sqrt(1. - sinc * sinc); /* in this range OK */ if (fabs(rh) <= EPS10) { lp.phi = P->phi0; lp.lam = 0.0; diff --git a/src/projections/patterson.cpp b/src/projections/patterson.cpp index 16e7b746..71099cdb 100644 --- a/src/projections/patterson.cpp +++ b/src/projections/patterson.cpp @@ -76,7 +76,7 @@ static PJ_XY patterson_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, fo static PJ_LP patterson_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double yc, tol, y2, f, fder; + double yc; int i; (void) P; @@ -90,10 +90,11 @@ static PJ_LP patterson_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, in } for (i = MAX_ITER; i ; --i) { /* Newton-Raphson */ - y2 = yc * yc; - f = (yc * (K1 + y2 * y2 * (K2 + y2 * (K3 + K4 * y2)))) - xy.y; - fder = C1 + y2 * y2 * (C2 + y2 * (C3 + C4 * y2)); - yc -= tol = f / fder; + const double y2 = yc * yc; + const double f = (yc * (K1 + y2 * y2 * (K2 + y2 * (K3 + K4 * y2)))) - xy.y; + const double fder = C1 + y2 * y2 * (C2 + y2 * (C3 + C4 * y2)); + const double tol = f / fder; + yc -= tol; if (fabs(tol) < EPS11) { break; } diff --git a/src/projections/poly.cpp b/src/projections/poly.cpp index fc7fc8bd..10d93ed2 100644 --- a/src/projections/poly.cpp +++ b/src/projections/poly.cpp @@ -33,8 +33,10 @@ static PJ_XY poly_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward xy.y = -Q->ml0; } else { sp = sin(lp.phi); - ms = fabs(cp = cos(lp.phi)) > TOL ? pj_msfn(sp, cp, P->es) / sp : 0.; - xy.x = ms * sin(lp.lam *= sp); + cp = cos(lp.phi); + ms = fabs(cp) > TOL ? pj_msfn(sp, cp, P->es) / sp : 0.; + lp.lam *= sp; + xy.x = ms * sin(lp.lam); xy.y = (pj_mlfn(lp.phi, sp, cp, Q->en) - Q->ml0) + ms * (1. - cos(lp.lam)); } @@ -45,14 +47,14 @@ static PJ_XY poly_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward static PJ_XY poly_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double cot, E; if (fabs(lp.phi) <= TOL) { xy.x = lp.lam; xy.y = Q->ml0; } else { - cot = 1. / tan(lp.phi); - xy.x = sin(E = lp.lam * sin(lp.phi)) * cot; + const double cot = 1. / tan(lp.phi); + const double E = lp.lam * sin(lp.phi); + xy.x = sin(E) * cot; xy.y = lp.phi - P->phi0 + cot * (1. - cos(E)); } @@ -69,26 +71,28 @@ static PJ_LP poly_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse lp.lam = xy.x; lp.phi = 0.; } else { - double r, c, sp, cp, s2ph, ml, mlb, mlp, dPhi; int i; - r = xy.y * xy.y + xy.x * xy.x; + const double r = xy.y * xy.y + xy.x * xy.x; lp.phi = xy.y; for (i = I_ITER; i ; --i) { - sp = sin(lp.phi); - s2ph = sp * ( cp = cos(lp.phi)); + const double sp = sin(lp.phi); + const double cp = cos(lp.phi); + const double s2ph = sp * cp; if (fabs(cp) < ITOL) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } - c = sp * (mlp = sqrt(1. - P->es * sp * sp)) / cp; - ml = pj_mlfn(lp.phi, sp, cp, Q->en); - mlb = ml * ml + r; + double mlp = sqrt(1. - P->es * sp * sp); + const double c = sp * mlp / cp; + const double ml = pj_mlfn(lp.phi, sp, cp, Q->en); + const double mlb = ml * ml + r; mlp = P->one_es / (mlp * mlp * mlp); - lp.phi += ( dPhi = + const double dPhi = ( ml + ml + c * mlb - 2. * xy.y * (c * ml + 1.) ) / ( P->es * s2ph * (mlb - 2. * xy.y * ml) / c + - 2.* (xy.y - ml) * (c * mlp - 1. / s2ph) - mlp - mlp )); + 2.* (xy.y - ml) * (c * mlp - 1. / s2ph) - mlp - mlp ); + lp.phi += dPhi; if (fabs(dPhi) <= ITOL) break; } @@ -96,7 +100,7 @@ static PJ_LP poly_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; } - c = sin(lp.phi); + const double c = sin(lp.phi); lp.lam = asin(xy.x * tan(lp.phi) * sqrt(1. - P->es * c * c)) / sin(lp.phi); } @@ -106,21 +110,20 @@ static PJ_LP poly_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, inverse static PJ_LP poly_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double B, dphi, tp; - int i; if (fabs(xy.y = P->phi0 + xy.y) <= TOL) { lp.lam = xy.x; lp.phi = 0.; } else { lp.phi = xy.y; - B = xy.x * xy.x + xy.y * xy.y; - i = N_ITER; + const double B = xy.x * xy.x + xy.y * xy.y; + int i = N_ITER; while(true) { - tp = tan(lp.phi); - lp.phi -= (dphi = (xy.y * (lp.phi * tp + 1.) - lp.phi - - .5 * ( lp.phi * lp.phi + B) * tp) / - ((lp.phi - xy.y) / tp - 1.)); + const double tp = tan(lp.phi); + const double dphi = (xy.y * (lp.phi * tp + 1.) - lp.phi - + .5 * ( lp.phi * lp.phi + B) * tp) / + ((lp.phi - xy.y) / tp - 1.); + lp.phi -= dphi; if( !(fabs(dphi) > CONV) ) break; --i; diff --git a/src/projections/putp2.cpp b/src/projections/putp2.cpp index a11e479e..5ff44000 100644 --- a/src/projections/putp2.cpp +++ b/src/projections/putp2.cpp @@ -17,39 +17,39 @@ PROJ_HEAD(putp2, "Putnins P2") "\n\tPCyl, Sph"; static PJ_XY putp2_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; - double p, c, s, V; - int i; - (void) P; - - p = C_p * sin(lp.phi); - s = lp.phi * lp.phi; - lp.phi *= 0.615709 + s * ( 0.00909953 + s * 0.0046292 ); - for (i = NITER; i ; --i) { - c = cos(lp.phi); - s = sin(lp.phi); - lp.phi -= V = (lp.phi + s * (c - 1.) - p) / - (1. + c * (c - 1.) - s * s); - if (fabs(V) < EPS) - break; - } - if (!i) - lp.phi = lp.phi < 0 ? - PI_DIV_3 : PI_DIV_3; - xy.x = C_x * lp.lam * (cos(lp.phi) - 0.5); - xy.y = C_y * sin(lp.phi); - - return xy; + int i; + (void) P; + + const double p = C_p * sin(lp.phi); + const double phi_pow_2 = lp.phi * lp.phi; + lp.phi *= 0.615709 + phi_pow_2 * ( 0.00909953 + phi_pow_2 * 0.0046292 ); + for (i = NITER; i ; --i) { + const double c = cos(lp.phi); + const double s = sin(lp.phi); + const double V = (lp.phi + s * (c - 1.) - p) / + (1. + c * (c - 1.) - s * s); + lp.phi -= V; + if (fabs(V) < EPS) + break; + } + if (!i) + lp.phi = lp.phi < 0 ? - PI_DIV_3 : PI_DIV_3; + xy.x = C_x * lp.lam * (cos(lp.phi) - 0.5); + xy.y = C_y * sin(lp.phi); + + return xy; } static PJ_LP putp2_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; - double c; - lp.phi = aasin(P->ctx,xy.y / C_y); - lp.lam = xy.x / (C_x * ((c = cos(lp.phi)) - 0.5)); - lp.phi = aasin(P->ctx,(lp.phi + sin(lp.phi) * (c - 1.)) / C_p); + lp.phi = aasin(P->ctx,xy.y / C_y); + const double c = cos(lp.phi); + lp.lam = xy.x / (C_x * (c - 0.5)); + lp.phi = aasin(P->ctx,(lp.phi + sin(lp.phi) * (c - 1.)) / C_p); - return lp; + return lp; } diff --git a/src/projections/putp4p.cpp b/src/projections/putp4p.cpp index a5728b74..365f7c1b 100644 --- a/src/projections/putp4p.cpp +++ b/src/projections/putp4p.cpp @@ -22,7 +22,8 @@ static PJ_XY putp4p_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwa lp.phi = aasin(P->ctx,0.883883476 * sin(lp.phi)); xy.x = Q->C_x * lp.lam * cos(lp.phi); - xy.x /= cos(lp.phi *= 0.333333333333333); + lp.phi *= 0.333333333333333; + xy.x /= cos(lp.phi); xy.y = Q->C_y * sin(lp.phi); return xy; diff --git a/src/projections/putp6.cpp b/src/projections/putp6.cpp index bcf9ad8e..db334ff9 100644 --- a/src/projections/putp6.cpp +++ b/src/projections/putp6.cpp @@ -23,15 +23,15 @@ PROJ_HEAD(putp6p, "Putnins P6'") "\n\tPCyl, Sph"; static PJ_XY putp6_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double p, r, V; int i; - p = Q->B * sin(lp.phi); + const double p = Q->B * sin(lp.phi); lp.phi *= 1.10265779; for (i = NITER; i ; --i) { - r = sqrt(1. + lp.phi * lp.phi); - lp.phi -= V = ( (Q->A - r) * lp.phi - log(lp.phi + r) - p ) / - (Q->A - 2. * r); + const double r = sqrt(1. + lp.phi * lp.phi); + const double V = ( (Q->A - r) * lp.phi - log(lp.phi + r) - p ) / + (Q->A - 2. * r); + lp.phi -= V; if (fabs(V) < EPS) break; } diff --git a/src/projections/robin.cpp b/src/projections/robin.cpp index 5a3b081c..3c4d9f07 100644 --- a/src/projections/robin.cpp +++ b/src/projections/robin.cpp @@ -101,8 +101,7 @@ static PJ_XY robin_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar 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; + double t; struct COEFS T; int iters; @@ -119,7 +118,7 @@ static PJ_LP robin_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers } } else { /* general problem */ /* in Y space, reduce to table interval */ - i = isnan(lp.phi) ? -1 : lround(floor(lp.phi * NODES)); + long i = isnan(lp.phi) ? -1 : lround(floor(lp.phi * NODES)); if( i < 0 || i >= NODES ) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return lp; @@ -133,7 +132,8 @@ static PJ_LP robin_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers /* first guess, linear interp */ t = 5. * (lp.phi - T.c0)/(Y[i+1].c0 - T.c0); for (iters = MAX_ITER; iters ; --iters) { /* Newton-Raphson */ - t -= t1 = (V(T,t) - lp.phi) / DV(T,t); + const double t1 = (V(T,t) - lp.phi) / DV(T,t); + t -= t1; if (fabs(t1) < EPS) break; } diff --git a/src/projections/rpoly.cpp b/src/projections/rpoly.cpp index 58decf66..b065861f 100644 --- a/src/projections/rpoly.cpp +++ b/src/projections/rpoly.cpp @@ -34,7 +34,8 @@ static PJ_XY rpoly_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar xy.y = - P->phi0; } else { xy.y = 1. / tan(lp.phi); - xy.x = sin(fa = 2. * atan(fa * sin(lp.phi))) * xy.y; + fa = 2. * atan(fa * sin(lp.phi)); + xy.x = sin(fa) * xy.y; xy.y = lp.phi - P->phi0 + (1. - cos(fa)) * xy.y; } return xy; @@ -48,7 +49,9 @@ PJ *PROJECTION(rpoly) { return pj_default_destructor(P, ENOMEM); P->opaque = Q; - if ((Q->mode = (Q->phi1 = fabs(pj_param(P->ctx, P->params, "rlat_ts").f)) > EPS)) { + Q->phi1 = fabs(pj_param(P->ctx, P->params, "rlat_ts").f); + Q->mode = Q->phi1 > EPS; + if (Q->mode) { Q->fxb = 0.5 * sin(Q->phi1); Q->fxa = 0.5 / Q->fxb; } diff --git a/src/projections/sconics.cpp b/src/projections/sconics.cpp index 60965c21..f305e291 100644 --- a/src/projections/sconics.cpp +++ b/src/projections/sconics.cpp @@ -90,7 +90,8 @@ static PJ_LP sconics_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, (and ellipsoi struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); double rho; - rho = hypot (xy.x, xy.y = Q->rho_0 - xy.y); + xy.y = Q->rho_0 - xy.y; + rho = hypot (xy.x, xy.y); if (Q->n < 0.) { rho = - rho; xy.x = - xy.x; @@ -164,14 +165,16 @@ static PJ *setup(PJ *P, enum Type type) { Q->n = sin (Q->sig); Q->c2 = cos (del); Q->c1 = 1./tan (Q->sig); - if (fabs (del = P->phi0 - Q->sig) - EPS10 >= M_HALFPI) + del = P->phi0 - Q->sig; + if (fabs (del) - EPS10 >= M_HALFPI) return pj_default_destructor(P, PJD_ERR_LAT_0_HALF_PI_FROM_MEAN); Q->rho_0 = Q->c2 * (Q->c1 - tan (del)); break; case VITK1: - Q->n = (cs = tan (del)) * sin (Q->sig) / del; + cs = tan (del); + Q->n = cs * sin (Q->sig) / del; Q->rho_c = del / (cs * tan (Q->sig)) + Q->sig; Q->rho_0 = Q->rho_c - P->phi0; break; diff --git a/src/projections/somerc.cpp b/src/projections/somerc.cpp index be1f660d..fe6477fa 100644 --- a/src/projections/somerc.cpp +++ b/src/projections/somerc.cpp @@ -82,7 +82,9 @@ PJ *PROJECTION(somerc) { cp *= cp; Q->c = sqrt (1 + P->es * cp * cp * P->rone_es); sp = sin (P->phi0); - Q->cosp0 = cos( phip0 = aasin (P->ctx, Q->sinp0 = sp / Q->c) ); + Q->sinp0 = sp / Q->c; + phip0 = aasin (P->ctx, Q->sinp0); + Q->cosp0 = cos(phip0); sp *= P->e; Q->K = log (tan (M_FORTPI + 0.5 * phip0)) - Q->c * ( log (tan (M_FORTPI + 0.5 * P->phi0)) - Q->hlf_e * diff --git a/src/projections/stere.cpp b/src/projections/stere.cpp index a1bd41e7..d95bb7fa 100644 --- a/src/projections/stere.cpp +++ b/src/projections/stere.cpp @@ -44,13 +44,14 @@ static double ssfn_ (double phit, double sinphi, double eccen) { static PJ_XY stere_e_forward (PJ_LP lp, PJ *P) { /* Ellipsoidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double coslam, sinlam, sinX = 0.0, cosX = 0.0, X, A = 0.0, sinphi; + double coslam, sinlam, sinX = 0.0, cosX = 0.0, A = 0.0, sinphi; coslam = cos (lp.lam); sinlam = sin (lp.lam); sinphi = sin (lp.phi); if (Q->mode == OBLIQ || Q->mode == EQUIT) { - sinX = sin (X = 2. * atan(ssfn_(lp.phi, sinphi, P->e)) - M_HALFPI); + const double X = 2. * atan(ssfn_(lp.phi, sinphi, P->e)) - M_HALFPI; + sinX = sin (X); cosX = cos (X); } @@ -116,7 +117,8 @@ oblcon: proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - xy.x = (xy.y = Q->akm1 / xy.y) * cosphi * sinlam; + xy.y = Q->akm1 / xy.y; + xy.x = xy.y * cosphi * sinlam; xy.y *= (Q->mode == EQUIT) ? sinphi : cosph0 * sinphi - sinph0 * cosphi * coslam; break; @@ -129,7 +131,8 @@ oblcon: proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } - xy.x = sinlam * ( xy.y = Q->akm1 * tan (M_FORTPI + .5 * lp.phi) ); + xy.y = Q->akm1 * tan (M_FORTPI + .5 * lp.phi); + xy.x = sinlam * xy.y; xy.y *= coslam; break; } @@ -147,11 +150,12 @@ static PJ_LP stere_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, invers switch (Q->mode) { case OBLIQ: case EQUIT: - cosphi = cos ( tp = 2. * atan2 (rho * Q->cosX1 , Q->akm1) ); + tp = 2. * atan2 (rho * Q->cosX1 , Q->akm1); + cosphi = cos (tp); sinphi = sin (tp); - if ( rho == 0.0 ) + if ( rho == 0.0 ) phi_l = asin (cosphi * Q->sinX1); - else + else phi_l = asin (cosphi * Q->sinX1 + (xy.y * sinphi * Q->cosX1 / rho)); tp = tan (.5 * (M_HALFPI + phi_l)); @@ -164,7 +168,8 @@ static PJ_LP stere_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, invers xy.y = -xy.y; /*-fallthrough*/ case S_POLE: - phi_l = M_HALFPI - 2. * atan (tp = - rho / Q->akm1); + tp = - rho / Q->akm1; + phi_l = M_HALFPI - 2. * atan (tp); halfpi = -M_HALFPI; halfe = -.5 * P->e; break; @@ -190,9 +195,11 @@ static PJ_LP stere_e_inverse (PJ_XY xy, PJ *P) { /* Ellipsoidal, invers static PJ_LP stere_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double c, rh, sinc, cosc; + double c, sinc, cosc; - sinc = sin (c = 2. * atan ((rh = hypot (xy.x, xy.y)) / Q->akm1)); + const double rh = hypot (xy.x, xy.y); + c = 2. * atan (rh / Q->akm1); + sinc = sin (c); cosc = cos (c); lp.lam = 0.; @@ -210,7 +217,8 @@ static PJ_LP stere_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers lp.phi = P->phi0; else lp.phi = asin (cosc * sinph0 + xy.y * sinc * cosph0 / rh); - if ((c = cosc - sinph0 * sin (lp.phi)) != 0. || xy.x != 0.) + c = cosc - sinph0 * sin (lp.phi); + if (c != 0. || xy.x != 0.) lp.lam = atan2 (xy.x * sinc * cosph0, c * rh); break; case N_POLE: @@ -248,8 +256,8 @@ static PJ *setup(PJ *P) { /* general initialization */ Q->akm1 = 2. * P->k0 / sqrt (pow (1+P->e,1+P->e) * pow (1-P->e,1-P->e)); else { - Q->akm1 = cos (Q->phits) / - pj_tsfn (Q->phits, t = sin (Q->phits), P->e); + t = sin (Q->phits); + Q->akm1 = cos (Q->phits) / pj_tsfn (Q->phits, t, P->e); t *= P->e; Q->akm1 /= sqrt(1. - t * t); } diff --git a/src/projections/sts.cpp b/src/projections/sts.cpp index 4d682a53..cbc36b7d 100644 --- a/src/projections/sts.cpp +++ b/src/projections/sts.cpp @@ -23,49 +23,48 @@ struct pj_opaque { static PJ_XY sts_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double c; - - xy.x = Q->C_x * lp.lam * cos(lp.phi); - xy.y = Q->C_y; - lp.phi *= Q->C_p; - c = cos(lp.phi); - if (Q->tan_mode) { - xy.x *= c * c; - xy.y *= tan (lp.phi); - } else { - xy.x /= c; - xy.y *= sin (lp.phi); - } - return xy; + + xy.x = Q->C_x * lp.lam * cos(lp.phi); + xy.y = Q->C_y; + lp.phi *= Q->C_p; + const double c = cos(lp.phi); + if (Q->tan_mode) { + xy.x *= c * c; + xy.y *= tan (lp.phi); + } else { + xy.x /= c; + xy.y *= sin (lp.phi); + } + return xy; } static PJ_LP sts_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, inverse */ PJ_LP lp = {0.0,0.0}; struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque); - double c; - - xy.y /= Q->C_y; - c = cos (lp.phi = Q->tan_mode ? atan (xy.y) : aasin (P->ctx, xy.y)); - lp.phi /= Q->C_p; - lp.lam = xy.x / (Q->C_x * cos(lp.phi)); - if (Q->tan_mode) - lp.lam /= c * c; - else - lp.lam *= c; - return lp; + + xy.y /= Q->C_y; + lp.phi = Q->tan_mode ? atan (xy.y) : aasin (P->ctx, xy.y); + const double c = cos (lp.phi); + lp.phi /= Q->C_p; + lp.lam = xy.x / (Q->C_x * cos(lp.phi)); + if (Q->tan_mode) + lp.lam /= c * c; + else + lp.lam *= c; + return lp; } static PJ *setup(PJ *P, double p, double q, int mode) { - P->es = 0.; - P->inv = sts_s_inverse; - P->fwd = sts_s_forward; - static_cast<struct pj_opaque*>(P->opaque)->C_x = q / p; - static_cast<struct pj_opaque*>(P->opaque)->C_y = p; - static_cast<struct pj_opaque*>(P->opaque)->C_p = 1/ q; - static_cast<struct pj_opaque*>(P->opaque)->tan_mode = mode; - return P; + P->es = 0.; + P->inv = sts_s_inverse; + P->fwd = sts_s_forward; + static_cast<struct pj_opaque*>(P->opaque)->C_x = q / p; + static_cast<struct pj_opaque*>(P->opaque)->C_y = p; + static_cast<struct pj_opaque*>(P->opaque)->C_p = 1/ q; + static_cast<struct pj_opaque*>(P->opaque)->tan_mode = mode; + return P; } diff --git a/src/projections/tcc.cpp b/src/projections/tcc.cpp index 3dd47940..9413b567 100644 --- a/src/projections/tcc.cpp +++ b/src/projections/tcc.cpp @@ -12,10 +12,10 @@ PROJ_HEAD(tcc, "Transverse Central Cylindrical") "\n\tCyl, Sph, no inv"; static PJ_XY tcc_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0, 0.0}; - double b, bt; - b = cos (lp.phi) * sin (lp.lam); - if ((bt = 1. - b * b) < EPS10) { + const double b = cos (lp.phi) * sin (lp.lam); + const double bt = 1. - b * b; + if (bt < EPS10) { proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION); return xy; } diff --git a/src/projections/tmerc.cpp b/src/projections/tmerc.cpp index 5c76b141..4b2a96f0 100644 --- a/src/projections/tmerc.cpp +++ b/src/projections/tmerc.cpp @@ -152,21 +152,21 @@ static PJ_XY approx_s_fwd (PJ_LP lp, PJ *P) { static PJ_LP approx_e_inv (PJ_XY xy, PJ *P) { PJ_LP lp = {0.0,0.0}; struct pj_opaque_approx *Q = static_cast<struct pj_opaque_approx*>(P->opaque); - double n, con, cosphi, d, ds, sinphi, t; lp.phi = pj_inv_mlfn(P->ctx, Q->ml0 + xy.y / P->k0, P->es, Q->en); if (fabs(lp.phi) >= M_HALFPI) { lp.phi = xy.y < 0. ? -M_HALFPI : M_HALFPI; lp.lam = 0.; } else { - sinphi = sin(lp.phi); - cosphi = cos(lp.phi); - t = fabs (cosphi) > 1e-10 ? sinphi/cosphi : 0.; - n = Q->esp * cosphi * cosphi; - d = xy.x * sqrt (con = 1. - P->es * sinphi * sinphi) / P->k0; + const double sinphi = sin(lp.phi); + const double cosphi = cos(lp.phi); + double t = fabs (cosphi) > 1e-10 ? sinphi/cosphi : 0.; + const double n = Q->esp * cosphi * cosphi; + double con = 1. - P->es * sinphi * sinphi; + const double d = xy.x * sqrt (con) / P->k0; con *= t; t *= t; - ds = d * d; + const double ds = d * d; lp.phi -= (con * ds / (1.-P->es)) * FC2 * (1. - ds * FC4 * (5. + t * (3. - 9. * n) + n * (1. - 4 * n) - ds * FC6 * (61. + t * (90. - 252. * n + diff --git a/src/projections/tpeqd.cpp b/src/projections/tpeqd.cpp index 0bc3be87..58aeb8e1 100644 --- a/src/projections/tpeqd.cpp +++ b/src/projections/tpeqd.cpp @@ -28,7 +28,8 @@ static PJ_XY tpeqd_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwar z1 *= z1; z2 *= z2; - xy.x = Q->r2z0 * (t = z1 - z2); + t = z1 - z2; + xy.x = Q->r2z0 * t; t = Q->z02 - t; xy.y = Q->r2z0 * asqrt (4. * Q->z02 * z2 - t * t); if ((Q->ccs * sp - cp * (Q->cs * sin(dl1) - Q->sc * sin(dl2))) < 0.) @@ -53,14 +54,16 @@ static PJ_LP tpeqd_s_inverse (PJ_XY xy, PJ *P) { /* Spheroidal, invers /* lam--phi now in system relative to P1--P2 base equator */ sp = sin (lp.phi); cp = cos (lp.phi); - lp.phi = aasin (P->ctx, Q->sa * sp + Q->ca * cp * (s = cos(lp.lam -= Q->lp))); + lp.lam -= Q->lp; + s = cos(lp.lam); + lp.phi = aasin (P->ctx, Q->sa * sp + Q->ca * cp * s); lp.lam = atan2 (cp * sin(lp.lam), Q->sa * cp * s - Q->ca * sp) + Q->lamc; return lp; } PJ *PROJECTION(tpeqd) { - double lam_1, lam_2, phi_1, phi_2, A12, pp; + double lam_1, lam_2, phi_1, phi_2, A12; struct pj_opaque *Q = static_cast<struct pj_opaque*>(pj_calloc (1, sizeof (struct pj_opaque))); if (nullptr==Q) return pj_default_destructor(P, ENOMEM); @@ -94,7 +97,8 @@ PJ *PROJECTION(tpeqd) { Q->hz0 = .5 * Q->z02; A12 = atan2(Q->cp2 * sin (Q->dlam2), Q->cp1 * Q->sp2 - Q->sp1 * Q->cp2 * cos (Q->dlam2)); - Q->ca = cos(pp = aasin(P->ctx, Q->cp1 * sin(A12))); + const double pp = aasin(P->ctx, Q->cp1 * sin(A12)); + Q->ca = cos(pp); Q->sa = sin(pp); Q->lp = adjlon ( atan2 (Q->cp1 * cos(A12), Q->sp1) - Q->hz0); Q->dlam2 *= .5; diff --git a/src/projections/vandg2.cpp b/src/projections/vandg2.cpp index 05314833..223620d6 100644 --- a/src/projections/vandg2.cpp +++ b/src/projections/vandg2.cpp @@ -24,7 +24,8 @@ static PJ_XY vandg2_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwa double x1, at, bt, ct; bt = fabs(M_TWO_D_PI * lp.phi); - if ((ct = 1. - bt * bt) < 0.) + ct = 1. - bt * bt; + if (ct < 0.) ct = 0.; else ct = sqrt(ct); diff --git a/src/projections/vandg4.cpp b/src/projections/vandg4.cpp index a5cfd8e6..5d3e1885 100644 --- a/src/projections/vandg4.cpp +++ b/src/projections/vandg4.cpp @@ -32,7 +32,8 @@ static PJ_XY vandg4_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forwa dt = sqrt(dt * dt - 4.); if ((fabs(lp.lam) - M_HALFPI) < 0.) dt = -dt; dt2 = dt * dt; - x1 = bt + ct; x1 *= x1; + x1 = bt + ct; + x1 *= x1; t = bt + 3.*ct; ft = x1 * (bt2 + ct2 * dt2 - 1.) + (1.-bt2) * ( bt2 * (t * t + 4. * ct2) + diff --git a/src/projections/wag7.cpp b/src/projections/wag7.cpp index 3afb5a03..ee891243 100644 --- a/src/projections/wag7.cpp +++ b/src/projections/wag7.cpp @@ -11,15 +11,18 @@ PROJ_HEAD(wag7, "Wagner VII") "\n\tMisc Sph, no inv"; static PJ_XY wag7_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0, 0.0}; - double theta, ct, D; - (void) P; /* Shut up compiler warnnings about unused P */ - - theta = asin (xy.y = 0.90630778703664996 * sin(lp.phi)); - xy.x = 2.66723 * (ct = cos (theta)) * sin (lp.lam /= 3.); - xy.y *= 1.24104 * (D = 1/(sqrt (0.5 * (1 + ct * cos (lp.lam))))); - xy.x *= D; - return (xy); + (void) P; /* Shut up compiler warnnings about unused P */ + + xy.y = 0.90630778703664996 * sin(lp.phi); + const double theta = asin (xy.y); + const double ct = cos (theta); + lp.lam /= 3.; + xy.x = 2.66723 * ct * sin (lp.lam); + const double D = 1/(sqrt (0.5 * (1 + ct * cos (lp.lam)))); + xy.y *= 1.24104 * D; + xy.x *= D; + return (xy); } diff --git a/src/projections/wink2.cpp b/src/projections/wink2.cpp index 4aaf1972..7fd08649 100644 --- a/src/projections/wink2.cpp +++ b/src/projections/wink2.cpp @@ -20,25 +20,25 @@ struct pj_opaque { static PJ_XY wink2_s_forward (PJ_LP lp, PJ *P) { /* Spheroidal, forward */ PJ_XY xy = {0.0, 0.0}; - double k, V; - int i; - - xy.y = lp.phi * M_TWO_D_PI; - k = M_PI * sin (lp.phi); - lp.phi *= 1.8; - for (i = MAX_ITER; i ; --i) { - lp.phi -= V = (lp.phi + sin (lp.phi) - k) / - (1. + cos (lp.phi)); - if (fabs (V) < LOOP_TOL) - break; - } - if (!i) - lp.phi = (lp.phi < 0.) ? -M_HALFPI : M_HALFPI; - else - lp.phi *= 0.5; - xy.x = 0.5 * lp.lam * (cos (lp.phi) + static_cast<struct pj_opaque*>(P->opaque)->cosphi1); - xy.y = M_FORTPI * (sin (lp.phi) + xy.y); - return xy; + int i; + + xy.y = lp.phi * M_TWO_D_PI; + const double k = M_PI * sin (lp.phi); + lp.phi *= 1.8; + for (i = MAX_ITER; i ; --i) { + const double V = (lp.phi + sin (lp.phi) - k) / + (1. + cos (lp.phi)); + lp.phi -= V; + if (fabs (V) < LOOP_TOL) + break; + } + if (!i) + lp.phi = (lp.phi < 0.) ? -M_HALFPI : M_HALFPI; + else + lp.phi *= 0.5; + xy.x = 0.5 * lp.lam * (cos (lp.phi) + static_cast<struct pj_opaque*>(P->opaque)->cosphi1); + xy.y = M_FORTPI * (sin (lp.phi) + xy.y); + return xy; } @@ -48,8 +48,8 @@ PJ *PROJECTION(wink2) { return pj_default_destructor(P, ENOMEM); P->opaque = Q; - static_cast<struct pj_opaque*>(P->opaque)->cosphi1 = cos(pj_param(P->ctx, P->params, "rlat_1").f); - P->es = 0.; + static_cast<struct pj_opaque*>(P->opaque)->cosphi1 = cos(pj_param(P->ctx, P->params, "rlat_1").f); + P->es = 0.; P->inv = nullptr; P->fwd = wink2_s_forward; |