1.空间计算
地理数据处理（geoprocessing）计算函数：
多边形（Polygon）：
1、交：poly3.Intersection(poly2)
2、并：poly3.Union(poly2)
3、差：poly3.Difference(poly2)
4、补：poly3.SymmetricDifference(poly2)
下面我们来做一段练习，首先我们有如下一份数据，
该图层内有两个几何体，我们就用这两个集合体来进行练习
from osgeo import ogr, osr
shp = r'D:/work/全国行政区划/江苏省.shp'
drivier = ogr.GetDriverByName('ESRI Shapefile')
ds = drivier.Open(shp)
layer = ds.GetLayer(0)
srs = layer.GetSpatialRef() # 获取投影
print(layer.GetFeatureCount())
feat_js = layer.GetNextFeature()
feat_rec = layer.GetNextFeature()
geom_js = feat_js.GetGeometryRef()
geom_rec = feat_rec.GetGeometryRef()
geom_In = geom_js.Intersection(geom_rec) # 相交
geom_un = geom_js.Union(geom_rec) # 并集
geom_diff = geom_js.Difference(geom_rec) # 差集
geom_sy = geom_js.SymmetricDifference(geom_rec) # 补充
# 写入shp文件
ds2 = drivier.CreateDataSource('test2.shp')
layer = ds2.CreateLayer('test', srs, geom_type=ogr.wkbPolygon)
#定义字段
fieldDefn = ogr.FieldDefn('id', ogr.OFTString)
fieldDefn.SetWidth(4)
layer.CreateField(fieldDefn)
#定义要素
featureDefn = layer.GetLayerDefn()
feature1 = ogr.Feature(featureDefn)
feature2 = ogr.Feature(featureDefn)
feature3 = ogr.Feature(featureDefn)
feature4 = ogr.Feature(featureDefn)
#设置几何形状
feature1.SetGeometry(geom_In)
feature2.SetGeometry(geom_un)
feature3.SetGeometry(geom_diff)
feature4.SetGeometry(geom_sy)
#设定某字段的数值
feature1.SetField('id', 1)
feature2.SetField('id', 2)
feature3.SetField('id', 3)
feature4.SetField('id', 4)
#将feature写入layer
layer.CreateFeature(feature1)
layer.CreateFeature(feature2)
layer.CreateFeature(feature3)
layer.CreateFeature(feature4)
ds2.Destroy()
得到的结果
1、两个几何体的交集
2、两个几何体的并集
3、差集（geom_js去除相交的部分后剩余的部分）
4、补集，并集减去相交的部分
几何形状(Geometry)：
1.2给geometry加buffer，效果就是把点、线变成多边形，线条加粗：
<geom>.Buffer(<distance>)；
1.3两个geometry是否相等：<geom1>.Equal(<geom2>) ；
1.4两个geometry之间的最短距离：<geom1>.Distance(<geom2>) ；
1.5用方框边界框住这个几何形状，并返回四个角的点坐标(minx,maxx, miny, maxy)：<geom>.GetEnvelope()。
geom_js.GetENvelope()
(116.35518309500014, 121.92747178000022, 30.760280365000085, 35.127189627500144)
2.空间分析的方法
2.1.获取面积
建议采用GetArea()获得整个多边形覆盖的面积。
geom_js.GetArea()`
>>> 9.896416947183278
不过此用法以“度”为量纲，已经失去它作为面积的意义了。若要较为准确的面积，可根据空间参考将度转换为米，再将面积计算为平方米。
判断两个对象的关系
2.2Intersect语句是判断两个要素是否相交：
poly2.Intersect(poly1)
返回0表示不相交，返回1表示相交。
2.3Disjoint是判断两个要素是否不相交：
poly2.Disjoint(poly1)
返回1表示不相交，返回0表示相交，跟Intersect正好相反。
2.4Touch表示相邻：
poly2.Touches(poly1)
返回0表示不相邻，返回1表示相邻。
2.5Crosses表示横穿，一般是指一条线穿过一个多边形：
poly2.Crosses(line)
返回0表示不穿过，返回1表示穿过。
2.6Within表示包含，指一个要素完全被另一个要素包含：
ptB.Within(poly1)
返回0表示点在多边形外，返回1表示点在多边形内。
2.7Contains也表示包含，跟Within正好相反：
poly1.Contains(ptB)
与Within的主调对象与参数对换。
2.8Overlaps重叠，只有两个多边形之间才能算作overlap
poly2.Overlaps(poly3)
返回0表示不重叠，返回1表示重叠。
