PostGIS:

Oh, the Places You’ll Index!

Jackson Voelkel

Kaiser Permanente | Portland State University

Overview

Who am I
What is GIS
What are Spatial Analytics
What is PostGIS, and where does it fit in?
PostGIS Examples

Who am I?

Jackson Voelkel

jacksonvoelkel.com

I’ve also been saying for years I’ll start a blog…

Health Data Analyst, Kaiser Permanente (National GIS team)

Spatial statistics and custom analysis
Workflow automation
Database design/compilation
R, Python, SQL

Adjunct Professor of Geography, Portland State University

GEOG 597: Advanced Spatial Quantitative Analysis in R
GEOG 575: Spatial Database Design and Compilation
GEOG 4/588: GIS 1: Introduction to GIS

Principal / Lead Analyst, Lone Fir Labs (Consulting)

Same as above, but far more specialized/specific tasks
Remote Sensing
Machine Learning
Web map application development

What is GIS?

That’s a hard question to answer

It depends on who you ask
It depends on the audience’s interests

Common Responses:

“Digital Cartography”
“Math with maps”
“Spatial data management”
“Geodatabases”

… or you can just watch a video.

Esri is the main GIS software company around today.

Here’s William Shatner talking about one of the first GIS tools!

Why is GIS Important?

Everything happens somewhere.

Geographic space is often the lowest common denomenator for compiling and analyzing disparate data.

Spatial Analysis

Four Types of Spatial Analysis

Spatial data manipulation
Spatial data analysis
Spatial statistical analysis
Spatial modeling

Spatial data manipulation

Typically done via GIS
E.g., merging, overlay, dissolve.

Spatial data analysis

First steps in analysis: exploratory and descriptive

Spatial statistical analysis

Hypothesis testing

Spatial Modeling

Forecasts and predictions
E.g., climate change modeling

Spatial Data Analysis

Spatial data analysis is concerned with that branch of data analysis where the geographical referencing of objects contains important information. In many areas of data collection, especially in some areas of experimental science, the indexes that distinguish different cases can be exchanged without any loss of information. All the information relevant to understanding the variation in the data set is contained in the observations and no relevant information is contained in the indexing. In the case of spatial data the indexing (by location and time) may contain crucial information. A definition of spatial analysis (of which spatial data analysis is one element) is that it represents a collection of techniques and models that explicitly use the spatial referencing of each data case.

Goodchild, Michael, and Haining. 2004. “GIS and spatial data analysis: converging perspectives”. Papers in Regional Science 83:363-385

Spatial data analysis is concerned with that branch of data analysis where the geographical referencing of objects contains important information. In many areas of data collection, especially in some areas of experimental science, the indexes that distinguish different cases can be exchanged without any loss of information. All the information relevant to understanding the variation in the data set is contained in the observations and no relevant information is contained in the indexing. In the case of spatial data the indexing (by location and time) may contain crucial information. A definition of spatial analysis (of which spatial data analysis is one element) is that it represents a collection of techniques and models that explicitly use the spatial referencing of each data case.

Goodchild, Michael, and Haining. 2004. “GIS and spatial data analysis: converging perspectives”. Papers in Regional Science 83:363-385

How is it different than GIS?

Different perspective on spatial data

GIS is built around the entity-attribute model. Spatial analysis uses these data, and sees data as patterns that are the outcomes of processes. This is similar… but different than spatial data manipulation.

Proprietary Tools

ArcGIS
ENVI
MicroStation

FOSS Tools

R
Python
QGIS
GRASS (still around, and amazing!)

Enterprise Spatial Analysis

Arcgis

They ‘own’ GIS
Desktop/Server software
Python API
Some DB development (though nothing like PostgreSQL)

Example 1: Enriching data

We wish to filter a specific cohort of people by distance to rivers, perform and analysis, and push the results back to our DB.

Create Oracle Query to grab cohort of interest (not filtered by distance to rivers!)
Run query and download as .csv
Open .csv in ArcGIS
Create points from XY columns
Run analysis (buffer, spatial join, etc)
Export results to csv
Upload to Oracle

Example 2: Joining data

We wish to join the newest Census Blocks (~12 million) to every point we have in a database (lets imagine ~12 million patients, trees, or aardvarks).

Select all XY from Oracle
Download a CSV
Import into ArcGIS
Join
Wait many hours.
Export results as a CSV
Upload to Oracle.

Aside: Oracle Spatial

Yes, I could do some of this in Oracle Spatial, but I’m not made of money … anyways, I think PostGIS is better!

Limits

Time: this takes forever!
Reproducability: too much hands-on work
Dependencies: Between Oracle and ArcGIS there are enough bugs to drive you mad… imagine when they update!
Versioning: Floating files of different vintages

Where Does PostgreSQL Fit In?

Everywhere!

Storage

Well, obviously!

(Spatial) Analysis

PostGIS!

PostGIS: Spatial Extention

As easy as…

$ sudo apt install postgresql-10-postgis-2.4

Or the equivalent for your inferior (just kidding) system/version

and

CREATE EXTENSION POSTGIS;

Adds in Massive Amounts of Power

Data types
Analytical functions
So big, it has more extensions for processes like volumetrics and street network analysis.

Geometry

Spatial data works in a standard entity-attribute format, with various geometry types.

Geometry

The Open Geographic Consortium sets standards, which PostGIS follows.

OGC Simple Feature Access - Part 1: Common Architecture

Geometry

Let’s look at a very simple data set: fire stations in the Portland Metro Area:

SELECT * FROM fire_sta;

 id  |                                 geom                                 | station |               address               |      city       |               district               
-----+----------------------------------------------------------------------+---------+-------------------------------------+-----------------+--------------------------------------
   1 | 01040000201E6900000100000001010000006848EEED96711F410D8911649E305341 | 50      | 12617 SW WALNUT ST                  | TIGARD          | TUALATIN VALLEY FIRE & RESCUE
   2 | 01040000201E6900000100000001010000003057F1DB66D61F41779684954A2D5341 | 211     | 4555 JEAN RD                        | LAKE OSWEGO     | LAKE OSWEGO F&R&LIFE SAFETY
   3 | 01040000201E6900000100000001010000000B36E71F0FDB1D412D5B860991315341 | 11      | 102 E MAIN ST                       | GASTON          | GASTON RFPD
   4 | 01040000201E6900000100000001010000001E8D7697A0D61F41D67BD6BA5B415341 | 26      | 5247 N LOMBARD ST                   | PORTLAND        | PORTLAND FIRE & RESCUE
   5 | 01040000201E6900000100000001010000006F3DEF04B1A51F4197EB934FA6265341 | 56      | 8455 SW ELLIGSEN RD                 | WILSONVILLE     | TUALATIN VALLEY FIRE & RESCUE
   6 | 01040000201E6900000100000001010000008FE7EC0CDF06204146A2177BF1195341 | 65      | 26815 HIGHWAY 170                   |                 | CANBY RFPD #62
   7 | 01040000201E6900000100000001010000002FC149D18B392241CBD9801B12245341 | 253     | 87600 E GOVERNMENT CAMP LOOP        | GOVERNMENT CAMP | HOODLAND RFPD #74
   8 | 01040000201E690000010000000101000000B046D46B25882041D7EC8E183C3C5341 | 74      | 1520 NE 192ND AVE                   | GRESHAM         | GRESHAM FIRE & EMERG SRVCS
   9 | 01040000201E690000010000000101000000C668D7930D901F41350FEBB7F2335341 | 53      | 8480 SW SCHOLLS FERRY RD            | BEAVERTON       | TUALATIN VALLEY FIRE & RESCUE
  10 | 01040000201E69000001000000010100000082EDAD44C9E91F417153C4A654345341 | 18      | 8720 SW 30TH AVE                    | PORTLAND        | PORTLAND FIRE & RESCUE

SELECT  id
        , station
        , address
        , city
        , st_astext(geom) geom 
FROM fire_sta;


 id | station |           address            |      city       |                     geom                      
----+---------+------------------------------+-----------------+-----------------------------------------------
  1 | 50      | 12617 SW WALNUT ST           | TIGARD          | MULTIPOINT(515173.732354289 5030521.56357027)
  2 | 211     | 4555 JEAN RD                 | LAKE OSWEGO     | MULTIPOINT(521625.714787829 5027114.33621751)
  3 | 11      | 102 E MAIN ST                | GASTON          | MULTIPOINT(489155.781155438 5031492.14882545)
  4 | 26      | 5247 N LOMBARD ST            | PORTLAND        | MULTIPOINT(521640.147913175 5047662.91934105)
  5 | 56      | 8455 SW ELLIGSEN RD          | WILSONVILLE     | MULTIPOINT(518508.254818878 5020313.24340334)
  6 | 65      | 26815 HIGHWAY 170            |                 | MULTIPOINT(525167.525244938 5007301.9233175)
  7 | 253     | 87600 E GOVERNMENT CAMP LOOP | GOVERNMENT CAMP | MULTIPOINT(597189.908765828 5017672.42973943)
  8 | 74      | 1520 NE 192ND AVE            | GRESHAM         | MULTIPOINT(541714.710603913 5042416.38372346)
  9 | 53      | 8480 SW SCHOLLS FERRY RD     | BEAVERTON       | MULTIPOINT(517123.394376409 5033930.87372189)
 10 | 18      | 8720 SW 30TH AVE             | PORTLAND        | MULTIPOINT(522866.317069732 5034322.60573278)

Notice the projection of the data gives hard-to-read coordinates!

SELECT  id
        , station
        , address
        , city
        , (ST_ASTEXT(ST_TRANSFORM(geom,4326))) geom 
FROM fire_sta;


 id | station |           address            |      city       |                      geom                      
----+---------+------------------------------+-----------------+------------------------------------------------
  1 | 50      | 12617 SW WALNUT ST           | TIGARD          | MULTIPOINT(-122.806027592858 45.428051325105)
  2 | 211     | 4555 JEAN RD                 | LAKE OSWEGO     | MULTIPOINT(-122.723699509766 45.3972124570865)
  3 | 11      | 102 E MAIN ST                | GASTON          | MULTIPOINT(-123.138647942155 45.4368683505389)
  4 | 26      | 5247 N LOMBARD ST            | PORTLAND        | MULTIPOINT(-122.722608693335 45.5821707668332)
  5 | 56      | 8455 SW ELLIGSEN RD          | WILSONVILLE     | MULTIPOINT(-122.763784232417 45.3360833873664)
  6 | 65      | 26815 HIGHWAY 170            |                 | MULTIPOINT(-122.679455050874 45.2187565962888)
  7 | 253     | 87600 E GOVERNMENT CAMP LOOP | GOVERNMENT CAMP | MULTIPOINT(-121.760251085664 45.3058277709189)
  8 | 74      | 1520 NE 192ND AVE            | GRESHAM         | MULTIPOINT(-122.46574190944 45.5340339071676)
  9 | 53      | 8480 SW SCHOLLS FERRY RD     | BEAVERTON       | MULTIPOINT(-122.780985680916 45.4586940268239)
 10 | 18      | 8720 SW 30TH AVE             | PORTLAND        | MULTIPOINT(-122.707514263889 45.462055584884)

Polygons

Here we have the Neighborhoods in the Portland Metro Area:

SELECT  nbo_hood_id
        , name
        , ST_ASTEXT(ST_TRANSFORM(geom,4326)) 
FROM nbo_hood;

 nbo_hood_id |          name           |                                                                                                                                                                                           
-------------+-------------------------+-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
           8 | LLOYD DISTRICT          | MULTIPOLYGON(((-122.650480796952 45.5333589655569,-122.650510431462 45.5314726993377,-122.654567304327 45.531497367206,-122.654583838372 45.5293616414892,-122.653311176239 45.52876729323
         211 | GRESHAM - NORTH GRESHAM | MULTIPOLYGON(((-122.460068956153 45.5591750869026,-122.460071836589 45.5581701759594,-122.460071474271 45.5581033059463,-122.459943494087 45.5581034443263,-122.459921124039 45.5581034656
          96 | PARKROSE HEIGHTS        | MULTIPOLYGON(((-122.537634551912 45.5337757045262,-122.537636878715 45.5335736333397,-122.537824367438 45.5335715422847,-122.545188064426 45.533490336749,-122.547000762669 45.53347027081
          28 | NORTHFIELD              | MULTIPOLYGON(((-122.497756840315 45.6083488130372,-122.495572604578 45.6083495710326,-122.495576108117 45.6101271754758,-122.493112204749 45.6101643722772,-122.49309425842 45.60917751321
         287 | NORTH SALMON CREEK      | MULTIPOLYGON(((-122.655676626356 45.7289193353581,-122.655477162054 45.7282245128809,-122.654758732645 45.7257234828009,-122.654380547849 45.7247388828086,-122.654560664795 45.7235452739
         154 | CASCADE HIGHLANDS       | MULTIPOLYGON(((-122.503939226142 45.6102082210104,-122.503923730676 45.6072320292005,-122.503908694318 45.6043453309102,-122.504661932669 45.6043443269167,-122.504804951482 45.6043411063
         162 | FISHERS LANDING EAST    | MULTIPOLYGON(((-122.495576108117 45.6101271754758,-122.495572604578 45.6083495710326,-122.497756840315 45.6083488130372,-122.497753173955 45.6066773378049,-122.492101178464 45.6066973944
         350 | CHILDS                  | MULTIPOLYGON(((-122.712797021563 45.3873736667375,-122.71278317843 45.3862934000188,-122.712778032241 45.3859530356456,-122.712773822213 45.3856748658333,-122.712769287943 45.38535516713
         147 | FOURTH PLAIN VILLAGE    | MULTIPOLYGON(((-122.639197560032 45.6458942572195,-122.638937368859 45.6458664162537,-122.638389870996 45.6463356847865,-122.63783388227 45.646320363385,-122.637376439233 45.646147518485
         316 | WEST HAZEL DELL         | MULTIPOLYGON(((-122.714306015297 45.6971360464785,-122.713822245802 45.6969213359745,-122.706154919092 45.6970914144172,-122.703045332396 45.6971602242795,-122.702985164086 45.6971605011

I’ve cut off the entire polygon description.

PostGIS Analysis

Ok, so we can store geographic information… let’s put it to work!

Spatial Functions

Standard GIS operates with “tools”/“widgets”:

ArcGIS Toolbox

QGIS Toolbox

In PostGIS, we have functions (lots of them):

Post GIS Functions

PostGIS has a massive amount of high-performance tools
Between PostGIS and R/Python, you can do anything
Analysis can take place where the data lives - no unnecessary transfers!

The next few examples are very basic! PostGIS is capable of highly advanced analytics, which we wont get into here.

Buffering

Here we have a cohort of schools:

SELECT * FROM schools WHERE grade = '9-12';

Buffering

Here we have a buffer for each (converting to different coordinate system):

SELECT ST_BUFFER(ST_TRANSFORM(geom,26910),2000) 
FROM schools WHERE grade = '9-12';

Buffer + Union

Here we have a single ‘unioned’ buffer:

SELECT ST_UNION(ST_BUFFER(ST_TRANSFORM(geom,26910),2000)) 
FROM schools WHERE grade = '9-12';

Spatial Join

One of the most common functions in GIS is a “spatial join”, wherein tables are joined based on geographic relation.

In most desktop software this can take hours on large data.
In PostGIS, it is blazingly fast.

Let’s join the neighborhood name onto each school in the region!

Spatial Join

As a reminder, here are the schools and neighborhoods:

Spatial Join

SELECT s.*, n.name nbo_name 
FROM schools s 
JOIN nbo_hood n 
  ON ST_Intersects(s.geom,st_transform(n.geom,4326)) 
WHERE s.grade = '9-12';

Spatial Join: Point-in-Polygon

Intersections

Let’s get all of the Census Block Groups (CBGs) within Downtown Portland. There are many options for this:

“Clip” CBGs to neighborhood
Get all CBGs with centroids in the neighborhood
Get all CBGs that touch or ‘intersect’ the neighborhood

Intersections: ‘Clipping’

SELECT  c.fips
        , c.pop10
        , ST_Intersection(n.geom, c.geom) geom
FROM  cbg c, 
      nbo_hood n 
WHERE n.name = 'DOWNTOWN'
  AND ST_Intersects(n.geom,c.geom);

Intersections: ‘Touching’

SELECT c.*
FROM  cbg c, 
      nbo_hood n 
WHERE ST_INTERSECTS(c.geom, n.geom) 
  AND n.name = 'DOWNTOWN';

Intersections: ‘Centroids’

SELECT c.*
FROM  cbg c, 
      nbo_hood n 
WHERE ST_INTERSECTS(ST_CENTROID(c.geom), n.geom) 
  AND n.name = 'DOWNTOWN';

Intersections: Density

While we’re at it, lets also use ST_AREA to create a population density map!

SELECT  fips
        , pop10/(ST_AREA(c.geom)/1000) pop_density
        , c.geom
FROM  cbg c, 
      nbo_hood n 
WHERE ST_INTERSECTS(ST_CENTROID(c.geom), n.geom) 
  AND n.name = 'DOWNTOWN';

Intersections

These aren’t the only options! There are far more options such as:

Completely within
Within a distance
Share (or don’t share) edges

Length

Trails

Let’s find the top 10 longest trails in the Portland Metro Area!

Trails

SELECT  systemname
        , sum(st_length(geom))/5280 length
        , ST_UNION(geom) -- merging them
FROM trails WHERE systemname IS NOT NULL
GROUP BY systemname ORDER BY length DESC LIMIT 10;

Distance

Firestations and Schools

What is the closest fire station to Riverdale High?

SELECT  f.*,
        ST_DISTANCE(f.geom,ST_TRANSFORM(s.geom,26910)) dist
FROM  fire_sta f, 
      (SELECT * FROM schools WHERE name = 'Riverdale High') s
ORDER BY dist
LIMIT 1;

 id | station |                   st_astext                   |       dist       
----+---------+-----------------------------------------------+------------------
 97 | 10      | MULTIPOINT(524763.806430033 5034486.74243309) | 1032.21244503901

Not bad, but we can do more!

Distance: Line

WITH s AS 
  (SELECT ST_Transform(geom,26910) geom 
  FROM schools 
  WHERE name = 'Riverdale High')
SELECT ST_Transform(ST_Makeline(st_centroid(f.geom), ST_Centroid(s.geom)),4326)
FROM s,
(SELECT  f.*, ST_DISTANCE(f.geom,s.geom) dist
FROM  fire_sta f, s
ORDER BY dist
LIMIT 1) f;

Distance: All of the Lines!

We can use a LATERAL JOIN to calculate all nearest neighbors!

Distance: All of the Lines!

WITH s AS (SELECT * FROM schools WHERE grade = '9-12')
SELECT ST_Transform(ST_Makeline(ST_Centroid(f.geom), ST_Centroid(ST_Transform(s.geom,26910))),4326)
FROM s CROSS JOIN LATERAL
  (SELECT address, geom
  FROM fire_sta ORDER BY ST_Transform(s.geom,26910) <#> geom LIMIT 1) f;

Distance: Filtering

We can use distance to features to filter our cohort based on other features:

SELECT  s.name
        , ROUND(ST_Distance(ST_Transform(s.geom,2913),t.geom)) as dist
        , s.geom
FROM schools s, trails t 
WHERE ST_Distance(ST_Transform(s.geom,2913),t.geom) < 100;

Street Networks

What??

Load in data from regional governments, vendors, or even Open Streetmaps!

However, we’ll need one more extension….

pgRouting

Like PostGIS, pgRouting is an extension that adds data types and functions for the handling of network analysis. Some cool functions are:

All Pairs Shortest Path, Floyd-Warshall Algorithm
Shortest Path Dijkstra
Driving Distance
K-Dijkstra, One to Many Shortest Path
Traveling Sales Person

Open Streetmaps (OSM)

Crowd sourced
Various qualities, though good in major American/European urban areas
Free as .osm.pbf files from GEOFabrik
Can be directly parsed into

OSM2PO

I primarily use OSM data because of this tool! It will parse raw OSM to networks with:

Speed Limits
Road Types
Trails
Stairwells
Ferry lines

Routing

Routing: A-to-B

First, lets look at the roads data, which consists of Nodes and Edges. We’re mostly interested in the for visualization and nodes for analysis:

 id | osm_id  |        osm_name         | osm_meta | osm_source_id | osm_target_id | clazz | flags | source | target |    km     | kmh |   cost    | reverse_cost |      x1      |     y1     |      x2      |     y2     |    geom    
----+---------+-------------------------+----------+---------------+---------------+-------+-------+--------+--------+-----------+-----+-----------+--------------+--------------+------------+--------------+------------+--------
  1 | 4440143 | Northeast Hodes Street  |          |      38756004 |      29936590 |    41 |     7 |      1 |     50 |  0.059913 |  40 | 0.0014978 |    0.0014978 | -122.8929526 | 45.5297334 | -122.8922158 | 45.5296458 | 0105000...
  2 | 4440143 | Northeast Hodes Street  |          |      29936590 |      27286355 |    41 |     7 |     50 |      6 | 0.0115988 |  40 |   0.00029 |      0.00029 | -122.8922158 | 45.5296458 | -122.8920669 | 45.5296452 | 0105000...
  3 | 4440143 | Northeast Hodes Street  |          |      27286355 |      29936591 |    41 |     7 |      6 |     54 | 0.0530484 |  40 | 0.0013262 |    0.0013262 | -122.8920669 | 45.5296452 | -122.8913859 | 45.5296425 | 0105000...
  4 | 4440143 | Northeast Hodes Street  |          |      29936591 |      29936595 |    41 |     7 |     54 |     58 | 0.0657458 |  40 | 0.0016436 |    0.0016436 | -122.8913859 | 45.5296425 | -122.8905419 | 45.5296391 | 0105000...
  5 | 4440143 | Northeast Hodes Street  |          |      29936595 |      27199860 |    41 |     7 |     58 |     13 | 0.0685188 |  40 |  0.001713 |     0.001713 | -122.8905419 | 45.5296391 | -122.8896623 | 45.5296356 | 0105000...

and target in the edges (above) links up with the node/intersection ID’s (below)!

 id | clazz |  osm_id  | osm_name | ref_count | restrictions | geom_vertex  
----+-------+----------+----------+-----------+--------------+---------------
  1 |     0 | 38756004 |          |         2 |              | 010100002...
 50 |     0 | 29936590 |          |         1 |              | 010100002...
  6 |     0 | 27286355 |          |         1 |              | 010100002...
 54 |     0 | 29936591 |          |         1 |              | 010100002...
 58 |     0 | 29936595 |          |         1 |              | 010100002...

Routing: A-to-B

For this first route, let’s route between these targets: 45667 and 34045

SELECT * FROM pdx_roads_node WHERE id IN (92015,251821);

Routing: A-to-B

SELECT   p.osm_name, x.cost, p.geom
FROM pdx_roads p join (
    SELECT * FROM pgr_dijkstra(
        'SELECT id, source, target, cost, reverse_cost FROM pdx_roads'
        ,92015
        ,251821
    )
) x ON p.id = x.edge;

              osm_name               |   cost    |      geom            
-------------------------------------+-----------+-------------------
 Southeast Helene Street             |  0.001856 | 0105000020E6100000...
 Southeast 58th Avenue               | 0.0027587 | 0105000020E6100000...
 Southeast Concorde Place            | 0.0011906 | 0105000020E6100000...
 Southeast Imlay Avenue              | 0.0037218 | 0105000020E6100000...
 Southeast Lois Street               | 0.0059255 | 0105000020E6100000...
 Southeast Lois Street               | 0.0007229 | 0105000020E6100000...
 Southeast Lois Street               | 0.0021492 | 0105000020E6100000...
 Southeast Lois Street               | 0.0014191 | 0105000020E6100000...
 Southeast Lois Street               | 0.0026148 | 0105000020E6100000...
 Southeast Lois Street               |  0.002863 | 0105000020E6100000...
 Southeast Lois Street               | 0.0007276 | 0105000020E6100000...
 Southeast Lois Street               | 0.0009849 | 0105000020E6100000...
 Southeast Lois Street               | 0.0024699 | 0105000020E6100000...
 Southeast Lois Street               | 0.0011073 | 0105000020E6100000...
 Southeast Lois Street               |  0.001348 | 0105000020E6100000...
 Southeast Lois Street               | 0.0025831 | 0105000020E6100000...
 Southeast Lois Street               |  0.001708 | 0105000020E6100000...
 Southeast Lois Street               | 0.0015847 | 0105000020E6100000...
 Southeast Lois Street               | 0.0074722 | 0105000020E6100000...
 Southeast Lois Street               | 0.0019735 | 0105000020E6100000...
 Southeast Lois Street               | 0.0025121 | 0105000020E6100000...
 Southwest Cornelius Pass Road       | 0.0027643 | 0105000020E6100000...
 Southwest Cornelius Pass Road       | 0.0008814 | 0105000020E6100000...
 Southwest Cornelius Pass Road       | 0.0018234 | 0105000020E6100000...
 Southeast Cornelius Pass Road       | 0.0007404 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0007444 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0028978 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0009953 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0009114 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0004179 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0004582 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0006979 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0005732 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0004901 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0027577 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0004096 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0004396 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0009592 | 0105000020E6100000...
 Southwest Baseline Road             |   0.00151 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0011425 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0028631 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0007222 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0028163 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0005794 | 0105000020E6100000...
 Southwest Baseline Road             |  0.001034 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0028915 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0010045 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0058289 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0022246 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0002601 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0012276 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0001367 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0001642 | 0105000020E6100000...
 Southwest Baseline Road             | 0.0006795 | 0105000020E6100000...

And on, and on, and on…

Routing: A-to-B

We can even get stats about the ride:

WITH route as (SELECT p.osm_name, x.cost, st_transform(p.geom,26910) geom FROM pdx_roads p join (SELECT * FROM pgr_dijkstra('SELECT id, source, target, cost, reverse_cost FROM pdx_roads' ,92015,251821)) x ON p.id = x.edge)

SELECT  SUM(ST_Length(geom)/1000) Total_km
        , SUM(cost)*60 Total_Min
FROM route;

 
     total_km     | total_min 
------------------+-----------
 32.6837767909505 | 27.968016

Routing: [A,B,C…]-to-Z

We can route many-to-many if we’d like!

The following app shows commute flows based off of the LEHD LODES data set and OSM data.

COMNET

Routing: Driving Distance

Driving Distance

What if we wanted to see everywhere we could drive within x minutes?

The “Driving Distance” tool gives us the ability to definite a starting point, and identify all points within a specified cost from it.

Driving Distance

SELECT  x.agg_cost as cost, x.node, x.edge, p.geom
FROM pdx_roads_node p join (
    SELECT * FROM pgr_drivingdistance(
        'SELECT id, source, target, cost, reverse_cost FROM pdx_roads'
        ,250376 -- start point, Downtown Portland
        ,(SELECT (SELECT 10.000) / (SELECT 60.000))) -- 10min
) x ON p.id = x.node;

Walking Distance

What if we wanted to see everywhere we could walk within x minutes?

Walking Distance

SELECT  x.agg_cost as cost, x.node, x.edge, p.geom
FROM pdx_roads_node p join (
    SELECT * FROM pgr_drivingdistance(
        'SELECT id, source, target, km/5.00 as cost, reverse_cost FROM pdx_roads --!!!
        WHERE clazz > 14' --!!! OSM has road types!
        ,250376 -- start point, Downtown Portland
        ,(SELECT (SELECT 10.000) / (SELECT 60.000))) -- 10min
) x ON p.id = x.node;

Well, we don’t know nodes…

We don’t know which node we live by… so we can use our ‘nearest neighbor’ function to find it!

Routing: Search

Using the following, it would be easy to create an application that created driving distances based on lat/lon (such as a click on a map):

SELECT  x.agg_cost as cost, x.node, x.edge, p.geom
FROM pdx_roads_node p join (
    SELECT * FROM pgr_drivingdistance(
        'SELECT id, source, target, cost, reverse_cost FROM pdx_roads'
        , (SELECT p.id 
          FROM pdx_roads_node p, 
            (SELECT ST_SetSRID(ST_MakePoint(-122.551,45.4668),4326) as geom) s
          ORDER BY ST_DISTANCE(s.geom,p.geom) 
          LIMIT 1)
        ,(SELECT (SELECT 6.500) / (SELECT 60.000)))
) x ON p.id = x.node;

Routing: Search

Routing

Things I didn’t create an example for:

Route from directly in front of a home, part way down edge
Create driving distance polygons part way down edges
Turn nodes into isochrones

Routing: Isochrones:

I’m not going over the code to do this, but it’s from Abel Vázquez

The GIS Stack

What Does a Spatial Database Do?

Esri Stack (Proprietary)

Open Source Stack

Here’s a simple example!

Indexes

Spatial Indexing

“Spatial Indexing”

Indexes and PostGIS functions

Many of PostGIS’s functions work naitively with spatial indexes, which add to it’s incredible performance when compared to other tools.

Big Take Aways

Fast!

Not just in terms of time saved by being an all-in-one-shop, but the actual performance vs. proprietary systems!

Simple!

No more moving data between applications, servers, and storage: just run analyses where the data lives!

Great in Enterprise!

Allows teams to use the same data.

This seems like such a ‘duh’ comment, but in the GIS world it’s common to have multiple versions floating around multiple computers…

Accessible!

Though I didn’t provide examples, in a large enterprise it is great to ‘drive’ PostGIS from anywhere!

In fact, this entire presentation (including calling the queries and creating the maps) was made in R!

Expansive

It is hard to give a talk on PostGIS and feel satisfied… I feel I just scratched the surface
It does everything
It connects to everything
It’s all in SQL
It’s shockingly fast
It’s free

Thank You!

Questions?

jacksonvoelkel.com

PostGIS:

Oh, the Places You’ll Index!

Jackson VoelkelKaiser Permanente | Portland State University

Overview

Who am I?

Who am I?

Health Data Analyst, Kaiser Permanente (National GIS team)

Adjunct Professor of Geography, Portland State University

Principal / Lead Analyst, Lone Fir Labs (Consulting)

What is GIS?

That’s a hard question to answer

Common Responses:

Why is GIS Important?

Everything happens somewhere.

Spatial Analysis

Four Types of Spatial Analysis

Spatial data manipulation

Spatial data analysis

Spatial statistical analysis

Spatial Modeling

Spatial Data Analysis

How is it different than GIS?

Proprietary Tools

FOSS Tools

Enterprise Spatial Analysis

Arcgis

Example 1: Enriching data

Example 2: Joining data

Aside: Oracle Spatial

Limits

Where Does PostgreSQL Fit In?

Everywhere!

Storage

(Spatial) Analysis

PostGIS: Spatial Extention

Adds in Massive Amounts of Power

Geometry

Geometry

Geometry

Polygons

PostGIS Analysis

Spatial Functions

Post GIS Functions

Buffering

Buffering

Buffering

Buffer + Union

Spatial Join

Spatial Join

Spatial Join

Spatial Join

Spatial Join: Point-in-Polygon

Intersections

Intersections

Intersections: ‘Clipping’

Intersections: ‘Touching’

Intersections: ‘Centroids’

Intersections: Density

Intersections

Length

Trails

Trails

Distance

Firestations and Schools

Distance: Line

Distance: All of the Lines!

Distance: All of the Lines!

Distance: Filtering

Street Networks

What??

pgRouting

Open Streetmaps (OSM)

OSM2PO

Routing

Routing: A-to-B

Routing: A-to-B

Routing: A-to-B

Routing: A-to-B

Routing: [A,B,C…]-to-Z

COMNET

Jackson Voelkel

Kaiser Permanente | Portland State University