无投影信息
自己设置
def intersection(shp1, shp2, outshp):
"""
:param shp1: 缓冲区矢量路径
:param shp2: 目标矢量路径
:param outshp: 输出矢量路径
:return:
"""
driver = ogr.GetDriverByName("ESRI Shapefile")
dataSource1 = driver.Open(shp1, 1)
layer1 = dataSource1.GetLayer()
dataSource2 = driver.Open(shp2, 1)
layer2 = dataSource2.GetLayer()
# 获取到shp的投影坐标系信息
prosrs = layer1.GetSpatialRef()
geosrs = osr.SpatialReference()
# 设置输出坐标系为WGS84
geosrs.SetProjection("GCS_WGS_1984")
# 新建DataSource，Layer
out_ds = driver.CreateDataSource(outshp)
out_lyr = out_ds.CreateLayer(outshp, prosrs, ogr.wkbPolygon)
def_feature = out_lyr.GetLayerDefn()
for feature in layer1:
geom = feature.GetGeometryRef()
for feature_ in layer2:
geom_ = feature_.GetGeometryRef()
if geom.Intersect(geom_) == 1:
out_feature = ogr.Feature(def_feature)
out_feature.SetGeometry(geom_)
out_lyr.CreateFeature(out_feature)
out_ds.FlushCache()
del dataSource1, dataSource2
# 写入投影文件
geosrs.MorphFromESRI()
prjfile = open(outshp.replace('.shp','.prj'),'w')
prjfile.write(geosrs.ExportToWkt())
prjfile.close()
有投影信息
投影信息来源于.prj文件中。
def intersection(shp1, shp2, outshp):
"""
:param shp1: 缓冲区矢量路径
:param shp2: 目标矢量路径
:param outshp: 输出矢量路径
:return:
"""
driver = ogr.GetDriverByName("ESRI Shapefile")
dataSource1 = driver.Open(shp1, 1)
layer1 = dataSource1.GetLayer()
spatialref = layer1.GetSpatialRef()
dataSource2 = driver.Open(shp2, 1)
layer2 = dataSource2.GetLayer()
# 新建DataSource，Layer
out_ds = driver.CreateDataSource(outshp)
out_lyr = out_ds.CreateLayer(outshp, spatialref, ogr.wkbPolygon)
def_feature = out_lyr.GetLayerDefn()
for feature in layer1:
geom = feature.GetGeometryRef()
for feature_ in layer2:
geom_ = feature_.GetGeometryRef()
if geom.Intersect(geom_) == 1:
out_feature = ogr.Feature(def_feature)
out_feature.SetGeometry(geom_)
out_lyr.CreateFeature(out_feature)
out_ds.FlushCache()
del dataSource1, dataSource2
