split a lines table using a points table

[Andrea Borruso]

Hi,

I would like to split a lines table using a points table. The points has been snapped to the lines table, but these points are not coincident with any of the vertex of the lines and than I cannot use LinesCutAtNodes.

Is there a SQL way to create new vertices in the line tables, for every point of the points table?

Or could I use directly Split SQL function?

Thank you,

Andrea

[a.fu...@lqt.it]

> Or could I use directly Split SQL function?
>

Hi Andrea,

ST_Split() seems to be the tool of the trade you are loocking for:

https://www.gaia-gis.it/fossil/libspatialite/wiki?name=liblwgeom-4.0

bye Sandro

[andy]

Ciao Sandro,

On 8 April 2015 at 20:26, <a.fu...@lqt.it> wrote:

ST_Split() seems to be the tool of the trade you are loocking for:

I have applied this SQL to the attached file:

SELECT ST_Split(input."Geometry",blade."Geometry") AS Geometry,

ST_SplitLeft(input.Geometry, blade.Geometry), 

ST_SplitRight(input.Geometry, blade.Geometry)

FROM "ContornoPerCosta_elem" as input,"coords_mooved" as blade 

I have no splitting. I have snapped again the points by hand in QGIS.  As I written these points are not coincident with any of the vertex of the lines, but these are coincident with lines segments. 

Does ST_Split work with this kind of input layers?

Thank you very much,

Andrea

--

Andrea Borruso
website: http://blog.spaziogis.it
GEO+ geomatica in Italia http://bit.ly/GEOplus 

38° 7' 48" N, 13° 21' 9" E, EPSG:4326
--

"cercare e saper riconoscere chi e cosa,

 in mezzo all’inferno, non è inferno, 

e farlo durare, e dargli spazio"

Italo Calvino

[a.fu...@lqt.it]

On Wed, 8 Apr 2015 21:32:42 +0200, andy wrote:
> I have no splitting. I have snapped again the points by hand in
> QGIS.  As I written these points are not coincident with any of the
> vertex of the lines, but these are coincident with lines segments. 
>

Andrea,

"snapping by hand" does not necessarily means creating a real
genuine geometric snap; and floating point arithmetic is always
subject to nastry rounding/truncation approximations.

SELECT pt.id_t, ln.pk_uid
FROM coords_mooved AS pt
LEFT JOIN ContornoPerCosta_elem AS ln
ON (ST_Intersects(ln.geometry, pt.geometry) = 1);
------------------------
1I NULL
8I NULL
8F NULL

SELECT ST_Intersection(ln.geometry, pt.geometry)
FROM ContornoPerCosta_elem AS ln,
coords_mooved AS pt;
-------------------------------
NULL
NULL
NULL

conclusion: in rigorous numeric terms there is no real
intersection between the linestring and the points;
and this fully explains why ST_Split() will then fail.


solution: creating a really "robust" snap
=========================================
CREATE TABLE contorno2 AS
SELECT ln.pk_elem AS pk_elem,
ST_Snap(ln.geometry, mpt.geom, 0.000001) AS geom
FROM ContornoPerCosta_elem AS ln,
(SELECT ST_Collect(geometry) AS geom
FROM coords_mooved) AS mpt;
SELECT RecoverGeometryColumn('contorno2', 'geom',
32633, 'LINESTRING', 'XY');

by calling ST_Snap() we'll introduce a slight change:
three new vertices will be inserted into the Linestring
so to exactly match each Point.
you can easily check this small detail:

SELECT ST_NumPoints(geometry)
FROM ContornoPerCosta_elem;
---------
40

SELECT ST_NumPoints(geom)
FROM contorno2;
---------
43


SELECT ST_Split(input.geom, blade.geometry) AS Geometry,
ST_SplitLeft(input.geom, blade.Geometry),
ST_SplitRight(input.geom, blade.Geometry)
FROM contorno2 as input, coords_mooved as blade;

and now ST_Split() finally works as expected ;-)

bye Sandro

[andy]

Caro Sandro,

I start with my violin: you're always very kind, clear and professional. You are unique. Not rhetoric, is what I think

On 8 April 2015 at 22:48, <a.fu...@lqt.it> wrote:

"snapping by hand" does not necessarily means creating a real
genuine geometric snap; and floating point arithmetic is always
subject to nastry rounding/truncation approximations.

Two things.

One is a question. With ST_Snap() we have three new vertices inserted into the Linestring that exactly match each Point. Is this required to use ST_Split, or I can use it also when my points simply are snappep properly along a line segment?

Another one is a perspective, and in particular that one of a software final user. I have made my "snapping by hand" using QGIS 2.6. Now I have discovered that the snap to line does not work and it's in some way only a visual snapping. Am I wrong? 

I know that everything we made with GIS tools has a degree of approximation, but if I use a snap tool, often it will be useful to me to make some network analysis. But if there is no snap, it will be impossible to make it properly or my newtwork tools must know the snapping approximation of the tools I have used.

But before using "snapping by hand", I have used SQL:

SELECT ST_Closestpoint(ST_Collect(a.Geometry), b.Geometry) AS Geometry
FROM my_lines a INNER JOIN my_points b
ON PtDistWithin(a.Geometry, b.Geometry, 100)
GROUP BY b.PK_UID, b.Geometry;

And I have had the some not real genuine geometric snap. But probably, it's normal and this query does not produce a real snap.

And once again thank you,

[a.fu...@lqt.it]

ciao Andrea,

> Two things.
>
> One is a question. With ST_Snap() we have three new vertices
> inserted into the Linestring that exactly match each Point. Is this
> required to use ST_Split, or I can use it also when my points simply
> are snappep properly along a line segment?
>

interpolating more "snapping point(s)" it's not a strict requirement
for ST_Split(); anyway it could become a practical requirement in
many real world cases.
the real problem is that ST_Split() requires a strict intersection
between the linestring (input) and the point (blade), and due to
the finite precision of DOUBLE this could easily become a painful
nightmare.
if the linestring is exactly horizontal or vertical this is an
absolutely trivial condition, and we can expect the same for few
others "simple" slope values (e.g. 45, 30 or 60 degrees); but
assuming an arbitrary slope you can easily encounter many failures.

just as a recall and to clarify better the general terms of the
question: you can amuse yourself by testing various intersections
between two straight segments, then checking if the intersection
point do really intersects both lines. something like this:

CREATE TABLE segm1 (id INTEGER);
SELECT AddGeometryColumn('segm1', 'geom', -1, 'LINESTRING', 'XY');
INSERT INTO segm1 VALUES (1,
ST_GeomFromText('LINESTRING(0 0, 100 100)', -1));
INSERT INTO segm1 VALUES (2,
ST_GeomFromText('LINESTRING(0.0001 0, 100 100)', -1));
INSERT INTO segm1 VALUES (3,
ST_GeomFromText('LINESTRING(0 0.0000001, 100 100)', -1));
INSERT INTO segm1 VALUES (4,
ST_GeomFromText('LINESTRING(0 0, 99.99999999 100)', -1));
INSERT INTO segm1 VALUES (5,
ST_GeomFromText('LINESTRING(0 0, 100 99.99999999999)', -1));

CREATE TABLE segm2 (id INTEGER);
SELECT AddGeometryColumn('segm2', 'geom', -1, 'LINESTRING', 'XY');
INSERT INTO segm2 VALUES (1,
ST_GeomFromText('LINESTRING(100 0, 0 100)', -1));
INSERT INTO segm2 VALUES (2,
ST_GeomFromText('LINESTRING(100 0.0001, 0 100)', -1));
INSERT INTO segm2 VALUES (3,
ST_GeomFromText('LINESTRING(100.0000001 0, 0 100)', -1));
INSERT INTO segm2 VALUES (4,
ST_GeomFromText('LINESTRING(100 0, 0.00000001 100)', -1));
INSERT INTO segm2 VALUES (5,
ST_GeomFromText('LINESTRING(100 0, 0 99.99999999999)', -1));

SELECT id_a, id_b, ST_AsText(pt) AS point,
ST_Intersects(ln_a, pt) AS ok_line1,
ST_Intersects(ln_b, pt) AS ok_line2
FROM (SELECT l1.id AS id_a, l2.id AS id_b,
ST_Intersection(l1.geom, l2.geom) AS pt,
l1.geom AS ln_a, l2.geom AS ln_b
FROM segm1 AS l1, segm2 AS l2);

as you can easily check by yourself, an apparently valid intersection
point is always returned; but not always such point do really
intersects both segments.
there is no "black magic" in all this; it simply is what you can
expect from using arithmetic computations based on finite precision
DOUBLE. there is no way to circvumvent such issues.

> Another one is a perspective, and in particular that one of a
> software
> final user. I have made my "snapping by hand" using QGIS 2.6. Now I
> have discovered that the snap to line does not work and it's in some
> way only a visual snapping. Am I wrong?
> I know that everything we made with GIS tools has a degree of
> approximation,
>

you have not to specifically blame GIS tools; this is a general
issue affecting any scientific/mathematical sw based on DOUBLE
(floating point numbers).
you always have to carefully consider that after all you are
usning finite precision computations, and you have to safely
shield your code (and data) against the many hidden pitfalls
you can easily encounter.

> but if I use a snap tool, often it will be useful to me
> to make some network analysis.
>

this is a typical topology task: and producing sane and sound
topology data surely isn't a trivial task.
some kind of appropriate "preliminary preparation/validation"
is often required, and using ST_Snap() if frequently required
so to create robust and affordable snap-points not subject to
the fancy amenities of finite arithmetics.

baciamo le mani a Voscenza,
Sandro

[andy]

Ciao Sandro,

I have applied what you have explained and ST_Split() seems to work. Thank you.

After running the queries I have opened SplitRight and SplitLeft results. 

First point is that I must run another process to combine results, remove overlapping segments and obtain the 4 segments of my goal. I hope to find a solution myself.

The second point is a question: is it normal that in the ST_Split() results I do not have the A-B segment?

If yes, I must add something more to my procedure.

Thank you

[a.fu...@lqt.it]

On Thu, 9 Apr 2015 13:03:22 +0200, andy wrote:
> The second point is a question: is it normal that in the ST_Split()
> results I do not have the A-B segment?
>

Andrea,

I'm not really sure to correctly understand your question.
in the following examples I'll adopt your A-B-C names;
"start" and "end" are the extreme points of the initial
linestring.

case #1
---------------------
- you'll split the same original Linestring three times
using each Point as a blade; so you'll get:
- start-A / A-end
- start-B / B-end
- start-C / C-end
- you'll never get A-B and B-C for rather obvious reasons
(I imagine you've followed this approach).

case #2
---------------------
- you'll split the original Linestring using the first
point (A), so you'll now get start-A and A-end.
discard the initial Linestring and insert the two news
then marking A as already processed.
- next iteration: B intersect A-end; so you'll now
split A-end using B as a blade. you'll get A-B
and B-end; discard A-end and insert A-B and B-end
then marking B as already processed.
- ... and so on until you've consumed all your "blades".
- in this second case at the end of the process you'll
surely get: start-A, A-B, B-C and C-end

bye sandro

[andy]

Sandro,

On 9 April 2015 at 17:19, <a.fu...@lqt.it> wrote:

I'm not really sure to correctly understand your question.
in the following examples I'll adopt your A-B-C names;

I'm a normal, near to basic, GIS user, and  I was looking for a function that starting from 3 connected in parallel blades gives me the option to cut a cilindric cake in 4 parts in one pass :)

Thank you again for your time