3

I'm facing a problem dealing with Australian census collection shapes provided in MapInfo format by the Australian Bureau of Statistics. I am loading these into a PostGIS database using the ogr2ogr tool, which works for a majority of the shapes, but not all of them.

A simple example of the problem I am facing is a query like this (requires loading the NSW data set):

SELECT st_union(wkb_geometry) FROM cd06answ WHERE cd_code_2006 LIKE '1291%'

The result of this query is not the expected shape but NULL.

There are no null values in the table, but there are invalid geometries. For example

SELECT cd_code_2006 FROM cd06answ 
WHERE cd_code_2006 LIKE '1291%' AND NOT st_isvalid(wkb_geometry)

retrieves the values '1291301' and '1291321'. If I exclude invalid geometries the st_union succeeds.

Connecting Quantum GIS to the database allows rendering both shapes in question. They should be part of the geometric union, so I need to fix the problem somehow.

Are there better ways to load the MapInfo data into PostGIS? Or some means of fixing the data inside PostGIS? Since the database data renders ok it should be possible to save it, shouldn't it?

EDIT: based on Christophe's feedback I experimented a bit more with st_buffer and st_snaptogrid. The result of this query:

SELECT 
    cd_code_2006, 
    st_isvalid(st_buffer(wkb_geometry,0)), 
    st_isvalid(st_snaptogrid(wkb_geometry, 0.00000001)),
    st_isvalid(st_snaptogrid(wkb_geometry, 0.0000001)) 
FROM
    cd06answ 
WHERE 
    cd_code_2006 LIKE '1291%' 
AND
    NOT st_isvalid(wkb_geometry)

Is that for both affected geometries the first and the last of the three st_isvalids is true, the middle one isn't.

Unfortunately neither approach fixes the union, only a

SELECT st_union(st_buffer(wkb_geometry,0.4)) FROM cd06answ 
WHERE cd_code_2006 LIKE '1291%'

results in a geometry, but

SELECT st_union(st_buffer(wkb_geometry,0.3)) FROM cd06answ 
WHERE cd_code_2006 LIKE '1291%'

does not (I had tried the small buffer trick earlier, but didn't push it up to this level).

This seems a bit too much for a fix.

Glorfindel
  • 21,988
  • 13
  • 81
  • 109
Peter Becker
  • 8,795
  • 7
  • 41
  • 64
  • The ABS also provides those boundaries in ESRI shapefile format. Not sure whether you would have got better results using them (after all, a utility for loading shapefiles comes bundled with PostGIS). – Todd Owen Jul 31 '13 at 03:32

2 Answers2

2

Try running st_buffer (with radius 0 first and then 0.000000001 etc) or st_snaptogrid on these invalid geometries to 'repair' them (link to the docs here and here).

I've seen these errors pop up when importing from Mapinfo or other sources with ogr2ogr2 (basically caused by higher precision and/or rounding issues). I think the Postgis developers were planning to include a specific precision reducer function but if I recall correctly there isn't one in 1.4.

If that does not help please post your current postgis version and a wkt version of your polygon and projection. There are other possible causes for invalidness of polygons.

ChristopheD
  • 112,638
  • 29
  • 165
  • 179
  • Thanks, Christophe. See edit on the entry. The WKT for just one shape is thousands of characters long, I don't think I can easily post that. PostGIS is 1.5.2, the current stable release. I've encountered the problem on a 32bit Windows stack using Postgres 9.0.x as well as a 64bit Ubuntu with a Postgres 8.4.x. – Peter Becker Jan 27 '11 at 03:31
  • @Peter Becker: I'll post some more suggestions later this evening (Europe time ;-) – ChristopheD Jan 27 '11 at 06:50
  • I tried it again and I thought I did the same thing, but this time it worked. I have no idea why it didn't work the first few times I tried, it might have been a case of PEBKAC. – Peter Becker Feb 02 '11 at 01:39
  • @Peter Becker: glad it worked out for you! I see I forgot to add other suggestions as well... – ChristopheD Feb 03 '11 at 18:49
1

Hallo

Have you tried ST_IsValidReason(geometry) to get a clue about what is wrong?

/Nicklas

Nicklas Avén
  • 4,706
  • 1
  • 18
  • 15
  • Being the noob I am I didn't even know about that function. Both geometries are self-intersecting, oddly enough in exactly the same point. – Peter Becker Jan 30 '11 at 23:37
  • This leads to the following question: http://stackoverflow.com/questions/4846391/how-do-i-clean-up-self-intersecting-polygons-in-a-spatial-database – Peter Becker Jan 30 '11 at 23:44