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 /src | |
| 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...
Diffstat (limited to 'src')
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; |