Cartopy 系列教程：安装并绘制地图

bugsuse

发布于 2020-04-21 17:45:53

10K0

发布于 2020-04-21 17:45:53

Cartopy 也是一个 python 地图绘制包，同样能完成很多 Basemap 能实现的功能，而且旨在使数据分析及可视化尽可能简单。其利用了强大的 PROJ.4，numpy，shapely库，提供了简单直观的绘图接口，可以创建满足出版质量的地图。

Cartopy具有如下特点：

面向对象的投影定义
point，line，vector，polygon 和 image 可在投影间转换
和 matplotlib 绘图包具有很高集成度，从而具有简单易用的接口
集成了 shapely 库，从而具有强大的 vector 数据处理能力

你可以从其 Github 主页获取源代码 [注1]。

安装

Cartopy 的安装相对于 Basemap 来说就要容易多了（windows 表示不服）。

使用 Conda 可在所有平台进行安装

conda install -c scitools cartopy

当然你也可以使用源码安装

下载源码包后，解压并从终端进入到解压路径下

python setup.py install

如果你是在 windows 上使用 cartopy 的话，建议使用 conda 安装，如果安装时出现冲突提示，可以创建一个虚拟环境然后再安装。

打开终端

conda create --name pycartopy python=2.7 matplotlib

这里仅先安装 matplotlib 绘图包，创建好环境之后再安装其它需要的包和 cartopy

安装好后进入环境继续安装需要的库（当然你也可以创建 python3 的环境）

activate pycartopy

conda install -c scitools cartopy

绘制地图

Cartopy 可以非常简单的创建地图，而不需要像 Basemap 一样先创建 map 实例：

import cartopy.crs as ccrs
import matplotlib.pyplot as plt

ax = plt.axes(projection=ccrs.PlateCarree())
ax.coastlines()

plt.show()

matplotlib 可用的投影在 Cartopy 投影列表页可以查看。

plt.axes(projection=ccrs.PlateCarree()) 设置一个 GeoAxes 类，其具有大量的地图相关的方法。在之前的例子中，我们使用 coastlines 方法添加海岸线到地图中。

下面使用另一种投影创建地图，然后使用 stock_img 方法添加背景图片到地图中。

import cartopy.crs as ccrs
import matplotlib.pyplot as plt

ax = plt.axes(projection=ccrs.Mollweide())
ax.stock_img()
plt.show()

添加数据到地图

一旦你创建好了地图，添加数据的方式就和常规的 matplotlib 绘图方法相同。默认情况下，添加到 GeoAxes 的任何数据的坐标系统和 GeoAxes 本身的坐标系统是相同的。为了控制给定数据的坐标系，你可以添加 transform 关键词，并传递一个合适的 cartopy.crs.CRS 实例。

import cartopy.crs as ccrs
import matplotlib.pyplot as plt

ax = plt.axes(projection=ccrs.PlateCarree())
ax.stock_img()

ny_lon, ny_lat = -75, 43
delhi_lon, delhi_lat = 77.23, 28.61

plt.plot([ny_lon, delhi_lon], [ny_lat, delhi_lat],
         color='blue', linewidth=2, marker='o',
         transform=ccrs.Geodetic(),
         )

plt.plot([ny_lon, delhi_lon], [ny_lat, delhi_lat],
         color='gray', linestyle='--',
         transform=ccrs.PlateCarree(),
         )

plt.text(ny_lon - 3, ny_lat - 12, 'New York',
         horizontalalignment='right',
         transform=ccrs.Geodetic())

plt.text(delhi_lon + 3, delhi_lat - 12, 'Delhi',
         horizontalalignment='left',
         transform=ccrs.Geodetic())

plt.show()

注意：在 PlateCarree 投影平面地图上， New York 和 Delhi 之间的蓝色线并不是直线，这是因为 Geodetic 坐标系是真正的球面坐标系，两点之间的线被定义为在球坐标，而不是2D笛卡尔坐标上的最短路径。

注意：默认情况下，matplotlib 会根据绘图数据自动确定坐标轴上下限。因为 cartopy 使用的 GeoAxes 类，这相当于询问地图的坐标轴范围。有时候会得到预期效果，但有时候并不会。

有以下几种方式可以设置 cartopy GeoAxes 的范围：

使用 set_global 方法可以绘制全球图
任何坐标系中都可以使用 set_extent 方法根据边界框来设置地图范围
也可以使用 set_xlim 或 set_ylim 方法设置坐标轴界限

关于其它的一些功能，后面再进行介绍。既然上一篇介绍了 Basemap 的白化方法，趁热打铁此篇介绍一下使用 Cartopy 进行白化的方法。

白化

Cartopy 同样也能像 Basemap 那样对地图进行白化，而且效果也是很好的。

下面就利用 Cartopy 来白化地图：

import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.io.shapereader as shpreader
from cartopy.examples.waves import sample_data

# subplot_kw 是一个用于设置子图属性的关键词参数：设置投影
fig, ax = plt.subplots(subplot_kw=dict(projection=ccrs.PlateCarree()))

plate_carre_data_transform = ccrs.PlateCarree()._as_mpl_transform(ax)

lons, lats, data = sample_data()

con = ax.contourf(lons, lats, data, transform=ccrs.PlateCarree())

shpfilename = shpreader.natural_earth(resolution='10m',
                                      category='cultural',
                                      name='admin_0_countries')
ax.set_xlim([60, 140])
ax.set_ylim([0, 60])                                      

# pathkw 用于设置 patch 的属性                                      
path = crsmask(shpfilename, 'China', ax, ccrs.PlateCarree(), con, 
       pathkw =  dict(facecolor='none',edgecolor='black', linewidth=1.5))

plt.show()

使用的 shp 文件由 Natural Earth 提供。 Cartopy 中提供了获取文件的接口。

上例中进行白化时使用了 crsmask 函数：

def crsmask(shpfilename, region, ax, crs, conf, pathkw = None):
    '''
    目的:
        对选定区域进行白化
    输入参数:
        @parameter:   shpfilename  
                   natural earth下载的shp文件
        @type:  string
        
        @parameter: region  shp文件中包含的区域名
        @type:  string
        
        @parameter: ax  坐标轴实例
        @type:  GeoAxes instance
        
        @parameter: crs  投影
        @type:  cartopy.crs
        
        @parameter:  conf  绘图返回值，比如contourf和contour返回值
        @type:
        
        @parameter: pathkw  PathPatch关键字参数
        @type:  dict
        
    输出参数:
        @parameter:   path  clipping path
        @type:  matplotlib.patches.PathPatch
    '''

    import cartopy.io.shapereader as shpreader
    from cartopy.mpl.patch import geos_to_path
    from matplotlib.patches import PathPatch
    
    reader = shpreader.Reader(shpfilename)
    countries = reader.records()
    
    try:
        multipoly, = [country.geometry for country in countries
                        if country.attributes['name'] == region]
    except KeyError:
        multipoly, = [country.geometry for country in countries
                        if country.attributes['NAME'] == region]
                        
    main_geom = sorted(multipoly.geoms, key=lambda geom: geom.area)[-1]

    path, = geos_to_path(main_geom)

    plate_carre_data_transform = crs._as_mpl_transform(ax)

    for collection in conf.collections:
        collection.set_clip_path(path, plate_carre_data_transform)
    
# Draw the path of interest.
    path = PathPatch(path, transform = crs, **pathkw)
    
    return path

程序中 try 和 except 中对文件中的区域进行了判断，这是用于白化某一块区域。如果你想白化多个区域的话，可以将