Как найти пересекающиеся географические регионы между двумя таблицами рекурсивно

Я запускаю Postgres 9.6.1 и PostGIS 2.3.0 r15146 и имею две таблицы.
geographies может иметь 150 000 000 строк, paths может иметь 10 000 000 строк:

CREATE TABLE paths (id uuid NOT NULL, path path NOT NULL, PRIMARY KEY (id))
CREATE TABLE geographies (id uuid NOT NULL, geography geography NOT NULL, PRIMARY KEY (id))

Учитывая массив/набор ids для таблицы geographies, каков наилучший способ нахождения всех пересекающихся путей и геометрий?

Другими словами, если начальный geography имеет соответствующий пересекающийся path, нам также нужно найти все остальные geographies, которые пересекает этот path. Оттуда нам нужно найти все остальные paths, которые только что нашли geographies, и так далее, пока мы не найдем все возможные пересечения.

Идентификаторы начальной географии (наш вход) могут быть от 0 до 700. В среднем около 40.
Минимальные пересечения будут равны 0, максимальное значение будет около 1000. Среднее значение, вероятно, около 20, обычно менее 100.

Я создал функцию, которая делает это, но я новичок в GIS в PostGIS и Postgres вообще. Я отправил мое решение в качестве ответа на этот вопрос.

Я чувствую, что должен быть более красноречивый и быстрый способ сделать это, чем то, что я придумал.

Ответы

Ответ 1

Ваша функция может быть радикально упрощена.

Настройка

Я предлагаю вам преобразовать столбец paths.path в тип данных geography (или, по крайней мере, geometry). path является родным типом Postgres и не очень хорошо работает с функциями PostGIS и пространственными индексами. Вам нужно было бы сделать path::geometry или path::geometry::geography (в результате LINESTRING внутренне), чтобы он работал с функциями PostGIS, такими как ST_Intersects().

Мой ответ основан на этих адаптированных таблицах:

CREATE TABLE paths (
   id uuid PRIMARY KEY
 , path geography NOT NULL
);

CREATE TABLE geographies (
   id uuid PRIMARY KEY
 , geography geography NOT NULL
 , fk_id text NOT NULL
);

Все работает с типом данных geometry для обоих столбцов. geography, как правило, более точный, но и более дорогой. Что использовать? Читайте здесь FAQ по PostGIS.

Решение 1: оптимизирована функция

CREATE OR REPLACE FUNCTION public.function_name(_fk_ids text[])
  RETURNS TABLE(id uuid, type text) AS
$func$
DECLARE
   _row_ct int;
   _loop_ct int := 0;

BEGIN
   CREATE TEMP TABLE _geo ON COMMIT DROP AS  -- dropped at end of transaction
   SELECT DISTINCT ON (g.id) g.id, g.geography, _loop_ct AS loop_ct -- dupes possible?
   FROM   geographies g
   WHERE  g.fk_id = ANY(_fk_ids);

   GET DIAGNOSTICS _row_ct = ROW_COUNT;

   IF _row_ct = 0 THEN  -- no rows found, return empty result immediately
      RETURN;           -- exit function
   END IF;

   CREATE TEMP TABLE _path ON COMMIT DROP AS
   SELECT DISTINCT ON (p.id) p.id, p.path, _loop_ct AS loop_ct
   FROM   _geo  g
   JOIN   paths p ON ST_Intersects(g.geography, p.path);  -- no dupes yet

   GET DIAGNOSTICS _row_ct = ROW_COUNT;

   IF _row_ct = 0 THEN  -- no rows found, return _geo immediately
      RETURN QUERY SELECT g.id, text 'geo' FROM _geo g;
      RETURN;   
   END IF;

   ALTER TABLE _geo  ADD CONSTRAINT g_uni UNIQUE (id);  -- required for UPSERT
   ALTER TABLE _path ADD CONSTRAINT p_uni UNIQUE (id);

   LOOP
      _loop_ct := _loop_ct + 1;

      INSERT INTO _geo(id, geography, loop_ct)
      SELECT DISTINCT ON (g.id) g.id, g.geography, _loop_ct
      FROM   _paths      p
      JOIN   geographies g ON ST_Intersects(g.geography, p.path)
      WHERE  p.loop_ct = _loop_ct - 1   -- only use last round!
      ON     CONFLICT ON CONSTRAINT g_uni DO NOTHING;  -- eliminate new dupes

      EXIT WHEN NOT FOUND;

      INSERT INTO _path(id, path, loop_ct)
      SELECT DISTINCT ON (p.id) p.id, p.path, _loop_ct
      FROM   _geo  g
      JOIN   paths p ON ST_Intersects(g.geography, p.path)
      WHERE  g.loop_ct = _loop_ct - 1
      ON     CONFLICT ON CONSTRAINT p_uni DO NOTHING;

      EXIT WHEN NOT FOUND;
   END LOOP;

   RETURN QUERY
   SELECT g.id, text 'geo'  FROM _geo g
   UNION ALL
   SELECT p.id, text 'path' FROM _path p;

END
$func$  LANGUAGE plpgsql;

Вызов:

SELECT * FROM public.function_name('{foo,bar}');

Много быстрее, чем у вас.

Основные моменты

  • Вы используете запросы на весь набор, а не только последние дополнения к набору. Это становится все медленнее с каждым циклом без необходимости. Я добавил счетчик циклов (loop_ct) в избежать избыточной работы.

  • Обязательно иметь пространственные значения GiST индексы на geographies.geography и paths.path:

    CREATE INDEX geo_geo_gix ON geographies USING GIST (geography);
    CREATE INDEX paths_path_gix ON paths USING GIST (path);
    

    Так как Postgres 9.5 индексирование только сканирование будет вариантом для индексов GiST. Вы можете добавить id в качестве второго столбца индекса. Выгода зависит от многих факторов, вам придется протестировать. Тем не менее, для типа uuid не существует подходящего оператора GiST-класса. Он будет работать с bigint после установки расширения btree_gist:

  • Иметь подходящий индекс на g.fk_id. Опять же, индекс многоколонок на (fk_id, id, geography) может заплатить, если вы можете получить от него только индексные проверки. Значение по умолчанию btree index, fk_id должно быть первым столбцом индекса. Особенно, если вы часто запускаете запрос и редко обновляете таблицу, а строки таблицы намного шире, чем индекс.

  • Вы можете инициализировать переменные во время объявления. Требуется только один раз после перезаписи.

  • ON COMMIT DROP автоматически отбрасывает временные таблицы в конце транзакции. Таким образом, я удалил явные таблицы. Но вы получаете исключение, если вы дважды вызываете функцию в одной транзакции. В этой функции я бы проверял существование таблицы temp и использовал TRUNCATE в этом случае. Связанный:

  • Используйте GET DIAGNOSTICS, чтобы получить счетчик строк вместо запуска другого запроса для счетчика.

    Вам не нужно рассчитывать вообще после перезаписи. Дешево проверить FOUND.
    На самом деле вам нужно GET DIAGNOSTICS. CREATE TABLE не устанавливает FOUND (как указано в руководстве). У меня была ваша INSERT в моей оригинальной (проверенной) функции, которая устанавливает FOUND, следовательно, надзор. Исправлено.

  • Быстрее добавить индекс или ограничение PK/UNIQUE после заполнения таблицы. И не раньше, чем мы действительно нуждаемся в этом.

  • ON CONFLICT ... DO ... - это более простой и дешевый способ для UPSERT с Postgres 9.5.

    Для простой формы команды вы просто указываете столбцы или выражения индекса (например, ON CONFLICT (id) DO ...) и пусть Postgres выполняют уникальный вывод индекса для определения ограничения или индекса арбитра. Позже я оптимизировал, предоставив ограничение напрямую. Но для этого нам нужно фактическое ограничение - уникального индекса недостаточно. Исправлено. Подробности в руководстве здесь.

  • Это может помочь временным таблицам ANALYZE вручную, чтобы помочь Postgres найти лучший план запроса. (Но я не думаю, что вам это нужно в вашем случае.)

  • _geo_ct - _geographyLength > 0 - неудобный и дорогой способ сказать _geo_ct > _geographyLength. Но это полностью исчезло.

  • Не указывайте имя языка. Просто LANGUAGE plpgsql.

  • Функциональный параметр varchar[] для массива fk_id, но вы позже прокомментировали:

    Это поле bigint, которое представляет географическую область (это фактически вычисленный s2cell id на уровне 15).

    Я не знаю s2cell id на уровне 15, но в идеале вы передаете массив соответствующего типа данных, или если это не опция по умолчанию text[].

    Также, поскольку вы прокомментировали:

    Всегда есть только 13 fk_id.

    Это кажется идеальным вариантом использования для параметра функции VARIADIC. Таким образом, ваше определение функции будет:

    CREATE OR REPLACE FUNCTION public.function_name(_fk_ids VARIADIC text[]) ...

    Подробнее:

Решение 2: простой SQL с рекурсивным CTE

Трудно обернуть rCTE вокруг двух чередующихся циклов, но возможно с некоторой тонкостью SQL:

WITH RECURSIVE cte AS (
   SELECT g.id, g.geography::text, NULL::text AS path, text 'geo' AS type
   FROM   geographies g
   WHERE  g.fk_id = ANY($kf_ids)  -- your input array here

   UNION
   SELECT p.id, g.geography::text, p.path::text
        , CASE WHEN p.path IS NULL THEN 'geo' ELSE 'path' END AS type
   FROM   cte              c
   LEFT   JOIN paths       p ON c.type = 'geo'
                            AND ST_Intersects(c.geography::geography, p.path)
   LEFT   JOIN geographies g ON c.type = 'path'
                            AND ST_Intersects(g.geography, c.path::geography)
   WHERE (p.path IS NOT NULL OR g.geography IS NOT NULL)
   )
SELECT id, type FROM cte;

Это все. Вам нужны те же самые индексы, что и выше. Вы можете перевести его в функцию SQL или PL/pgSQL для повторного использования.

Дополнительные дополнительные точки

  • Приведение в text необходимо, потому что тип geography не является "хешируемым" (тот же для geometry). (Подробнее см. эту открытую проблему PostGIS.) Обходите ее, выбрав text. Строки уникальны только в силу (id, type), мы можем игнорировать столбцы geography для этого. Вернитесь к geography для соединения. Не стоит дорожать слишком много.

  • Нам нужно два LEFT JOIN, чтобы не исключать строки, потому что на каждой итерации только одна из двух таблиц может вносить больше строк.
    Последнее условие гарантирует, что мы еще не закончили:

    WHERE (p.path IS NOT NULL OR g.geography IS NOT NULL)
    

    Это работает, потому что дублирующие выводы исключаются из временного промежуточная таблица. Руководство:

    Для UNION (но не UNION ALL) отбросить повторяющиеся строки и строки, которые дублируйте любую предыдущую строку результатов. Включите все остальные строки в результат рекурсивного запроса, а также поместить их во временный промежуточная таблица.

Итак, что быстрее?

rCTE, вероятно, быстрее, чем функция для небольших наборов результатов. Временные таблицы и индексы в функции означают значительно больше накладных расходов. Однако для больших наборов результатов функция может быть быстрее. Только тестирование с вашей фактической настройкой может дать вам окончательный ответ. *
 * См. обратную связь OP в комментарии.

Ответ 2

Я подумал, что было бы неплохо разместить мое собственное решение здесь, даже если оно не оптимально.

Вот что я придумал (используя совет Стива Чамберса):

CREATE OR REPLACE FUNCTION public.function_name(
    _fk_ids character varying[])
    RETURNS TABLE(id uuid, type character varying)
    LANGUAGE 'plpgsql'
    COST 100.0
    VOLATILE
    ROWS 1000.0
AS $function$

    DECLARE
        _pathLength bigint;
        _geographyLength bigint;

        _currentPathLength bigint;
        _currentGeographyLength bigint;
    BEGIN
        DROP TABLE IF EXISTS _pathIds;
        DROP TABLE IF EXISTS _geographyIds;
        CREATE TEMPORARY TABLE _pathIds (id UUID PRIMARY KEY);
        CREATE TEMPORARY TABLE _geographyIds (id UUID PRIMARY KEY);

        -- get all geographies in the specified _fk_ids
        INSERT INTO _geographyIds
            SELECT g.id
            FROM geographies g
            WHERE g.fk_id= ANY(_fk_ids);

        _pathLength := 0;
        _geographyLength := 0;
        _currentPathLength := 0;
        _currentGeographyLength := (SELECT COUNT(_geographyIds.id) FROM _geographyIds);
        -- _pathIds := ARRAY[]::uuid[];

        WHILE (_currentPathLength - _pathLength > 0) OR (_currentGeographyLength - _geographyLength > 0) LOOP
            _pathLength := (SELECT COUNT(_pathIds.id) FROM _pathIds);
            _geographyLength := (SELECT COUNT(_geographyIds.id) FROM _geographyIds);

            -- gets all paths that have paths that intersect the geographies that aren't in the current list of path ids

            INSERT INTO _pathIds 
                SELECT DISTINCT p.id
                    FROM paths p
                    JOIN geographies g ON ST_Intersects(g.geography, p.path)
                    WHERE
                        g.id IN (SELECT _geographyIds.id FROM _geographyIds) AND
                        p.id NOT IN (SELECT _pathIds.id from _pathIds);

            -- gets all geographies that have paths that intersect the paths that aren't in the current list of geography ids
            INSERT INTO _geographyIds 
                SELECT DISTINCT g.id
                    FROM geographies g
                    JOIN paths p ON ST_Intersects(g.geography, p.path)
                    WHERE
                        p.id IN (SELECT _pathIds.id FROM _pathIds) AND
                        g.id NOT IN (SELECT _geographyIds.id FROM _geographyIds);

            _currentPathLength := (SELECT COUNT(_pathIds.id) FROM _pathIds);
            _currentGeographyLength := (SELECT COUNT(_geographyIds.id) FROM _geographyIds);
        END LOOP;

        RETURN QUERY
            SELECT _geographyIds.id, 'geography' AS type FROM _geographyIds
            UNION ALL
            SELECT _pathIds.id, 'path' AS type FROM _pathIds;
    END;

$function$;

Ответ 3

Пример графика и данных из этого script Пример графика

Он может быть чистым реляционным с агрегатной функцией. Эта реализация использует одну таблицу path и одну таблицу point. Оба являются геометриями. Точка проще создавать тестовые данные и тестировать, чем общую географию, но ее нужно легко адаптировать.

create table path (
    path_text text primary key,
    path geometry(linestring) not null
);
create table point (
   point_text text primary key,
   point geometry(point) not null
);

Тип, чтобы сохранить состояние агрегатной функции:

create type mpath_mpoint as (
    mpath geometry(multilinestring),
    mpoint geometry(multipoint)
);

Функция построения штата:

create or replace function path_point_intersect (
    _i mpath_mpoint[], _e mpath_mpoint
) returns mpath_mpoint[] as $$

    with e as (select (e).mpath, (e).mpoint from (values (_e)) e (e)),
    i as (select mpath, mpoint from unnest(_i) i (mpath, mpoint))
    select array_agg((mpath, mpoint)::mpath_mpoint)
    from (
        select
            st_multi(st_union(i.mpoint, e.mpoint)) as mpoint,
            (
                select st_collect(gd)
                from (
                    select gd from st_dump(i.mpath) a (a, gd)
                    union all
                    select gd from st_dump(e.mpath) b (a, gd)
                ) s
            ) as mpath
        from i inner join e on st_intersects(i.mpoint, e.mpoint)

        union all
        select i.mpoint, i.mpath
        from i inner join e on not st_intersects(i.mpoint, e.mpoint)

        union all
        select e.mpoint, e.mpath
        from e
        where not exists (
            select 1 from i
            where st_intersects(i.mpoint, e.mpoint)
        )
    ) s;
$$ language sql;

Агрегат:

create aggregate path_point_agg (mpath_mpoint) (
    sfunc = path_point_intersect,
    stype = mpath_mpoint[]
);

Этот запрос вернет набор строк multilinestring, multipoint, содержащих согласованные пути/точки:

select st_astext(mpath), st_astext(mpoint)
from unnest((
    select path_point_agg((st_multi(path), st_multi(mpoint))::mpath_mpoint)
    from (
        select path, st_union(point) as mpoint
        from
            path 
            inner join
            point on st_intersects(path, point)
        group by path
    ) s
)) m (mpath, mpoint)
;
                         st_astext                         |          st_astext          
-----------------------------------------------------------+-----------------------------
 MULTILINESTRING((-10 0,10 0,8 3),(0 -10,0 10),(2 1,4 -1)) | MULTIPOINT(0 0,0 5,3 0,5 0)
 MULTILINESTRING((-9 -8,4 -8),(-8 -9,-8 6))                | MULTIPOINT(-8 -8,2 -8)
 MULTILINESTRING((-7 -4,-3 4,-5 6))                        | MULTIPOINT(-6 -2)