Python与开源GIS:根据属性条件选择要素

根据属性来对数据库中的记录进行选择是数据库应用中必不可少的一项功能,对于地理数据,也是如此,譬如人口在百万以上的城市,面积在百万平方千米以上的国家。

根据属性条件选择要素

Layer对象有一个方法叫SetAttributeFilter(),可以将根据属性,将Layer中符合某一条件的Feature过滤出来。设定了Filter之后就可以用GetNextFeature()方法依次取出符合条件的Feature了。

    from osgeo import ogr
    import os
    shpfile = '/gdata/world_borders.shp'
    ds = ogr.Open(shpfile)
    layer = ds.GetLayer(0)
    lyr_count = layer.GetFeatureCount()
    print(lyr_count)
    layer.SetAttributeFilter("AREA > 800000")
    lyr_count = layer.GetFeatureCount()
    print(lyr_count)

上面的代码首先将图层中的要素的数目打印出来。 然后,使用过滤器,来选择面积在80万平方公里以上的国家。这里使用了“AREA”字段,这个字段是早已经在Shapefile中生成好的。 然后再打印图层中要素的数目,就会发现数据变少了。

根据属性条件生成要素

大多数时候,选择要素的目的并不是为了简单地查看一下数目, 而是要将数据保存了,以方便以后的使用。请仔细查看下面的代码。注意,生成矢量数据可以在不同的层次上。

    driver = ogr.GetDriverByName("ESRI Shapefile")
    extfile = '/v/filter_demo.shp'
    if os.access( extfile, os.F_OK ):
    driver.DeleteDataSource( extfile )
    newds = driver.CreateDataSource(extfile)
    layernew = newds.CreateLayer('rect',None,ogr.wkbPolygon)
    feat = layer.GetNextFeature()
    while feat is not None:
    layernew.CreateFeature(feat)
    feat = layer.GetNextFeature()
    newds.Destroy()

看一下生成的结果吧。这是一幅不完整的世界地图,因为它只包含了面积在80万平方公里以上的国家。

对于这个结果,首先要注意的就是有一些面积比较小的图斑也被选了出来,这些图斑是属于某些国家的,其面积定义是与某个大图斑一致的,也就是这此处的面积是国家的面积,而并非是图斑的面积;另外要注意的是台湾并没有被选择出来:这是个政治问题,而非技术问题。