SPLIT LINES WITH POINTS, THE SPATIALITE WAY

 

PROBLEMA

Dato uno strato di linee e uno di punti (quest’ultimo non necessariamente sopra le linee), dividere le linee usando i punti.

Workflow e note:

• Poiché i punti di input NON sono sulle linee, per ottenere il punto più vicino sulla linea (probabilmente su un segmento) occorre usare la funzione st_closetpoint;
• I vertici delle linee di input sono densificati con st_segmentize: questo permetterà in seguito di scattare i punti di input – ora sui segmenti – su un nodo/vertice, la precisione dipende dalle opzioni addensamento;
• Estrai i nodi/vertici delle linee densificate con la funzione ST_DissolvePoints;
• Punti di snap (quelli sovrapposti) sul nodo/vertice più vicino delle linee densificate usando la funzione st_snap;
• Dividi le linee con i punti usando la funzione LinesCutAtNodes:
• Il processo sopra riportato è necessario in quanto i punti possono essere usati per dividere le linee solo se si trovano esattamente sopra un nodo/vertice di linea (LinesCutAtNodes non funzionerà nei punti su un segmento di linea);
• dump geometria utilizzando la VirtualElementary è un nuovo driver per la tabella virtuale introdotto a partire dalla versione 4.2.1 di SpatiaLite.

Step1: inserire punti ‘attaccati’ ai segmenti di linea (calamita)

CREATE TABLE points_over_lines
AS SELECT a.id,ST_ClosestPoint(ST_Union(b.geom), a.geom) AS geom
FROM civici a, strade b
GROUP BY a.geom,a.id;

SELECT RecoverGeometryColumn('points_over_lines','geom',3004,'POINT','XY');

Step2: densificare i vertici delle linee ed estrarre i nodi/vertici come unica geometria multipoint

CREATE TABLE lines_nodes_densified AS
SELECT pk AS id, ST_Union(ST_DissolvePoints(st_segmentize(geom,1))) AS geom 
FROM strade;

SELECT RecoverGeometryColumn('lines_nodes_densified','geom',3004,'MULTIPOINT','XY');

Step3: snap point su linee nodi/vertici

CREATE TABLE points_snapped AS
SELECT b.id, ST_snap(ST_Union(b.geom),a.geom, ST_Distance(a.geom,b.geom)*1.01) AS geom 
FROM lines_nodes_densified a, points_over_lines b
GROUP BY a.geom, b.geom, b.id;

SELECT RecoverGeometryColumn('points_snapped','geom',3004,'POINT','XY');

Step4: dividere le linee

CREATE TABLE lines_split AS
SELECT a.id, ST_LinesCutAtNodes(st_segmentize(a.geom,1),ST_Union(b.geom)) AS geom
FROM strade a, points_snapped b
GROUP BY a.id,a.geom;

SELECT RecoverGeometryColumn('lines_split','geom',3004,'MULTILINESTRING','XY');

Step5: dump geometrie – tabella con geometria finale

create table lines_split_dump as
SELECT t.id as multi_id, e.item_no,e.geometry
FROM lines_split AS t
JOIN ElementaryGeometries AS e ON (e.f_table_name = 'lines_split' AND e.origin_rowid = t.id);

SELECT RecoverGeometryColumn('lines_split_dump','geometry',3004,'LINESTRING','XY');

  NB: il risultato della divisione è TOPOLOGICAMENTE corretto!!!  

prima: linee e punti

dopo: linee divise tramite i punti

 

Versione di Spatialite utilizzata:

Spatialite_gui

 

EDIT:  unico script SQL :

	—
	— split lines with point, the Spatialite way
	— di Salvatore Fiandaca
	— e-mail: pigrecoinfinito@gmail.com
	—
	CREATE TABLE "points_over_lines" AS
	SELECT a.id AS id,ST_ClosestPoint(ST_Union(b.geom), a.geom) AS geom
	FROM civici a, strade b
	GROUP BY a.geom,a.id;
	SELECT RecoverGeometryColumn('points_over_lines','geom',3004,'POINT','XY');
	CREATE TABLE "lines_nodes_densified" AS
	SELECT pk AS id, ST_Union(ST_DissolvePoints(st_segmentize(geom,1))) AS geom
	FROM strade;
	SELECT RecoverGeometryColumn('lines_nodes_densified','geom',3004,'MULTIPOINT','XY');
	CREATE TABLE "points_snapped" AS
	SELECT b.id AS id, ST_snap(ST_Union(b.geom),a.geom, ST_Distance(a.geom,b.geom)*1.01) AS geom
	FROM "lines_nodes_densified" a, "points_over_lines" b
	GROUP BY a.geom, b.geom, b.id;
	SELECT RecoverGeometryColumn('points_snapped','geom',3004,'POINT','XY');
	CREATE TABLE "lines_split" AS
	SELECT a.id AS id, ST_LinesCutAtNodes(st_segmentize(a.geom,1),ST_Union(b.geom)) AS geom
	FROM strade a, "points_snapped" b
	GROUP BY a.id,a.geom;
	SELECT RecoverGeometryColumn('lines_split','geom',3004,'MULTILINESTRING','XY');
	CREATE TABLE "lines_split_dump" AS
	SELECT t.id AS multi_id, e.item_no,e.geometry
	FROM "lines_split" t
	JOIN ElementaryGeometries e ON (e.f_table_name = 'lines_split' AND e.origin_rowid = t.id);
	SELECT RecoverGeometryColumn('lines_split_dump','geometry',3004,'LINESTRING','XY');
	—
	— cancello le tabelle inutili, i primi 4 step
	—
	drop table points_over_lines;
	drop table lines_nodes_densified;
	drop table points_snapped;
	drop table lines_split;
	—
	— aggiorno statistiche e VACUUM
	—
	UPDATE geometry_columns_statistics set last_verified = 0;
	SELECT UpdateLayerStatistics('geometry_table_name');
	VACUUM;

view raw

scrip_SQL_split_lines_with_points.sql

hosted with ❤ by GitHub

 

database per test – articolo

database usato nel video

spatialite_gui-4.3.0a-win-amd64

 

sitografia:

idea presa da qui – stesso argomento da con postgis

forum spatialite user – qui trovate il papà di spatialite, A. Furieri

gaia gis

spatialite

SQLite

 

video demo – usando DBmanager di QGIS 3.0:

 

———-

donate