ESRIシェイプファイル内の特定の緯度/経度で情報を検索するためにpythonを使用するにはどうすればよいですか？

Question

私はESRIシェイプファイルを持っています（ここから：http://pubs.usgs.gov/ds/425/）。私は、指定された緯度/経度でシェイプファイル（この場合は表面的なマテリアル）から情報を検索するためにpythonを使用したいと考えています。

この問題を解決する最善の方法は何ですか？

ありがとうございました。

最終的な解決策：

#!/usr/bin/python 

from osgeo import ogr, osr 

dataset = ogr.Open('./USGS_DS_425_SHAPES/Surficial_materials.shp') 
layer = dataset.GetLayerByIndex(0) 
layer.ResetReading() 

# Location for New Orleans: 29.98 N, -90.25 E 
point = ogr.CreateGeometryFromWkt("POINT(-90.25 29.98)") 

# Transform the point into the specified coordinate system from WGS84 
spatialRef = osr.SpatialReference() 
spatialRef.ImportFromEPSG(4326) 
coordTransform = osr.CoordinateTransformation(
     spatialRef, layer.GetSpatialRef()) 

point.Transform(coordTransform) 

for feature in layer: 
    if feature.GetGeometryRef().Contains(point): 
     break 

for i in range(feature.GetFieldCount()): 
    print feature.GetField(i) 

Answer 1

gdal/ogrツールキットに対するPythonバインディングを使用できます。ここでは例です：

from osgeo import ogr 

ds = ogr.Open("somelayer.shp") 
lyr = ds.GetLayerByName("somelayer") 
lyr.ResetReading() 

point = ogr.CreateGeometryFromWkt("POINT(4 5)") 

for feat in lyr: 
    geom = feat.GetGeometryRef() 
    if geom.Contains(point): 
     sm = feat.GetField(feat.GetFieldIndex("surface_material")) 
     # do stuff... 

Answer 2

アウトPython Shapefile Library

これは、ジオメトリと異なる情報を与える必要があります。

Answer 3

別のオプションは、格好の良い（GEOS、PostGISのためのエンジンをベースにPythonライブラリ）と（ファイルを書き込む/読み取るため、基本的です）フィオナを使用することです：

import fiona 
import shapely 

with fiona.open("path/to/shapefile.shp") as fiona_collection: 

    # In this case, we'll assume the shapefile only has one record/layer (e.g., the shapefile 
    # is just for the borders of a single country, etc.). 
    shapefile_record = fiona_collection.next() 

    # Use Shapely to create the polygon 
    shape = shapely.geometry.asShape(shapefile_record['geometry']) 

    point = shapely.geometry.Point(32.398516, -39.754028) # longitude, latitude 

    # Alternative: if point.within(shape) 
    if shape.contains(point): 
     print "Found shape for point." 

ポリゴンが大きく/複雑な場合（ポリゴンが極端に不規則な海岸線を持ついくつかの国のシェイプファイルなど）、ポイントポリゴンテストを行うことはコストがかかることに注意してください。それはすぐに、より集中的なテストを行う前に、物事を除外するために境界ボックスを使用するように助けることができるいくつかのケースでは：

minx, miny, maxx, maxy = shape.bounds 
bounding_box = shapely.geometry.box(minx, miny, maxx, maxy) 

if bounding_box.contains(point): 
    ... 

最後に、残念ながら（それは大/不規則なシェープファイルをロードし、解析に時間がかかることに注意してくださいこれらのタイプのポリゴンは、しばしばメモリに保存するのに高価です）。