Python地学分析 — 矢量数据读取 01

欢迎关注博主的微信公众号：“智能遥感”。

该公众号将为您奉上Python地学分析、爬虫、数据分析、Web开发、机器学习、深度学习等热门源代码。

代码和数据可以直接在我的GitHub上下载：

https://github.com/xbr2017/PyGdal_vector_01

********************************************************************

在上一节矢量数据集的介绍后，本节开始Python实战！！！

说到读取矢量数据或者遥感图像，我们不得不介绍今天的主角儿 — GDAL库。是的，你没看错，就是C++里面经常用到的那个库，在这里，Python也可以有哦（突然想到网上比较流行的一句话：Python除了不能生小孩，其他什么都能）。

那么，我们先来看一下如何安装GDAL。这里以Windows环境为例，有两种常用方法：

1. https://www.lfd.uci.edu/~gohlke/pythonlibs/上下载GDAL对应的*.whl格式文件。如Python2.7版本32位对应的GDAL文件：GDAL-2.2.4-cp27-cp27m-win32.whl。进入cmd命令行，cd到该GDAL文件所在路径下，输入命令pip install  GDAL-2.2.4-cp27-cp27m-win32.whl即可安装好GDAL库。

2. 如果觉得上面操作太繁琐（其实蛮锻炼人，你可能在方法1中遇到坑，找度娘或翻墙谷歌去吧），就直接打开cmd命令行在线安装，输入命令pip install gdal。

下面终于上代码啦，主要实现打印shp矢量文件的前10个要素。

# _*_ coding: utf-8 _*_
__author__ = 'xbr'
__date__ = '2018/10/30 22:44'

import sys

from osgeo import ogr

# 矢量图shp文件
fn = r'D:\osgeopy-data\global\ne_50m_populated_places.shp'
# 打开文件
ds = ogr.Open(fn, 0)
# 如果文件为空，则报错
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))
# 在确保打开文件后，从数据源中检索第一个图层
lyr = ds.GetLayer(0)
# 下面for循环是获取图层中前10个要素，
# 每个要素所含字段：城市名、最大人口数、经度、纬度
i = 0
for feat in lyr:
    pt = feat.geometry()
    x = pt.GetX()   # 经度
    y = pt.GetY()   # 纬度
    name = feat.GetField('NAME')   # 城市名
    pop = feat.GetField('POP_MAX') # 最大人口数
    print(name, pop, x, y)
    i += 1
    if i == 10:
        break
# 删除ds变量以强制文件关闭
del ds

该代码数据结果如下：

('Bombo', 75000, 32.533299524864844, 0.5832991056146284)
('Fort Portal', 42670, 30.27500161597942, 0.671004121125236)
('Potenza', 69060, 15.798996495640267, 40.642002130098206)
('Campobasso', 50762, 14.655996558921856, 41.56299911864397)
('Aosta', 34062, 7.315002595706176, 45.7370010670723)
('Mariehamn', 10682, 19.949004471869102, 60.09699618489543)
('Ramallah', 24599, 35.206209378189556, 31.90294475142406)
('Vatican City', 832, 12.453386544971766, 41.903282179960115)
('Poitier', 85960, 0.3332765285345545, 46.58329225573658)
('Clermont-Ferrand', 233050, 3.080008095928406, 45.779982115759424)