From 4de3ddf71b743d190b1ff2c0e12f3c08a961e9a8 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Martin=20T=C5=AFma?= Date: Mon, 17 Apr 2017 18:03:04 +0200 Subject: [PATCH] Switched from hardcoded values to CSV files for datum definitions. Non-WGS84 maps now fit correctly when defined using map grid. --- gpxsee.pro | 6 ++- pkg/datums.csv | 123 +++++++++++++++++++++++++++++++++++++++++++++ pkg/ellipsoids.csv | 30 +++++++++++ src/config.h | 20 +++++--- src/datum.cpp | 83 ++++++++++++++++++++++++++++++ src/datum.h | 38 ++++++++++++++ src/ellipsoid.cpp | 88 +++++++++++++++++++------------- src/ellipsoid.h | 33 +++++++----- src/gui.cpp | 42 ++++++++++++++++ src/gui.h | 1 + src/offlinemap.cpp | 75 ++++++++++++++++++--------- 11 files changed, 459 insertions(+), 80 deletions(-) create mode 100644 pkg/datums.csv create mode 100644 pkg/ellipsoids.csv create mode 100644 src/datum.cpp create mode 100644 src/datum.h diff --git a/gpxsee.pro b/gpxsee.pro index c2e8fb53..d8cdba29 100644 --- a/gpxsee.pro +++ b/gpxsee.pro @@ -92,7 +92,8 @@ HEADERS += src/config.h \ src/utm.h \ src/lambertconic.h \ src/ellipsoid.h \ - src/ozf.h + src/ozf.h \ + src/datum.h SOURCES += src/main.cpp \ src/gui.cpp \ src/poi.cpp \ @@ -158,7 +159,8 @@ SOURCES += src/main.cpp \ src/utm.cpp \ src/lambertconic.cpp \ src/ellipsoid.cpp \ - src/ozf.cpp + src/ozf.cpp \ + src/datum.cpp RESOURCES += gpxsee.qrc TRANSLATIONS = lang/gpxsee_cs.ts \ lang/gpxsee_sv.ts diff --git a/pkg/datums.csv b/pkg/datums.csv new file mode 100644 index 00000000..664af133 --- /dev/null +++ b/pkg/datums.csv @@ -0,0 +1,123 @@ +Adindan,4201,5,-162,-12,206 +Afgooye,4205,15,-43,-163,45 +Ain el Abd 1970,4204,14,-150,-251,-2 +Anna 1 Astro 1965,4708,2,-491,-22,435 +Arc 1950,4209,5,-143,-90,-294 +Arc 1960,4210,5,-160,-8,-300 +Ascension Island 1958,4712,14,-207,107,52 +Astro B4 Sorol Atoll,4707,14,114,-116,-333 +Astro Beacon 1945,4709,14,145,75,-272 +Astro DOS 71/4,4710,14,-320,550,-494 +Astronomic Stn 1952,4711,14,124,-234,-25 +Australian Geodetic 1966,4202,2,-133,-48,148 +Australian Geodetic 1984,4203,2,-134,-48,149 +Australian Geocentric 1994 (GDA94),4283,11,0,0,0 +Austrian,4312,3,594,84,471 +Bellevue (IGN),4714,14,-127,-769,472 +Bermuda 1957,4216,4,-73,213,296 +Bogota Observatory,4218,14,307,304,-318 +Campo Inchauspe,4221,14,-148,136,90 +Canton Astro 1966,4716,14,298,-304,-375 +Cape,4222,5,-136,-108,-292 +Cape Canaveral,4717,4,-2,150,181 +Carthage,4223,5,-263,6,431 +CH-1903,4149,3,674,15,405 +Chatham 1971,4672,14,175,-38,113 +Chua Astro,4224,14,-134,229,-29 +Corrego Alegre,4225,14,-206,172,-6 +Djakarta (Batavia),4211,3,-377,681,-50 +DOS 1968,,14,230,-199,-752 +Easter Island 1967,4719,14,211,147,111 +Egypt,,14,-130,-117,-151 +European 1950,4230,14,-87,-98,-121 +European 1950 (Mean France),,14,-87,-96,-120 +European 1950 (Spain and Portugal),,14,-84,-107,-120 +European 1979,4668,14,-86,-98,-119 +Finland Hayford,4123,14,-78,-231,-97 +Gandajika Base,4233,14,-133,-321,50 +Geodetic Datum 1949,4272,14,84,-22,209 +GGRS 87,4121,11,-199.87,74.79,246.62 +Guam 1963,4675,4,-100,-248,259 +GUX 1 Astro,4718,14,252,-209,-751 +Hartebeeshoek94,4148,20,0,0,0 +Hermannskogel,3906,3,653,-212,449 +Hjorsey 1955,4658,14,-73,46,-86 +Hong Kong 1963,4739,14,-156,-271,-189 +Hu-Tzu-Shan,4236,14,-634,-549,-201 +Indian Bangladesh,4682,6,289,734,257 +Indian Thailand,4240,6,214,836,303 +Israeli,4281,23,-235,-85,264 +Ireland 1965,4299,1,506,-122,611 +ISTS 073 Astro 1969,4724,14,208,-435,-229 +Johnston Island,4725,14,191,-77,-204 +Kandawala,4244,6,-97,787,86 +Kerguelen Island,4698,14,145,-187,103 +Kertau 1948,4245,7,-11,851,5 +L.C. 5 Astro,4726,4,42,124,147 +Liberia 1964,4251,5,-90,40,88 +Luzon Mindanao,,4,-133,-79,-72 +Luzon Philippines,4253,4,-133,-77,-51 +Mahe 1971,4256,5,41,-220,-134 +Marco Astro,4616,14,-289,-124,60 +Massawa,4262,3,639,405,60 +Merchich,4261,5,31,146,47 +Midway Astro 1961,4727,14,912,-58,1227 +Minna,4263,5,-92,-93,122 +NAD27 Alaska,,4,-5,135,172 +NAD27 Bahamas,,4,-4,154,178 +NAD27 Canada,,4,-10,158,187 +NAD27 Canal Zone,,4,0,125,201 +NAD27 Caribbean,,4,-7,152,178 +NAD27 Central,,4,0,125,194 +NAD27 CONUS,,4,-8,160,176 +NAD27 Cuba,,4,-9,152,178 +NAD27 Greenland,,4,11,114,195 +NAD27 Mexico,,4,-12,130,190 +NAD27 San Salvador,,4,1,140,165 +NAD83,4269,11,0,0,0 +Nahrwn Masirah Ilnd,,5,-247,-148,369 +Nahrwn Saudi Arbia,,5,-231,-196,482 +Nahrwn United Arab,,5,-249,-156,381 +Naparima BWI,4271,14,-2,374,172 +NGO1948,4273,27,315,-217,528 +NTF France,4275,24,-168,-60,320 +Norsk,4817,27,278,93,474 +NZGD1949,4272,14,84,-22,209 +NZGD2000,4167,20,0,0,0 +Observatorio 1966,4182,14,-425,-169,81 +Old Egyptian,4229,12,-130,110,-13 +Old Hawaiian,4135,4,61,-285,-181 +Oman,4232,5,-346,-1,224 +Ord Srvy Grt Britn,4277,0,375,-111,431 +Pico De Las Nieves,4728,14,-307,-92,127 +Pitcairn Astro 1967,4729,14,185,165,42 +Potsdam Rauenberg DHDN,4314,3,606,23,413 +Prov So Amrican 1956,4248,14,-288,175,-376 +Prov So Chilean 1963,4254,14,16,196,93 +Puerto Rico,4139,4,11,72,-101 +Pulkovo 1942 (1),4284,15,28,-130,-95 +Pulkovo 1942 (2),4284,15,28,-130,-95 +Qatar National,4285,14,-128,-283,22 +Qornoq,4287,14,164,138,-189 +Reunion,4626,14,94,-948,-1262 +Rijksdriehoeksmeting,4289,3,593,26,478 +Rome 1940,4806,14,-225,-65,9 +RT 90,4124,3,498,-36,568 +S42,4179,15,28,-121,-77 +Santo (DOS),4730,14,170,42,84 +Sao Braz,4184,14,-203,141,53 +Sapper Hill 1943,4292,14,-355,16,74 +Schwarzeck,4293,21,616,97,-251 +South American 1969,4291,16,-57,1,-41 +South Asia,,8,7,-10,-26 +Southeast Base,4615,14,-499,-249,314 +Southwest Base,4183,14,-104,167,-38 +Timbalai 1948,4298,6,-689,691,-46 +Tokyo,4301,3,-128,481,664 +Tristan Astro 1968,4734,14,-632,438,-609 +Viti Levu 1916,4731,5,51,391,-36 +Wake-Eniwetok 1960,4732,13,101,52,-39 +WGS 72,4322,19,0,0,5 +WGS 84,4326,20,0,0,0 +Yacare,4309,14,-155,171,37 +Zanderij,4311,14,-265,120,-358 diff --git a/pkg/ellipsoids.csv b/pkg/ellipsoids.csv new file mode 100644 index 00000000..b761a2e7 --- /dev/null +++ b/pkg/ellipsoids.csv @@ -0,0 +1,30 @@ +0,Airy 1830,6377563.396,299.3249646 +1,Modified Airy,6377340.189,299.3249646 +2,Australian National,6378160.0,298.25 +3,Bessel 1841,6377397.155,299.1528128 +4,Clarke 1866,6378206.4,294.9786982 +5,Clarke 1880,6378249.145,293.465 +6,Everest (India 1830),6377276.345,300.8017 +7,Everest (1948),6377304.063,300.8017 +8,Modified Fischer 1960,6378155.0,298.3 +9,Everest (Pakistan),6377309.613,300.8017 +10,Indonesian 1974,6378160.0,298.247 +11,GRS 80,6378137.0,298.257222101 +12,Helmert 1906,6378200.0,298.3 +13,Hough 1960,6378270.0,297.0 +14,International 1924,6378388.0,297.0 +15,Krassovsky 1940,6378245.0,298.3 +16,South American 1969,6378160.0,298.25 +17,Everest (Malaysia 1969),6377295.664,300.8017 +18,Everest (Sabah Sarawak),6377298.556,300.8017 +19,WGS 72,6378135.0,298.26 +20,WGS 84,6378137.0,298.257223563 +21,Bessel 1841 (Namibia),6377483.865,299.1528128 +22,Everest (India 1956),6377301.243,300.8017 +23,Clarke 1880 Palestine,6378300.789,293.466 +24,Clarke 1880 IGN,6378249.2,293.466021 +25,Hayford 1909,6378388.0,296.959263 +26,Clarke 1858,6378350.87,294.26 +27,Bessel 1841 (Norway),6377492.0176,299.1528 +28,Plessis 1817 (France),6376523.0,308.6409971 +29,Hayford 1924,6378388.0,297.0 diff --git a/src/config.h b/src/config.h index 3f4334ba..63db33b3 100644 --- a/src/config.h +++ b/src/config.h @@ -13,6 +13,8 @@ #define FONT_SIZE 12 #define SCREEN_DPI 96.0 +#define ELLIPSOID_FILE QString("ellipsoids.csv") +#define DATUM_FILE QString("datums.csv") #define MAP_FILE QString("maps.txt") #define MAP_DIR QString("maps") #define POI_DIR QString("POI") @@ -29,12 +31,16 @@ #define GLOBAL_DIR QString("/usr/share/gpxsee") #endif -#define USER_MAP_DIR USER_DIR + QString("/") + MAP_DIR -#define USER_MAP_FILE USER_DIR + QString("/") + MAP_FILE -#define USER_POI_DIR USER_DIR + QString("/") + POI_DIR -#define GLOBAL_MAP_DIR GLOBAL_DIR + QString("/") + MAP_DIR -#define GLOBAL_MAP_FILE GLOBAL_DIR + QString("/") + MAP_FILE -#define GLOBAL_POI_DIR GLOBAL_DIR + QString("/") + POI_DIR -#define TILES_DIR USER_DIR + QString("/tiles") +#define USER_ELLIPSOID_FILE USER_DIR + QString("/") + ELLIPSOID_FILE +#define USER_DATUM_FILE USER_DIR + QString("/") + DATUM_FILE +#define USER_MAP_DIR USER_DIR + QString("/") + MAP_DIR +#define USER_MAP_FILE USER_DIR + QString("/") + MAP_FILE +#define USER_POI_DIR USER_DIR + QString("/") + POI_DIR +#define GLOBAL_ELLIPSOID_FILE GLOBAL_DIR + QString("/") + ELLIPSOID_FILE +#define GLOBAL_DATUM_FILE GLOBAL_DIR + QString("/") + DATUM_FILE +#define GLOBAL_MAP_DIR GLOBAL_DIR + QString("/") + MAP_DIR +#define GLOBAL_MAP_FILE GLOBAL_DIR + QString("/") + MAP_FILE +#define GLOBAL_POI_DIR GLOBAL_DIR + QString("/") + POI_DIR +#define TILES_DIR USER_DIR + QString("/tiles") #endif /* CONFIG_H */ diff --git a/src/datum.cpp b/src/datum.cpp new file mode 100644 index 00000000..c1e934a1 --- /dev/null +++ b/src/datum.cpp @@ -0,0 +1,83 @@ +#include +#include "wgs84.h" +#include "datum.h" + + +static QMap WGS84() +{ + QMap map; + map.insert("WGS 84", Datum(Ellipsoid(WGS84_RADIUS, WGS84_FLATTENING), + 0, 0, 0)); + return map; +} + +QMap Datum::_datums = WGS84(); +QString Datum::_errorString; +int Datum::_errorLine = 0; + +Datum Datum::datum(const QString &name) +{ + QMap::const_iterator it = _datums.find(name); + + if (it == _datums.end()) + return Datum(); + + return it.value(); +} + +bool Datum::loadList(const QString &path) +{ + QFile file(path); + bool res; + + if (!file.open(QFile::ReadOnly)) { + _errorString = qPrintable(file.errorString()); + return false; + } + + _errorLine = 1; + _errorString.clear(); + + while (!file.atEnd()) { + QByteArray line = file.readLine(); + QList list = line.split(','); + if (list.size() != 6) { + _errorString = "Format error"; + return false; + } + + int eid = list[2].trimmed().toInt(&res); + if (!res) { + _errorString = "Invalid ellipsoid id"; + return false; + } + double dx = list[3].trimmed().toDouble(&res); + if (!res) { + _errorString = "Invalid dx"; + return false; + } + double dy = list[4].trimmed().toDouble(&res); + if (!res) { + _errorString = "Invalid dy"; + return false; + } + double dz = list[5].trimmed().toDouble(&res); + if (!res) { + _errorString = "Invalid dz"; + return false; + } + + Ellipsoid e = Ellipsoid::ellipsoid(eid); + if (e.isNull()) { + _errorString = "Unknown ellipsoid ID"; + return false; + } + + Datum d(e, dx, dy, dz); + _datums.insert(list[0].trimmed(), d); + + _errorLine++; + } + + return true; +} diff --git a/src/datum.h b/src/datum.h new file mode 100644 index 00000000..128b003d --- /dev/null +++ b/src/datum.h @@ -0,0 +1,38 @@ +#ifndef DATUM_H +#define DATUM_H + +#include +#include "ellipsoid.h" + +class Datum +{ +public: + Datum() : _ellipsoid(Ellipsoid()), _dx(0), _dy(0), _dz(0) {} + Datum(const Ellipsoid &ellipsoid, double dx, double dy, double dz) + : _ellipsoid(ellipsoid), _dx(dx), _dy(dy), _dz(dz) {} + + const Ellipsoid &ellipsoid() const {return _ellipsoid;} + double dx() const {return _dx;} + double dy() const {return _dy;} + double dz() const {return _dz;} + + bool isNull() const {return _ellipsoid.isNull();} + bool isWGS84() const + {return _ellipsoid.isWGS84() && _dx == 0 && _dy == 0 && _dz == 0;} + + static bool loadList(const QString &path); + static const QString &errorString() {return _errorString;} + static int errorLine() {return _errorLine;} + + static Datum datum(const QString &name); + +private: + Ellipsoid _ellipsoid; + double _dx, _dy, _dz; + + static QMap _datums; + static QString _errorString; + static int _errorLine; +}; + +#endif // DATUM_H diff --git a/src/ellipsoid.cpp b/src/ellipsoid.cpp index b28c299a..844803d7 100644 --- a/src/ellipsoid.cpp +++ b/src/ellipsoid.cpp @@ -1,42 +1,62 @@ -#include -#include -#include "wgs84.h" +#include #include "ellipsoid.h" +QMap Ellipsoid::_ellipsoids; +QString Ellipsoid::_errorString; +int Ellipsoid::_errorLine = 0; -#define INTERNATIONAL_RADIUS 6378388.0 -#define INTERNATIONAL_FLATTENING (1.0/297.0) -#define KRASSOVSKY_RADIUS 6378245.0 -#define KRASSOVSKY_FLATTENING (1.0/298.3) -#define CLARKE1866_RADIUS 6378206.4 -#define CLARKE1866_FLATTENING (1.0/294.9786982) -#define GRS80_RADIUS 6378137.0 -#define GRS80_FLATTENING (1.0/298.257222101) -#define WGS72_RADIUS 6378135.0 -#define WGS72_FLATTENING (1.0/298.26) -#define GRS67_RADIUS 6378160.0 -#define GRS67_FLATTENING (1.0/298.25) - - -// 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() +Ellipsoid Ellipsoid::ellipsoid(int id) { - _radius = WGS84_RADIUS; - _flattening = WGS84_FLATTENING; + QMap::const_iterator it = _ellipsoids.find(id); + + if (it == _ellipsoids.end()) + return Ellipsoid(); + + return it.value(); } -Ellipsoid::Ellipsoid(Name name) +bool Ellipsoid::loadList(const QString &path) { - _radius = ellipsoids[name].radius; - _flattening = ellipsoids[name].flattening; + QFile file(path); + bool res; + + if (!file.open(QFile::ReadOnly)) { + _errorString = qPrintable(file.errorString()); + return false; + } + + _errorLine = 1; + _errorString.clear(); + + while (!file.atEnd()) { + QByteArray line = file.readLine(); + QList list = line.split(','); + if (list.size() != 4) { + _errorString = "Format error"; + return false; + } + + int id = list[0].trimmed().toInt(&res); + if (!res) { + _errorString = "Invalid ellipsoid id"; + return false; + } + double radius = list[2].trimmed().toDouble(&res); + if (!res) { + _errorString = "Invalid ellipsoid radius"; + return false; + } + double flattening = list[3].trimmed().toDouble(&res); + if (!res) { + _errorString = "Invalid ellipsoid flattening"; + return false; + } + + Ellipsoid e(radius, 1.0/flattening); + _ellipsoids.insert(id, e); + + _errorLine++; + } + + return true; } diff --git a/src/ellipsoid.h b/src/ellipsoid.h index 75a55474..c128ede6 100644 --- a/src/ellipsoid.h +++ b/src/ellipsoid.h @@ -1,28 +1,37 @@ #ifndef ELLIPSOID_H #define ELLIPSOID_H +#include +#include +#include "wgs84.h" + class Ellipsoid { public: - enum Name { - Clarke1866, - GRS80, - International1924, - Krassowsky1940, - WGS84, - WGS72, - GRS67 - }; - - Ellipsoid(); - Ellipsoid(Name name); + Ellipsoid() : _radius(-1.0), _flattening(-1.0) {} + Ellipsoid(double radius, double flattening) + : _radius(radius), _flattening(flattening) {} double radius() const {return _radius;} double flattening() const {return _flattening;} + bool isNull() const {return _radius < 0 || _flattening < 0;} + bool isWGS84() const + {return _radius == WGS84_RADIUS && _flattening == WGS84_FLATTENING;} + + static bool loadList(const QString &path); + static const QString &errorString() {return _errorString;} + static int errorLine() {return _errorLine;} + + static Ellipsoid ellipsoid(int id); + private: double _radius; double _flattening; + + static QMap _ellipsoids; + static QString _errorString; + static int _errorLine; }; #endif // ELLIPSOID_H diff --git a/src/gui.cpp b/src/gui.cpp index 13680252..353e16fd 100644 --- a/src/gui.cpp +++ b/src/gui.cpp @@ -26,6 +26,8 @@ #include "keys.h" #include "settings.h" #include "data.h" +#include "ellipsoid.h" +#include "datum.h" #include "map.h" #include "maplist.h" #include "mapdir.h" @@ -47,6 +49,7 @@ GUI::GUI() { + loadDatums(); loadMaps(); loadPOIs(); @@ -117,6 +120,45 @@ void GUI::createBrowser() _browser->setFilter(filter); } +void GUI::loadDatums() +{ + QString ef, df; + + if (QFile::exists(USER_ELLIPSOID_FILE)) + ef = USER_ELLIPSOID_FILE; + else if (QFile::exists(GLOBAL_ELLIPSOID_FILE)) + ef = GLOBAL_ELLIPSOID_FILE; + else + qWarning("No ellipsoids file found."); + + if (QFile::exists(USER_DATUM_FILE)) + df = USER_DATUM_FILE; + else if (QFile::exists(GLOBAL_DATUM_FILE)) + df = GLOBAL_DATUM_FILE; + else + qWarning("No datums file found."); + + if (!ef.isNull() && !df.isNull()) { + if (!Ellipsoid::loadList(ef)) { + if (Ellipsoid::errorLine()) + qWarning("%s: parse error on line %d: %s", qPrintable(ef), + Ellipsoid::errorLine(), qPrintable(Ellipsoid::errorString())); + else + qWarning("%s: %s", qPrintable(ef), qPrintable( + Ellipsoid::errorString())); + } else { + if (!Datum::loadList(df)) { + if (Datum::errorLine()) + qWarning("%s: parse error on line %d: %s", qPrintable(ef), + Datum::errorLine(), qPrintable(Datum::errorString())); + else + qWarning("%s: %s", qPrintable(ef), qPrintable( + Datum::errorString())); + } + } + } +} + void GUI::loadMaps() { QList online, offline; diff --git a/src/gui.h b/src/gui.h index efe8ed05..0898d51f 100644 --- a/src/gui.h +++ b/src/gui.h @@ -79,6 +79,7 @@ private slots: private: typedef QPair DateRange; + void loadDatums(); void loadMaps(); void loadPOIs(); void closeFiles(); diff --git a/src/offlinemap.cpp b/src/offlinemap.cpp index f4235bf3..6152c108 100644 --- a/src/offlinemap.cpp +++ b/src/offlinemap.cpp @@ -12,7 +12,7 @@ #include "wgs84.h" #include "coordinates.h" #include "matrix.h" -#include "ellipsoid.h" +#include "datum.h" #include "latlon.h" #include "mercator.h" #include "transversemercator.h" @@ -22,17 +22,39 @@ #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} -}; +// Abridged Molodensky transformation +static Coordinates toWGS84(Coordinates c, const Datum &datum) +{ + Ellipsoid WGS84(WGS84_RADIUS, WGS84_FLATTENING); + + double dX = datum.dx(); + double dY = datum.dy(); + double dZ = datum.dz(); + + double slat = sin(deg2rad(c.lat())); + double clat = cos(deg2rad(c.lat())); + double slon = sin(deg2rad(c.lon())); + double clon = cos(deg2rad(c.lon())); + double ssqlat = slat * slat; + + double from_f = datum.ellipsoid().flattening(); + double df = WGS84.flattening() - from_f; + double from_a = datum.ellipsoid().radius(); + double da = WGS84.radius() - from_a; + double from_esq = datum.ellipsoid().flattening() + * (2.0 - datum.ellipsoid().flattening()); + double adb = 1.0 / (1.0 - from_f); + double rn = from_a / sqrt(1 - from_esq * ssqlat); + double rm = from_a * (1 - from_esq) / pow((1 - from_esq * ssqlat), 1.5); + double from_h = 0.0; // we're flat! + + double dlat = (-dX * slat * clon - dY * slat * slon + dZ * clat + da * rn + * from_esq * slat * clat / from_a + +df * (rm * adb + rn / adb) * slat + * clat) / (rm + from_h); + double dlon = (-dX * slon + dY * clon) / ((rn + from_h) * clat); + + return Coordinates(c.lon() + rad2deg(dlon), c.lat() + rad2deg(dlat)); +} int OfflineMap::parseMapFile(QIODevice &device, QList &points, QString &projection, ProjectionSetup &setup, QString &datum) @@ -161,32 +183,30 @@ bool OfflineMap::createProjection(const QString &datum, return false; } - Ellipsoid::Name ellipsoid = Ellipsoid::WGS84; - for (size_t i = 0; i < ARRAY_SIZE(datums); i++) { - if (datum.startsWith(datums[i].name)) { - ellipsoid = datums[i].ellipsoid; - break; - } + Datum d = Datum::datum(datum); + if (d.isNull()) { + qWarning("%s: %s: unknown datum", qPrintable(_name), qPrintable(datum)); + return false; } if (projection == "Mercator") _projection = new Mercator(); else if (projection == "Transverse Mercator") - _projection = new TransverseMercator(Ellipsoid(ellipsoid), + _projection = new TransverseMercator(d.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(ellipsoid), + _projection = new LambertConic(d.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(Ellipsoid(ellipsoid), setup.zone); + _projection = new UTM(d.ellipsoid(), setup.zone); else if (!points.first().ll.isNull()) - _projection = new UTM(Ellipsoid(ellipsoid), points.first().ll); + _projection = new UTM(d.ellipsoid(), points.first().ll); else { qWarning("%s: Can not determine UTM zone", qPrintable(_name)); return false; @@ -197,9 +217,14 @@ bool OfflineMap::createProjection(const QString &datum, return false; } - for (int i = 0; i < points.size(); i++) - if (points.at(i).ll.isNull()) - points[i].ll = _projection->xy2ll(points.at(i).pp); + for (int i = 0; i < points.size(); i++) { + if (points.at(i).ll.isNull()) { + if (d.isWGS84()) + points[i].ll = _projection->xy2ll(points.at(i).pp); + else + points[i].ll = toWGS84(_projection->xy2ll(points.at(i).pp), d); + } + } return true; }