User:Maxbe/Kartenversuch
Contents |
Vorbemerkung und Warnung
Ich hab die Karte nur zum Spielen gemacht, um halt auch mal eine Karte zu rendern. Und nur für eigene Zwecke, deshalb reicht mir auch eine kleine Gegend (10°Ost bis 13°Ost und 47°Nord bis 48.655°Nord, erzeugt aus dem alps-Extrakt der Geofabrik) südlich und östlich von München sowie die Alpen südlich von München.
Ich kann kaum SQL und ich mag eigentlich keine Datenbanken. Viele Abfragen im Mapfile und im preprocessing sind echte quick'n dirty-Programmierung, wo ganze Tabellen durchsucht werden statt eines Indexes. Insbesondere wird häufig der spatial index von PostGIS ignoriert. Sowas ist pfui, das funktioniert mit winzigen Gebieten, skaliert aber nicht mit größeren Datenmengen. Wenn diese Programme auf ein größeres Gebiet in der Datenbank losgelassen werden, wird der Server stehen bleiben. Egal wie viel CPU, RAM und Plattenplatz Du hast, Dein Server wird stehen bleiben, wenn Du diese Skripte ohne Änderung übernimmst und auf Europa loslässt!
Trotzdem erzähle ich aber gern, wie ich die Karte mache. Vielleicht kann ja ein besserer Datenbanker was damit anfangen...
Worum geht es?
Ziel war eine Karte, die als Hintergrund für die typische "Das ist mein Laden"- oder "das ist der Wanderweg auf meinen Lieblingsberg"-Karte im Web dienen kann. Dazu werden relativ wenige POIs gebraucht (z.B. nicht der Konkurrenzladen in der Nachbarschaft), aber alle wichtigen Straße, Sehenswürdigkeiten, Bahnhöfe und Ortsnamen sollten drauf sein. Für die Berge sollten Höhenlinien und Schummerung drauf sein. Zu besichtigen ist die Karte hier: http://geo.dianacht.de/topo/
Das ganze sollte mit dem mapserver (Version 6) von mapserver.org laufen. Weil damit kenne ich mich aus und kann recht einfach Grafikexporte in bliebiger Größe
erzeugen. Betriebssystem ist Debian Lenny 32 Bit, die Datenbank sollte nicht mehr als 3GB auf der Platte brauchen und auf einen vServer der 10-Euro-Klasse 1-2 gleichzeitige User bedienen können.
Vorbereitung
OSM-Datenbank
Man braucht eine laufende PostGIS-DB mit OSM-Daten. Da habe ich mich an http://trac.osgeo.org/mapserver/wiki/RenderingOsmData gehalten und fürs die gelegentlichen Updates an http://wiki.openstreetmap.org/wiki/DE:HowTo_minutely_hstore
Das Einspielen und Updaten der Datenbank ist dort auch gut beschrieben. Letztendlich hatte ich eine Datenbank der gewüschten Gegend namens "osm" mit Koordinaten in EPSG 900913.
Höhendaten
Für den Import der Höhendaten braucht man die Tools von http://www.gdal.org/ und für die Schummerung die perrygeo-Tools von http://www.perrygeo.net/wordpress/ (zumindest das Programm "hillshade").
Die Gebrauchsanweisung für das Einspielen der Höhendaten steht bei http://wiki.openstreetmap.org/wiki/Shaded_relief_maps_using_mapnik
Dummerweise ist dort das Einspielen für OSM-Datenbanken mit EPSG 4326-Daten beschreiben, weswegen ich die Daten noch nachträglich konvertiere. Die eigentlichen Daten nehme ich in Form eines geoTIFFs von http://www.opendem.info die sind schön geglättet und geflickt an den Löchern der originalen NASA-Daten.
Das Einspielen sieht dann bei mir so aus:
Runterladen und auspacken
wget "http://koenigstuhl.geog.uni-heidelberg.de/~mover/srtm_germany_dtm.zip" unzip srtm_germany_dtm.zip mv srtm_germany_dtm.tif srtm.tif rm srtm_germany_dtm.zip
Höhenlinien im 25m-Abstand erzeugen und in die DB namens "srtm" einspielen:
psql srtm -c "drop table contours;" gdal_translate -projwin 10 48.655 13 47 -of GTiff -co "TILED=YES" -a_srs "+proj=latlong" srtm.tif input_adapted.tif gdal_contour -i 25 -snodata 32767 -a height input_adapted.tif input_adapted.shp shp2pgsql -p -I -g way input_adapted.shp contours | psql -q srtm shp2pgsql -a -g way input_adapted.shp contours | psql -q srtm
Konvertieren in spherical mercator
psql srtm -c "UPDATE geometry_columns SET SRID=900913 WHERE f_table_name='contours';" psql srtm -c "alter table contours drop constraint enforce_srid_way;" psql srtm -c "UPDATE contours set way=Transform(ST_SetSRID(way,4326),900913);"
In ein Zusätzliches Feld den "Typ" der Höhenlinie eintragen. "25", "50" oder "100", das erleichtert später das rendern von z.B. nur den 100er-Linien
psql srtm -c "alter table contours add column typ text;" psql srtm -c "update contours set typ='25';" psql srtm -c "update contours set typ='50' where CAST(height as int)%50=0;" psql srtm -c "update contours set typ='100' where CAST(height as int)%100=0;"
Die Schummerung wird ebenfalls aus der oben erzeugten Datei input_adapted.tif erzeugt, nach
http://wiki.openstreetmap.org/wiki/Hillshading_with_Mapnik
gdalwarp -of GTiff -co "TILED=YES" -srcnodata 32767 -t_srs "+proj=merc +ellps=sphere +R=6378137 +a=6378137 +units=m" \ -rcs -order 3 -tr 24 24 -wt Float32 -ot Float32 -wo SAMPLE_STEPS=100 -multi input_adapted.tif warped.tif hillshade warped.tif hillshade.tif -z 3 cp hillshade.tif /home/www/wms/schummerung/hillshade.tif
Jetzt hab ich also eine DB "osm mit den Openstreetmapdaten, eine DB "srtm" mit den Höhen und ein TIFF mit den Graustufen der Schummerung.
Preprocessing
Für meine Zwecke reichen erstmal zwei größere Aktionen im preprocessing. Die werden 1x täglich durchgeführt. Beide hab ich aus Diskussionen im Forum entnommen: http://forum.openstreetmap.org/viewtopic.php?id=13907 und http://forum.openstreetmap.org/viewtopic.php?id=14210
S- und U-Bahnhöfe
Zum einen wollte ich alle U-Bahnstationen mit einem "U" und alle S-Bahnen mit "S" beschriftet haben. Keine ganz leichte Aufgabe, weil sowas nur unvollständig in den Daten steht. Aber in kleinen Schritten kommt man hin:
Verwendet wird das Feld "service". Das braucht bei Bahnhöfen keiner und ist gut geeignet, die Begriffe "station", "ubahn" und "sbahn aufzunehmen".
S-Bahnhöfe sind Stationen, die zu Bahnstrecken gehören, die der MVV betreibt:
update osm_point set service='station' where railway='station'; update osm_point set service='sbahn' where railway='station' and exists \ (select tags from osm_rels where 'n'||osm_id=ANY(members) and \ ARRAY['route','train']<@tags and ARRAY['network','MVV']<@tags);
U-Bahnhöfe sind Stationen, die zu U-Bahnstrecken gehören, die der MVV betreibt:
update osm_point set service='ubahn' where railway='station' and exists \ (select tags from osm_rels where 'n'||osm_id=ANY(members) and \ ARRAY['route','subway']<@tags and ARRAY['network','MVV']<@tags);
Manche sind aber nicht Teil der Relation, sondern Teil der Schiene:
for i in `psql osm -c "select nodes from osm_ways where exists (select osm_id \
from osm_line where id=osm_id and railway='subway');" | \
sed 's/{//g; s/}/,/g; s/,/ /g; s/[^0-9 ]//g'` ; do
echo $i
psql osm -c "update osm_point set service='ubahn' where osm_id='"$i"' and railway='station' ;"
done
Und jetzt der Holzhammer für alles was übrigbleibt. Ich habe eine Tabelle "bahnhof" mit allen S- und U-Bahnhaltestellen:
update osm_point as o1 set service='sbahn' where railway='station' and service='station' \ and exists(select name from bahnhof where bahn='s-bahn' and name=o1.name) and not exists \ (select * from osm_point as o2 where railway='station' and service='sbahn' and o2.name=o1.name);
Sowas funktioniert natürlich nur, wenn man in seinem Gebiet nur einen Betreiber für S-Bahnen hat und nur eine Stadt mit U- und S-Bahn (an dieser Stelle nochmal der Hinweis auf die Vorbemerkung...).
Dominanz von Gipfeln
Um mit der Angabe in "ele" rechnen zu können, muss man dafür sorgen, dass nur Zahlen drinstehen. Sofern "300m" oder "2055 Meter" drin steht, kann man das "m" einfach wegwerfen, "566,6"kann man durch 566.6" ersetzen, aber den Rest muss man einfach rauswerfen und auf "222 estimated" verzichten:
psql osm -c "update osm_point set ele=replace(ele,' ',) where ele like '% %';" psql osm -c "update osm_point set ele=replace(ele,'meter',) where ele like '%m%';" psql osm -c "update osm_point set ele=replace(ele,'Meter',) where ele like '%M%';" psql osm -c "update osm_point set ele=replace(ele,'m',) where ele like '%m%';" psql osm -c "update osm_point set ele=replace(ele,',','.') where ele like '%,%';" psql osm -c "update osm_point set ele=NULL where ele similar to '%[a-zA-Z]%';" psql osm -c "update osm_point set ele=NULL where ele not similar to '[0-9.+-]+';"
1-Stellige Genauigkeit reicht auch für Höhen, sieht besser aus als ein übernommener Waypoint mit höchster Präzision "1234.44444563":
psql osm -c "update osm_point set ele=round(cast(ele as float)*10)/10 where ele is not null and round(cast(ele as float))!=cast(ele as float);"
zur besseren Übersichtlichkeit von Bergregionen in kleineren Maßstäben wird bei Gipfeln noch die Entfernung zum nächstgelegenen höheren Gipfel bestimmt und im unbenutzten Feld "population" abgelegt:
php /home/osmupdate/diffs/populatepeaks.php | psql osm
mit diesem "populatepeaks.php"
<?php
$db=pg_pconnect('host=localhost port=5432 dbname=osm user=osm password=xxx') or die("Keine DB, sonst alles in Ordnung");
$query=" select st_astext(st_transform(st_PointOnSurface(way),4326)) as p, osm_id,name,ele from osm_point where \"natural\"='peak';";
$result=pg_query($db,$query);
while($gipfel=pg_fetch_array($result)) {
if($gipfel["ele"]>100){
$query2="select st_distance(way,(select way from osm_point where osm_id=".$gipfel["osm_id"].")) as dominanz,ele,name \
from osm_point where osm_id!=".$gipfel["osm_id"]." and cast(ele as float)>".$gipfel["ele"]." and \"natural\"='peak' order by \
dominanz limit 1;";
$result2=pg_query($db,$query2);
$gipfel2=pg_fetch_array($result2);
$p=preg_split('/ /',preg_replace('/[A-Z\(\)]+/',,$gipfel["p"]));
$dom=cos($p[1]/180*3.14)*$gipfel2["dominanz"];
$d=3+$dom/100;if($d>30){$d=30;}
if($d=0){$d=30000;}
$query2="update osm_point set population=".(int)$dom." where osm_id=".$gipfel["osm_id"]." and \"natural\"='peak';";
echo $query2."\n";
}
}
pg_close($db);
echo "update osm_point set population=0 where \"natural\"='peak' and population is null;\n";
echo "update osm_point set population=100000 where \"natural\"='peak' and name like 'Gro%glockner' and population='0';\n";
?>
Dabei werden übrigens sämtliche Gipfel der Datenbank durchsucht und berechnet. Bei größeren Gebieten wäre das eher belastend. Wer einen höheren Berg als den Glockner in seiner DB hat, muss die letzte Zeile ändern...
Views
Für die ganzen Extrawürste, die beim Postprocessing gebraten werden, brauchen wir natürlich auch wieder views. Ausserdem brauchen wir die schon alleine deshalb, weil OSM das Feld "natural" angelegt hat, was in Anführungszeichen geschrieben werden muss. Anführungszeichen darf man aber in mapserver-mapfiles nicht verwenden. Und wir brauchen sie für einige Felder, die im Standard-osm2pgsql nicht importiert werden, sondern im hstore abgelegt werden.
Views haben keinen index und wenn man sie als Tabelle anlegt, brauchen sie Plattenplatz. Hier liegt der beste Ansatzpunkt für Optimierung. Weil im Gegensatz zu den bisherigen Dingen, wird das hier beim Rendern abgefragt, und nicht nur einmal beim Updaten oder Füllen der DB.
create or replace view osm_polygon_view as select osm_id, "natural" as nature, tags->'public_transport' as public_transport, \ landuse, waterway, sport, highway,railway, name, leisure, amenity, aeroway, wood, tags->'way_area' as area, way from osm_polygon; create or replace view osm_line_view as select way,osm_id,waterway,highway,ref,name,bridge,tunnel,railway,aeroway,z_order,man_made, \ layer,surface,tracktype,aerialway, tags->'public_transport' as public_transport,tags->'sac_scale' as sac_scale,"natural" as nature \ from osm_line; create or replace view osm_addr_view as select osm_id,way,"addr:housenumber" as housenumber from osm_point where "addr:housenumber" \ is not null union select osm_id,centroid(way) as way,"addr:housenumber" as housenumber from osm_polygon where \ "addr:housenumber" is not null; create or replace view osm_peak_view as select *,coalesce(name||'|'||ele,name,ele) as elename,population as isolation, \ tags->'mountain_pass' as mountain_pass,"natural" as nature from osm_point where "natural"='peak' or tags->'mountain_pass'='yes'; create or replace view osm_locality_view as select *,coalesce(name||'|'||ele,name,ele) as elename from osm_point \ where place is not null;
Rendern
Die Icons für Sportplätze, Kirchen... sind von http://www.sjjb.co.uk/mapicons/ Die sind CC-0, was praktisch ist, wenn man die Karte beliebig verwenden möchte. Ein paar sind selbstgemalt. Der Dank für die 2,5D-Darstellung der Häser gebührt User:SunCobalt.
Das Mapfile gibt es hier zum runterladen: User:Maxbe/Kartenversuch/mapfile. Allerdings sei auch hier nochmal zur Warnung auf die völlig fehlende Optimierung hingewiesen! Um zu vermeiden, dass doppelt getagte Kirchen auch doppelt auf der Karte erscheinen, wird für jede Kirche bei jedem rendern abgefragt, ob im Umriss des Gebäudes bereits ein Kirchenpunkt existiert.
select way,osm_id ,place,name,amenity,religion from osm_point where amenity='place_of_worship' and not exists (select way from osm_polygon where amenity='place_of_worship' and ST_Contains(osm_polygon.way,osm_point.way)) union select centroid(way) as way,osm_id ,place,name,amenity,religion from osm_polygon where amenity='place_of_worship';
Sowas gehört eigentlich ins preprocessing...
Wers brauche kann, kann sich bedienen. Der Stil steht unter CC-0-Lizenz, wer ihn findet, darf ihn behalten, verändern, verbessern, nutzen und verbreiten. Der Stand entspricht der Stilkritik im Forum http://forum.openstreetmap.org/viewtopic.php?id=14436 , die dort vorgebrachten Vorschläge wären also noch zu berücksichtigen.
Viel Spaß damit!
