From 29efa84075311a079eefd79f05011af91fff14a8 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Martin=20T=C5=AFma?= Date: Mon, 19 Mar 2018 19:13:48 +0100 Subject: [PATCH] Do not depend on the initializing ellipsoid --- src/map/albersequal.cpp | 9 +++++---- src/map/albersequal.h | 3 +-- src/map/transversemercator.cpp | 23 ++++++++++++----------- src/map/transversemercator.h | 2 +- 4 files changed, 19 insertions(+), 18 deletions(-) diff --git a/src/map/albersequal.cpp b/src/map/albersequal.cpp index cedde0d4..d4447c68 100644 --- a/src/map/albersequal.cpp +++ b/src/map/albersequal.cpp @@ -70,7 +70,6 @@ AlbersEqual::AlbersEqual(const Ellipsoid *ellipsoid, double standardParallel1, double sp1, sp2; - _e = ellipsoid; _latitudeOrigin = deg2rad(latitudeOrigin); _longitudeOrigin = deg2rad(longitudeOrigin); _falseEasting = falseEasting; @@ -79,7 +78,9 @@ AlbersEqual::AlbersEqual(const Ellipsoid *ellipsoid, double standardParallel1, sp1 = deg2rad(standardParallel1); sp2 = deg2rad(standardParallel2); - _es2 = 2 * _e->flattening() - _e->flattening() * _e->flattening(); + _a2 = ellipsoid->radius() * ellipsoid->radius(); + _es2 = 2 * ellipsoid->flattening() - ellipsoid->flattening() + * ellipsoid->flattening(); _es = sqrt(_es2); _one_minus_es2 = 1 - _es2; _two_es = 2 * _es; @@ -109,7 +110,7 @@ AlbersEqual::AlbersEqual(const Ellipsoid *ellipsoid, double standardParallel1, _n = sin_lat1; _C = sqr_m1 + _n * q1; - _a_over_n = _e->radius() / _n; + _a_over_n = ellipsoid->radius() / _n; nq0 = _n * q0; _rho0 = (_C < nq0) ? 0 : _a_over_n * sqrt(_C - nq0); } @@ -172,7 +173,7 @@ Coordinates AlbersEqual::xy2ll(const QPointF &p) const if (rho != 0.0) theta = atan2(dx, rho0_minus_dy); rho_n = rho * _n; - q = (_C - (rho_n * rho_n) / (_e->radius() * _e->radius())) / _n; + q = (_C - (rho_n * rho_n) / _a2) / _n; qc = 1 - ((_one_minus_es2) / (_two_es)) * log((1.0 - _es) / (1.0 + _es)); if (fabs(fabs(qc) - fabs(q)) > 1.0e-6) { q_over_2 = q / 2.0; diff --git a/src/map/albersequal.h b/src/map/albersequal.h index fb6af6b0..3842f05b 100644 --- a/src/map/albersequal.h +++ b/src/map/albersequal.h @@ -18,13 +18,12 @@ public: virtual Coordinates xy2ll(const QPointF &p) const; private: - const Ellipsoid *_e; - double _latitudeOrigin; double _longitudeOrigin; double _falseEasting; double _falseNorthing; + double _a2; double _rho0; double _C; double _n; diff --git a/src/map/transversemercator.cpp b/src/map/transversemercator.cpp index a6ae85a0..b276d3a4 100644 --- a/src/map/transversemercator.cpp +++ b/src/map/transversemercator.cpp @@ -47,14 +47,14 @@ Defense. #define SPHSN(lat) \ - ((double)(_e->radius() / sqrt(1.e0 - _es * pow(sin(lat), 2)))) + ((double)(_a / sqrt(1.e0 - _es * pow(sin(lat), 2)))) #define SPHTMD(lat) \ ((double)(_ap * lat - _bp * sin(2.e0 * lat) + _cp * sin(4.e0 * lat) \ - _dp * sin(6.e0 * lat) + _ep * sin(8.e0 * lat))) #define DENOM(lat) \ ((double)(sqrt(1.e0 - _es * pow(sin(lat),2)))) #define SPHSR(lat) \ - ((double)(_e->radius() * (1.e0 - _es) / pow(DENOM(lat), 3))) + ((double)(_a * (1.e0 - _es) / pow(DENOM(lat), 3))) TransverseMercator::TransverseMercator(const Ellipsoid *ellipsoid, @@ -65,31 +65,32 @@ TransverseMercator::TransverseMercator(const Ellipsoid *ellipsoid, double b; - _e = ellipsoid; + _a = ellipsoid->radius(); _longitudeOrigin = deg2rad(longitudeOrigin); _latitudeOrigin = deg2rad(latitudeOrigin); _scale = scale; _falseEasting = falseEasting; _falseNorthing = falseNorthing; - _es = 2 * _e->flattening() - _e->flattening() * _e->flattening(); + _es = 2 * ellipsoid->flattening() - ellipsoid->flattening() + * ellipsoid->flattening(); _ebs = (1 / (1 - _es)) - 1; - b = _e->radius() * (1 - _e->flattening()); + b = _a * (1 - ellipsoid->flattening()); - tn = (_e->radius() - b) / (_e->radius() + b); + tn = (_a - b) / (_a + b); tn2 = tn * tn; tn3 = tn2 * tn; tn4 = tn3 * tn; tn5 = tn4 * tn; - _ap = _e->radius() * (1.e0 - tn + 5.e0 * (tn2 - tn3) / 4.e0 + 81.e0 + _ap = _a * (1.e0 - tn + 5.e0 * (tn2 - tn3) / 4.e0 + 81.e0 * (tn4 - tn5) / 64.e0); - _bp = 3.e0 * _e->radius() * (tn - tn2 + 7.e0 * (tn3 - tn4) / 8.e0 + 55.e0 + _bp = 3.e0 * _a * (tn - tn2 + 7.e0 * (tn3 - tn4) / 8.e0 + 55.e0 * tn5 / 64.e0 ) / 2.e0; - _cp = 15.e0 * _e->radius() * (tn2 - tn3 + 3.e0 * (tn4 - tn5 ) / 4.e0) / 16.0; - _dp = 35.e0 * _e->radius() * (tn3 - tn4 + 11.e0 * tn5 / 16.e0) / 48.e0; - _ep = 315.e0 * _e->radius() * (tn4 - tn5) / 512.e0; + _cp = 15.e0 * _a * (tn2 - tn3 + 3.e0 * (tn4 - tn5 ) / 4.e0) / 16.0; + _dp = 35.e0 * _a * (tn3 - tn4 + 11.e0 * tn5 / 16.e0) / 48.e0; + _ep = 315.e0 * _a * (tn4 - tn5) / 512.e0; } QPointF TransverseMercator::ll2xy(const Coordinates &c) const diff --git a/src/map/transversemercator.h b/src/map/transversemercator.h index 1c11f319..35ce0a45 100644 --- a/src/map/transversemercator.h +++ b/src/map/transversemercator.h @@ -18,12 +18,12 @@ public: virtual Coordinates xy2ll(const QPointF &p) const; private: - const Ellipsoid *_e; double _longitudeOrigin; double _latitudeOrigin; double _scale; double _falseEasting; double _falseNorthing; + double _a; double _es; double _ebs;