1 #!C:\Program Files\pythonxy\python\python.exe
2 #-*- coding:gb2312 -*-
3
4 from osgeo importogr,osr,gdal5 importos6
7 """
8 Understanding OGR Data Type:9 Geometry - wkbPoint,wkbLineString,wkbPolygon,wkbMultiPoint,wkbMultiLineString,wkbMultiPolygon10 Attribute - OFTInteger,OFTReal,OFTString,OFTDateTime11 """
12
13 classARCVIEW_SHAPE:14 #------------------------------
15 #read shape file
16 #------------------------------
17 defread_shp(self,file):18 #open
19 ds = ogr.Open(file,False) #False - read only, True - read/write
20 layer =ds.GetLayer(0)21 #layer = ds.GetLayerByName(file[:-4])
22 #fields
23 lydefn =layer.GetLayerDefn()24 spatialref =layer.GetSpatialRef()25 #spatialref.ExportToProj4()
26 #spatialref.ExportToWkt()
27 geomtype =lydefn.GetGeomType()28 fieldlist =[]29 for i inrange(lydefn.GetFieldCount()):30 fddefn =lydefn.GetFieldDefn(i)31 fddict = {'name':fddefn.GetName(),'type':fddefn.GetType(),32 'width':fddefn.GetWidth(),'decimal':fddefn.GetPrecision()}33 fieldlist +=[fddict]34 #records
35 geomlist =[]36 reclist =[]37 feature =layer.GetNextFeature()38 while feature is notNone:39 geom =feature.GetGeometryRef()40 geomlist +=[geom.ExportToWkt()]41 rec ={}42 for fd infieldlist:43 rec[fd['name']] = feature.GetField(fd['name'])44 reclist +=[rec]45 feature =layer.GetNextFeature()46 #close
47 ds.Destroy()48 return(spatialref,geomtype,geomlist,fieldlist,reclist)49
50 #------------------------------
51 #write shape file
52 #------------------------------
53 defwrite_shp(self,file,data):54 gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8","YES");55 gdal.SetConfigOption("SHAPE_ENCODING","UTF-8");56 spatialref,geomtype,geomlist,fieldlist,reclist =data57 #create
58 driver = ogr.GetDriverByName("ESRI Shapefile")59 ifos.access(file, os.F_OK ):60 driver.DeleteDataSource(file)61 ds =driver.CreateDataSource(file)62 #spatialref = osr.SpatialReference( 'LOCAL_CS["arbitrary"]' )
63 #spatialref = osr.SpatialReference().ImportFromProj4('+proj=tmerc ...')
64 layer = ds.CreateLayer(file[:-4],srs=spatialref,geom_type=geomtype)65 #print type(layer)
66 #fields
67 for fd infieldlist:68 field = ogr.FieldDefn(fd['name'],fd['type'])69 if fd.has_key('width'):70 field.SetWidth(fd['width'])71 if fd.has_key('decimal'):72 field.SetPrecision(fd['decimal'])73 layer.CreateField(field)74 #records
75 for i inrange(len(reclist)):76 geom =ogr.CreateGeometryFromWkt(geomlist[i])77 feat =ogr.Feature(layer.GetLayerDefn())78 feat.SetGeometry(geom)79 for fd infieldlist:80 #print(fd['name'],reclist[i][fd['name']])
81 feat.SetField(fd['name'],reclist[i][fd['name']])82 layer.CreateFeature(feat)83 #close
84 ds.Destroy()85
86 #--------------------------------------
87 #main function
88 #--------------------------------------
89 if __name__ == "__main__":90 test =ARCVIEW_SHAPE()91 data = test.read_shp(r'../data/chn_adm2.shp')92 spatialref,geomtype,geomlist,fieldlist,reclist =data93 test.write_shp(r'../data/chn_adm2_bak.shp',[spatialref,geomtype,geomlist,fieldlist,reclist])