Erfassen der ArcGIS-ähnlichen Geschwindigkeit in Postgis

62

Ich benutze Postgis 2.0 jetzt seit einem Vierteljahr und obwohl ich es sehr gerne benutze, hat eine übermäßige Abfrageverarbeitungszeit es für meinen Anwendungsfall im Grunde unbrauchbar gemacht.

Ich neige dazu, starke Geoverarbeitung für kommunale Datensätze durchzuführen, die häufig Hunderttausende von Multipolygonen enthalten. Diese Multipolygone sind manchmal sehr unregelmäßig geformt und können von 4 bis 78.000 Punkten pro Multipolygon variieren.

Wenn ich beispielsweise einen Paketdatensatz mit 329.152 Multipolygonen mit einem Gerichtsbarkeitsdatensatz schneide, der 525 Multipolygone enthält, erhalte ich die folgenden Statistiken für die Gesamtverbrauchszeit:

ArcGIS 10.0 (on same host with windows 7 OS): 3 minutes
Postgis:56 minutes (not including geometry pre-processing queries)

Mit anderen Worten, es dauert in Postgis 1500% länger als in ArcGIS - und dies ist eine meiner einfacheren Abfragen!

Einer der Gründe, warum ArcGIS angeblich schneller ausgeführt wird, sind bessere Indizes. Einige Programmierer haben kürzlich herausgefunden, wie diese Indizes funktionieren, und ich frage mich, ob jemand weiß, wie man diese Indizes in Postgis erstellt (oder Tabellen erstellt, die die Indizes imitieren). Vielleicht würde dies die meisten Geschwindigkeitsprobleme in Postgis lösen. Ich kann nur hoffen, dass es eine Möglichkeit gibt, zumal ArcGIS nur 4 GB RAM verwenden kann, während ich für meinen Postgis-Server das bis zu 4-fache verwenden konnte!

Natürlich gibt es viele Gründe, warum Postgis langsam laufen kann. Deshalb werde ich eine detaillierte Version meiner Systemspezifikationen bereitstellen:

Machine: Dell XPS 8300 
Processor: i7-2600 CPU @ 3.40 GHz 3.40 GHz 
Memory: Total Memory 16.0 GB (10.0 GB on virtual machine)

Platform: Ubuntu Server 12.04 Virtual Box VM

Potgres Version: 9.1.4
Postgis Version: POSTGIS="2.0.1 r9979" GEOS="3.3.5-CAPI-1.7.5" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.9.1, released 2012/05/15" LIBXML="2.7.8" LIBJSON="UNKNOWN" TOPOLOGY RASTER

Ich erläutere auch den gesamten Installationsprozess, den ich zum Einrichten von Postgis verwendet habe, einschließlich der Erstellung der VM selbst .

Ich habe außerdem den gemeinsamen Speicher von 24 MB auf 6 GB in der conf-Datei erhöht und die folgenden Befehle ausgeführt, um die Ausführung von postgres zu ermöglichen:

sudo sysctl -w kernel.shmmax=7516192768 (I know this setting is deleted every time you restart the OS)
sudo /etc/init.d/postgresql restart

Soweit ich das beurteilen kann, macht sich an der Leistung absolut nichts bemerkbar.

Hier sind Links zu den Daten, die ich für diesen Test verwendet habe:

  1. Pakete: tcad_parcels_06142012.shp.zip aus City of Austin, TX
  2. Gerichtsstände: Gerichtsstandsgrenzen von City of Austin, TX

Hier sind die Schritte, die ich unternommen habe, um die Daten zu verarbeiten:

ArcGIS

  1. Hinzufügen von Datensätzen zu ArcMap
  2. Stellen Sie das Koordinatensystem auf die zentralen Texas-Füße ein (srid 2277)
  3. Verwenden Sie das Schnittwerkzeug aus dem Dropdown-Menü

Postgis

Pakete importieren mit:

shp2pgsql -c -s 2277 -D -i -I -W UTF-8 "tcad_parcels_06142012.shp" "public"."tcad_parcels_06142012" |psql -d postgis_testing -U postgres -h local_ip -p 5432

Gerichtsbarkeiten importieren mit:

shp2pgsql -c -s 2277 -D -i -I -W UTF-8 "jurisdictions.shp" "public"."jurisdictions" |psql -d postgis_testing -U postgres -h local_ip -p 5432

Ungültige Geometrie in Paketen bereinigen:

DROP TABLE IF EXISTS valid_parcels;
CREATE TABLE valid_parcels(
  gid serial PRIMARY KEY,
  orig_gid integer,
  geom geometry(multipolygon,2277)
);
CREATE INDEX ON valid_parcels USING gist (geom);
INSERT INTO valid_parcels(orig_gid,geom)
  SELECT 
    gid 
    orig_gid,
    st_multi(st_makevalid(geom)) 
  FROM 
    tcad_parcels_06142012;
CLUSTER valid_parcels USING valid_parcels_geom_idx;

Ungültige Geometrie in Gerichtsbarkeiten bereinigen:

DROP TABLE IF EXISTS valid_jurisdictions;
CREATE TABLE valid_jurisdictions(
  gid serial PRIMARY KEY,
  orig_gid integer,
  geom geometry(multipolygon,2277)
);
CREATE INDEX ON valid_jurisdictions USING gist (geom);
INSERT INTO valid_jurisdictions(orig_gid,geom)
  SELECT 
    gid 
    orig_gid,
    st_multi(st_makevalid(geom)) 
  FROM 
    jurisdictions;
CLUSTER valid_jurisdictions USING valid_jurisdictions_geom_idx;

Cluster ausführen:

cluster;

Führen Sie eine Vakuumanalyse durch:

vacuum analyze;

Führen Sie eine Überschneidung für bereinigte Tabellen durch:

CREATE TABLE parcel_jurisdictions(
  gid serial primary key,
  parcel_gid integer,
  jurisdiction_gid integer,
  isect_geom geometry(multipolygon,2277)
);
CREATE INDEX ON parcel_jurisdictions using gist (isect_geom);

INSERT INTO parcel_jurisdictions(parcel_gid,jurisdiction_gid,isect_geom)
  SELECT
    a.orig_gid parcel_gid,
    b.orig_gid jurisdiction_gid,
    st_multi(st_intersection(a.geom,b.geom))
  FROM
    valid_parcels a, valid_jurisdictions b
  WHERE
    st_intersects(a.geom,b.geom);

Erläutern Sie die Abfrage zum Analysieren von Kreuzungen:

Total runtime: 3446860.731 ms
        Index Cond: (geom && b.geom)
  ->  Index Scan using valid_parcels_geom_idx on valid_parcels a  (cost=0.00..11.66 rows=2 width=1592) (actual time=0.030..4.596 rows=1366 loops=525)
  ->  Seq Scan on valid_jurisdictions b  (cost=0.00..113.25 rows=525 width=22621) (actual time=0.009..0.755 rows=525 loops=1)
Nested Loop  (cost=0.00..61428.74 rows=217501 width=24213) (actual time=2.625..3445946.889 rows=329152 loops=1)
  Join Filter: _st_intersects(a.geom, b.geom)

Nach allem, was ich gelesen habe, ist meine Kreuzungsabfrage effizient und ich habe absolut keine Ahnung, was ich falsch mache, damit die Abfrage bei sauberer Geometrie 56 Minuten dauert!

THX1138
quelle
2
In PostGIS ist es eine gebräuchliche Redewendung, eine Begrenzungsrahmen-Schnittpunktprüfung hinzuzufügen, um die Dinge zu beschleunigen. Fügen Sie Ihrer WHERE-Klausel 'AND a.geom && b.geom' hinzu und sehen Sie, welchen Unterschied dies macht.
Sean
2
st_intersects () enthält eine Bounding-Box-Abfrage, bevor Kreuzungstests in postgis 2.x durchgeführt werden. Das spart also leider keine Zeit.
THX1138
1
Können Sie Ihre Abfrage mit EXPLAIN ANALYZE ausführen und die Ergebnisse veröffentlichen
Nathan W
1
Sie sollten sich auch darüber im Klaren sein, dass Sie unterschiedliche Datensätze für Postgis und Arcgis ausführen, da Sie sagen, dass Sie die Gültigkeit dieser Datensätze für Postgis erhöhen möchten.
Nicklas Avén
2
Ist es möglich, die Datensätze zu erhalten, um einen Blick darauf zu werfen?
Nicklas Avén

Antworten:

87

Anderer Ansatz. Wenn Sie wissen, dass sich der Schmerz in ST_Intersection befindet und dass True / False-Tests schnell sind, kann der Versuch, die Menge an Geometrie zu minimieren, die durch die Kreuzung verläuft, die Dinge beschleunigen. Beispielsweise müssen Parzellen, die vollständig in einer Jurisdiktion enthalten sind, nicht abgeschnitten werden, aber ST_Intersection wird wahrscheinlich immer noch die Mühe machen, einen Teil der Schnittmengenüberlagerung zu erstellen, bevor erkannt wird, dass keine neue Geometrie generiert werden muss. Also das

INSERT INTO parcel_jurisdictions(parcel_gid,jurisdiction_gid,isect_geom)
SELECT
  a.orig_gid AS parcel_gid,
  b.orig_gid AS jurisdiction_gid,

  st_multi(st_intersection(a.geom,b.geom)) AS geom
FROM
  valid_parcels a, valid_jurisdictions b
WHERE
  st_intersects(a.geom, b.geom) and not st_within(a.geom, b.geom)
UNION ALL
SELECT
  a.orig_gid AS parcel_gid,
  b.orig_gid AS jurisdiction_gid,
  a.geom AS geom
FROM
  valid_parcels a, valid_jurisdictions b
WHERE
  st_within(a.geom, b.geom);

Oder noch knapper

INSERT INTO parcel_jurisdictions(parcel_gid,jurisdiction_gid,isect_geom)
SELECT
  a.orig_gid AS parcel_gid,
  b.orig_gid AS jurisdiction_gid,
  CASE 
     WHEN ST_Within(a.geom,b.geom) 
     THEN a.geom
     ELSE ST_Multi(ST_Intersection(a.geom,b.geom)) 
  END AS geom
FROM valid_parcels a
JOIN valid_jurisdictions b
ON ST_Intersects(a.geom, b.geom)

Könnte ohne die UNION sogar noch schneller sein.

Paul Ramsey
quelle
13
Danke, das bringt mich auf 3,63 Minuten! Ich hätte nie gedacht, dass eine Gewerkschaft schneller sein würde. Diese Antwort wird mich von nun an dazu bringen, meine Abfragen zu überdenken.
THX1138,
2
Das ist sehr cool. Ich hatte einen Fall bei der Arbeit, in dem meine st_intersection-Abfrage 30 Minuten + dauerte und jetzt weiß ich, wie ich das vermeiden kann :)
Nathan W
1
Diese Frage brachte mich dazu, Postgis zu lernen! Ich werde heute gut schlafen und sehe, wie Postgis Schulter an Schulter mit Arcgis rennt :-)
vinayan
2
Eine weitere Verbesserung von Martin Davis, könnten Sie das "rein oder raus?" Frage in die SELECT mit einer CASE-Anweisung und vermeiden Sie die UNION auf diese Weise.
Paul Ramsey
2
UNIONEntfernt doppelte Zeilen aus den beiden Abfragen, diese beiden Abfragen dürfen jedoch nicht dieselbe Zeile in ihrer Ergebnismenge enthalten. Also UNION ALL, was die Duplikatprüfung überspringt, wäre hier angebracht. (Ich benutze nicht UNIONso viel, aber ich finde im Allgemeinen, dass aus der Zeit, die ich benutze, ich viel häufiger UNION ALL.)
jpmc26
4

Was würde passieren, wenn Sie das "st_multi(st_intersection(a.geom,b.geom))"Teil weglassen ?

Bedeutet die folgende Abfrage nicht dasselbe ohne sie? Ich habe es anhand der von Ihnen angegebenen Daten ausgeführt.

INSERT INTO parcel_jurisdictions(parcel_gid,jurisdiction_gid,isect_geom)
  SELECT
    a.orig_gid parcel_gid,
    b.orig_gid jurisdiction_gid,
    a.geom
  FROM
    valid_parcels a, valid_jurisdictions b
  WHERE
    st_intersects(a.geom,b.geom);

Aufbau

Processor: AMD Athlon II X4 635 2.9 GHz 
Memory: 4 GB
Platform: Windows 7 Professional
Potgres Version: 8.4
Postgis Version: "POSTGIS="2.0.1 r9979" GEOS="3.3.5-CAPI-1.7.5" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.9.1, released 2012/05/15" LIBXML="2.7.8" LIBJSON="UNKNOWN" TOPOLOGY RASTER"

Ergebnisse analysieren

"Nested Loop  (cost=0.00..7505.18 rows=217489 width=1580) (actual time=1.994..248405.616 rows=329150 loops=1)"
"  Join Filter: _st_intersects(a.geom, b.geom)"
"  ->  Seq Scan on valid_jurisdictions b  (cost=0.00..37.25 rows=525 width=22621) (actual time=0.054..1.732 rows=525 loops=1)"
"  ->  Index Scan using valid_parcels_index on valid_parcels a  (cost=0.00..11.63 rows=2 width=1576) (actual time=0.068..6.423 rows=1366 loops=525)"
"        Index Cond: (a.geom && b.geom)"
"Total runtime: 280087.497 ms"
Vinayan
quelle
Nein, er möchte die resultierenden Kreuzungspolygone, aber Ihre Abfrage zeigt sehr gut, dass der gesamte Schmerz in der Kreuzungsgeneration liegt, nicht im binären True / False-Test-Teil der Abfrage. Und das ist durchaus zu erwarten, da der True / False-Code stark optimiert ist, während dies bei der Schnittmengenerzeugung nicht der Fall ist.
Paul Ramsey