Sem PostGIS Juliane
Wiki-Seite zum Selbststudium im GISpunkt-Seminar PostGIS. Zurück zu Sem PostGIS Selbststudium.
Inhaltsverzeichnis
Daten: OpenStreetMap von Zürich und Umgebung
Als Übungsdaten werden "LINK OpenStreetMap"-Daten der Schweiz verwendet. Diese stehen kostenfrei unter LINK geofabrik zum Download zur Verfügung. Der Datensatz switzerland.shp.zip enthält Gebäude (buildings.shp), Bodenbedeckung (natural.shp), Bahnen (railways.shp), Strassen (roads.shp), Gewässernetz (waterways.shp) und Points of Interest (points.shp). Aus dem Datensatz points.shp werden zur beispielhaften Bearbeitung der Aufgabe die Bushaltestellen extrahiert (vgl. Query 1).
Zudem werden die Schweizer Landesgrenze (ponet.e00) und die Grenzen der Siedlungsgebiete (pppoly.e00) des "LINK Digital Chart of the World" (DCW)-Datensatzes über LINK maproom bezogen. Diese Daten können ebenfalls kostenfrei heruntergeladen werden.
Metadaten-Records
OpenStreetMap
dc:title = OpenStreetMap
Buildings, Natural, Railways, Roads, Waterways, Points
dct:abstract = Die Dateien werden regelmässig mit dem Tool Osmosis aus dem so genannten "Planet File" ausgeschnitten, das
einen Gesamt-Abzug der "OpenStreetMap"-Datenbank darstellt. Solche Auszüge werden von der Geofabrik GmbH
erstellt, aktualisiert und als Shapefiles bereitgestellt. Folgende Elemente sind in den einzelnen Shapefiles
gespeichert und teilweise bereits vollständig attributiert:
buildings.shp: Gebäude, Öffentliche Gebäude, etc.
natural.shp: Wälder, Parkanlagen, Seen, etc.
railways.shp: Eisenbahnlinien, Tramlinien, Bahnhöfe, U-Bahnen, Schmalspurbahnen, etc.
roads.shp: Fusswege, Strassen (1. Klasse, 2. Klasse usw.), Radwege, etc.
waterways.shp: Kanäle, Bäche, Flüsse, Wehre, etc.
points.shp: Busstationen, Almhütten, Museen, Cafés, Hotels, Kindergärten, Universitäten, Tierparks, etc.
dc:format = Shapefile
dct:spatial = buildings.shp: northlimit=47.734204; southlimit=46.006147; eastlimit=10.409632; westlimit=6.037668
natural.shp: northlimit=47.808295; southlimit=45.905356; eastlimit=10.481553; westlimit=5.957567
railways.shp: northlimit=47.745569; southlimit=45.832066; eastlimit=10.287146; westlimit=5.993640
roads.shp: northlimit=47.808603; southlimit=45.821394; eastlimit=10.488137; westlimit=5.965783
waterways.shp: northlimit=47.783465; southlimit=45.904293; eastlimit=10.487938; westlimit=5.978676
points.shp: northlimit=47.792457; southlimit=45.832330; eastlimit=10.453054; westlimit=6.000755
dct:modified = 02-Feb-2009 07:41
dc:publisher = Geofabrik GmbH Karlsruhe
dc:language = en
dc:rights = Alle Daten aus dem "OpenStreetMap"-Projekt unterliegen der LINK Creative Commons Attribution Share-Alike-
Lizenz.
Digital Chart of the World
Bemerkung: Die Angabe beim DC-Element modified bezieht sich auf das Datum des Datenbezugs. Die Original DCW-Daten stammen vermutlich aus dem Jahr 1993.
dc:title = Digital Chart of the World
Ponet, Pppoly
dct:abstract = Verschiedene Kartenebenen zahlreicher Länder können im Arc/INFO Exportformat von ESRI's "Digital Chart of the
World"-Datensatz heruntergeladen werden.
ponet.e00: Landesgrenze der Schweiz
pppoly.e00: Siedlungsgebiete der Schweiz
dc:format = Arc/INFO coverage
dct:spatial = ponet.e00: northlimit=47.810918; southlimit=45.796475; eastlimit=10.552661; westlimit=5.912654
pppoly.e00: northlimit=47.721973; southlimit=45.982600; eastlimit=9.668165; westlimit=6.040542
dct:modified = 04-Feb-2009 11:28
dc:publisher = maproom, Pennsylvania State University
dc:language = en
dc:rights = Die ursprünglichen Daten der DCW wurden von der US Defence Mapping Agency erstellt und sind über die LINK
National Imagery and Mapping Agency erhältlich.
Laden der Datensätze
Die Coverages ponet.e00 und pppoly.e00 werden zuerst in LINK QuantumGIS geladen und als Shapefiles exportiert (rechte Maustaste auf entsprechenden Layer > Save as shapefile). Zu Übungszwecken für eine anschliessende Umprojektion (vgl. Query 2) wird als Coordinate reference System CH1903/LV03 (EPSG: 21781) ausgewählt.
Alle Shapefiles werden mit shp2pgsql.exe in sql-Files konvertiert.
Statement:
shp2pgsql -I -W codierung shapefilename tablename > outputsql.sql
Konvertieren der "OpenStreetMap"-Daten:
shp2pgsql –s 4326 –I –W ISO-8859-15 buildings gebaeude > buildings.sql shp2pgsql –s 4326 –I –W ISO-8859-15 natural flaechen > natural.sql shp2pgsql –s 4326 –I –W ISO-8859-15 points points > points.sql shp2pgsql –s 4326 –I –W ISO-8859-15 railways bahnen > railways.sql shp2pgsql –s 4326 –I –W ISO-8859-15 roads strassen > roads.sql shp2pgsql –s 4326 –I –W ISO-8859-15 waterways gewaesser > waterways.sql
Konvertieren der Daten der "Digital Chart of the World":
shp2pgsql –s 21781 –I –W ISO-8859-15 ponet grenze > ponet.sql shp2pgsql –s 21781 –I –W ISO-8859-15 pppoly siedlungen > pppoly.sql
Anschliessend werden die sql-Files über Terminal mit psql.exe in die PostGIS-Datenbank auf dem Server geladen.
Statement:
psql –U rolename –d datenbankname –h host –p port –f filename.sql
Queries
Query 1: Datensatz mit mindestens einem Thema
Der Datensatz points.shp beinhaltet zahlreiche Points of Interest (siehe Metadaten-Records). Für die weitere Bearbeitung werden alle Bushaltestellen selektiert und in eine neue Tabelle gespeichert.
-- Neue Tabelle busstationen anlegen CREATE TABLE busstationen ( gid serial NOT NULL, osm_id bigint, "name" character varying(32), "type" character varying(16), the_geom geometry, CONSTRAINT busstationen_pkey PRIMARY KEY (gid), CONSTRAINT enforce_dims_the_geom CHECK (ndims(the_geom) = 2), CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text OR the_geom IS NULL), CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 4326) ) WITHOUT OIDS; ALTER TABLE busstationen OWNER TO swai; -- Index kreieren CREATE INDEX busstationen_the_geom_gist ON busstationen USING gist (the_geom); -- Busstationen aus Tabelle points selektieren und in die zuvor angelegte Tabelle busstationen einfügen INSERT INTO busstationen ( SELECT gid, osm_id, "name", "type", the_geom FROM points WHERE "type" = 'bus_station' );
Query 2: Umprojizieren (CH03->WGS84)
Die "Digital Chart of the World"-Daten liegen im Koordinatensystem CH03/LV03 (EPSG: 21781) vor. Für den anschliessenden KML-Export müssen sie jedoch in das geografische Koordinatensystem WGS84 (EPSG: 4326) transformiert werden. Die Umprojektion wird mit Hilfe einer zweiten Geometriespalte (wgs84) durchgeführt.
-- Geometriespalte hinzufügen
SELECT AddGeometryColumn('public','ponet','wgs84',4326,'MULTILINESTRING', 2);
SELECT AddGeometryColumn('public','pppoly','wgs84',4326,'MULTILINESTRING', 2);
-- Daten umprojizieren
UPDATE ponet SET wgs84 = transform(the_geom, 4326);
UPDATE pppoly SET wgs84 = transform(the_geom, 4326);
Die "OpenStreetMap"-Daten müssen nicht umprojiziert werden. Sie liegen bereits im WGS84 Koordinatensystem vor.
Query 3: KML speichern
Die KML-Dateien werden mit der Funktion AsKML(text, text, geometrie) LINK google function erstellt. Mit dieser Funktion wird die Ausgabe KML-konform formatiert.
-- KML ponet SELECT askml( 'schweiz', 'landesgrenze', geomunion(transform(wgs84,4326))) FROM ponet;
Die Konvertierung alle anderen Datensätze dauert sehr lange. Aus diesem Grund werden alle weiteren KML’s aus den selektierten BBox Dateien (vgl. Query 4) erzeugt. Im Folgenden beispielhaft die Speicherung als KML für die Tabellen pppoly_zh und buildings_zh:
-- KML aus BBox pppoly_zh SELECT askml( 'zuerich', 'siedlungsflaeche', geomunion(transform(wgs84,4326))) FROM pppoly_zh; -- KML aus BBox buildings_zh SELECT askml( 'zuerich', 'gebaeude', geomunion(transform(the_geom,4326))) FROM buildings_zh;
Query 4: Ausschnitt anfragen (Perimeter, definiert mittels BBox)
Als Bounding Box-Umgebung wird die Region Zürich mit folgenden Eckkoordinaten definiert: 8.45 47.45, 8.6 47.35. Die in der Bounding Box enthaltenen Objekte der verschiedenen Tabellen Tankanlagen werden jeweils in eine separate Tabelle tabellenname_zh abgelegt. Anschliessend wird noch ein Primärschlüssel und ein räumlicher Index hinzugefügt (hier beispielhaft erläutert für die Tabellen pppoly und buildings).
-- Tabelle pppoly_zh aus BBox-Abfrage erstellen und füllen CREATE TABLE pppoly_zh AS ( SELECT * FROM pppoly WHERE wgs84 && SetSRID( 'BOX3D( 8.45 47.45, 8.6 47.35 )'::box3d, 4326 ) ); -- Primärschlüssel hinzufügen ALTER TABLE pppoly_zh ADD PRIMARY KEY( gid ); -- Räumlichen Index kreieren CREATE INDEX pppoly_zh_wgs84_gist ON pppoly_zh using gist (wgs84); -- Tabelle buildings_zh aus BBox-Abfrage erstellen und füllen CREATE TABLE buildings_zh AS ( SELECT * FROM buildings WHERE the_geom && SetSRID( 'BOX3D( 8.45 47.45, 8.6 47.35 )'::box3d, 4326 ) ); -- Primärschlüssel hinzufügen ALTER TABLE buildings_zh ADD PRIMARY KEY( gid ); -- Räumlichen Index kreieren CREATE INDEX buildings_zh_the_geom_gist ON buildings_zh using gist (the_geom);
Query 5: Buffer rund um ein Linien-Element
Es wird ein Buffer von ca. 400 m um die Limmat gelegt und anschliessend ebenfalls als KML-Datei gespeichert.
-- ca. 0.4 km Buffer um den Fluss Limmat SELECT Buffer(the_geom,0.005) as buffer, * FROM waterways_zh WHERE "name"='Limmat'; -- Buffer als KML speichern SELECT askml( 'zuerich', 'buffer limmat', geomunion(transform(Buffer(the_geom,0.005),4326))) FROM waterways_zh WHERE "name"='Limmat';
Bemerkung: Umrechnung ° in km
Längengradabstand a = cos(φ) × 2 π r / 360 = xx km Für Zürich φ = 47°40’ und r = 6370 km -> 0.1° ca. 75 km
