5.13. 关于获取要素：

== 获取要素最好用GetNextFeature而不要用索引值 == 因为shapefile在实现的时候是以0为起始索引。mapinfo的tab是以1为起始值，这样在获取的时候需要判断读取数据源类型。但是用!GetNextFeature就会自动从正确的起始值开始读取。如果需要重新开始读取，就用!ResetReading重新设置读取位置到开头。

== 获取要素中的坐标值最好用ExportToWkb而不要用GetX,GetY == !ExportToWkb可以批量获取所有的坐标值。而GetX，GetY只能获取单独的一个值。在绘制的时候还要费神去把它们合成一个数组，而获取了wkb，就可以进行批量解析，一下生成数组，这样，我们就可以利用到了强大的Numeric模块，进行神奇的矩阵魔术。至于wkb的格式，[[http://publib.boulder.ibm.com/infocenter/db2luw/v8/index.jsp?topic=/com.ibm.db2.udb.doc/opt/rsbp4121.htm|这里]]有简单的说明。 我写了一个小模块，用于解析wkb

{{{ #! python import struct import ogr from array import array

import sys endian_name = sys.byteorder

wkbXDR = ‘>’ # Big Endian wkbNDR = ‘<’ # Little Endian

def choose(bool, a, b):

return (bool and [a] or [b])[0]

BTOR = choose(endian_name == ‘little’,wkbNDR,wkbNDR)

def up_endian_type(wkb):

endian_t = struct.unpack(‘b’,wkb[0])[0] endian = choose(endian_t,’<’,’>’) wkbtype = struct.unpack(endian+’I’,wkb[1:5])[0] return endian,wkbtype,endian_t

def up_len(wkb,beg,endian):

return struct.unpack(endian+’I’,wkb[beg:beg+4])[0]

def up_point(wkb):

endian,wkbtype,et = up_endian_type(wkb) points = struct.unpack(endian+”2d”,wkb[5:]) return points

def up_linestring(wkb):

endian,wkbtype,et = up_endian_type(wkb) lenght = up_len(wkb,5,endian) points = array(‘d’,wkb[9:9+lenght*16]) if endian != BTOR : points.byteswap() return points

def up_linearring(wkb,ringcount,endian):

#endian,wkbtype,et = up_endian_type(wkb) points = [] ptr = 0 for i in range(ringcount):

length = up_len(wkb,ptr,endian) ps = array(‘d’,wkb[ptr+4:ptr+4+length*16]) if endian != BTOR : ps.byteswap() points.append(ps) ptr += 4+length*16

return points,ptr

def up_polygon(wkb,sub=-1):

endian,wkbtype,et = up_endian_type(wkb) if sub == -1:

ringcount = up_len(wkb,5,endian) points = up_linearring(wkb[9:],ringcount,endian)[0] return points

else:

points = [] ptr = 5 ringcount = up_len(wkb,ptr,endian) ps,ringlen = up_linearring(wkb[ptr+4:],ringcount,endian) points.append(ps) ptr += 4+ringlen return points,ptr

def up_mpoint(wkb):

endian,wkbtype,et = up_endian_type(wkb) subcount = up_len(wkb,5,endian) points = [] ptr = 9 for i in range(subcount):

subps = up_point(wkb[ptr:]) points.append(subps) ptr += 9+len(subps)*8

return points

def up_mlinestring(wkb):

endian,wkbtype,et = up_endian_type(wkb) subcount = up_len(wkb,5,endian) points = [] ptr = 9 for i in range(subcount):

subps = up_linestring(wkb[ptr:]) points.append(subps) ptr += 9+len(subps)*8

return points

def up_mpolygon(wkb):

endian,wkbtype,et = up_endian_type(wkb) subcount = up_len(wkb,5,endian) points = [] ptr = 9 for i in range(subcount):

subps,size = up_polygon(wkb[ptr:],i) points.append(subps) ptr += size

return points

fun_map = {

ogr.wkbPoint : up_point, ogr.wkbLineString : up_linestring, ogr.wkbPolygon : up_polygon, ogr.wkbMultiPoint : up_mpoint, ogr.wkbMultiLineString : up_mlinestring, ogr.wkbMultiPolygon : up_mpolygon }

def WkbUnPacker(wkb):

endian,wkbtype,endian_t = up_endian_type(wkb) foo = fun_map[wkbtype] points = foo(wkb) return [endian_t,wkbtype,points]

if __name__ == “__main__”:

import time ds = ogr.Open(“e:/gisdata/data/country.shp”) layer = ds.GetLayer() begt = time.time() #count = layer.GetFeatureCount() #count = 2 feature = layer.GetNextFeature() #for i in [1]:#range(count): while feature is not None:

#feature = layer.GetFeature(i) geom = feature.GetGeometryRef() #print geom.ExportToWkt() wkb = geom.ExportToWkb() wkbarr = WkbUnPacker(wkb) feature = layer.GetNextFeature() #print wkbarr

print time.time()-begt

}}} 因为手上数据不全，无法测试一些特殊情况。如果谁有数据，可以测试一下，如果有错，希望能反馈给我！

== 加快读取大量feature时的循环速度 == ogr在python绑定中读取大量feature时速度总是很慢。16000条记录单单只是获取Feature竟然运行了130秒，简直令人无法忍受。后来打开了ogr.py察看了一下绑定代码，发现在Feature类中有“__setarr__”和”__getarr__”两个自省方法在创建Feature后不断运行。试着注释掉它们，速度一下就提高到2秒。简直是飞跃啊！不过注释掉自省方法会不会有后遗症，我现在还不敢说（如果大家发现注释掉这两个自省方法有问题，或者不需要注释掉调代码而有更好的方法，希望能反馈给我）。

= 反馈 = 如果您发现我写的东西中有问题，或者您对我写的东西有意见，请登陆[[http://groups.google.com/group/geosings?hl=zh-CN|这个论坛]]。