经常看文献的小伙伴们一定会发现,在顶级期刊中总有一些全球地图,将矢量/栅格数据进行了精美的可视化。每当看到这些图片,都不禁感叹如果自己也会画就好了,本期小编就为大家带来此类全球要素分布图的绘制方法。
01
图片展示
在正式教学之前,我们先欣赏一下Nature中的案例。下图展示了最长年度干旱期(LAD)变化的全球分布:配色上彰显主题但不喧宾夺主,结构上巧妙地搭配纬向LAD变化的平均值和标准差,各司其职又浑然一体,简直是科技与美的结合。
图源:https://www.nature.com/articles/s41586-024-07887-y
02
绘制方法
我们以ERA5下载的2020年1月2 m高度的温度数据为例。首先,我们需要准备全球矢量数据和全球边界数据(其中全球边界数据可以自己绘制矩形)。
其次,根据我们需要表达的信息设置填充色并调整符号标注。为了体现示例数据的空间分布,我们选择纯色梯度渐变色卡。
第三步是为数据框设置投影坐标系,常用的坐标系一般有两种选择,分别为“Sphere_Equal_Earth_Greenwich(Equal_Earth)”与“Sphere_Aitoff(Aitoff)”。
最后,以“Sphere_Equal_Earth_Greenwich”为例,我们设置了经纬线并增加图例。
Python代码如下,代码中数据同样下载自ERA5,需要数据及相关文件的小伙伴可联系作者索要。
#加载包
import rasterio
import geopandas as gpd
import matplotlib.pyplot as plt
from rasterio.warp import calculate_default_transform, reproject, Resampling
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import numpy as np
from matplotlib.font_manager import FontProperties
# 文件路径
raster_path = "C:/Users/ASUS/Desktop/t2m/t2m_202001.tif"
shapefile_path = "C:/Users/ASUS/Desktop/t2m/World_countries.shp"
# 加载 shapefile
world = gpd.read_file(shapefile_path)
# 设置 Equal Earth 投影
equal_earth_crs = ccrs.EqualEarth()
# 打开栅格文件
with rasterio.open(raster_path) as src:
# 计算变换并重投影
transform, width, height = calculate_default_transform(
src.crs, equal_earth_crs, src.width, src.height, *src.bounds)
kwargs = src.meta.copy()
kwargs.update({
'crs': equal_earth_crs,
'transform': transform,
'width': width,
'height': height
})
with rasterio.open('reprojected.tif', 'w', **kwargs) as dst:
for i in range(1, src.count + 1):
reproject(
source=rasterio.band(src, i),
destination=rasterio.band(dst, i),
src_transform=src.transform,
src_crs=src.crs,
dst_transform=transform,
dst_crs=equal_earth_crs,
resampling=Resampling.nearest)
# 绘制
fig, ax = plt.subplots(figsize=(12, 8), subplot_kw={'projection': equal_earth_crs})
ax.set_global()
# 打开重投影后的栅格文件并使用 imshow 绘制
with rasterio.open('reprojected.tif') as src:
data = src.read(1)
extent = [src.bounds.left, src.bounds.right, src.bounds.bottom, src.bounds.top]
img = ax.imshow(data, cmap='RdPu', extent=extent, transform=equal_earth_crs, vmin=220, vmax=320)
# 添加 shapefile
world.boundary.plot(ax=ax, transform=ccrs.PlateCarree(), linewidth=0.5, edgecolor='black')
# 设置字体以支持中文
font = FontProperties(fname=r'C:\Windows\Fonts\simsun.ttc', size=14)
# 样式
ax.coastlines()
ax.gridlines(draw_labels=True, linewidth=0.5, color='gray', alpha=0.5)
ax.set_title('2020年1月2m高度温度变化', fontproperties=font, fontsize=14)
# 添加自定义特征
ax.add_feature(cfeature.BORDERS, linestyle=':', edgecolor='black')
# 添加色带
cbar = plt.colorbar(img, ax=ax, orientation='vertical', shrink=0.6, pad=0.05)
cbar.set_label('温度 (K)', fontproperties=font, fontsize=12)
cbar.set_ticks(np.linspace(220, 320, 11))
cbar.set_ticklabels([str(int(tick)) for tick in np.linspace(220, 320, 11)])
plt.show()
以上便是本期图绘的诚意分享啦,Python偏运动,Arcgis偏商务,希望大家能找到适合自己的方法画出自己满意的科研绘图,欢迎大家持续关注我们!
一图胜千言!水文图绘改版后致力于分享水文相关的精美图表,为读者提供作图思路和经验,帮助大家制作更漂亮丰富的图表。同时欢迎留言咨询绘图难点,我们会针对性地分享相关绘制经验。另外也期待读者踊跃来稿,分享更好的构图思维和技巧,稿件可发送至邮箱hydro90@126.com, 或者联系微信17339888901投稿。
编辑: 张恒飞 于明业 |校稿:Hydro90编委团