| 166 | == Function wrapper == |
| 167 | |
| 168 | The following code could be wrapped into a function like this. (It uses unlogged and temp tables, but not sure if it really helps to get speed): |
| 169 | |
| 170 | {{{ |
| 171 | -- takes a table name, gid and geom column and simplify the layer based on given tolerance. |
| 172 | -- returns a setof (gid int, geom geometry). where geom is the simplified geometry and gid is the tablename gid. |
| 173 | create or replace function simplifyLayerPreserveTopology (tablename text, id_col text, geom_col text, tolerance float) returns setof record as $$ |
| 174 | DECLARE |
| 175 | tabname alias for $1; |
| 176 | tid alias for $2; |
| 177 | geo alias for $3; |
| 178 | tol alias for $4; |
| 179 | numpoints int:=0; |
| 180 | |
| 181 | BEGIN |
| 182 | |
| 183 | EXECUTE 'select sum(st_npoints(geom)) from '||quote_ident(tabname) into numpoints; |
| 184 | raise notice 'Num points in %: %', tabname, numpoints; |
| 185 | |
| 186 | EXECUTE 'create unlogged table poly as (' |
| 187 | ||'select '||quote_ident(tid)||', (st_dump('||quote_ident(geo)||')).* from '||quote_ident(tabname)||')'; |
| 188 | |
| 189 | -- extract rings out of polygons |
| 190 | create unlogged table rings as |
| 191 | select st_exteriorRing((st_dumpRings(geom)).geom) as g from poly; |
| 192 | |
| 193 | raise notice 'rings created...'; |
| 194 | |
| 195 | drop table poly; |
| 196 | |
| 197 | -- Simplify the rings. Here, no points further than 10km: |
| 198 | create unlogged table simplerings as select st_simplifyPreserveTopology(st_linemerge(st_union(g)), tol) as g from rings; |
| 199 | |
| 200 | raise notice 'rings simplified...'; |
| 201 | |
| 202 | drop table rings; |
| 203 | |
| 204 | -- extract lines as individual objects, in order to rebuild polygons from these |
| 205 | -- simplified lines |
| 206 | create unlogged table simplelines as select (st_dump(g)).geom as g from simplerings; |
| 207 | |
| 208 | drop table simplerings; |
| 209 | |
| 210 | -- Rebuild the polygons, first by polygonizing the lines, with a |
| 211 | -- distinct clause to eliminate overlaping segments that may prevent polygon to be created, |
| 212 | -- then dump the collection of polygons into individual parts, in order to rebuild our layer. |
| 213 | create temp table simplepolys as |
| 214 | select (st_dump(st_polygonize(distinct g))).geom as g |
| 215 | from simplelines; |
| 216 | |
| 217 | select count(*) from simplepolys into numpoints; |
| 218 | |
| 219 | raise notice 'rings polygonized. num rings: %', numpoints; |
| 220 | |
| 221 | drop table simplelines; |
| 222 | |
| 223 | -- some spatial indexes |
| 224 | create index simplepolys_geom_gist on simplepolys using gist(g); |
| 225 | |
| 226 | raise notice 'spatial index created...'; |
| 227 | |
| 228 | -- comparing percentage of overlaping area gives good results to identify simplified objects. |
| 229 | -- (the 50% of overlapping area is empirical. Should fail with some odd-shaped polygons) |
| 230 | RETURN QUERY EXECUTE 'select '||quote_ident(tid)||', st_collect(geom) as geom ' |
| 231 | ||'from (' |
| 232 | --||' select distinct on (d.'||quote_ident(tid)||') d.'||quote_ident(tid)||', s.g as geom ' |
| 233 | ||' select d.'||quote_ident(tid)||', s.g as geom ' |
| 234 | ||' from '||quote_ident(tabname)||' d, simplepolys s ' |
| 235 | ||' where st_intersects(d.'||quote_ident(geo)||', s.g) ' |
| 236 | ||' and st_area(st_intersection(s.g, d.'||quote_ident(geo)||'))/st_area(s.g) > 0.5 ' |
| 237 | ||' order by d.'||quote_ident(tid)||', st_area(st_intersection(s.g, d.'||quote_ident(geo)||'))/st_area(s.g) desc ' |
| 238 | ||' ) as foo ' |
| 239 | ||'group by '||quote_ident(tid); |
| 240 | |
| 241 | drop table simplepolys; |
| 242 | |
| 243 | RETURN; |
| 244 | |
| 245 | END; |
| 246 | $$ language plpgsql strict; |
| 247 | }}} |
| 248 | |
| 249 | === Example usage === |
| 250 | |
| 251 | {{{ |
| 252 | select * from simplifyLayerPreserveTopology('departement', 'gid', 'geom', 10000) as (gid int, geom geometry); |
| 253 | }}} |
| 254 | |
| 255 | == ToDo == |
| 256 | |
| 257 | Handle enclosed polygons cases. |