Optimisation d'un très grand point dans une requête polygonale

J'ai un ensemble de données national de points d'adresse (37 millions) et un ensemble de données de polygones de contours d'inondation (2 millions) de type MultiPolygonZ, certains des polygones sont très complexes, max ST_NPoints est d'environ 200 000. J'essaie d'identifier à l'aide de PostGIS (2.18) quels points d'adresse se trouvent dans un polygone d'inondation et de les écrire dans une nouvelle table avec l'identifiant de l'adresse et les détails du risque d'inondation. J'ai essayé du point de vue de l'adresse (ST_Within) mais j'ai ensuite échangé cela à partir du point de vue de la zone inondable (ST_Contains), la raison étant qu'il y a de grandes zones sans risque d'inondation du tout. Les deux jeux de données ont été reprojetés en 4326 et les deux tableaux ont un index spatial. Ma requête ci-dessous fonctionne depuis 3 jours maintenant et ne montre aucun signe de fin dans un avenir proche!

select a.id, f.risk_factor_1, f.risk_factor_2, f.risk_factor_3
into gb.addresses_with_flood_risk
from gb.flood_risk_areas f, gb.addresses a
where ST_Contains(f.the_geom, a.the_geom);

Existe-t-il un moyen plus optimal de gérer cela? De plus, pour les requêtes à exécution longue de ce type, quelle est la meilleure façon de surveiller les progrès autres que l'examen de l'utilisation des ressources et de pg_stat_activity?

Ma requête d'origine s'est terminée OK bien que pendant 3 jours et j'ai été détourné avec d'autres travaux, donc je n'ai jamais pu consacrer du temps à essayer la solution. Cependant, je viens de revisiter cela et de travailler sur les recommandations, jusqu'ici tout va bien. J'ai utilisé ce qui suit:

1. Création d'une grille de 50 km sur le Royaume-Uni à l'aide de la solution ST_FishNet suggérée ici
2. Définissez le SRID de la grille générée sur British National Grid et construisez un index spatial sur celle-ci
3. Coupé mes données d'inondation (MultiPolygon) en utilisant ST_Intersection et ST_Intersects (seulement obtenu ici, je devais utiliser ST_Force_2D sur le geom car shape2pgsql a ajouté un index Z
4. Coupé mes données de point en utilisant la même grille
5. Index créés sur la ligne et col et index spatial sur chacune des tables

Je suis prêt à exécuter mon script maintenant, je vais parcourir les lignes et les colonnes en remplissant les résultats dans un nouveau tableau jusqu'à ce que j'aie couvert tout le pays. Mais je viens de vérifier mes données d'inondation et certains des plus grands polygones semblent avoir été perdus en translation! Voici ma requête:

SELECT g.row, g.col, f.gid, f.objectid, f.prob_4band, ST_Intersection(ST_Force_2D(f.geom), g.geom) AS geom 
INTO rofrse.tmp_flood_risk_grid 
FROM rofrse.raw_flood_risk f, rofrse.gb_grid g
WHERE (ST_Intersects(ST_Force_2D(f.geom), g.geom));

Mes données d'origine ressemblent à ceci:

Cependant, après coupure, cela ressemble à ceci:

Voici un exemple de polygone "manquant":

Pour répondre à votre dernière question en premier, consultez cet articlesur l'opportunité de pouvoir suivre l'avancement des requêtes. Le problème est difficile et serait aggravé par une requête spatiale, car le fait de savoir que 99% des adresses avaient déjà été analysées pour le confinement dans un polygone d'inondation, que vous pouviez obtenir du compteur de boucle dans l'implémentation de l'analyse de table sous-jacente, ne serait pas nécessairement aider si le dernier 1% des adresses croisent un polygone d'inondation avec le plus de points, tandis que les 99% précédents croisent une petite zone. C'est l'une des raisons pour lesquelles EXPLAIN peut parfois être inutile avec l'espace, car il donne une indication des lignes qui seront balayées, mais, pour des raisons évidentes, ne prend pas en compte la complexité des polygones (et donc une grande proportion de l'exécution) de toutes les requêtes de type intersections / intersections.

Un deuxième problème est que si vous regardez quelque chose comme

EXPLAIN 
SELECT COUNT(a.id) 
FROM sometable a, someothertable b
WHERE ST_Intersects (a.geom, b.geom)

vous verrez quelque chose comme, après avoir raté beaucoup de détails:

_st_intersects(a.geom, b.geom)
   ->  Bitmap Index Scan on ix_spatial_index_name  (cost...rows...width...))
   Index Cond: (a.geom && geom)

La condition finale, &&, signifie faire une case à cocher englobante, avant de faire une intersection plus précise des géométries réelles. C'est évidemment sensé et au cœur du fonctionnement des R-Trees. Cependant, et j'ai également travaillé sur des données sur les inondations au Royaume-Uni dans le passé, je connais donc bien la structure des données, si les (multi) polygones sont très étendus - ce problème est particulièrement aigu si une rivière coule à, disons, 45 degrés - vous obtenez d'énormes boîtes de délimitation, ce qui pourrait forcer la vérification d'un grand nombre d'intersections potentielles sur des polygones très complexes.

La seule solution que j'ai pu trouver pour le problème "ma requête fonctionne depuis 3 jours et je ne sais pas si nous sommes à 1% ou 99%" est d'utiliser une sorte de division et de conquête pour les nuls approche, par laquelle je veux dire, divisez votre zone en petits morceaux et exécutez-les séparément, soit en boucle dans plpgsql, soit explicitement dans la console. Cela a l'avantage de couper des polygones complexes en parties, ce qui signifie que les points suivants dans les vérifications de polygones fonctionnent sur des polygones plus petits et que les boîtes de délimitation des polygones sont beaucoup plus petites.

J'ai réussi à exécuter des requêtes en une journée en divisant le Royaume-Uni en blocs de 50 km sur 50 km, après avoir tué une requête qui fonctionnait depuis plus d'une semaine dans tout le Royaume-Uni. En passant, j'espère que votre requête ci-dessus est CREATE TABLE ou UPDATE et pas seulement un SELECT. Lorsque vous mettez à jour une table, les adresses, basées sur le fait d'être dans un polygone d'inondation, vous devrez de toute façon analyser la table en cours de mise à jour, les adresses, donc avoir un index spatial dessus n'est d'aucune utilité.

EDIT: Sur la base d'une image qui vaut mille mots, voici une image de certaines données sur les inondations au Royaume-Uni. Il y a un très grand multipolygone, dont la boîte englobante couvre toute cette zone, il est donc facile de voir comment, par exemple, en coupant d'abord le polygone d'inondation avec la grille rouge, le carré dans le coin sud-ouest ne serait soudainement testé que contre un minuscule sous-ensemble du polygone.