2016-04-04 15 views
2

Angesichts Land Polygone als ShapelyMultiPolygon, ich möchte das (Multi-) Polygon, das die z. 12 Seemeilen Puffer entlang der Küste.Geodätische Pufferung in Python

Die Verwendung der formschönen Methode buffer funktioniert nicht, da euklidische Berechnungen verwendet werden.

Kann mir jemand sagen, wie man geodätische Puffer in Python berechnet?

+0

Ich schlage vor, diese Frage auf gis.SE zu verschieben/stellen –

+0

@MikeT Ich bin in Ordnung mit der Frage auf gis.SE zu bewegen. Bitte fahre fort. Vielen Dank. – ARF

+0

Es könnte schneller sein, diese Frage selbst zu kopieren/einzufügen/zu löschen, da die Administratoren auf SO mit anderen Abfragen überlastet sind. –

Antwort

3

Dies ist kein formschönes Problem, da shapely in seiner Dokumentation ausdrücklich angibt, dass die Bibliothek nur für die planare Berechnung gedacht ist. Um Ihre Frage zu beantworten, sollten Sie jedoch die Koordinatensysteme angeben, die Sie für Ihre Multipolygone verwenden. Angenommen, Sie verwenden WGS84-Projektion (lat, lon), dies ist ein Rezept, das ich in einer anderen Frage SO gefunden habe (fix-up-shapely-polygon-object-when-discontinuous-after-map-projection). Sie benötigen pyproj Bibliothek.

import pyproj 
from shapely.geometry import MultiPolygon, Polygon 
from shapely.ops import transform as sh_transform 
from functools import partial 

wgs84_globe = pyproj.Proj(proj='latlong', ellps='WGS84') 

def pol_buff_on_globe(pol, radius): 
    _lon, _lat = pol.centroid.coords[0] 
    aeqd = pyproj.Proj(proj='aeqd', ellps='WGS84', datum='WGS84', 
         lat_0=_lat, lon_0=_lon) 
    project_pol = sh_transform(partial(pyproj.transform, wgs84_globe, aeqd), pol) 
    return sh_transform(partial(pyproj.transform, aeqd, wgs84_globe), 
          project_pol.buffer(radius)) 

def multipol_buff_on_globe(multipol, radius): 
    return MultiPolygon([pol_buff_on_globe(g, radius) for g in multipol]) 

pol_buff_on_globe Funktion hat die folgenden. Erstellen Sie zunächst eine äquidistante Azimut-Projektion, die im Mittelpunkt des Polygons zentriert ist. Ändern Sie dann das Koordinatensystem des Polygons zu dieser Projektion. Danach baut er den Puffer auf und ändert dann das Koordinatensystem des gepufferten Polygons in das WGS84-Koordinatensystem.

Einige besondere Sorgfalt erforderlich:

  • Sie müssen herausfinden, wie die Strecke, die Sie auf den Abstand in aeqd Projektion verwendeten übersetzen möchten.
  • Achten Sie darauf, die Pole nicht zu puffern (siehe die genannte SO-Frage).
  • Die Tatsache, dass wir den Mittelpunkt des Polygons verwenden, um die Projektion zu zentrieren, sollte garantieren, dass die Antwort gut genug ist, aber wenn Sie spezifische Genauigkeitsanforderungen haben, sollten Sie diese Lösung NICHT verwenden oder zumindest eine Charakterisierung des Fehlers vornehmen das typische Polygon, das Sie verwenden.