读者答疑：如何简单绘制全球土地覆盖图
个人信息
公众号：气 python 风雨
关注我获取更多学习资料，第一时间收到我的 Python 学习资料，也可获取我的联系方式沟通合作
温馨提示
由于可视化代码过长隐藏，可点击运行 Fork 查看
若没有成功加载可视化图，点击运行可以查看
ps:隐藏代码在【代码已被隐藏】所在行，点击所在行，可以看到该行的最右角，会出现个三角形，点击查看即可
前言
读者来信想看 python 读取 tif 数据绘制全球土地覆盖
数据为长时序动态土地覆盖数据（GLASS-GLC）
镜像为气象分析 3.7
安装库
!pip install rioxarray -i https://pypi.mirrors.ustc.edu.cn/simple/Looking in indexes: https://pypi.mirrors.ustc.edu.cn/simple/
Collecting rioxarray
Downloading https://mirrors.ustc.edu.cn/pypi/packages/e0/e1/c3c95e9f39485539a9dbcc2272ca98d7558a976b3476699aa6a614c50503/rioxarray-0.9.1.tar.gz (47 kB)
[K |████████████████████████████████| 47 kB 900 kB/s eta 0:00:011
[?25h Installing build dependencies ... [?25ldone
[?25h Getting requirements to build wheel ... [?25ldone
[?25h Preparing wheel metadata ... [?25ldone
[?25hRequirement already satisfied: packaging in /opt/conda/lib/python3.7/site-packages (from rioxarray) (20.4)
Requirement already satisfied: pyproj>=2.2 in /opt/conda/lib/python3.7/site-packages (from rioxarray) (3.2.1)
Requirement already satisfied: xarray>=0.17 in /opt/conda/lib/python3.7/site-packages (from rioxarray) (0.19.0)
Requirement already satisfied: rasterio in /opt/conda/lib/python3.7/site-packages (from rioxarray) (1.1.0)
Requirement already satisfied: pyparsing>=2.0.2 in /opt/conda/lib/python3.7/site-packages (from packaging->rioxarray) (2.4.7)
Requirement already satisfied: six in /opt/conda/lib/python3.7/site-packages (from packaging->rioxarray) (1.15.0)
Requirement already satisfied: certifi in /opt/conda/lib/python3.7/site-packages (from pyproj>=2.2->rioxarray) (2021.10.8)
Requirement already satisfied: setuptools>=40.4 in /opt/conda/lib/python3.7/site-packages (from xarray>=0.17->rioxarray) (49.3.0.post20200809)
Requirement already satisfied: pandas>=1.0 in /opt/conda/lib/python3.7/site-packages (from xarray>=0.17->rioxarray) (1.3.4)
Requirement already satisfied: numpy>=1.17 in /opt/conda/lib/python3.7/site-packages (from xarray>=0.17->rioxarray) (1.21.2)
Requirement already satisfied: affine in /opt/conda/lib/python3.7/site-packages (from rasterio->rioxarray) (2.3.0)
Requirement already satisfied: attrs in /opt/conda/lib/python3.7/site-packages (from rasterio->rioxarray) (19.3.0)
Requirement already satisfied: click<8,>=4.0 in /opt/conda/lib/python3.7/site-packages (from rasterio->rioxarray) (7.1.2)
Requirement already satisfied: cligj>=0.5 in /opt/conda/lib/python3.7/site-packages (from rasterio->rioxarray) (0.5.0)
Requirement already satisfied: snuggs>=1.4.1 in /opt/conda/lib/python3.7/site-packages (from rasterio->rioxarray) (1.4.7)
Requirement already satisfied: click-plugins in /opt/conda/lib/python3.7/site-packages (from rasterio->rioxarray) (1.1.1)
Requirement already satisfied: pytz>=2017.3 in /opt/conda/lib/python3.7/site-packages (from pandas>=1.0->xarray>=0.17->rioxarray) (2020.1)
Requirement already satisfied: python-dateutil>=2.7.3 in /opt/conda/lib/python3.7/site-packages (from pandas>=1.0->xarray>=0.17->rioxarray) (2.8.1)
Building wheels for collected packages: rioxarray
Building wheel for rioxarray (PEP 517) ... [?25ldone
[?25h Created wheel for rioxarray: filename=rioxarray-0.9.1-py3-none-any.whl size=54579 sha256=173d9d1df43e98c0008f8dfea7f21e1f9311ec3a008d842140d2293f646df976
Stored in directory: /home/mw/.cache/pip/wheels/41/08/3b/933ebe31dd607b64cbca07b137f8e666f553200fbf7219f1f2
Successfully built rioxarray
Installing collected packages: rioxarray
Successfully installed rioxarray-0.9.1简易代码
import rioxarray
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
# 读取GeoTIFF文件
ds = rioxarray.open_rasterio('/home/mw/input/landuse2841/GLASS-GLC/GLASS-GLC_7classes_2015.tif')
# 创建一个使用PlateCarree坐标系统的图形
plt.figure(figsize=(10, 5))
ax = plt.axes(projection=ccrs.PlateCarree())
# 绘制数据
ds.plot(ax=ax, transform=ccrs.PlateCarree(), cmap='tab20c')
# 添加地理特征
ax.add_feature(cfeature.COASTLINE)
ax.gridlines(draw_labels=True)
# 显示图表
plt.show()Warning: ecCodes 2.21.0 or higher is recommended. You are running version 2.14.1
/opt/conda/lib/python3.7/site-packages/cartopy/io/__init__.py:260: DownloadWarning: Downloading: https://naciscdn.org/naturalearth/110m/physical/ne_110m_coastline.zip
warnings.warn('Downloading: {}'.format(url), DownloadWarning)细化绘图
增加指北针，图例等
# from csdn https://blog.csdn.net/qq_32832803/article/details/110910540
#-----------函数：添加指北针--------------
def add_north(ax, labelsize=18, loc_x=0.88, loc_y=0.89, width=0.04, height=0.13, pad=0.14):
"""
画一个比例尺带'N'文字注释
主要参数如下
:param ax: 要画的坐标区域 Axes实例 plt.gca()获取即可
:param labelsize: 显示'N'文字的大小
:param loc_x: 以文字下部为中心的占整个ax横向比例
:param loc_y: 以文字下部为中心的占整个ax纵向比例
:param width: 指南针占ax比例宽度
:param height: 指南针占ax比例高度
:param pad: 文字符号占ax比例间隙
:return: None
"""
minx, maxx = ax.get_xlim()
miny, maxy = ax.get_ylim()
ylen = maxy - miny
xlen = maxx - minx
left = [minx + xlen*(loc_x - width*.5), miny + ylen*(loc_y - pad)]
right = [minx + xlen*(loc_x + width*.5), miny + ylen*(loc_y - pad)]
top = [minx + xlen*loc_x, miny + ylen*(loc_y - pad + height)]
center = [minx + xlen*loc_x, left[1] + (top[1] - left[1])*.4]
triangle = mpatches.Polygon([left, top, right, center], color='k')
ax.text(s='N',
x=minx + xlen*loc_x,
y=miny + ylen*(loc_y - pad + height),
fontsize=labelsize,
horizontalalignment='center',
verticalalignment='bottom')
ax.add_patch(triangle)import rioxarray
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
from cartopy.mpl.geoaxes import GeoAxes
from cartopy.mpl.gridliner import LongitudeFormatter, LatitudeFormatter
from cartopy.io.shapereader import Reader
from cartopy.feature import ShapelyFeature
import numpy as np
from meteva.base.tool.plot_tools import add_china_map_2basemap
import matplotlib.patches as mpatches
# 读取GeoTIFF文件
ds = rioxarray.open_rasterio('/home/mw/input/landuse2841/GLASS-GLC/GLASS-GLC_7classes_2015.tif')
# 创建一个使用PlateCarree坐标系统的图形
plt.figure(figsize=(16, 9))
ax = plt.axes(projection=ccrs.PlateCarree())
# 绘制数据
im = ds.plot(ax=ax, transform=ccrs.PlateCarree(), cmap='tab20c', add_colorbar=False)
# 添加地理特征，如国家边界、海岸线等
ax.add_feature(cfeature.COASTLINE)
add_china_map_2basemap(ax, name="river", edgecolor='k', lw=0.5, encoding='gbk', grid0=None) # "河流"
add_china_map_2basemap(ax, name="nation",edgecolor='k', lw=0.5, encoding='gbk', grid0=None) # "国界"
# 设置网格线样式
gl = ax.gridlines(draw_labels=True, linewidth=1, color='gray', alpha=0.5, linestyle='--')
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
# 添加颜色条
cbar = plt.colorbar(im, ax=ax, orientation='vertical', shrink=0.8)
cbar.set_label('Land Cover Class')
# 设置图表标题
plt.title('GLASS Land Cover Map for 2015')
# 添加指北针
add_north(ax)
# 创建图例
classes = ds.values.flatten()
unique_classes, counts = np.unique(classes, return_counts=True)
class_names = ['Class {}'.format(i) for i in unique_classes]
colors = [im.cmap(im.norm(value)) for value in unique_classes]
# 创建图例
for i, (class_name, color) in enumerate(zip(class_names, colors)):
ax.scatter([], [], s=100, color=color, label=class_name, marker='s')
# 添加图例到图表
plt.legend(scatterpoints=1, frameon=False, loc='lower left', bbox_to_anchor=(0, 0))
# 显示图表
plt.show()补充
请注意，图例中的类别名称默认为 "Class X" 形式，小编不太清楚土地覆盖怎么个分类法
如果你有实际的土地覆盖类别名称，你可以替换 class_names 列表中的字符串。
想要增加比例尺的话可参考
https://stackoverflow.com/questions/32333870/how-can-i-show-a-km-ruler-on-a-cartopy-matplotlib-plot
