Reducing the processing time when handling millions of 3D points to 3D lines in PostGIS

Question

I am trying to calculate the minimum distance between 16 millions of 3D points and 3D lines.

The query is working for more or less 1 million data, but after that, it took so long and never ends. All the tables have spatial index.

The query is:

SELECT
    public.points.orden,
    st_length(st_shortestline(public.points.point_geometry, public.line01.geom)) as line01_length,
    st_length(st_shortestline(public.points.point_geometry, public.line02.geom)) as line02_length,
    st_length(st_shortestline(public.points.point_geometry, public.line03.geom)) as line03_length,
    st_length(st_shortestline(public.points.point_geometry, public.line04.geom)) as line04_length
ST_3DLENGTH(ST_3DSHORTESTLINE(public.points.point_geometry, public.line01.geom)) as line01_3dlength,
ST_3DLENGTH(ST_3DSHORTESTLINE(public.points.point_geometry, public.line02.geom)) as line02_3dlength,
ST_3DLENGTH(ST_3DSHORTESTLINE(public.points.point_geometry, public.line03.geom)) as line03_3dlength,
ST_3DLENGTH(ST_3DSHORTESTLINE(public.points.point_geometry, public.line04.geom)) as line04_3dlength
into public.lengthdist
FROM public.points, public.line01, public.line02, public.line03, public.line04
LIMIT (SELECT COUNT(*) FROM public.points);

After, I select the min distance:

SELECT public.lengthdist.orden,
    min(line01_length),
    min(line02_length),
    min(line03_length),
    min(line04_length),
    min(line01_3dlength),
    min(line02_3dlength),
    min(line03_3dlength),
    min(line04_3dlength)
FROM public.lengthdist
GROUP BY public.lengthdist.orden

And also y set the endpoint for each min line with ST_EndPoint, and I have the same problem.

I tried with CREATE TABLE - INSERT INTO statement as well, same results. Also with an UPDATE. Also, I tried with just ST_SHORTESTLINE without ST_LENGTH.

How can I do to reduce the never-ending processing time?

I am working in Windows 10 x64, 16GB RAM.

DB Version: 11
OS Type: windows
DB Type: desktop
Total Memory (RAM):16 GB 
CPUs num: 4
Connections num: 20 Data Storage: ssd
max_connections = 20
shared_buffers = 512MB
effective_cache_size = 4GB
maintenance_work_mem = 1GB
checkpoint_completion_target = 0.5
wal_buffers = 16MB
default_statistics_target = 100
random_page_cost = 1.1
work_mem = 22573kB
min_wal_size = 100MB
max_wal_size = 2GB
max_worker_processes = 4
max_parallel_workers_per_gather = 2
max_parallel_workers = 4
max_parallel_maintenance_workers = 2

P.D.: I've VACUUM ANALYZE all tables too.

Answer 1

Updated methodology:

A more versatile, and a lot more performant way is to use an index driven (K)NN approach:

DROP TABLE IF EXISTS public.lengthdist;

CREATE TABLE public.lengthdist AS (
    SELECT  ln.id AS ln_id,
            cat.orden AS orden,
            ST_Distance(pts.geom, ln.geom),
            ST_3DDistance(pts.geom, ln.geom)
    FROM    (
        SELECT  1 AS id, geom
        FROM    public.line01
        UNION ALL
        SELECT  2 AS id, geom
        FROM    public.line02
        UNION ALL
        SELECT  3 AS id, geom
        FROM    public.line03
        UNION ALL
        SELECT  4 AS id, geom
        FROM    public.line04
    ) AS ln
    CROSS JOIN (
        SELECT orden
        FROM   public.points
        GROUP BY
               1
    ) AS cats
    CROSS JOIN LATERAL (
        SELECT  id, geom
        FROM    public.points
        WHERE   orden = cat.orden
        ORDER BY
                geom <-> ln.geom
        LIMIT   1
    ) AS pts
    ORDER BY
            1, 2, 3
);

This assumes

• a spatial index on all geom columns
• a BTREE index on public.points.orden

and runs a (K) Nearest Neighbor search on each set of points defined by orden via a double CROSS JOIN.

Execution time on the setup as described in comments is 2.5 seconds.

For more info on the overall concept:

• Efficient way to find nearest feature between huge postgres tables
• Get barrier edge id
• Find closest polygon from point and get its attributes?
• How to UPDATE with LATERAL Nearest-Neighbour query?

---

Old answer:

If it's the minimum distance per group of points to each line that you need, better ST_Collect the points by the desired category and get the ST_Distance/ST_3DDistance:

DROP TABLE IF EXISTS public.lengthdist;

CREATE TABLE public.lengthdist AS (
    SELECT  ln.id AS line_layer_id,
            pts.orden,
            ST_Distance(ln.geom, pts.geom) AS dist2d,
            ST_3DDistance(ln.geom, pts.geom) AS dist3d
    FROM    (
        SELECT  orden,
                ST_Collect(geom) AS geom
        FROM    public.points
        GROUP BY
                orden
    ) AS pts
    CROSS JOIN (
        SELECT  1 AS id, geom
        FROM    public.line01
        UNION ALL
        SELECT  2 AS id, geom
        FROM    public.line02
        UNION ALL
        SELECT  3 AS id, geom
        FROM    public.line03
        UNION ALL
        SELECT  4 AS id, geom
        FROM    public.line04
    ) AS ln
    ORDER BY
            ln.id, pts.orden
);

An example query with 16 million points and 4 lines finished in under 30 seconds, on a mid tech setup.

Note that

• spatial indexes are irrelevant here
• having a table for each individual line is somewhat pointless within a RDBMS; also, cross joining multiple tables (with >1 rows) has some serious side effects! I UNION ALL'ed them into a single table and added their table suffix as id
• if your geometries are in any geographical reference system, output distance units will be degrees, which are useless as a measure; you may want to CAST(geom AS GEOGRAPHY) (or ::GEOGRAPHY) both geom columns in ST_Distance, or use ST_DistanceSphere/ST_DistanceSpheroid. ST_3DDistance does not support the GEOGRAPHY type. Adding spheroidal distance calculation adds significant overhead to the overall query performance!

Answer 2

Since you are computing distance between every pair, spatial indices don't help, the computation has to be done regardless of the index.

But since you only care about the minimum distance, you don't need every single pair. If you do want to use index, you need to use it in filter. How to do it? If you can guess the approximate minimum distance, add a filter to only return the points within (ST_DWithin function) some k * guessed-min-distance, with k being something between 2 and 10, depending on how confident you are about your guess. The idea is to get as little results as possible (less is faster), but not an empty set. If you got some results - great, select shortest one. If you got no results - increase search radius, until you find points.

See also couple resources - this one from PostGIS how to use box for faster search (that is not guaranteed to be precise):

https://postgis.net/workshops/postgis-intro/knn.html

and one from BigQuery how to use script to do this search

https://medium.com/@mentin/nearest-neighbor-using-bq-scripting-373241f5b2f5.

Answer 3

Restructure the program so that the minimum distance is always known as it goes along so you don't need to waste time calculating it every time, i.e. just update a variable whenever a new minimum is reached, so you don't need to calculate it, then it should be known instantaneously at any given point in time.