[Gfoss] clip con spatialite

Maurizio Trevisani maurizio.trevisani a gmail.com
Lun 5 Gen 2015 12:36:35 CET


Ti riporto due esempi che uso io, per poligoni o linestring.

CREATE TABLE layerpoligonale_clippato(
  lpoly_id TEXT,
  lpoly_att1 TEXT,
  lpoly_att2 TEXT,
  lpoly_att3 TEXT,
  lpoly_att4 TEXT);
select addgeometrycolumn('layerpoligonale_clippato', 'geom', 25832,
'MULTIPOLYGON');
.timer ON
insert into  layerpoligonale_clippato(
  lpoly_id,
  lpoly_att1,
  lpoly_att2,
  lpoly_att3,
  lpoly_att4, geom)
   select
  a.lpoly_id,
  a.lpoly_att1,
  a.lpoly_att2,
  a.lpoly_att3,
  a.lpoly_att4,
CastToMulti(SanitizeGeometry(ST_SnapToGrid(ST_CollectionExtract(CastToGeometryCollection(ST_Intersection(a.geom,ww.geomclip)),3),
0.0, 0.0, 0.01, 0.01)))
 from layerpoligonaleoriginale as a, layer_geometriadiclip as ww
where ((ST_Intersects(a.geom, ww.geomclip)=1)
	AND ww.ROWID IN (
		SELECT ROWID FROM SpatialIndex
		WHERE f_table_name='layer_geometriadiclip' AND search_frame=a.geom));
.timer OFF
update  layerpoligonale_clippato set geom = st_makevalid(geom) where
st_isvalid(geom) = 0;
select * from layerpoligonale_clippato where geom is null;
delete from layerpoligonale_clippato where geom is null;
.dumpshp layerpoligonale_clippato geom layerpoligonale_clippato CP1252





CREATE TABLE layerlineare_clippato(
  lline_id TEXT,
  lline_att1 TEXT,
  lline_att2 TEXT,
  lline_att3 TEXT,
  lline_att4 TEXT);
select addgeometrycolumn('layerlineare_clippato', 'geom', 25832,
'MULTILINESTRING');
.timer ON
insert into  layerlineare_clippato(
  lline_id,
  lline_att1,
  lline_att2,
  lline_att3,
  lline_att4, geom)
select
  a.lline_id,
  a.lline_att1,
  a.lline_att2,
  a.lline_att3,
  a.lline_att4,
CastToMulti(SanitizeGeometry(ST_SnapToGrid(ST_CollectionExtract(CastToGeometryCollection(ST_Intersection(a.geom,ww.geomclip)),2),
0.0, 0.0, 0.01, 0.01)))
 from layerlineare_originale as a, layer_geometriadiclip as ww
where ((ST_Intersects(a.geom, ww.geomclip)=1)
	AND ww.ROWID IN (
		SELECT ROWID FROM SpatialIndex
		WHERE f_table_name='layer_geometriadiclip' AND search_frame=a.geom));
.timer OFF
update  layerlineare_clippato set geom = st_makevalid(geom) where
st_isvalid(geom) = 0;
select * from layerlineare_clippato where geom is null;
delete from layerlineare_clippato where geom is null;
.dumpshp layerlineare_clippato geom layerlineare_clippato CP1252

Dentro, per mi esigenze, trovi anche un arrotondamento delle
coordinate al centimetro, alcune ripuliture delle geometrie (sanitize
e poi makevalid, e l'eliminazione di eventuali geometrie nulle), e
soprattutto trovi come usare gli spatial index (indispensabili su
tabelle con molte geometrie).

Considera che un campo geometria va effettivamente dichiarato ed aggiunto
select addgeometrycolumn

Nella tua query tu selezioni tutte le LINES che intersecano la
EXTENSION.GEOM, ma non effettui alcun taglio (la st_intersects
verifica che vi sia una intersezione, non la esegue: per intersecare
la funzione è st_intersection).

Prova a ragionare sui due sample che ti ho proposto, cercando tutti i
comandi usati in
http://www.gaia-gis.it/gaia-sins/spatialite-sql-4.2.0.html
e prova anche qualche esempio che trovi nel cookbook in
https://www.gaia-gis.it/fossil/libspatialite/wiki?name=misc-docs

E' importante lo studio del SQL spaziale (Spatialite o Postgis) perchè
insegna ad affrontare le questioni geografiche con la necessaria
consapevolezza e non, come di solito succede con i GIS desktop, in
maniera empirica e senza una completa comprensione di ciò che succede
dietro le quinte!

Ciao,
Maurizio

Il 05/01/15, matteo<matteo.ghetta a gmail.com> ha scritto:
> Ciao a tutti,
> ho un problema (sono alle prime armi, abbiate pietà..) con spatialite.
> Molto in breve, ho due layer un DB di tutta l'italia e un layer su cui
> vorrei clippare l'italia.
>
> Ora, la query che ho lanciato è:
>
> create table clipped_lines as
> select * from lines, estensione
> where st_intersects(lines.GEOMETRY, extension.geom)
>
> dove *lines* sono le polilinee dell'italia e estensione è il layer di clip.
> La clip viene effettivamente fatta, ma il risultato è una tabella e non
> un layer spaziale che posso aggiungere alla mappa.
>
> Qualche esperto sa darmi una mano?
>
> Grazie mille
>
> Matteo
> _______________________________________________
> Gfoss a lists.gfoss.it
> http://lists.gfoss.it/cgi-bin/mailman/listinfo/gfoss
> Questa e' una lista di discussione pubblica aperta a tutti.
> I messaggi di questa lista non hanno relazione diretta con le posizioni
> dell'Associazione GFOSS.it.
> 666+40 iscritti al 5.6.2014


Maggiori informazioni sulla lista Gfoss