PostGIS Terminal Examples

Aus Geoinformation HSR
Version vom 17. April 2017, 18:56 Uhr von Stefan (Diskussion | Beiträge) (Alle OeV-Haltestellen)

Wechseln zu: Navigation, Suche

Here are some examples queries for the PostGIS Terminal. All SQL queries work in the terminal (like this) but mostly also on any PostGIS database created with osm2pgsql.

See also PostGIS - Tipps und Tricks.

Fun with SQL Queries

All Roundabouts

-- This query works but all attributes from joined lines in first query are "lost".
-- NOTE: Add "ST_Within(way, (SELECT way FROM osm_polygon WHERE osm_id=-51701))" 
-- for exact results within Switzerland.
WITH lines AS (
    SELECT 
        osm_id as id, 
        way as geom 
    FROM osm_line
    WHERE 
      (tags @> hstore('junction','roundabout') 
      OR tags @> hstore('highway','mini_roundabout'))
      AND NOT (tags ? 'proposed')
),
lines2 AS (
    SELECT 
        osm_id as id, 
        way as geom 
    FROM osm_polygon
    WHERE 
      (tags @> hstore('junction','roundabout') 
      OR tags @> hstore('highway','mini_roundabout'))
      AND NOT (tags ? 'proposed') 
)
SELECT
    'Roundabout' as label, 
    ST_AsText((ST_Dump(ST_LineMerge(ST_Multi(ST_Collect(geom))))).geom ) as geom
FROM lines
UNION
SELECT
    id::text as label, 
    ST_AsText(ST_ExteriorRing(geom))
FROM lines2
LIMIT 9999

/*

Other solution attempts:

A. Most promising but complex and not yet working: WITH RECURSIVE http://workshops.boundlessgeo.com/postgis-intro/advanced_geometry_construction.html
http://blog.cleverelephant.ca/2010/07/network-walking-in-postgis.html

B. This could be still an approach to follow as part of WITH RECURSIVE :

(
    SELECT *, COUNT(1) FROM lines l1 
    JOIN lines l2 ON ST_Intersects(l1.geom, l2.geom) WHERE l1.id <> l2.id
    GROUP BY l1.id, l1.geom
) as tmp
WHERE count = 2

Source: http://gis.stackexchange.com/questions/147316/merging-linestrings-in-qgis-using-postgresql-and-postgis


C. This query works and is from OpenMapTiles and assumes lines to join have common "name, highway, ref" - which is not the optimal condition here: 

SELECT
    -- NOTE: Make sure no MULII... are generated
    (ST_Dump(geometry)).geom AS geometry,
    -- NOTE: The osm_id is no longer the original one which can make it difficult
    -- to lookup road names by OSM ID
    member_osm_ids[0] AS osm_id,
    member_osm_ids,
    name,
    ref,
    highway,
    z_order
FROM (
    SELECT
        ST_LineMerge(ST_Collect(way)) AS geometry,
        name,
        ref,
        highway,
        min(z_order) AS z_order,
        array_agg(DISTINCT osm_id) AS member_osm_ids
    FROM osm_line
    -- NOTE: We only care about highways (not railways) for labeling
    WHERE (name <> '' OR ref <> '') AND NULLIF(highway, '') IS NOT NULL
    GROUP BY name, highway, ref
) AS highway_union

*/ 

Alle Geonamen, die bereits in Schweizerdeutscher Mundart ("name:gsw") erfasst sind

 SELECT 
   ST_AsText(way) geom,
   tags->'name:gsw' AS label
 FROM osm_poi
 WHERE tags ? 'name:gsw'

Finding a string part of a tag key or a tag value

 SELECT osm_id, name, key
 FROM (SELECT osm_id, name, skeys(tags) AS key 
       from osm_point) tmp
 WHERE key LIKE 'generator%';
 SELECT osm_id, name, value
 FROM (SELECT osm_id, name, svals(tags) AS value 
       FROM osm_polygon) tmp
 WHERE value LIKE 'nuclear%';

Alle Zoos der Schweiz

(Tipp: http://labs.geometa.info/postgisterminal/?xapi=*[tourism=zoo]):

  SELECT ST_AsText(way) AS geom, name||' '||osm_id AS label
  FROM osm_all
  WHERE tags @> hstore('tourism','zoo')  

Schweizer Kernkraftwerke mit 40 Km-Puffer

  SELECT
   ST_AsText(ST_Buffer(ST_Centroid(way),40000)) AS geom,
   name AS label
  FROM osm_all
  WHERE tags @> hstore('generator:source','nuclear')  

Alle Restaurants mit Namen 'Rössli' der Schweiz

  SELECT ST_AsText(way) AS geom, name AS label
  FROM osm_point
  WHERE amenity = 'restaurant'
  AND name ILIKE '%rössli%'  

Alle 4000er Berggipfel der Schweiz

  SELECT ST_AsText(way) AS geom, name||','||ele AS label
  FROM osm_point
  WHERE "natural" = 'peak'
  AND  to_number('0'||ele, '99999999999.000')::int >= 4000
  
  oder:
  
  SELECT ST_AsText(way) AS geom, COALESCE(name, '')||' '||osm_id AS label
  FROM osm_poi
  WHERE tags @> hstore('natural', 'peak')
  AND CAST(regexp_replace(hstore("tags")->'ele', '[^0-9\.]', '', 'g') AS real) >= 4000  

Alle Aussichtspunkte im Kanton Zürich, die höher als 500 m ü.M. sind

  SELECT ST_AsText(a.way) geom, COALESCE(name, '')||' '||ele||' m ü.M' AS label
  FROM 
    osm_point AS a, 
    (SELECT way FROM osm_polygon 
     WHERE name = 'Zürich'
     AND tags @> hstore('admin_level','4')) AS b
  WHERE ST_Contains(b.way,a.way)
  AND tags @> hstore('tourism','viewpoint')
  AND  to_number('0'||ele, '99999999999.000')::int >= 500   

Alle Picnic-Plätze und Aussichtspunkte im aktuellen Kartenausschnitt

  SELECT ST_AsText(way) AS geom, name AS label
  FROM osm_point
  WHERE tourism IN ('picnic_site','viewpoint')
  AND ST_Contains(mapextent(), way)  

Alle Schulhäuser im Umkreis von 40 km aller Kernkraftwerke

  SELECT ST_AsText(a.way) AS geom, '' AS label
  FROM 
    osm_poi AS a, 
    (SELECT ST_Buffer(ST_Centroid(way),40000) AS way
     FROM osm_poi
     WHERE tags @> hstore('generator:source','nuclear')) AS b
  WHERE ST_Within(a.way,b.way)
    AND a.tags @> hstore('amenity', 'school')
  UNION
  SELECT ST_AsText(ST_Buffer(ST_Centroid(way),40000)) AS geom, COALESCE(name, '') AS label
  FROM osm_poi
  WHERE tags @> hstore('generator:source','nuclear')
  LIMIT 2000  

Die nächsten 10 Bars in der Nähe von 'mylocation' (ungeachtet der Distanz)

  SELECT ST_AsText(osm_poi.way) AS geom, COALESCE(name, '') AS label
  FROM 
    osm_poi, 
    (SELECT ST_Transform(ST_GeomFromText('POINT(8.81638 47.22666)', 4326), 900913) AS way) AS mylocation
  WHERE osm_poi.tags @> hstore('amenity', 'bar')
  ORDER BY osm_poi.way <-> mylocation.way
  LIMIT 10  

Die nächsten 100 Restaurants in der Nähe von 'mylocation' im Umkreis von max. 20 Km (Luftlinie)

  SELECT ST_AsText(osm_poi.way) AS geom, COALESCE(name, '') AS label
  FROM 
    osm_poi, 
    (SELECT ST_Transform(ST_GeomFromText('POINT(8.81638 47.22666)', 4326), 900913) AS way) AS mylocation
  WHERE ST_DWithin(osm_poi.way, mylocation.way, 20000)
  AND osm_poi.tags @> hstore('amenity', 'restaurant')
  ORDER BY osm_poi.way <-> mylocation.way
  LIMIT 100  

Alle Strassen mit Namen "Bahnhofstrasse" im Kt.ZH

 
  SELECT ST_AsText(ST_LineMerge(w.way))
  FROM osm_line w, 
       (SELECT way FROM osm_polygon 
        WHERE name = 'Zürich'
        AND tags @> hstore('admin_level','4')) AS p
  WHERE w.name = 'Bahnhofstrasse'
  AND (hstore("tags")->'highway') IS NOT NULL
  AND ST_Within(w.way, p.way)   

Bounding Box gegeben ein geografischer Name (Beispiel "Entlebuch")

  • Enter Entlebuch in http://boundingbox.klokantech.com/
  • Select DublinCore (and reduce to 4 digits) => westlimit=7.9891; southlimit=46.9419; eastlimit=8.3543; northlimit=47.1263
  • Place coords. in Polygon:
 
  SELECT ST_AsText( ST_Transform( ST_GeomFromText(
    'MULTIPOLYGON(((7.9891 47.1263, 8.3543 47.1263, 
                    8.3543 46.9419, 7.9891 46.9419, 
                    7.9891 47.1263)))', 4326), 900913) )   

Alle OeV-Haltestellen

  • Es hat 4000(!) mehr als in OSM als Didok/uic_name eingetragene Nodes?
  • Tags-Zusammenstellung aus DIDOK: https://github.com/datendelphin/didok/blob/master/db-import/import_stops.py
  • "SELECT way FROM osm_polygon WHERE osm_id=-51701" gibt die Schweizer Grenze zurück (EOSMDBOne enthält auch Teile des angrenzenden Auslands!).
  • Performance Note: The presence of a subquery is blocking the inlining of st_within, so putting the subquery into the from-clause instead.
 -- Just run this e.g. in PostGIS Terminal:
  SELECT count(*) 
  FROM 
    osm_point AS osm,
    (SELECT way FROM osm_polygon WHERE osm_id=-51701) AS ch
  WHERE ST_Within(osm.way, ch.way)
  AND NOT (tags ? 'uic_name')
  AND (tags @> 'railway=>station' OR tags @> 'railway=>halt' OR tags @> 'railway=>tram_stop' 
       OR tags @> 'highway=>bus_stop' OR tags @> 'amenity=>bus_station' 
       OR tags @> 'amenity=>ferry_terminal' OR tags @> 'railway=>funicular' OR tags @> 'aerialway=>station')
  -- 5291, query runtime: 1913ms.

  -- Use PostGIS Terminal to visualize the stops not in Didok (uic_name):
  SELECT ST_AsText(way) as geom, coalesce(name,'') as label
  FROM 
    osm_point AS osm,
    (SELECT way FROM osm_polygon WHERE osm_id=-51701) AS ch
  WHERE ST_Within(osm.way, ch.way)
  AND NOT (tags ? 'uic_name')
  AND (tags @> 'railway=>station' OR tags @> 'railway=>halt' OR tags @> 'railway=>tram_stop' 
       OR tags @> 'highway=>bus_stop' OR tags @> 'amenity=>bus_station' 
       OR tags @> 'amenity=>ferry_terminal' OR tags @> 'railway=>funicular' OR tags @> 'aerialway=>station')

  -- Get result as GeoCSV/WKT - to be run e.g. in psql or pgAdmin and export result to file:
  -- Round to mm precision (0.00000001 fractional part) from float using SnapToGrid. 
  -- See https://en.wikipedia.org/wiki/Decimal_degrees 
  FROM 
    osm_point AS osm,
    (SELECT way FROM osm_polygon WHERE osm_id=-51701) AS ch
  WHERE ST_Within(osm.way, ch.way)
  AND NOT (tags ? 'uic_name')
  AND (tags @> 'railway=>station' OR tags @> 'railway=>halt' OR tags @> 'railway=>tram_stop' 
       OR tags @> 'highway=>bus_stop' OR tags @> 'amenity=>bus_station' 
       OR tags @> 'amenity=>ferry_terminal' OR tags @> 'railway=>funicular' OR tags @> 'aerialway=>station')
  -- Resultat in https://gist.github.com/sfkeller/a18675d2190c4e121549e745474dc80f    

Die nächsen 15 Berge rund um Zürich

  • Next 15 mountains near Zurich
  • Kriterien: Höher als 1111m und mit Wikipedia-Tag (und natürlich mit Namen).
 SELECT 
    COALESCE(name, '') AS name,
    to_number('0'||ele, '99999999999.000')::int AS height,
    round(ST_Distance(start.geom, way))::int AS distance,
    CASE 
      WHEN round(degrees(ST_Azimuth(start.geom, way))) BETWEEN 315 AND 359 THEN 'north'
      WHEN round(degrees(ST_Azimuth(start.geom, way))) BETWEEN   0 AND  45 THEN 'north'
      WHEN round(degrees(ST_Azimuth(start.geom, way))) BETWEEN  45 AND 135 THEN 'east'
      WHEN round(degrees(ST_Azimuth(start.geom, way))) BETWEEN 135 AND 225 THEN 'south'
      WHEN round(degrees(ST_Azimuth(start.geom, way))) BETWEEN 225 AND 315 THEN 'west'
    END as direction,
    'https://de.wikipedia.org/wiki/'||(tags->'wikipedia') AS wikipedia,
    osm_id,
    ST_AsText(ST_SnapToGrid(ST_Transform(way,4326),0.00000001)) AS wkt 
  FROM 
    osm_point peak, 
    (SELECT ST_Transform(ST_GeomFromText('POINT(8.5393635 47.3781008)', 4326), 900913) AS geom) AS start
  WHERE tags @> hstore('natural', 'peak')
  AND to_number('0'||ele, '99999999999.000')::int >= 1111
  AND NULLIF(name, '') > ''
  AND tags ? 'wikipedia'
  ORDER BY way <-> start.geom
  LIMIT 15 

Special SQL Queries

Extracting numbers out of OSM values

See also regex above...

 with tmp(txt) as (
   values ('1'),('2 '),('9.999999'),('12345678901.999'),('5.5'),('0.111')
     ,('-1.1'),('-0.00000001'),('-9999')
     ,(' 5.65 '),('xxx')
 )
 select 'Numeric', txt, coalesce(substring(txt, '([-]?0\.\d*[1-9]\d*|[-]?[1-9]\d*(\.\d+)?)')::numeric, 0)
 from tmp
 union
 select 'Integer', txt, coalesce(substring(txt, '([-]?[0-9]{1,19})')::bigint, 0)
 from tmp
 union
 select 'Number ', txt, to_number('0'||txt, 'S99999999999.000')::numeric
 from tmp
 order by 1,2

Regex: ....

Name string handling

Select names of all OSM-objects containing 'zoo' at start middle or end (using wildcard '%' in String):

 SELECT ST_AsText(way) geom,
 COALESCE(name, ' ')||' '||osm_id label
 FROM osm_all
 WHERE hstore(tags)->'name' ILIKE '%zoo%'

Select all parking lots for disabled persons (Rollstuhlparkplatz / wheel parking) within visible map area (without FIXME). Note: Tag capacity:disabled currently occurs always together with amenity={parking, parking_space, parking_entrance}):

 SELECT ST_AsText(way) geom, 
        COALESCE(name,' ')||' ('||COALESCE(SUBSTRING((tags->'capacity:disabled') FROM E'[0-9]+'),' ')||')' label
 FROM osm_poi
 WHERE (tags @> '"capacity:disabled"=>"yes"' 
        OR SUBSTRING((tags->'capacity:disabled') FROM E'[0-9]+')::int > 0)
 AND ST_Contains(mapextent(), way)

Select all motorways (higways) with speed limit greater or equal than 100:

 SELECT ST_AsText(way) geom
 FROM osm_line
 WHERE tags @> '"highway"=>"motorway"'
 AND COALESCE(SUBSTRING((tags->'maxspeed') FROM E'[0-9]+')::int,0) >= 100
 -- 1806!

Count tags except xxx...:

 SELECT array_upper(%# tags, 1) AS "#tags", count(*) AS "count"
 FROM osm_point
 WHERE array_upper(%# tags, 1) > 0
 AND NOT EXISTS (SELECT skeys FROM skeys(tags) WHERE skeys LIKE 'xxx:%')
 GROUP BY 1
 ORDER BY 1


Point Clustering

Point Clustering with round/truncate to grid and centroid. Parameter 9000 defines the grid width (meter) and is dependent from the units of the coordinate reference systems. Fast - reduces 13358 Restaurants to 1163 - but with a "visual grid effect".

 SELECT
   ST_AsText( ST_Centroid(ST_Collect(position)) ) AS geom,
   CASE COUNT(position)
     WHEN 1 THEN COALESCE(max(name), '1') ELSE COUNT(position)::text
   END AS label,
   array_agg(osm_id) AS list_of_ids 
 FROM (
   SELECT way AS position, name, osm_id 
   FROM osm_poi
   WHERE tags @> hstore('amenity','restaurant')
 ) tmp
 GROUP BY ST_SnapToGrid(position, 9000);

Marker Queries

marker query (must have exactly field names lon, lat, title and description):

 SELECT X(p2.way) AS lon, Y(p2.way) AS lat, 'Briefkasten' AS title, p2.ref AS description
 FROM planet_osm_polygon p1
 JOIN planet_osm_point p2 ON CONTAINS(p1.way, p2.way)
 WHERE p1.name = 'Uster'
 AND p2.amenity = 'post_box'

Extension with user defined marker icon (attribute 'icon'):

 SELECT X(p2.way) AS lon, Y(p2.way) AS lat, 'Briefkasten' AS title, p2.ref AS description, 'http://myserver/marker.png' as icon 
 FROM planet_osm_polygon p1
 JOIN planet_osm_point p2 ON CONTAINS(p1.way, p2.way)
 WHERE p1.name = 'Uster'
 AND p2.amenity = 'post_box'

Mouse Queries

Finding Nearest Point to Linestrings (Snapping)

  • Lösung mit Distanzen innerhalb 50m (siehe ST_DWithin).
  • Anwendung: PostGIS-Terminal aufrufen, Copy&Paste dieser Query, dann in Karte klicken.
  • Bemerkungen:
    • (_mouse_x, _mouse_y) liefern die Mauskoordinaten (anstelle mypos).
    • Die Lösung verwendet WITH-Klausel, da (_mouse_x, _mouse_y) vom Terminal nur einmal geparst werden (leider).
    • (tags->'highway') liefert die Strassenart, z.B. unclassified, pedestrian, etc.
    • (SELECT geom FROM tmp1) liefert die Mauskoordinaten
    • COALESCE(name,'NN') liefert NN wenn kein Name vorhanden
 WITH tmp1 AS ( 
   SELECT ST_SetSRID(ST_MakePoint(_mouse_x, _mouse_y), 900913) AS geom
 ), 
 tmp2 AS ( 
   SELECT osm_id, way, name, tags, ST_Distance((SELECT min(geom)::geometry FROM tmp1), way) as distance
   FROM osm_line
   WHERE (tags->'highway') IS NOT NULL
   AND ST_DWithin((SELECT geom FROM tmp1), way, 50)
 )
 SELECT 
   ST_AsText(way) AS geom, 
   COALESCE(name,'NN')||' ('||(tags->'highway')||')' AS label
 FROM tmp2
 ORDER BY distance
 LIMIT 1

(Fast) dasselbe wie oben aber mit KNN Index (<->) nun aber mit den 10 nächsten Linien (k-Nearest Neighbor (KNN) Index)

  • geordnet nach dem Zentrum der Bounding Boxes ('<->' Operator).
  • Hinweis: der <#> ordnet nach den Bounding Boxen selber (ab PostGIS 2.0 und PostgreSQL 9.1)
 WITH tmp1 AS ( 
   SELECT ST_SetSRID(ST_MakePoint(_mouse_x, _mouse_y), 900913) AS geom
 ), 
 tmp2 AS ( 
   SELECT osm_id, way, name, tags, ST_Distance((SELECT min(geom)::geometry FROM tmp1), way) as distance
   FROM osm_line
   WHERE (tags->'highway') IS NOT NULL
   ORDER BY way <-> (SELECT min(geom)::geometry FROM tmp1)
   LIMIT 10
 )
 SELECT 
   ST_AsText(way) AS geom, 
   COALESCE(name,'NN')||' ('||(tags->'highway')||')' AS label
 FROM tmp2
 ORDER BY distance
 LIMIT 1

Statistics

All tupels in all tables

 SELECT '1. '||to_char(count(*), '999G999G999')||' osm_point(s)'
 FROM osm_point
 UNION
 SELECT '2. '||to_char(count(*), '999G999G999')||' osm_poi(s)'
 FROM osm_poi
 UNION
 SELECT '3. '||to_char(count(*), '999G999G999')||' osm_line(s)'
 FROM osm_line
 UNION
 SELECT '4. '||to_char(count(*), '999G999G999')||' osm_polygon(s)'
 FROM osm_polygon
 UNION
 SELECT '5. '||to_char(count(*), '999G999G999')||' osm_nodes'
 FROM osm_nodes
 UNION
 SELECT '6. '||to_char(count(*), '999G999G999')||' osm_ways'
 FROM osm_ways
 UNION
 SELECT '7. '||to_char(count(*), '999G999G999')||' osm_rels'
 FROM osm_rels
 ORDER BY 1

All Tag-Value-Pairs of OSM data

Kann anstelle mit osm_point auch mit osm_all durchgeführt werden.

 --
 -- Key-Value Statistics of OSM Data
 -- Return all key-value-pairs of type 'enum' without some 
 -- types numeric, date/time etc. chosen by hand:
 -- 
 -- Alternative with separate tag-value: 
 -- SELECT tmp.key as tag, tmp.value as value, count(*)::integer as freq
 SELECT tmp.key||'='||tmp.value as kvp, count(*)::integer as freq 
 FROM (
   SELECT (each(tags)).key as key, (each(tags)).value as value 
   FROM osm_point) as tmp
 WHERE (trim(value) !~ '^[-]*[0-9,.:\ ]+[m]*$')
 AND NOT (value ILIKE '%fixme%' OR key ILIKE '%fixme%')
 AND key NOT LIKE '%:%'
 AND key NOT LIKE '%description%'
 AND key NOT LIKE '%comment%'
 AND key NOT LIKE '%name'
 AND key NOT LIKE 'uic_%'
 AND key NOT LIKE '%ref'
 AND key NOT ILIKE '%fixme%'
 AND key NOT ILIKE '%todo%'
 AND key NOT IN ('osm_timestamp', 'name','operator','_picture_','_waypoint_','address','alt','is_in','url','website','wikipedia','email',
                 'converted_by','phone','information','opening_hours','date','time','collection_times','colour','fee',
                 'population','access','noexit','towards','bus_routes','busline','lines','type','denotation',
                 'CONTINUE','continue','copyright','stop')
 GROUP BY tmp.key, tmp.value
 HAVING COUNT(*) > 1
 ORDER by freq DESC; -- ca. 500 sec.!