PostGIS-Beispiele
Aus Geoinformation HSR
Version vom 4. Dezember 2007, 10:15 Uhr von Stefan (Diskussion | Beiträge) (Die Seite wurde neu angelegt: Für weitere PostGIS-Beispiele siehe auch GISpunkt-Seminar PostGIS. == British Columbia-Daten == Quelle: Paul Ramsey, Refractions Research, 11. Jul. 2007, [ht...)
Für weitere PostGIS-Beispiele siehe auch GISpunkt-Seminar PostGIS.
British Columbia-Daten
Quelle: Paul Ramsey, Refractions Research, 11. Jul. 2007, [1]
SQL Commands for the Demonstration
Simple Spatial SQL
create table points ( pt geometry, name varchar ); insert into points values ( 'POINT(0 0)', 'Origin' ); insert into points values ( 'POINT(5 0)', 'X Axis' ); insert into points values ( 'POINT(0 5)', 'Y Axis' ); select name, AsText(pt), ST_Distance(pt, 'POINT(5 5)') from points; drop table points;
Loading Shape Files
cd \postgis-workshop\data pg_setenv shp2pgsql -i -s 3005 bc_pubs.shp bc_pubs > bc_pubs.sql notepad bc_pubs.sql pg_shpsql psql -U postgres -f bc_data.sql -d postgis
Creating Spatial Indexes
create index bc_roads_gidx on bc_roads using gist ( the_geom ); create index bc_pubs_gidx on bc_pubs using gist ( the_geom ); create index bc_voting_areas_gidx on bc_voting_areas using gist ( the_geom ); create index bc_municipality_gidx on bc_municipality using gist ( the_geom ); create index bc_hospitals_gidx on bc_hospitals using gist ( the_geom );
Using Spatial Indexes
-- Unindexed query select gid, name from bc_roads where _ST_Crosses(the_geom, ST_GeomFromText('LINESTRING(1220446 477473,1220417 477559)', 3005)); -- Indexed query select gid, name from bc_roads where ST_Crosses(the_geom, ST_GeomFromText('LINESTRING(1220446 477473,1220417 477559)', 3005));
Indexes and Query Plans
-- Unindexed query plan explain select gid, name from bc_roads where _ST_Crosses(the_geom, ST_GeomFromText('LINESTRING(1220446 477473,1220417 477559)', 3005)); -- Indexed query plan explain select gid, name from bc_roads where ST_Crosses(the_geom, ST_GeomFromText('LINESTRING(1220446 477473,1220417 477559)', 3005));
PostgreSQL Optimization
open postgresql.conf
Spatial Analysis in SQL
-- What is the total length of all roads in the province, in kilometers? select sum(ST_Length(the_geom))/1000 as km_roads from bc_roads;
-- How large is the city of Prince George, in hectares? select ST_Area(the_geom)/10000 as hectares from bc_municipality where name = 'PRINCE GEORGE';
-- What is the largest municipality in the province, by area? select name, ST_Area(the_geom)/10000 as hectares from bc_municipality order by hectares desc limit 1;
Data Integrity
-- How many invalid voting area polygons? select count(*) from bc_voting_areas where not ST_IsValid(the_geom);
Distance Queries
-- Where is the Tabor Arms pub? select ST_AsText(the_geom) from bc_pubs where name ilike 'Tabor Arms%';
-- How many Unity Party supporters live within 2km of Tabor Arms pub? select sum(upbc) as unity_voters from bc_voting_areas where ST_DWithin(the_geom, ST_GeomFromText('POINT(1209385 996204)', 3005), 2000);
Spatial Joins
-- Find all pubs located within 250 meters of a hospital. select h.name, p.name from bc_hospitals h, bc_pubs p where ST_DWithin(h.the_geom, p.the_geom, 250);
-- Summarize the 2000 provincial election results by municipality. select m.name, sum(v.ndp) as ndp, sum(v.lib) as liberal, sum(v.gp) as green, sum(v.upbc) as unity, sum(v.vtotal) as total from bc_voting_areas v, bc_municipality m where ST_Intersects(v.the_geom, m.the_geom) group by m.name order by m.name;
Overlays
-- Clip the voting areas by the Prince George boundary. create table pg_voting_areas as select ST_Intersection(v.the_geom, m.the_geom) as intersection_geom, ST_Area(v.the_geom) as va_area, v.*, m.name from bc_voting_areas v, bc_municipality m where ST_Intersects(v.the_geom, m.the_geom) and m.name = 'PRINCE GEORGE';
-- What is the area of the clipped features? select sum(ST_Area(intersection_geom)) from pg_voting_areas;
-- How does that compare to the clipping polygon? (Should be the same.) select ST_Area(the_geom) from bc_municipality where name = 'PRINCE GEORGE';
Coordinate Projection
-- What is the proj4 definition of srid 3005? select proj4text from spatial_ref_sys where srid=3005;
-- What are the native coordinate of the first road in the roads table? select ST_AsText(the_geom) from bc_roads limit 1;
-- What do those coordinates look like in lon/lat? select ST_AsText(ST_Transform(the_geom,4326)) from bc_roads limit 1;
Exercises
-- What is the perimeter of the municipality of Vancouver? select ST_Perimeter(the_geom) from bc_municipality where name = 'VANCOUVER'; -- What is the total area of all voting areas in hectares? select sum(ST_Area(the_geom))/10000 as hectares from bc_voting_areas; -- What is the total area of all voting areas with more than 100 voters in them? select sum(ST_Area(the_geom))/10000 as hectares from bc_voting_areas where vtotal > 100; -- What is the length in kilometers of all roads named Douglas St? select sum(ST_Length(the_geom))/1000 as kilometers from bc_roads where name = 'Douglas St';
Advanced Exercises
-- What is the length in kilometers of Douglas St in Victoria? select sum(ST_Length(r.the_geom))/1000 as kilometers from bc_roads r, bc_municipality m where ST_Contains(m.the_geom, r.the_geom) and r.name = 'Douglas St' and m.name = 'VICTORIA'; -- What two pubs have the most Green Party supporters within 500 meters of them? select p.name, p.city, sum(v.gp) as greens from bc_pubs p, bc_voting_areas v where ST_DWithin(p.the_geom, v.the_geom, 500) group by p.name, p.city order by greens desc limit 2;
-- What is the latitude of the most southerly hospital in the province? select ST_Y(ST_Transform(the_geom,4326)) as latitude from bc_hospitals order by latitude asc limit 1;
-- What were the percentage NDP and Liberal vote within the city limits of Prince George in the 2000 provincial election? select 100*sum(v.ndp)/sum(v.vtotal) as ndp, 100*sum(v.lib)/sum(v.vtotal) as liberal from bc_voting_areas v, bc_municipality m where ST_Contains(m.the_geom, v.the_geom) and m.name = 'PRINCE GEORGE';
-- What is the largest voting area polygon that has a hole? select gid, id, ST_Area(the_geom) as area from bc_voting_areas where ST_NRings(the_geom) > 1 order by area desc limit 1;
-- How many NDP voters live within 50 meters of Simcoe St in Victoria? select sum(v.ndp) as ndp from bc_voting_areas v, bc_municipality m, bc_roads r where ST_Contains(m.the_geom, r.the_geom) and r.name = 'Simcoe St' and m.name = 'VICTORIA' and ST_DWithin(r.the_geom, v.the_geom, 50);
Weitere
...