Spring naar bijdragen

Aanbevolen berichten

Geplaatst (bewerkt)

De volgende hyperlink beschrijft de free-ware plug-in ttMaps voor het gebruik van rasterkaarten (dus bijvoorbeeld ook de 1:25.000 topokaarten van TOP25 to move van het Kadaster) op een TomTom:

 

ttMaps

 

ttMaps is een ontwikkeling binnen de Open Tom Community die al een tijd aan de gang is (zie: http://www.opentom.org/OpenTom:Community_Portal ).

 

Recent is de ontwikkelaar zover, dat de software goed genoeg wordt geacht om breder bekend te stellen (versie 0.0.1 d.d. 23.06.2008; versie 0.1.3 d.d. 13.08.2008, sinds enkele dagen versie 0.1.4). De applicatie bevat een logging functie en kan in een "voetganger modus" worden gebruikt (schakelen tussen static en non-static mode van de SIRFIII GPS chip).

 

In Nederland verwacht ik belangstelling voor ttMaps vooral om topografische kaartbladen te tonen. De DVD van het kadaster met alle kaarten van Nederland in tif + tfw formaat kost slechts € 99,95.

 

Ik heb de ontwikkelaar een bestandje gestuurd met een nederlandse vertaling van de user interface, dat al in versie 0.1.4 is opgenomen.

 

Heb (volgens de beschrijving op de site van ttMaps) al een aantal topo kaarten geconverteerd van .tif naar .ecw, dat ging zonder problemen.

 

De projectiedatabase in ttMaps is nu ook uitgebreid met de nederlandse geodetische gegevens. Er is nog wel een probleem met de juiste vertaling van de parameters (zoals door het Kadaster verstrekt voor RDNAPTRANS2004) naar PROJ.4 code die voor de berekening wordt gebruikt. Daardoor is op dit moment voor de nederlandse RD situatie een positie shift van zo'n 150 meter aanwezig.

bewerkt door NAV
ReizenReizen

Geplaatst

De ontwikkelaar van ttMaps stuit op het probleem dat met de gebruikte PROJ.4 code de transformatie van RD coördinaten naar WGS84 tot nu toe tot afwijkingen leidt van zo’n 150 meter. Via fora en overige internet informatie komen hij en ik niet veel verder (eerder wordt de verwarring groter).

 

Ik hoop dat op deze site een geodeet of GIS deskundige aan kan geven wat er misschien fout gaat.

 

Als basis heb ik de officiële parameters uit de RDNAPTRANS2004 procedure van het Kadaster aangeleverd. De ontwikkelaar heeft deze al op verschillende manieren naar PROJ.4 code vertaald, met tot dusver geen goede resultaten.

 

De laatste pogingen zijn door de ontwikkelaar als volgt gerapporteerd :

 

QUOTE

I have tested both parameters set in the package on Texel coordinates :

 

cs2cs +init=epsg:4326 +to +proj=sterea +lat_0=52.15616055555555

+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000

+ellps=bessel

+towgs84=565.2369,50.0087,465.658,1.9725,-1.7004,9.0677,4.0812 +units=m

+no_defs

 

4.8247619722 53.1607530501

117230.55 574972.34 -41.75

 

 

cs2cs +init=epsg:4326 +to +proj=sterea +lat_0=52.15616055555555

+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000

+ellps=bessel

+towgs84=593.0297,26.0038,478.7534,1.9725,-1.7004,9.0677,4.0812 +units=m

+no_defs

 

4.8247619722 53.1607530501

117256.91 574984.83 -67.62

 

 

The expected result is :

117380.1200 575040.3400

 

You can see that none of them give a correct result.

UNQUOTE

 

Wie heeft er kennis van RD naar WGS84 transformaties én van PROJ.4, en kan ons uit de brand helpen?

Geplaatst

Het ziet er allemaal veelbelovend uit.

Zelf ga ik het waarschijnlijk niet gebruiken (tenzij je weer terug naar de tomtom ROM kunt flashen).

 

Ik ga het zeker in de gaten houden

 

Bedankt voor de tip

Geplaatst
Zelf ga ik het waarschijnlijk niet gebruiken (tenzij je weer terug naar de tomtom ROM kunt flashen).

 

Er wordt niet naar ROM geflashed. Je hebt alleen voldoende vrij geheugen nodig. Dat betekent dat ttMaps met name geschikt is voor TomTom systemen met losse SD kaarten. De TomTom applicatie hoeft niet verwijderd te worden.

Geplaatst

Zelf ga ik het waarschijnlijk niet gebruiken (tenzij je weer terug naar de tomtom ROM kunt flashen).

 

Er wordt niet naar ROM geflashed. Je hebt alleen voldoende vrij geheugen nodig. Dat betekent dat ttMaps met name geschikt is voor TomTom systemen met losse SD kaarten. De TomTom applicatie hoeft niet verwijderd te worden.

 

Ok, dat is dan wel voor mij van toepassing.

Maar goed, de GPSmap60 bied voldoende mogelijkheden, ik zou echter ttmaps wel willen uitproberen.

 

Ik hou het in de gaten....

Geplaatst (bewerkt)

Kun je dan op je tomtom gewoon een map toevoegen zoals bv benelux of europa ook werkt of moet er ook nog software of iets dergelijks extra op?

bewerkt door Tornado Bram
Geplaatst (bewerkt)
Kun je dan op je tomtom gewoon een map toevoegen zoals bv benelux of europa ook werkt of moet er ook nog software of iets dergelijks extra op?

 

ttMaps is de software die geïnstalleerd moet worden en waarmee rasterkaarten kunnen worden gebruikt. Die moet je ook zelf toevoegen, en wel in .ecw formaat. Topo kaarten die als .tif bestand geleverd worden, kunnen redelijk eenvoudig naar .ecw worden omgezet. Je kunt ook zelf papieren kaarten scannen en als .tif opslaan, en vervolgens die omzetting doen. Als .tif kaarten niet vergezeld zijn zijn van een .tfw world bestand of in het formaat Geotiff zijn, dan zul je eerst de kaarten nog moeten calibreren (dit kan bijvoorbeeld door aan te geven wat de coördinaten van de linkeronder en rechterbovenhoek zijn.

bewerkt door NAV
Geplaatst
Werkt dit ook met de PDA versie of alleen met een 'echte' TomTom?

 

De software draait onder Linux, dus werkt alleen op "echte" TomToms en niet op PDA's of GSM's.

Geplaatst

Het onderwerp van deze draad was mij even ontgaan. Maar ik ben vandaag opgeport om mijn professorale licht over deze zaak te laten schijnen.

Eerst over de input files.

twf bevat linker bovenhoek en pixelgrootte (plus eventuele rotatie) in metrische coordinaten.

Hetzelfde met een GeoTIF: ook metrisch. Maar met de complicatie dat de notatie niet standaard is. Divers GeoTIF bij mij laten Ozi Exploreer klagen bij import: geen geodata in de GeoTIF terwijl vele andere progjes het wel kunnen vinden.

Bijvoorbeeld GeoTIF Examiner (Mentor software) laat zien wat er is.

Met ECW loopt het beter: dat is een fabrikant-eigen format. De geoinfo in een ECW bekijk ik met ECWSpy of met ECW Header Editor van Earth Resource Mapping. De info is wat een twf ook bevat plus de datum (naam) en de projectie (naam).

Maar als een importerend programma beide niet kent, dan moet je externe info aandragen.

 

De projectie

Onze kaarten zijn met een De DubbelProjectie volgens Schreiber gemaakt.

Dat is een heel obscure die voor ons prima werkt omdat ons landje nagenoeg in een vierkant past.

Gelukkig kunnen we deze projectie benaderen met iets anders. De fout die optreed is niet merkbaar op kaarten.

Het lukt prima met de wel bekende Transverse Mercator. Wel moet je dan de juiste parameters hebben. Dat is uitgezocht en die zijn:

Lat origin= 0.000000000 N

Central meridian 5.387633333

Scale Factor 0.999907900

False Easting 155000.00,

False Northing -5316592.00

 

Als software niet met de 0.00 N overweg kan (zero divide!)dan kan je de lat org. zetten op het nulpunt van ons raster in Amersfoort. Dan krijg je wel een andere False Northing. Welke? dat staat in de TRANSNAP documentatie.

 

Als het programma de Bessel ellipsoïde kent dan is dit voldoende om terug te kunnen rekenen van meters naar graden (op de Bessel ellips).

 

Ik heb hier de zaak in MatLab geprogrammerd. Dat ziet er dan zo uit:

 % rd2bessel
%
% Omrekenen geo coordinaten (tov Bessel) naar RD-grid coordinaten
% Formules van Geodetic Datum of Australia website
% 1 december 2001

function geo = rd2bessel(rd);

RDx = rd(1);
RDy = rd(2);

deg = pi/180;		% omrekenfactor naar radialen

% Bessel 1841 datum:
a = 6377397.155;
f = 1/299.1528128;

% RD grid parameters:
lab0   = 5.3876333*deg;		% origin
Efalse = 155000;	   		% false easting
K0	 = 0.9999079;			% scale factor
Nfalse = -5316592;			% false northing

% omrekenen:

e2 = f*(2-f);
e  = sqrt(e2);

a0 = 1 - 1/4*e^2 - 3/64*e^4 - 5/256*e^6;
a2 = 3/8*(e^2 + 1/4*e^4 + 15/128*e^6);
a4 = 15/256*(e^4 + 3/4*e^6);
a6 = 35/3072*e^6;

m = (RDy-Nfalse)/K0;

n   = f/(2-f);
G   = a*(1-n)*(1-n^2)*(1+9/4*n^2 + 225/64*n^4)*pi/180;
sig = m/G*pi/180;

% phic = phi-accent.
phic = sig +((3*n/2)-(27*n^3/32))*sin(2*sig) + ((21*n^2/16)-(55*n^4/32))*sin(4*sig)...
+(151*n^3/96)*sin(6*sig) + (1097*n^4/512)*sin(8*sig);

rho = a*(1-e^2)/(1-e^2*(sin(phic))^2)^(3/2);
nu  = a/sqrt(1-e^2*(sin(phic))^2);
psi = nu/rho;

% Grid to Geographical:

t = tan(phic);
Ec = RDx - Efalse;
x = Ec/(K0*nu);

term1 = (t/(K0*rho))*(x*Ec/2);
term2 = (t/(K0*rho))*(Ec*x^3/24)*(-4*psi^2+9*psi*(1-t^2)+12*t^2);
term3 = (t/(K0*rho))*((Ec*x^5)/720)*(8*psi^4*(11-24*t^2)-12*psi^3*(21-71*t^2)+15*(psi^2*...
15-98*t^2+15*t^4)+180*psi*(5*t^2-3*t^4)+360*t^4);
term4 = (t/(K0*rho))*(Ec*x^7/40320)*(1385+3633*t^2+4095*t^4+1575*t^6);

phi = phic - term1 + term2 - term3 + term4;

term1 = x*sec(phi);
term2 = (x^3/6)*sec(phi)*(psi+2*t^2);
term3 = (x^5/120)*sec(phi)*(-4*psi^3*(1-6*t^2)+psi^2*(9-68*t^2)+72*psi*t^2+24*t^4);
term4 = (x^7/5040)*sec(phi)*(61+662*t^2+1320*t^4+720*t^6);

omega = term1-term2+term3-term4;

lab = lab0 + omega;

geo(1) = phi/deg;		% lengte [graden]
geo(2) = lab/deg;		% breedte [graden]

 

Heb je belangstelling vor het omgekeerde (laat ik het voor de volledigheid maar geven):

% bessel2rd
%
% Omrekenen geo coordinaten (tov Bessel) naar RD-grid coordinaten
% Formules van Geodetic Datum of Australia website
% Gecontroleerd met GPS
% 7 juni 2000

function coordinaten = bessel2rd(geo);

deg = pi/180;		% omrekenfactor naar radialen

phi = geo(1)*deg;		% lengte [graden]
lab = geo(2)*deg;		% breedte [graden]

% Bessel 1841 datum:
a = 6377397.155;
f = 1/299.1528128;

% RD grid parameters:
lab0   = 5.3876333*deg;		% origin
Efalse = 155000;	   		% false easting
K0	 = 0.9999079;			% scale factor
Nfalse = -5316592;			% false northing

% omrekenen:

e2 = f*(2-f);
e  = sqrt(e2);

a0 = 1 - 1/4*e^2 - 3/64*e^4 - 5/256*e^6;
a2 = 3/8*(e^2 + 1/4*e^4 + 15/128*e^6);
a4 = 15/256*(e^4 + 3/4*e^6);
a6 = 35/3072*e^6;

m = a*(a0*phi - a2*sin(2*phi) + a4*sin(4*phi) - a6*sin(6*phi));

n   = f/(2-f);
G   = a*(1-n)*(1-n^2)*(1+9/4*n^2 + 225/64*n^4)*pi/180;
sig = m/G*pi/180;

rho = a*(1-e^2)/(1-e^2*(sin(phi))^2)^(3/2);
nu  = a/sqrt(1-e^2*(sin(phi))^2);
psi = nu/rho;

% Geographical to grid:

t = tan(phi);
omega = lab - lab0;

term1 = 1/6*omega^2*(cos(phi))^2*(psi-t^2);
term2 = 1/120*omega^4*(cos(phi))^4*(4*psi^3*(1-6*t^2)+psi^2*(1+8*t^2)-psi*2*t^2+t^4);
term3 = 1/5040*omega^6*(cos(phi))^6*(61-479*t^2+179*t^4-t^6);

E = (K0*nu*omega*cos(phi))*(1+term1+term2+term3) + Efalse;

term1 = 1/2*omega^2*nu*sin(phi)*cos(phi);
term2 = 1/24*omega^4*nu*sin(phi)*(cos(phi))^3*(4*psi^2+psi-t^2);
term3 = 1/720*omega^6*nu*sin(phi)*(cos(phi))^5*(8*psi^4*(11-24*t^2)-28*psi^3*(1-6*t^2)+psi^2*(1-32*t^2)-psi*(2*t^2)+t^4);
term4 = 1/40320*omega^8*nu*sin(phi)*(cos(phi))^7*(1385-3111*t^2+543*t^4-t^6);

N = K0*(m + term1 + term2 + term3 + term4) + Nfalse;

coordinaten = [E,N];

 

 

Heb je phi en lambda (in graden) op de Bessel Ellips dan zijn er standaard manieren om die positie over te zetten naar een andere ellipsoïde. In ons geval die van WGS84.

Voor ons is de 5-parameter transformatie voldoende nauwkeurig.

Zijn de parameters van beide ellipsen bekend dan dit door te rekenen. Hier gebeurt dit zo:

% molodensky datum omzetting van Bessel1841 naar WGS84
% Formules volgens website Peter Dana.
% Gecontroleerd met rekenvoorbeeld en GPS.
% In/uitvoer in graden en meters

% verbeterd nav check met EPS: fout in berekening Rm gehaald.
% 8 feb 2001
% 1 dec 2001: hoogte niet *deg!

function [to_geo]=bessel2wgs(from_geo);

deg = pi/180;					% omrekenfactor naar radialen

from_phi = from_geo(1)*deg;		% from_hoogte boven evenaar
from_lab = from_geo(2)*deg;		% from_oosterlengte
from_h =   from_geo(3); 		% from_hoogte boven zeespiegel

% parameters voor WGS 84:
to_a = 6378137.000;	  % semi-major axis of ellipsoid
to_f = 1/298.257223563;  % 1/flattening

% parameters voor Bessel 1841:
from_a = 6377397.155;
from_f = 1/299.1528128;

% correctiefactoren uitrekenen:
da  = to_a - from_a;
df  = to_f - from_f;
bda = 1-from_f;
from_es = 2*from_f - from_f^2;

dx = 593;
dy = 26;
dz = 478;

% Nu coordinaten omrekenen van Bessel 1841 naar WGS84

sphi = sin(from_phi);
cphi = cos(from_phi);
slab = sin(from_lab);
clab = cos(from_lab);

Rn = from_a/sqrt(1.0 - from_es*(sin(from_phi))^2);

% oud/fout: Rm = from_a*(1-from_es)/(1.0-from_es*(sin(from_phi)^2)^(3/2));
Rm = from_a*(1-from_es)/((1.0-from_es*(sin(from_phi)^2))^(3/2));

term1 = (((-dx*sphi*clab-dy*sphi*slab) + dz*cphi) + da*Rn*from_es*sphi*cphi/from_a);
term2 = df*(Rm/bda + Rn*bda)*sphi*cphi;

dphi = (term1 + term2)/(Rm+from_h);
dlab = (-dx*slab + dy*clab)/((Rn+from_h)*cphi);
dh = dx*cphi*clab + dy*cphi*slab + dz*sphi - da*from_a/Rn + df*bda*Rn*sphi*sphi;

to_phi = from_phi+dphi;
to_lab = from_lab + dlab;
to_h   = from_h + dh;

to_phi = to_phi/deg;
to_lab = to_lab/deg;
to_h   = from_h + dh;

to_geo = [to_phi,to_lab,to_h];

Je ziet hoe gauw je een foutje maakt in de programmering van deze formules!

In omgekeerde richting is het:

% molodensky datum omzetting van WGS84 naar Bessel 1841
% Formules volgens website Peter Dana.
% Gecontroleerd met rekenvoorbeeld en GPS (echter later nog fout in berekening Rm gehaald)
% In/uitvoer in graden en meters

% 1 dec 2001: hoogte niet *deg omrekenen!

function [to_geo]=wgs2bessel(from_geo);

deg = pi/180;					% omrekenfactor naar radialen

from_phi = from_geo(1)*deg;		% from_hoogte boven evenaar
from_lab = from_geo(2)*deg;		% from_oosterlengte
from_h =   from_geo(3); 		% from_hoogte boven zeespiegel

% parameters voor WGS 84:
from_a = 6378137.000;	  % semi-major axis of ellipsoid
from_f = 1/298.257223563;  % 1/flattening

% parameters voor Bessel 1841:
to_a = 6377397.155;
to_f = 1/299.1528128;

% correctiefactoren uitrekenen:
da  = to_a - from_a;
df  = to_f - from_f;
bda = 1-from_f;
from_es = 2*from_f - from_f^2;

dx = -593;
dy = -26;
dz = -478;

% Nu coordinaten omrekenen van WGS 84 naar Bessel 1841:

sphi = sin(from_phi);
cphi = cos(from_phi);
slab = sin(from_lab);
clab = cos(from_lab);

Rn = from_a/sqrt(1.0 - from_es*(sin(from_phi))^2);
% fout en oud: Rm = from_a*(1-from_es)/(1.0-from_es*(sin(from_phi)^2)^(3/2));
Rm = from_a*(1-from_es)/((1.0-from_es*(sin(from_phi)^2))^(3/2));

term1 = (((-dx*sphi*clab-dy*sphi*slab) + dz*cphi) + da*Rn*from_es*sphi*cphi/from_a);
term2 = df*(Rm/bda + Rn*bda)*sphi*cphi;

dphi = (term1 + term2)/(Rm+from_h);
dlab = (-dx*slab + dy*clab)/((Rn+from_h)*cphi);
dh = dx*cphi*clab + dy*cphi*slab + dz*sphi - da*from_a/Rn + df*bda*Rn*sphi*sphi;

to_phi = from_phi+dphi;
to_lab = from_lab + dlab;
to_phi = to_phi/deg;
to_lab = to_lab/deg;
to_h   = from_h + dh;

to_geo = [to_phi,to_lab,to_h];

 

Zodat mijn transfer functie simpel is:

% omzetten RD x-y coordinaten (in meters) naar graden (tov WGS84)
% 1 december 2001

function geo = rd2wgs(rd);

bessel = rd2bessel(rd);
bessel(3) = 0;	  %hoogte nul meter
geo = bessel2wgs(bessel);

 

Al deze bovenstaande formules zijn in de praktijk getest en tot op heden correct bevonden.

 

Er is een andere manier natuurlijk om van meters RD naar graden WGS84 te gaan door middel van een directe polynome aanpassing. Het voordeel hiervan is dat er minder rekenkracht nodig is voor een polynoom.

Mocht de omrekening te veel rekenkracht slurpen in de omzetting RD-> Bessel-graden (veel sinussen en zo) dan is er een reeksontwikkeleing mogelijk.

Die staat beschreven op pag 31 en verder van het document Coördinaattransformaties en kaartprojecties van de meetkundige dienst.

Geplaatst
ik ben vandaag opgeport om mijn professorale licht over deze zaak te laten schijnen.

 

Mijn dank is groot. We zullen de gegeven informatie uitgebreid bestuderen. Als er vragen blijven kom ik er op terug.

Geplaatst (bewerkt)
De ontwikkelaar van ttMaps stuit op het probleem dat met de gebruikte PROJ.4 code de transformatie van RD coördinaten naar WGS84 tot nu toe tot afwijkingen leidt van zo’n 150 meter. Via fora en overige internet informatie komen hij en ik niet veel verder (eerder wordt de verwarring groter).

.....

Wie heeft er kennis van RD naar WGS84 transformaties én van PROJ.4, en kan ons uit de brand helpen?

 

Een eerste test vandaag met opnieuw nieuwe parameters geeft goede resultaten :beerchug: . Ik zal deze week een en ander verder in de praktijk testen.

 

Dit is de Proj.4 procedure die nu gebruikt gaat worden:

cs2cs +init=epsg:4326 +to +proj=sterea +lat_0=52.15616055555555

+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000

+ellps=bessel

+towgs84=565.2369,50.0087,465.658,-0.406857330322398,0.350732676542563,-

+1.8703473836068,4.0812

+no_defs

bewerkt door NAV
Geplaatst (bewerkt)

De ontwikkelaar van ttMaps stuit op het probleem dat met de gebruikte PROJ.4 code de transformatie van RD coördinaten naar WGS84 tot nu toe tot afwijkingen leidt van zo’n 150 meter. Via fora en overige internet informatie komen hij en ik niet veel verder (eerder wordt de verwarring groter).

.....

Wie heeft er kennis van RD naar WGS84 transformaties én van PROJ.4, en kan ons uit de brand helpen?

 

Een eerste test vandaag met opnieuw nieuwe parameters geeft goede resultaten :cry: . Ik zal deze week een en ander verder in de praktijk testen.

 

Dit is de Proj.4 procedure die nu gebruikt gaat worden:

cs2cs +init=epsg:4326 +to +proj=sterea +lat_0=52.15616055555555

+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000

+ellps=bessel

+towgs84=565.2369,50.0087,465.658,-0.406857330322398,0.350732676542563,-

+1.8703473836068,4.0812

+no_defs

 

Eventueel heb ik ook nog wel de code beschikbaar die ik gebruikt heb voor de conversietool op http://www.speurvuur.nl/geostats/tools.php

Het is gebaseerd op de tool van gpsgek.

Mocht het nodig zijn stuur maar een PM

bewerkt door speurvuur
Geplaatst
Eventueel heb ik ook nog wel de code beschikbaar .....

Mocht het nodig zijn stuur maar een PM

 

Hartelijk dank voor het aanbod, maar het probleem is naar alle waarschijnlijkheid geheel opgelost met de laatste Proj.4 procedure.

Maak een account aan of meld je aan om een opmerking te plaatsen

Je moet lid zijn om een opmerking achter te kunnen laten

Account aanmaken

Maak een account aan in onze gemeenschap. Het is makkelijk!

Registreer een nieuw account

Aanmelden

Ben je al lid? Meld je hier aan.

Nu aanmelden
  • Onlangs hier   0 leden

    • Er kijken geen geregistreerde gebruikers naar deze pagina.
×
×
  • Nieuwe aanmaken...