不使用arcpy地理处理shapefile的Python脚本示例

#1 楼

这很奇怪，好像人们突然发现Python的功能一样（没有ArcPy就是其中的一个Python模块），例如参见问题Visualize shapefile in Python：


geospatial processing in Python具有悠久的历史，比Arcpy（或arcgisscripting）要古老得多->保罗说，这里没有“模仿” ArcPy的功能，大多数已经在ArcPy之前存在。
Python模块的参考是Python软件包索引（Pypi）并有专门的部分：主题::科学/工程:: GIS

您可以使用这些模块执行任何操作，并且它比ArcPy更加容易和快捷，因为它是纯净的Python（无光标...）。

Shapely是用于处理地理空间几何的这些模块之一->计算多边形的面积并将多边形转换为点。.
如果要处理矢量图层，有osgeo / ogr，Fiona或Pyshp（以及其他较少使用的图层）->通过属性查询shapefile，创建新图层fr进行选择，计算多边形的面积，将多边形转换为点
用于处理栅格，标准是osgeo / gdal

用于空间分析，有Pysal

对于3D，您可以使用numpy或scipy等其他Scientific模块（3D算法，网格，也可以使用统计数据，地统计，2D或3D）
，我不会谈论mapnik，matplotlib / basemap，Geodjango和.. 。

您可以将所有（Pysal与shape，...）组合在一起，并将它们与其他Scientific模块混合。
因此获取Python脚本示例，搜索Pyshp Fiona，ogr，gdal或可以在gis.stackexchange或Internet中正常运行（很多示例，不仅是英语）。
其中一个是法语（脚本和图形是通用的！）：


Python：在地质角度上使用矢量和栅格图层，而没有GIS软件
其他英语版本：

具有Python，Shapely和Fiona的GIS和西班牙语

/>
使用顶点的坐标确定不规则多边形的区域
gis.stackexchange
高程剖面线的每侧10 km
使用Pyshp更新属性
如何创建3D从栅格获取shapefile？
用于获取两点之间的高程差异的Python脚本
等

仅使用Python字典的Fiona可以更简单地编写Aaron演示的脚本：

import fiona
with fiona.open('sites.shp', 'r') as input:
    with open('hw1a.txt', 'w') as output:
       for pt in input:
           id = pt['properties']['id']
           cover = pt['properties']['cover']
           x = str(point['geometry']['coordinates'][0])
           y = str(point['geometry']['coordinates'][21])
           output.write(id + ' ' + x + ' ' + y+ ' ' + cover + '\n')

，如果还要匀称使用：

from shapely.geometry import shape
with fiona.open('sites.shp', 'r') as input:
    with open('hw1a.txt', 'w') as output:
       for pt in input:
           id = pt['properties']['id']
           cover = pt['properties']['cover']
           x = str(shape(pt['geometry']).x)
           y = str(shape(pt['geometry']).y)
           output.write(id + ' ' + x + ' ' + y+ ' ' + cover + '\n')

还有两本书：
Eric Westra的Python地理空间开发。

使用Joel Lawhead的Python学习地理空间分析
（来源：cloudfront.net）
Python还被QGIS（Quantum GIS），GRASS GIS，gvSIG或OpenJump等其他GIS应用程序用作脚本语言。或Paraview等3D建模器（还有Blender！）。您可以在所有这些应用程序中使用大多数地理空间模块（请参阅使用Blender可视化QGIS数据）

#2 楼

我强烈推荐USU网站使用开源GIS使用Python进行地理处理，以帮助您入门。在整个练习中，他们主要使用GDAL / OGR库。安装GDAL / OGR可能会有些困难，因此此博客文章可能对您有所帮助：在Windows上为Python安装GDAL（和OGR）。另请查看在GIS.SE上使用Arcpy的替代方法。

下面的开源地理处理脚本示例（来自USU站点）用于提取属性数据并将其写入文本文件：

# import modules
import ogr, os, sys

# set the working directory
os.chdir('f:/data/classes/python/data')

# open the output text file for writing
file = open('hw1a.txt', 'w')

# get the shapefile driver
driver = ogr.GetDriverByName('ESRI Shapefile')

# open the data source
datasource = driver.Open('sites.shp', 0)
if datasource is None:
  print 'Could not open file'
  sys.exit(1)

# get the data layer
layer = datasource.GetLayer()

# loop through the features in the layer
feature = layer.GetNextFeature()
while feature:

  # get the attributes
  id = feature.GetFieldAsString('id')
  cover = feature.GetFieldAsString('cover')

  # get the x,y coordinates for the point
  geom = feature.GetGeometryRef()
  x = str(geom.GetX())
  y = str(geom.GetY())

  # write info out to the text file
  file.write(id + ' ' + x + ' ' + y + ' ' + cover + '\n')

  # destroy the feature and get a new one
  feature.Destroy()
  feature = layer.GetNextFeature()

# close the data source and text file
datasource.Destroy()
file.close()

#3 楼

您可能对GDAL / OGR感兴趣。

GDAL用于处理栅格，而OGR用于矢量。两者都是开源库。

如果要消除对ArcPy的依赖，可以通过将信息读取到数组并在纯Python中运行自己的计算来模仿某些功能。

我最近通过选择多边形中的点来做到这一点，如此处所示。给定多边形顶点的坐标，它利用射线投射算法确定点是否位于多边形内。

#4 楼

我从未亲自使用过，但是办公室中的其他人喜欢匀称使用：https://pypi.python.org/pypi/Shapely