[Talk-de] Tilekoordinaten
Dimitri Junker
OSM at dimitri-junker.de
Sa Mai 8 21:09:12 UTC 2010
Hallo,
bei dem Versuch Worldfiles zu erzeugen kam die Frage nach der richtigen
Berechnung der OSM-Tilekoordinaten auf. Lt.
<http://wiki.openstreetmap.org/wiki/Slippy_map_tilenames#tile_numbers_to_lon
..2Flat>
lat_rad = arctan(sinh(Pi * (1 - 2 * ytile / n)))
mit n = 2 ^ zoom
In taho wird fast die gleiche Formel verwendet, aber statt Pi
ProjectF(85.0511), dies ergibt: 3.1415868112787813, also fast genau Pi, aber
eben nicht genau. Bei Zoomlevel 12 ergibt sich da z.B. ein Unterschied von
etwa 8m für Aachen. Ich hänge mal die entsprechenden Funktionen aus taho.pl
und der c++ Version an. Letztere habe ich aus der perl Version abgeleitet.
Wer denkt sich soetwas kompliziertes aus wenn es falsch ist? Also was ist
richtig?
Dimitri
C++ Version aus taho.exe
//##############################################################################
//# transform the x,y tilenumbers to the S,W,N,E edge ccordinates
//##############################################################################
void CPixmap::Project(int x, int y, int zoom, CGeoRect &gr)
{
double unit = 1. / (1<<zoom);
double relY1 = y * unit;
double relY2 = relY1 + unit;
double limitY = ProjectF(85.0511);
limitY=3.14159265;
double rangeY = 2 * limitY;
relY1 = limitY - rangeY * relY1;
relY2 = limitY - rangeY * relY2;
gr.m_n = ProjectMercToLat(relY1);
gr.m_s = ProjectMercToLat(relY2);
unit = 360. / (1<<zoom);
gr.m_w = -180 + x * unit;
gr.m_e=gr.m_w+unit;
}
//##############################################################################
//# helper sub for Project()
//##############################################################################
double CPixmap::ProjectF(double winkel)
{
double lat = 3.14159265/180 * winkel;
return(log(tan(lat) + 1/cos(lat)));
}
double CPixmap::ProjectMercToLat(double mercY)
{
return( 180/3.14159265* atan(sinh(mercY)));
}
Perl Version von taho.pl:
sub Project{
my ($X,$Y, $Zoom) = @_;
my $Unit = 1 / (2 ** $Zoom);
my $relY1 = $Y * $Unit;
my $relY2 = $relY1 + $Unit;
my $LimitY = ProjectF(85.0511);
my $RangeY = 2 * $LimitY;
$relY1 = $LimitY - $RangeY * $relY1;
$relY2 = $LimitY - $RangeY * $relY2;
my $Lat1 = ProjectMercToLat($relY1);
my $Lat2 = ProjectMercToLat($relY2);
$Unit = 360 / (2 ** $Zoom);
my $Long1 = -180 + $X * $Unit;
return(($Lat2, $Long1, $Lat1, $Long1 + $Unit)); # S,W,N,E
}
##############################################################################
# helper sub for Project()
##############################################################################
sub ProjectMercToLat($){
my $MercY = shift();
return( 180/pi* atan(sinh($MercY)));
}
##############################################################################
# helper sub for Project()
##############################################################################
sub ProjectF($){
my $Lat = pi/180 * shift();
my $Y = log(tan($Lat) + sec($Lat));
return($Y);
}
Mehr Informationen über die Mailingliste Talk-de