6

I'm trying to replicate the result of the QGIS line intersections tool in PostGIS.

The result of this tool gives me 133 intersections represented as points as depicted in the map below.

enter image description here

Using the same data I've tried to replicate this in PostGIS using this code

WITH select_trips AS 
    (
    SELECT t.* FROM trips t
    WHERE t.id IN (25758,24730,23875,25886,24442,26723)
        )
    SELECT 
        a.id,
        st_intersection(a.the_geom_webmercator, b.the_geom_webmercator) as buffer
    FROM 
        select_trips a
    JOIN 
    select_trips b ON st_intersects(a.the_geom_webmercator, b.the_geom_webmercator)

My problems are two fold. 1: this only returns 36 rows, 2: when I try to add the layer using the QGIS DB manager, nothing gets added.

From what I understand st_intersection should return the shared geometry of the two lines represented as a point, so why isn't this working?

I know I can run the tool in QGIS and then import to my postgresql DB, but I'd like to be able to run the whole process using PostGIS.

UPDATE: Using this code returns some very strange results a line string that looks somewhat like the trip lines below and 1127 intersection points. (the light blue is the correct intersection points generated with QGIS line intersections tool)

SELECT      
    ST_Intersection(a.the_geom, b.the_geom),
    a.id
FROM
    trips as a,
    trips as b
WHERE
    st_intersects(a.the_geom, b.the_geom)
    and a.id in (25758,24730,23875,25886,24442,26723);

enter image description here

1

3 Answers 3

3

A start is to exclude lines from matching themselves. Add AND a.id!=b.id to your WHERE clause in your second block of code, as follows:

SELECT ST_Intersection(a.the_geom, b.the_geom), a.id FROM trips as a, trips as b WHERE st_intersects(a.the_geom, b.the_geom) and a.id in (25758,24730,23875,25886,24442,26723) AND a.id!=b.id;

How does that change things?

2
  • Thanks for the help, the a.id != to b.id only removed 6 geometries. The other odd thing is that if I create a table this query gives me lines and points, rather than just points. Totally confused as to what is going on.
    – hselbie
    Sep 2, 2016 at 16:46
  • 1
    It should only have removed 6 geometries because table a was restricted to 6 geometries. And lines not points will be generated if two of the lines in your table overlap, i.e. share a segment. Try restricting your query further to see what is happening. E.g. add AND b.id in (25758,24730,23875,25886,24442,26723)
    – amball
    Sep 2, 2016 at 17:13
2

I think that you want to do the intersection as a whole, therefore something like line merge will suit you problem better. This may not perform the best, but it should give you the results that you want.

WITH select_trips AS 
(
    SELECT t.* FROM trips t
    WHERE t.id IN (25758,24730,23875,25886,24442,26723)
)
SELECT int_geom
FROM (
    SELECT (
        ST_Dumppoints(                           -- 3) dump out all the points
            ST_LineMerge(                        -- 2) Merge the lines, causes lines to be noded
                ST_Union(a.the_geom_webmercator) -- 1) Union all the features into one
                    )
                )
            )
        ).geom AS int_geom -- 4) Keep the geometry part
    FROM select_trips a
    ) b
GROUP BY int_geom   -- 5) Group on the geometry
HAVING COUNT(*) > 1 -- 6) Keep where there is more than one instance
1

Though by definition ST_Intersects states

Returns TRUE if the Geometries/Geography "spatially intersect in 2D"

I don't see any spatially intersecting points in your demo which is leading to more number of intersecting points.

You can use ST_Touches which states

Returns TRUE if the geometries have at least one point in common, but their interiors do not intersect.

Which exactly match your requirement.

SELECT
    ST_Intersection(a.the_geom, b.the_geom),
    a.id
FROM
    trips as a,
    trips as b
WHERE
    st_touches(a.the_geom, b.the_geom)
    and a.id in (25758,24730,23875,25886,24442,26723);

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.