root/GPXSee
1
0
Fork 0
mirror of https://github.com/tumic0/GPXSee.git synced 2024-11-24 19:55:53 +01:00
Code Issues Projects Releases Activity

Do not depend on the initializing ellipsoid

Browse Source
This commit is contained in:
Martin Tůma 2018-03-19 19:13:48 +01:00
parent 3d06fe8831
commit 29efa84075
4 changed files with 19 additions and 18 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
src/map/albersequal.cpp
View File

@ -70,7 +70,6 @@ AlbersEqual::AlbersEqual(const Ellipsoid *ellipsoid, double standardParallel1,
double sp1, sp2; double sp1, sp2;
_e = ellipsoid;
_latitudeOrigin = deg2rad(latitudeOrigin); _latitudeOrigin = deg2rad(latitudeOrigin);
_longitudeOrigin = deg2rad(longitudeOrigin); _longitudeOrigin = deg2rad(longitudeOrigin);
_falseEasting = falseEasting; _falseEasting = falseEasting;
@ -79,7 +78,9 @@ AlbersEqual::AlbersEqual(const Ellipsoid *ellipsoid, double standardParallel1,
sp1 = deg2rad(standardParallel1); sp1 = deg2rad(standardParallel1);
sp2 = deg2rad(standardParallel2); 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); _es = sqrt(_es2);
_one_minus_es2 = 1 - _es2; _one_minus_es2 = 1 - _es2;
_two_es = 2 * _es; _two_es = 2 * _es;
@ -109,7 +110,7 @@ AlbersEqual::AlbersEqual(const Ellipsoid *ellipsoid, double standardParallel1,
_n = sin_lat1; _n = sin_lat1;
_C = sqr_m1 + _n * q1; _C = sqr_m1 + _n * q1;
_a_over_n = _e->radius() / _n; _a_over_n = ellipsoid->radius() / _n;
nq0 = _n * q0; nq0 = _n * q0;
_rho0 = (_C < nq0) ? 0 : _a_over_n * sqrt(_C - nq0); _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) if (rho != 0.0)
theta = atan2(dx, rho0_minus_dy); theta = atan2(dx, rho0_minus_dy);
rho_n = rho * _n; 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)); qc = 1 - ((_one_minus_es2) / (_two_es)) * log((1.0 - _es) / (1.0 + _es));
if (fabs(fabs(qc) - fabs(q)) > 1.0e-6) { if (fabs(fabs(qc) - fabs(q)) > 1.0e-6) {
q_over_2 = q / 2.0; q_over_2 = q / 2.0;

3
src/map/albersequal.h
View File

@ -18,13 +18,12 @@ public:
virtual Coordinates xy2ll(const QPointF &p) const; virtual Coordinates xy2ll(const QPointF &p) const;
private: private:
const Ellipsoid *_e;
double _latitudeOrigin; double _latitudeOrigin;
double _longitudeOrigin; double _longitudeOrigin;
double _falseEasting; double _falseEasting;
double _falseNorthing; double _falseNorthing;
double _a2;
double _rho0; double _rho0;
double _C; double _C;
double _n; double _n;

23
src/map/transversemercator.cpp
View File

@ -47,14 +47,14 @@ Defense.
#define SPHSN(lat) \ #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) \ #define SPHTMD(lat) \
((double)(_ap * lat - _bp * sin(2.e0 * lat) + _cp * sin(4.e0 * lat) \ ((double)(_ap * lat - _bp * sin(2.e0 * lat) + _cp * sin(4.e0 * lat) \
- _dp * sin(6.e0 * lat) + _ep * sin(8.e0 * lat))) - _dp * sin(6.e0 * lat) + _ep * sin(8.e0 * lat)))
#define DENOM(lat) \ #define DENOM(lat) \
((double)(sqrt(1.e0 - _es * pow(sin(lat),2)))) ((double)(sqrt(1.e0 - _es * pow(sin(lat),2))))
#define SPHSR(lat) \ #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, TransverseMercator::TransverseMercator(const Ellipsoid *ellipsoid,
@ -65,31 +65,32 @@ TransverseMercator::TransverseMercator(const Ellipsoid *ellipsoid,
double b; double b;
_e = ellipsoid; _a = ellipsoid->radius();
_longitudeOrigin = deg2rad(longitudeOrigin); _longitudeOrigin = deg2rad(longitudeOrigin);
_latitudeOrigin = deg2rad(latitudeOrigin); _latitudeOrigin = deg2rad(latitudeOrigin);
_scale = scale; _scale = scale;
_falseEasting = falseEasting; _falseEasting = falseEasting;
_falseNorthing = falseNorthing; _falseNorthing = falseNorthing;
_es = 2 * _e->flattening() - _e->flattening() * _e->flattening(); _es = 2 * ellipsoid->flattening() - ellipsoid->flattening()
* ellipsoid->flattening();
_ebs = (1 / (1 - _es)) - 1; _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; tn2 = tn * tn;
tn3 = tn2 * tn; tn3 = tn2 * tn;
tn4 = tn3 * tn; tn4 = tn3 * tn;
tn5 = tn4 * 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); * (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; * tn5 / 64.e0 ) / 2.e0;
_cp = 15.e0 * _e->radius() * (tn2 - tn3 + 3.e0 * (tn4 - tn5 ) / 4.e0) / 16.0; _cp = 15.e0 * _a * (tn2 - tn3 + 3.e0 * (tn4 - tn5 ) / 4.e0) / 16.0;
_dp = 35.e0 * _e->radius() * (tn3 - tn4 + 11.e0 * tn5 / 16.e0) / 48.e0; _dp = 35.e0 * _a * (tn3 - tn4 + 11.e0 * tn5 / 16.e0) / 48.e0;
_ep = 315.e0 * _e->radius() * (tn4 - tn5) / 512.e0; _ep = 315.e0 * _a * (tn4 - tn5) / 512.e0;
} }
QPointF TransverseMercator::ll2xy(const Coordinates &c) const QPointF TransverseMercator::ll2xy(const Coordinates &c) const

2
src/map/transversemercator.h
View File

@ -18,12 +18,12 @@ public:
virtual Coordinates xy2ll(const QPointF &p) const; virtual Coordinates xy2ll(const QPointF &p) const;
private: private:
const Ellipsoid *_e;
double _longitudeOrigin; double _longitudeOrigin;
double _latitudeOrigin; double _latitudeOrigin;
double _scale; double _scale;
double _falseEasting; double _falseEasting;
double _falseNorthing; double _falseNorthing;
double _a;
double _es; double _es;
double _ebs; double _ebs;