GDALpython教程（1）——⽤OGR读写⽮量数据本教程的讲义和源码都是取⾃Utah State University的openGIS课程
相关资料，包括讲义、源码、数据样例，请从此处下载
本⼈只是做点翻译，写写学习体会⽽已，版权属于原作者。
欢迎转载，不过别忘了上⾯这段话。
==================================================
为什么⽤open source？
优点
1. 免费，适合个⼈和⼩公司
2. 强⼤的开发⼯具，bug更容易
3. 跨平台，windows和linux都能⽤
4. 拉风！
缺点
1. 没有内嵌地理处理器
2. ⽤的⼈少
Open source RS/GIS模块
1. OGR⽮量库：简单的⽮量数据读写，是GDAL的⼀部分
2. GDAL地理空间数据抽象库：
a)        读写栅格数据
b)        ArcGIS也是基于GDAL开发的
c)        C++库，但是可以⽤python调⽤
相关模块
1. Numeric：⾼速的数组处理，对栅格数据尤其重要
2. NumPy：下⼀代的Numeric
3. 更强⼤的gis库
导⼊库：
import ogr
或者：
python怎么读取json文件
from osgeo import ogr
万能的⽅法是：
try:
from osgeo import ogr
except:
import ogr
要读取某种类型的数据，必须要先载⼊数据驱动，也就是初始化⼀个对象，让它“知道”某种数据结构。
import ogr
driver = ogr.GetDriverByName(‘ESRI Shapefile’)
数据驱动driver的open()⽅法返回⼀个数据源对象
open(<filename>, <update>)
其中update为0是只读，为1是可写
例如：
from osgeo import ogr
driver = ogr.GetDriverByName('ESRI Shapefile')
filename = 'C:/Users/gongwei/Documents/My eBooks/python_and_sage/GDAL python/test/ospy_data1/sites.shp' dataSource = driver.Open(filename,0)
if dataSource is None:
print 'could not open'
print 'done!'
注意filename⼀定要写绝对路径！
因为⼀定要⽤绝对路径，为了简化代码，经常会使⽤到os.chdir()
读取数据层
layer = dataSource.GetLayer(0)
⼀般ESRI的shapefile都是填0的，如果不填的话默认也是0.
再看看这个数据层⾥⾯有⼏个点呢？
n = layer.GetFeatureCount()
print 'feature count:', n
读出上下左右边界
extent = layer.GetExtent()
print 'extent:', extent
print 'ul:', extent[0], extent[3]
print 'lr:', extent[1], extent[2]
读取某⼀要素feature（总算切⼊正题了），这⾥读取的是⼀个点
feat = layer.GetFeature(41)
fid = feat.GetField('id')
print fid
feat = layer.GetFeature(0)
fid = feat.GetField('id') #should be a different id
print fid
另外还有按顺序读取feature，循环遍历所有的feature
feat = layer.GetNextFeature() #读取下⼀个
while feat:
feat = layer.GetNextFeature()
later.ResetReading() #复位
提取feature的⼏何形状
geom = feat.GetGeometryRef()
geom.GetX()
geom.GetY()
print geom.
释放内存
feature.Destroy()
关闭数据源，相当于⽂件系统操作中的关闭⽂件
dataSource.Destroy()
读完了再说怎么写
创建新⽂件
driver.CreateDataSource(<filename>)
但是这个⽂件不能已经存在了，否则会出错
创建新的layer
dataSource.CreateLayer(<name>,CreateLayer(<name>, geom_type=<OGRwkbGeometryType>, [srs])举个例⼦：
ds2 = driver.CreateDataSource('test.shp')
layer2 = ds2.CreateLayer('test', geom_type=ogr.wkbPoint)
要删除⼀个shp⽂件
driver.DeleteDataSource('test.shp')
要添加⼀个新字段，只能在layer⾥⾯加，⽽且还不能有数据
添加的字段如果是字符串，还要设定宽度
fieldDefn = ogr.FieldDefn('id', ogr.OFTString)
fieldDefn.SetWidth(4)
layer.CreateField(fieldDefn)
添加⼀个新的feature，⾸先得完成上⼀步，把字段field都添加齐了
然后从layer中读取相应的feature类型，并创建feature
featureDefn = layer.GetLayerDefn()
feature = ogr.Feature(featureDefn)
设定⼏何形状
feature.SetGeometry(point)
设定某字段的数值feature.SetField('id', 23)
将feature写⼊layer layer.CreateFeature(feature)

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。