8

I am usually an ArcGIS desktop user, but keen to start using PostGIS more and I have a really big bit of processing to do. Not sure what functions to use, hopefully someone can help.

I have a polygon dataset (several million features), based on a type of landuse/ landcover classfication (20 categories). I have a number of regions in another dataset.

For each of the regions, I would like to know the area of each landcover classfication.

In ArcGIS (if it was a smaller dataset) I would imagine first adding the region to each of the polygons in the attribute table using a join. Then using "summarize" on the table by region and by landcover classification.

Not sure where to start doing this in PostGIS / SQL.

Update:

Wow thanks that has been a huge help.

It has been running a long time (44 hours!) and now I get:

NOTICE:  TopologyException: found non-noded intersection between LINESTRING (coords
edited) and LINESTRING (coords edited) at (position edited)
ERROR:  GEOS Intersection() threw an error!

********** Error **********

ERROR: GEOS Intersection() threw an error!
SQL state: XX000

I assume this is a problem in the original data - just a case of reviewing the original data or can I first check the topology some how for the whole data? Is there something about accepting certain errors / processing tolerances?

Improve this question
3
  • 1
    Quoting PostGIS docs (postgis.refractions.net/docs/ch04.html): "If you encounter any strange error messages such as "GEOS Intersection() threw an error!" or "JTS Intersection() threw an error!" when calling PostGIS functions with valid input geometries, you likely found an error in either PostGIS or one of the libraries it uses, and you should contact the PostGIS developers.
    – underdark
    Commented Oct 1, 2011 at 21:19
  • ok, is there a way to evaluate the original data, should I post a new question?
    – newpostgisuser
    Commented Oct 1, 2011 at 21:21
  • Yes, you could open a new question. Have you checked this related question: gis.stackexchange.com/questions/12451/…
    – underdark
    Commented Oct 1, 2011 at 21:34

1 Answer 1

Reset to default
9

Assuming you have the following table layout

landcover(id,type,the_geom)
region(id,name,the_geom)

Area values per landcover type per region can be calculated using

SELECT r.id, r.name, r.the_geom, l.type, 
       Sum(ST_Area(ST_Intersection(l.the_geom,r.the_geom)))
FROM landcover l,
     region r
WHERE l.the_geom && r.the_geom
GROUP by r.id, r.name, r.the_geom, l.type

ST_Intersection is used to account for landcover polygons that are only partially within a region.

Improve this answer
1
  • 2
    You may want to add "WHERE l.the_geom && r.the_geom" to the query to make use of spatial indexes.
    – diciu
    Commented Sep 29, 2011 at 9:14

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.