Spatial Join in PostGIS according to Line Attribute

Question

I am in the process of building a data warehouse in PostGIS, that will store the length of the road network of a country (Lines) per Administrative Area (i.e. Province or State, referred to as AA hereafter and is a Polygon). This is my first time of building a large scale Relational data warehouse so I am not really familiar with this kind of operations.

The "id" field of the spatial join resulting table is used to relate to the unique "id" field of the AA polygon layer. The length of the road segments is already calculated and is derived from the "len" attribute of the Road Network table.
So far the operation yields the desired results by performed the following query:

SELECT  
  m.id,  
  m.order08,  
  sum(r.len) as length,  
FROM
  public.and_lines_2012_03 AS r, 
  public.and_a8_2010_12 AS m
WHERE
  ST_Intersects(m.geom,r.geom) 
GROUP BY 
  m.id, m.order08
ORDER BY 
  id;

And the resulting table looks something like this:

gid |   id   | length  
 1  | 000002 | 74118
 2  | 000012 | 128131
 3  | 000013 | 174296
 4  | 000014 | 82240
 ...|  ...   |  ...

The Road Network is divided in 3 different line types that are attributed in the Road Network table by the name "line_type" and can have values 0,1,2 My intention is to create a table that will have a column for the length of each line type for each Administrative Area. That is the table would ideally look something like this:

gid |   id   | length_type_00 | length_type_01 | length_type_02
 1  | 000002 |     30028      |     11564      |     32526
 2  | 000012 |     64270      |     15602      |     48259
 3  | 000013 |     104492     |     51981      |     17823
 4  | 000014 |     38153      |     9524       |     34563
 ...|   ...  |      ...       |      ...       |      ...

I have experimented with different approaches in SQL statements, but so far have not managed to grasp the approach of this kind of queries. The most relative resource I have found on the Internet is this one from OpenGeo:

Update 2012/05/23

@Nicklas

I have managed to produce what seems to be very accurate results with the way you have suggested. The query used is the following:

SELECT  
  b.id AS id,
  b.order08 AS order08,
  round(sum(b.line_type_00)) AS line_type_00,
  round(sum(b.line_type_01)) AS line_type_01,
  round(sum(b.line_type_02)) AS line_type_02
FROM
(SELECT  
  m.id,
  m.order08,
  (line_type=0)::int * sum(ST_Length(ST_Intersection(r.geom::geography, m.geom::geography))) line_type_00,  
  (line_type=1)::int * sum(ST_Length(ST_Intersection(r.geom::geography, m.geom::geography))) line_type_01,  
  (line_type=2)::int * sum(ST_Length(ST_Intersection(r.geom::geography, m.geom::geography))) line_type_02
FROM
  public.and_lines_2012_03 AS r, 
  public.and_a8_2010_12 AS m
WHERE
  ST_Intersects(m.geom,r.geom)
GROUP BY 
  m.id, m.order08,line_type
ORDER BY
  m.id
) 
AS
  b
GROUP BY
  b.id,b.order08
;  

The query is quite slow though (4846 tuples take around 41 seconds), since from what I understand it re-projects on the fly to produce the the results in meters. Have I understood/am I doing something wrong?

Answer 1

What you want is a pivot table. You can create that with the pivot functionality in PostgreSQL, but I think it is easier to create it like this in your case:

SELECT  
  m.id,  
  (length_type='00')::int * sum(ST_Length(r.geom)) length_type_00,  
  (length_type='01')::int * sum(ST_Length(r.geom)) length_type_01,  
  (length_type='02')::int * sum(ST_Length(r.geom)) length_type_02
FROM
  public.and_lines_2012_03 AS r, 
  public.and_a8_2010_12 AS m
WHERE
  ST_Intersects(m.geom,r.geom) 
GROUP BY 
  m.id, m.order08
;

But be aware here that you will get to many km of roads. The roads intersecting two polygons will be represented in both with their full length added to both polygons. Depending on your data you have two options to solve this.

I your roads stops at the border to the polygon, it should, in it's full length be added to one of the polygons. Then you can use ST_Line_Interpolate_Point(geom, 0.5); to get the middle of the road. Then use that point in your ST_Intersection. Then the road will only be joined to the polygon with the roads center point.

The second option is if your roads crosses the border. Then, to get the right result you will have to use ST_Intersection to get the portion of the road that is in that particular polygon. That will be slower if you have a lot of roads and or the roads or polygons have many vertex-points.

The query could then look:

SELECT  
  m.id,  
  (length_type='00')::int * sum(ST_Length(ST_Intersection(r.geom, m.geom))) length_type_00,  
  (length_type='01')::int * sum(ST_Length(ST_Intersection(r.geom, m.geom))) length_type_01,  
  (length_type='02')::int * sum(ST_Length(ST_Intersection(r.geom, m.geom))) length_type_02 
FROM
  public.and_lines_2012_03 AS r, 
  public.and_a8_2010_12 AS m
WHERE
  ST_Intersects(m.geom,r.geom) 
GROUP BY 
  m.id, m.order08
;

HTH

Nicklas