1
0
mirror of https://github.com/tumic0/GPXSee.git synced 2024-11-28 05:34:47 +01:00

Code cleanup

This commit is contained in:
Martin Tůma 2017-04-09 10:26:09 +02:00
parent 051a3ed303
commit 8822bf7e5f
9 changed files with 123 additions and 141 deletions

View File

@ -12,85 +12,31 @@
#define CLARKE1866_FLATTENING (1.0/294.9786982)
#define GRS80_RADIUS 6378137.0
#define GRS80_FLATTENING (1.0/298.257222101)
#define WGS70_RADIUS 6378135.0
#define WGS70_FLATTENING (1.0/298.26)
#define SAD69_RADIUS 6378160.0
#define SAD69_FLATTENING (1.0/298.25)
#define WGS72_RADIUS 6378135.0
#define WGS72_FLATTENING (1.0/298.26)
#define GRS67_RADIUS 6378160.0
#define GRS67_FLATTENING (1.0/298.25)
#ifndef M_PI_2
#define M_PI_2 1.57079632679489661923
#endif // M_PI_2
#ifndef M_PI_4
#define M_PI_4 0.78539816339744830962
#endif // M_PI_4
#define ARRAY_SIZE(array) (sizeof(array) / sizeof((array)[0]))
typedef struct {
const char *name;
double radius;
double flattening;
} Entry;
static Entry list[] = {
{"S42", KRASSOVSKY_RADIUS, KRASSOVSKY_FLATTENING},
{"Pulkovo 1942", KRASSOVSKY_RADIUS, KRASSOVSKY_FLATTENING},
{"European 1950", INTERNATIONAL_RADIUS, INTERNATIONAL_FLATTENING},
{"European 1979", INTERNATIONAL_RADIUS, INTERNATIONAL_FLATTENING},
{"NZGD1949", INTERNATIONAL_RADIUS, INTERNATIONAL_FLATTENING},
{"NAD27", CLARKE1866_RADIUS, CLARKE1866_FLATTENING},
{"NAD83", GRS80_RADIUS, GRS80_FLATTENING},
{"WGS 72", WGS70_RADIUS, WGS70_FLATTENING},
{"South American 1969", SAD69_RADIUS, SAD69_FLATTENING}
// Must be in Ellipsoid::Name order!
struct {double radius; double flattening;} static ellipsoids[] = {
{CLARKE1866_RADIUS, CLARKE1866_FLATTENING},
{GRS80_RADIUS, GRS80_FLATTENING},
{INTERNATIONAL_RADIUS, INTERNATIONAL_FLATTENING},
{KRASSOVSKY_RADIUS, KRASSOVSKY_FLATTENING},
{WGS84_RADIUS, WGS84_FLATTENING},
{WGS72_RADIUS, WGS72_FLATTENING},
{GRS67_RADIUS, GRS67_FLATTENING}
};
Ellipsoid::Ellipsoid()
{
_radius = WGS84_RADIUS;
_flattening = WGS84_FLATTENING;
}
Ellipsoid::Ellipsoid(const QString &datum)
Ellipsoid::Ellipsoid(Name name)
{
for (size_t i = 0; i < ARRAY_SIZE(list); i++) {
if (datum.startsWith(list[i].name)) {
_radius = list[i].radius;
_flattening = list[i].flattening;
return;
}
}
_radius = WGS84_RADIUS;
_flattening = WGS84_FLATTENING;
}
double Ellipsoid::q(double b) const
{
double e = sqrt(_flattening * (2. - _flattening));
double esb = e * sin(b);
return log(tan(M_PI_4 + b / 2.) * pow((1. - esb) / (1. + esb), e / 2.));
}
double Ellipsoid::iq(double q) const
{
double e = sqrt(_flattening * (2. - _flattening));
double b0 = 0.;
double b = 2. * atan(exp(q)) - M_PI_2;
do {
b0 = b;
double esb = e * sin(b);
b = 2. * atan(exp(q) * pow((1. - esb) / (1. + esb), -e / 2.)) - M_PI_2;
} while (fabs(b - b0) > 1e-10);
return b;
}
double Ellipsoid::nradius(double phi) const
{
double sin_phi = sin(phi);
return (_radius / sqrt(1. - (_flattening * (2. - _flattening)) * sin_phi
* sin_phi));
_radius = ellipsoids[name].radius;
_flattening = ellipsoids[name].flattening;
}

View File

@ -1,21 +1,25 @@
#ifndef ELLIPSOID_H
#define ELLIPSOID_H
class QString;
class Ellipsoid
{
public:
enum Name {
Clarke1866,
GRS80,
International1924,
Krassowsky1940,
WGS84,
WGS72,
GRS67
};
Ellipsoid();
Ellipsoid(const QString &datum);
Ellipsoid(Name name);
double radius() const {return _radius;}
double flattening() const {return _flattening;}
double q(double b) const;
double iq(double q) const;
double nradius(double phi) const;
private:
double _radius;
double _flattening;

View File

@ -2,24 +2,59 @@
#include "rd.h"
#include "lambertconic.h"
#ifndef M_PI_2
#define M_PI_2 1.57079632679489661923
#endif // M_PI_2
#ifndef M_PI_4
#define M_PI_4 0.78539816339744830962
#endif // M_PI_4
static double q(const Ellipsoid &el, double b)
{
double e = sqrt(el.flattening() * (2. - el.flattening()));
double esb = e * sin(b);
return log(tan(M_PI_4 + b / 2.) * pow((1. - esb) / (1. + esb), e / 2.));
}
static double iq(const Ellipsoid &el, double q)
{
double e = sqrt(el.flattening() * (2. - el.flattening()));
double b0 = 0.;
double b = 2. * atan(exp(q)) - M_PI_2;
do {
b0 = b;
double esb = e * sin(b);
b = 2. * atan(exp(q) * pow((1. - esb) / (1. + esb), -e / 2.)) - M_PI_2;
} while (fabs(b - b0) > 1e-10);
return b;
}
static double nradius(const Ellipsoid &el, double phi)
{
double sin_phi = sin(phi);
return (el.radius() / sqrt(1. - (el.flattening() * (2. - el.flattening()))
* sin_phi * sin_phi));
}
LambertConic::LambertConic(const Ellipsoid &ellipsoid, double standardParallel1,
double standardParallel2, double centralParallel, double centralMeridian,
double standardParallel2, double latitudeOrigin, double longitudeOrigin,
double scale, double falseEasting, double falseNorthing) : _e(ellipsoid)
{
_cm = centralMeridian;
_cm = longitudeOrigin;
_fe = falseEasting;
_fn = falseNorthing;
double sp1 = deg2rad(standardParallel1);
double sp2 = deg2rad(standardParallel2);
double N1 = _e.nradius(sp1);
double N2 = _e.nradius(sp2);
double N1 = nradius(_e, sp1);
double N2 = nradius(_e, sp2);
_q0 = ellipsoid.q(deg2rad(centralParallel));
double q1 = _e.q(sp1);
double q2 = _e.q(sp2);
_q0 = q(_e, deg2rad(latitudeOrigin));
double q1 = q(_e, sp1);
double q2 = q(_e, sp2);
_n = log((N1 * cos(sp1)) / (N2 * cos(sp2))) / (q2 - q1);
double R1 = N1 * cos(sp1) / _n;
@ -29,7 +64,7 @@ LambertConic::LambertConic(const Ellipsoid &ellipsoid, double standardParallel1,
QPointF LambertConic::ll2xy(const Coordinates &c) const
{
double dl = _n * (deg2rad(c.lon()) - deg2rad(_cm));
double R = _R0 * exp(_n * (_q0 - _e.q(deg2rad(c.lat()))));
double R = _R0 * exp(_n * (_q0 - q(_e, deg2rad(c.lat()))));
return QPointF(_fe + R * sin(dl), _fn + _R0 - R * cos(dl));
}
@ -42,5 +77,5 @@ Coordinates LambertConic::xy2ll(const QPointF &p) const
double R = sqrt(dx * dx + dy * dy);
double q = _q0 - log(R / _R0) / _n;
return Coordinates(rad2deg(deg2rad(_cm) + dl / _n), rad2deg(_e.iq(q)));
return Coordinates(rad2deg(deg2rad(_cm) + dl / _n), rad2deg(iq(_e, q)));
}

View File

@ -8,7 +8,7 @@ class LambertConic : public Projection
{
public:
LambertConic(const Ellipsoid &ellipsoid, double standardParallel1,
double standardParallel2, double centralParallel, double centralMeridian,
double standardParallel2, double latitudeOrigin, double longitudeOrigin,
double scale, double falseEasting, double falseNorthing);
virtual QPointF ll2xy(const Coordinates &c) const;

View File

@ -3,6 +3,9 @@
#include <QRectF>
#define ARRAY_SIZE(array) \
(sizeof(array) / sizeof((array)[0]))
double niceNum(double x, int round);
int str2int(const char *str, int len);
QRectF scaled(const QRectF &rect, qreal factor);

View File

@ -21,6 +21,18 @@
#include "offlinemap.h"
struct {const char *name; Ellipsoid::Name ellipsoid;} static datums[] = {
{"S42", Ellipsoid::Krassowsky1940},
{"Pulkovo 1942", Ellipsoid::Krassowsky1940},
{"European 1950", Ellipsoid::International1924},
{"European 1979", Ellipsoid::International1924},
{"NZGD1949", Ellipsoid::International1924},
{"NAD27", Ellipsoid::Clarke1866},
{"NAD83", Ellipsoid::GRS80},
{"WGS 72", Ellipsoid::WGS72},
{"South American 1969", Ellipsoid::GRS67}
};
int OfflineMap::parseMapFile(QIODevice &device, QList<ReferencePoint> &points,
QString &projection, ProjectionSetup &setup, QString &datum)
{
@ -108,12 +120,12 @@ int OfflineMap::parseMapFile(QIODevice &device, QList<ReferencePoint> &points,
} else if (key == "Projection Setup") {
if (list.count() < 8)
return ln;
setup.centralParallel = list.at(1).trimmed().toFloat(&res);
setup.latitudeOrigin = list.at(1).trimmed().toFloat(&res);
if (!res)
setup.centralParallel = 0;
setup.centralMeridian = list.at(2).trimmed().toFloat(&res);
setup.latitudeOrigin = 0;
setup.longitudeOrigin = list.at(2).trimmed().toFloat(&res);
if (!res)
setup.centralMeridian = 0;
setup.longitudeOrigin = 0;
setup.scale = list.at(3).trimmed().toFloat(&res);
if (!res)
setup.scale = 1.0;
@ -142,29 +154,39 @@ bool OfflineMap::createProjection(const QString &datum,
const QString &projection, const ProjectionSetup &setup,
QList<ReferencePoint> &points)
{
Ellipsoid::Name ellipsoid = Ellipsoid::WGS84;
if (points.count() < 2) {
qWarning("%s: insufficient number of reference points",
qPrintable(_name));
return false;
}
for (size_t i = 0; i < ARRAY_SIZE(datums); i++) {
if (datum.startsWith(datums[i].name)) {
ellipsoid = datums[i].ellipsoid;
break;
}
}
if (projection == "Mercator")
_projection = new Mercator();
else if (projection == "Transverse Mercator")
_projection = new TransverseMercator(Ellipsoid(datum),
setup.centralMeridian, setup.scale, setup.falseEasting,
_projection = new TransverseMercator(Ellipsoid(ellipsoid),
setup.longitudeOrigin, setup.scale, setup.falseEasting,
setup.falseNorthing);
else if (projection == "Latitude/Longitude")
_projection = new LatLon();
else if (projection == "Lambert Conformal Conic")
_projection = new LambertConic(Ellipsoid(datum), setup.standardParallel1,
setup.standardParallel2, setup.centralParallel, setup.centralMeridian,
setup.scale, setup.falseEasting, setup.falseNorthing);
_projection = new LambertConic(Ellipsoid(ellipsoid),
setup.standardParallel1, setup.standardParallel2,
setup.latitudeOrigin, setup.longitudeOrigin, setup.scale,
setup.falseEasting, setup.falseNorthing);
else if (projection == "(UTM) Universal Transverse Mercator") {
if (setup.zone)
_projection = new UTM(setup.zone);
_projection = new UTM(Ellipsoid(ellipsoid), setup.zone);
else if (!points.first().ll.isNull())
_projection = new UTM(points.first().ll);
_projection = new UTM(Ellipsoid(ellipsoid), points.first().ll);
else {
qWarning("%s: Can not determine UTM zone", qPrintable(_name));
return false;
@ -496,36 +518,6 @@ void OfflineMap::unload()
}
}
QRectF OfflineMap::bounds() const
{
return QRectF(QPointF(0, 0), _size);
}
qreal OfflineMap::zoomFit(const QSize &size, const QRectF &br)
{
Q_UNUSED(size);
Q_UNUSED(br);
return 1.0;
}
qreal OfflineMap::resolution(const QPointF &p) const
{
Q_UNUSED(p);
return _resolution;
}
qreal OfflineMap::zoomIn()
{
return 1.0;
}
qreal OfflineMap::zoomOut()
{
return 1.0;
}
void OfflineMap::draw(QPainter *painter, const QRectF &rect)
{
if (_tileSize.isValid()) {

View File

@ -22,13 +22,13 @@ public:
const QString &name() const {return _name;}
QRectF bounds() const;
qreal resolution(const QPointF &p) const;
QRectF bounds() const {return QRectF(QPointF(0, 0), _size);}
qreal resolution(const QPointF &) const {return _resolution;}
qreal zoom() const {return 1.0;}
qreal zoomFit(const QSize &size, const QRectF &br);
qreal zoomIn();
qreal zoomOut();
qreal zoomFit(const QSize &, const QRectF &) {return 1.0;}
qreal zoomIn() {return 1.0;}
qreal zoomOut() {return 1.0;}
QPointF ll2xy(const Coordinates &c) const
{return _transform.map(_projection->ll2xy(c));}
@ -57,8 +57,8 @@ private:
} ReferencePoint;
typedef struct {
double centralParallel;
double centralMeridian;
double latitudeOrigin;
double longitudeOrigin;
double scale;
double falseEasting;
double falseNorthing;

View File

@ -1,12 +1,14 @@
#include "ellipsoid.h"
#include "utm.h"
UTM::UTM(int zone) : _tm(Ellipsoid(), (qAbs(zone) - 1)*6 - 180 + 3, 0.9996,
500000, zone < 0 ? 10000000 : 0)
UTM::UTM(const Ellipsoid &ellipsoid, int zone)
: _tm(ellipsoid, (qAbs(zone) - 1)*6 - 180 + 3, 0.9996, 500000,
zone < 0 ? 10000000 : 0)
{
}
UTM::UTM(const Coordinates &c) : _tm(Ellipsoid(), 0, 0, 0, 0)
UTM::UTM(const Ellipsoid &ellipsoid, const Coordinates &c)
: _tm(ellipsoid, 0, 0, 0, 0)
{
int zone = int((c.lon() + 180)/6) + 1;
@ -24,6 +26,6 @@ UTM::UTM(const Coordinates &c) : _tm(Ellipsoid(), 0, 0, 0, 0)
}
double cm = (zone - 1)*6 - 180 + 3;
_tm = TransverseMercator(Ellipsoid(), cm, 0.9996, 500000,
_tm = TransverseMercator(ellipsoid, cm, 0.9996, 500000,
(c.lat() < 0) ? 10000000 : 0);
}

View File

@ -7,8 +7,8 @@
class UTM : public Projection
{
public:
UTM(int zone);
UTM(const Coordinates &c);
UTM(const Ellipsoid &ellipsoid, int zone);
UTM(const Ellipsoid &ellipsoid, const Coordinates &c);
virtual QPointF ll2xy(const Coordinates &c) const
{return _tm.ll2xy(c);}