Matplotlib地形图绘制：山体阴影与光照渲染实战指南

matplotlib绘制地形图——山体阴影和光照我一直觉得在数据可视化这个领域里matplotlib是被低估得太厉害的一个库。很多人一提到地形图、山体阴影、光照渲染这些词第一反应就是GIS软件或者那些重型的可视化框架。但实际工作里你总会遇到这样的场景手里只有一份高程数据想快速出图看看地形趋势或者要在论文/报告里放一张干净、可控、风格统一的局部地形示意图这时候再跑去折腾专业GIS软件未免太笨重而且还容易在配色和风格上跟文章整体脱节。matplotlib其实完全能接住这类需求而且效果不差关键看你懂不懂得用它的三维曲面和光照着色机制。这篇文章我想完整聊聊怎么用matplotlib把一份普通的高程栅格数据变成一张带有山体阴影和光照方向感的地形图。核心会放在光照模拟的原理、Normalize与LightSource的配合方式、以及如何避开那些会让图废掉的坑。所涉及的方法不依赖任何额外地理库纯用matplotlib加numpy就能跑通非常适合需要快速出图或者做深度定制渲染的场景。1. 为什么用matplotlib画地形图场景优势与固有边界先说清楚它适合什么不适合什么免得你做到一半发现方向错了回来骂我。matplotlib在处理地形图这块有两个天然优势一是它跟数据分析生态无缝衔接你从GeoTIFF或者NetCDF里读出来的数据经过numpy处理后可以直接喂给绘图接口不用做任何格式转换二是它的每一个视觉元素都是可编程控制的这点等同于给了你一张白纸不像某些封装过深的GIS库想调整一个光照角度还要翻半天文档。现实里我有几次用matplotlib画地形图的经历基本都是为了快速验证数据质量。比如拿到一个DEM数字高程模型文件如果只打印数值范围你很难判断这个区域是山区还是平原过渡带或者有没有明显的断裂线。这时候直接画一张带山体阴影的3D地形图所有的地形趋势一下子就在眼前展开了。这种目视检查的效率是任何统计学指标都比不了的。当然它也有明显的边界必须提前说清楚。matplotlib的三维渲染本质上是OpenGL式的简易投影极度复杂的模型面片多了会卡顿光线追踪级别的真实光影效果它是做不了的超大范围的高精度地形渲染性能也不够好。但凡是大于几百米、小于几十公里的区域性地形示意它的出图质量完全够用。换句话说如果你需要的是那种带有科学论文质感的地形示意图而不是实景级游戏画面那matplotlib在实用性和可控性上是非常能打的。我个人的评判标准很简单数据量在千万级像素以下、对光影效果要求是看起来有立体感就行、图面需要跟论文其他子图保持同一种风格——这三条满足任意两条用matplotlib就是最优解。反之如果你需要渲染整条山脉的高精度模型或者想做那种可以360度旋转查看的交互式地形场景那就该去用专业地理可视化工具了。2. 环境准备与数据获取从高程栅格到numpy数组这个环节看起来基础但其实80%的初学者都倒在第一步。很多人问为什么自己plot_surface出来的图是平的为什么纹理跟预期差那么多追根溯源都是数据本身没处理好。所以这里我仔细捋一遍从原始数据到绘图数组的完整链路。2.1 处理高程数据的常用途径如果你的高程数据是GeoTIFF格式那可以用rasterio或者GDAL来读取但出于轻量化的考虑我一般更倾向于用rasterio转换成numpy数组后直接交给matplotlib。如果手头没有现成的DEM数据也有公开的数据集可以下载比如SRTM或者ASTER GDEM它们的分辨率和精度对于区域性地形示意来说完全够用。不过为了让这篇内容专注在绘图和光照本身我建议用一段可以任意生成测试高程数组的代码来做示例。这样你不用跑出去下数据也能完整复现所有效果。等代码流程跑通了再换真实数据自然水到渠成。生成模拟地形的原理并不复杂多个不同频率、不同振幅的正弦波叠加在一起就可以模拟出类似自然山脊线和山谷走向的起伏特征。import numpy as np import matplotlib.pyplot as plt from matplotlib import cm from matplotlib.colors import LightSource from matplotlib.ticker import LinearLocator # 第一步生成模拟地形数据 np.random.seed(42) x np.linspace(0, 10, 200) y np.linspace(0, 10, 200) X, Y np.meshgrid(x, y) Z (3 * np.sin(X * 0.8) * np.cos(Y * 0.6) 1.5 * np.sin(X * 2.5 1.2) * np.cos(Y * 1.8) 0.8 * np.sin(X * 5.0) * np.sin(Y * 4.0) 0.4 * np.random.randn(*X.shape))真正重要的是下一步——对Z值做预处理。我这里生成了多组不同频率的正弦波叠加高频成分模拟了局部地形的粗糙感低频成分控制了山脊走向。实际用真实DEM数据的时候你的Z数组可能是几千乘几千的分辨率这个流程同样成立。2.2 高程数据必须是二维规则网格这一点必须反复强调。matplotlib的plot_surface要求X、Y、Z三个数组必须具有相同的二维形状分别表示网格每个节点的x坐标、y坐标和高程值。你从文本CSV文件读数据的时候经常遇到的坑是文件里只存了经度、纬度、高程三列长度是N行看着好像没问题但plot_surface会直接报错或者画出完全扭曲的图。解决办法是把这些散点数据栅格化。最简单的方式是判断经纬度是否是规则的间隔如果是直接通过reshape就能恢复成二维数组如果是不规则采样就必须先插值了。插值的工具有很多pandas的pivot_table是一个便捷选择scipy的griddata函数更通用。# 从散点数据转规则的二维网格示意 import pandas as pd # 假如df有三列: lon, lat, elevation # df_pivot df.pivot_table(indexlat, columnslon, valueselevation) # Z df_pivot.values # X df_pivot.columns.values # Y df_pivot.index.values做这一步的时候有两个坑我踩过说给你听第一pivot_table之后要检查有没有NaN产生这个通常是原始数据的采样点分布不均匀导致的可以用邻近值填充或者插值解决第二要确认行和列的排序是否正确如果经纬度坐标没有从小到大排列画出来的地形图往往是倒置或者扭曲的检查方式很粗暴但有效——print一下X和Y的前几个值看看是不是递增的。2.3 归一化的玄机直接画与着色之间隔着一层Normalize先把Z值直接画出来看看效果你会发现颜色分布可能全部集中在某一个色段。这其实不是数据的问题而是颜色映射的默认区间没有匹配上数据的实际范围。Normalize对象的存在就是解决这个问题的它的作用是把数据值线性映射到[0,1]区间内再传递给colormap。# 不指定normmatplotlib会自动使用数据的最小值和最大值 # 但当你需要多张子图对比的时候必须手动指定统一的norm from matplotlib.colors import Normalize norm Normalize(vminnp.percentile(Z, 2), vmaxnp.percentile(Z, 98))这里故意用了百分位数而不是min和max原因是真实地形数据里往往有异常值比如某个像素的高程因为传感器误差出现一个离谱的数字如果用min和max作为上下界整个色带的分布就会失衡颜色层次全部被压缩到一个非常窄的区间里等于白画。按2%到98%截断是个稳妥的默认选择视觉上损失的信息很少但能有效抵抗离群值的干扰。顺便说一句这招不仅用于地形图画任何带有长尾分布的数据图都通用。3. 三维地形曲面的基座plot_surface核心参数逐项拆解三维曲面是地形图呈现的骨架这次我们先从最基础的plot_surface入手把它的参数彻底吃透后面再叠加上光照和阴影效果。这里面有几个参数的作用我花了很多时间才真正理解而不是只看官方文档的描述。3.1 rstride与cstride控制网格密度的前后台逻辑在matplotlib的三维曲面中X和Y的网格点数量直接决定了曲面的精细度。如果你的数据是200x200的网格全部画出来意味着4万个点matplotlib需要在这个基础上构建三角剖分计算量不小而且对于很多地形图应用视觉上看不出太大差别。这时rstride和cstride参数就派上用场了——它们的含义是每隔多少个点取一个采样点进行绘制。比如rstride2, cstride2实际绘制网格就从原来的200x200变成了100x100计算量降低了四倍画出来的地形轮廓基本不受影响。但如果你把这两个值调得太大比如大于等于5地形细节就会丢失山脊线变得模糊阴影效果也会跟着失真。我这里建议一个经验值绘图数据的长和宽不超过300的时候直接用默认参数等大于0时自动优化完全没问题如果超过500就显式设置rstride和cstride控制在300以内可以保证交互式操作的流畅度。官方文档里说rstride和cstride在未来的版本可能会弃用但直到我写这篇文章的当前版本还在正常工作所以你完全可以放心用只是知道有这么回事就好等哪天真被移除了再迁移也不迟。3.2 cmap与facecolors颜色的平滑过渡与逐点着色cmap参数控制的是曲面的整体颜色映射方案而facecolors参数可以让你对每个面片单独指定颜色。在地形渲染里这两者结合可以产生一个非常实用的效果用facecolors叠加光照阴影卷积同时保留cmap对应的高程颜色信息。这个做法的实质是先根据高程值生成一份基础颜色然后把光照强度矩阵作为透明度混合到面片上。如果用colormap直接上色你会得到一个虽然立体感尚可但缺少阴影层次的图。但如果把光照强度单独乘进去整张图的立体感会瞬间提升一个档次山体看起来是鼓起来的而不是一张平平的彩色贴图。fig, ax plt.subplots(subplot_kwdict(projection3d)) surf ax.plot_surface(X, Y, Z, cmapcm.terrain, linewidth0, antialiasedTrue)这是最基本的一张3D地形图。单看这张图你可能会觉得它的立体感还不错山脊山谷的形状大致能辨认出来但跟加上山体阴影之后的效果一对比差距就很明显了。原因在于matplotlib默认的光照方向是固定的而且它并不会根据地形起伏产生局部阴影你看到的高程颜色只是数值映射的结果并没有真实的光影信息。这也是为什么专业的地形图渲染都要单独做山体阴影叠加我们下面就来解决这个关键问题。3.3 视角设置ax.view_init的方位角与仰角如何影响地形辨识度视角的选择对地形图的判读影响非常大这一点经常被忽略。同样的数据如果你从正上方俯视看到的可能只是一团模糊的色块换一个斜视角山地的层次感立刻出来了。ax.view_init(elev, azim)控制的就是这个视角参数——elev是视线与z轴的夹角仰角azim是视线在xy平面上的方位角。在实践里我总结出了两个比较可靠的推荐值组合。第一种是区域总览图适合体现整体地形大势仰角30度左右方位角在-60到-120度之间根据地势走向微调这种视角下山脊线和谷地的层次最明显第二种是局部细节图适合突出特定山体结构仰角要降到15-25度方位角要顺着主山脊的走向来设置这样阳光在侧面的投影会拉长阴影增强地形特征的可见度。ax.view_init(elev30, azim-60)一个很反直觉但真实的经验是把仰角设成90度即完全俯视得到的并不是一张好用的地形图因为在高程色带强烈对比下你会看到大量杂色的马赛克根本分不清山体结构。真正有经验的人都是把仰角降低让观察角度接近真实地形扫描的感觉立体感瞬间就出来了。4. 山体阴影的底层原理从坡度坡向到光照强度矩阵这里介绍的内容是整个流程的核心。很多人直接调用LightSource然后加shadeTrue就算完事但不知道为什么参数是那样设置的碰到效果不好时也不知道从哪里调。真正理解了光照强度的计算方式你才能举一反三做出自己想要的光影风格。4.1 什么是山体阴影从太阳方位角到局部坡度坡向的数学关系山体阴影Hillshade是一个专业的GIS术语它模拟的是太阳光从某个方向照射到地形表面时每个像素受到的光照强度。这个值取决于三个因素太阳的高度角简称太阳高度、太阳的方位角、以及地面该点的坡度和坡向。用一个简单日常生活化的类比来说你拿着手电筒从侧面照一个揉皱的纸团纸团的每个褶皱因为朝向不同亮度也不同。褶皱朝向光源的部分就亮背向光源的部分就暗。山体阴影就是这样只不过光源变成了设定好的太阳位置褶皱变成了地形网格上的无数个小面片。具体到数学上光照强度跟两个矢量有关一个是光线的方向向量一个是地形面片的法向量。当光照方向和面片法向量方向一致即太阳直射坡面时亮度最高当两者垂直时亮度为零当面片朝向完全背对光源时就处于阴影中。计算每个面片的法向量的方法是对地形曲面在x和y方向分别求偏导数然后做叉积。在matplotlib中这个计算过程已经封装在LightSource类里了你不用自己手动去算偏导数和叉积但理解了底层逻辑对你调参非常有帮助。因为光照方向的选择、光照强度的缩放因子、颜色混合的模式这些参数都是围绕这个数学关系展开的。4.2 matplotlib的LightSource核心机制从Hillshade到blendmatplotlib的LightSource类在较新版本里做了一次比较重要的升级增加了对真实光照模型的支持。在旧版本里shade函数只是根据高程值计算光照强度然后用灰度图直接替换原图效果虽然能用但总感觉生硬。新版本的LightSource提供了两种模式。默认的hsv模式是把光照信息编码到HSV颜色空间的V分量明度上然后跟原始颜色混合这样既保持了原先的色调又增加了明暗变化。absolute模式则是一种简化的光照计算适合不需要考虑光照方向的分析场景。对于绝大多数地形图来说hsv模式是最合适的因为它能够在保持高程颜色的同时叠加具有方向性的明暗变化。如果要手动控制混合程度shade函数里还有一个blend_mode参数它有三个可选值soft、hard、overlay。soft模式是平滑混合适合保留颜色层次但阴影感会弱一点hard模式是直接替换明度阴影感最强但颜色层次会损失比较明显overlay模式借鉴了图像处理中的叠加混合既能保留颜色细节又能体现出较强的立体感。我的习惯是如果想突出地形解剖学特征就用overlay如果想做一张柔和的背景底图就用soft。4.3 自己动手计算坡度和坡向一份可复制的numpy实现虽然matplotlib封装了LightSource理解它的内部计算逻辑仍然能帮助你解决很多实际问题。坡度的数学定义是地形曲面在水平方向上的最大变化率坡向则是这个最大变化率指向的方向。计算坡度和坡向需要用到数值差分在规则网格上可以用Central difference方法来实现用numpy的gradient函数即可。# 手动计算光照强度矩阵的示例hillshade核心逻辑 dx, dy np.gradient(Z, x[1]-x[0], y[1]-y[0]) slope np.pi/2 - np.arctan(np.sqrt(dx*dx dy*dy)) aspect np.arctan2(-dx, dy) azimuth np.radians(315) # 西北方向光源示例 altitude np.radians(45) shaded np.sin(altitude) * np.sin(slope) np.cos(altitude) * np.cos(slope) * np.cos(azimuth - aspect) shaded np.clip(shaded, 0, 1)这段代码的输出是一个范围在[0,1]之间的光照强度矩阵。0表示完全在阴影中1表示受到完全直射。拿到这个矩阵之后你可以自己决定怎么跟颜色混合而不必拘泥于matplotlib预设的方式。我经常把这个矩阵直接作为透明度通道或乘性遮罩用matplotlib的alpha或者cmap技巧来实现自定义渲染。这个自由度是直接用一句shadeTrue得不到的。5. 完整绘图流程实战把光照叠加到三维地形上原理讲清楚了下面进入实战环节。我会一步步给出从数据读取到最终出图的完整代码并对每处关键参数做解释。这份代码你可以直接复制去跑也可以根据自己的数据调整参数。5.1 使用LightSource进行山体阴影渲染的完整示例# 初始化光照模型设定太阳方位角和高度角 ls LightSource(azdeg315, altdeg45) # 方式一直接绘制带光照的3D曲面 fig, ax plt.subplots(figsize(10, 8), subplot_kwdict(projection3d)) rgb ls.shade(Z, cmapcm.terrain, blend_modesoft, vert_exag2) surf ax.plot_surface(X, Y, Z, facecolorsrgb, linewidth0, antialiasedTrue, rstride1, cstride1) ax.view_init(elev25, azim-60) plt.show()这里的关键逻辑是plot_surface函数里的facecolors参数接收一个数组数组的每个元素是一个RGBA颜色值对应网格中的每个面片。我们先用ls.shade计算出带有山体阴影效果的RGB颜色然后把它传给facecolors。值得注意的是vert_exag参数它的作用是垂直方向的夸张系数。真实地形的高程振幅往往远小于水平范围如果直接绘图山体看起来会非常平。提高vert_exag能让山体在视觉上更加显著。对于模拟数据vert_exag一般2到3就足够了对于真实DEM数据需要根据你的水平和垂直单位比例来估算。如果单位不一致这个值可能需要调到几十甚至上百才能看出地形起伏。5.2 两种出图风格对比直观高程色与山体阴影叠色为了让你更直观地看到山体阴影到底带来了什么提升我把两种出图风格放在一起做一次对比。第一种是直接plot_surface加cmap第二种是叠加了LightSource的shade效果。同时从四个不同视角各拍一张这样你能充分感受到视角和光照之间微妙的互动关系。风格立体感颜色层次地形结构辨识度推荐场景纯cmap着色中等丰富中等快速查看高程分布叠加光照shade强烈较丰富高正式汇报、论文插图在实际项目里如果是自己调试数据我建议先用纯cmap模式快速看一遍等找到感兴趣的局部区域了再切换到叠加光照模式进行调整和出图。这能节省不少时间因为shade模式下参数稍微一变就得重新渲染效率并不高。5.3 调整光照方位和强度做出清晨与正午的地形质感光照方向的改变不仅影响美观还直接影响地形信息的呈现。在LightSource中azdeg控制方位角altdeg控制高度角。如果把方位角设在315度西北高度角设在30度得到的光影效果会让人联想到清晨或傍晚的阳光阴影很长能突出微小的地形起伏如果把altdeg提高到70度就像正午的太阳阴影很短地形看起来比较平。这里给出我实验过的几个参数组合你直接抄就行。视觉效果azdegaltdeg适用场景清晨/黄昏长阴影28025突出微小地貌标准制图照明31545通用平衡效果高亮正午18070弱化阴影强调色带我的个人倾向是除非有特殊需求否则用315/45这个组合基本不会出错。这个方向的光线在绝大多数地形上都能产生不错的山体层次既不会因为阴影太重掩盖颜色信息也不会因为阴影太弱显得平淡。5.4 垂直夸张系数对最终视觉的重大影响前面提到vert_exag参数但它值得单独拿出来说。因为在实际渲染中它是决定图上地形够不够立体的最关键因素比光照方位的影响还要大。如果vert_exag设得太小山体跟平地差不多无论光照怎么打都白搭如果设得太大峰谷之间的过渡显得非常突兀像是夸张的微型沙盘模型。我给的调试建议是先用一个较小的值比如1渲染如果看起来太平坦就逐步翻倍。每次改完值以后最好固定视角和光照参数只调这一个参数这样你能清楚知道这个参数对画面效果的单独贡献。当阴影边缘开始出现明显的锐利转折时说明已经接近上限了往回调一点就是合适的值。6. 把二维Hillshade叠加到3D曲面轻松实现高级地形特效如果上面这些步骤你都跑通了我们还可以再进一步。把山体阴影从三维曲面的贴图变成二维平面的叠加底图再贴回三维曲面这种做法的效果比直接贴要细腻得多因为你可以分别控制光照纹理的清晰度和曲面的平滑度两者之间不会互相干扰。6.1 在3D曲面上叠加纹理的底层操作思路是这样的先用LightSource生成一张带光照效果的RGB图像二维然后把这个图像作为facecolors传回plot_surface而plot_surface本身并不关心facecolors里的信息具体是怎么算出来的——它只负责把颜色贴上去。这就给了你一个非常灵活的接口你可以在把图像传给plot_surface之前对图像做任意处理比如锐化、去噪、对比度增强甚至叠加其他地理图层。ls LightSource(azdeg315, altdeg45) cmap cm.terrain # 生成二维光照图 rgb2d ls.shade(Z, cmapcmap, blend_modeoverlay, vert_exag3) # 直接作为3D曲面的facecolors fig plt.figure(figsize(12, 8)) ax fig.add_subplot(111, projection3d) ax.plot_surface(X, Y, Z, facecolorsrgb2d, linewidth0, antialiasedTrue) ax.view_init(elev20, azim-70) plt.show()这种方法的自由度非常高。比如你可以把光照图像的分辨率做得比曲面网格高很多再使用高分辨率的光照图配合较稀疏的曲面网格这样曲面的渲染速度快但地形明暗变化依然细腻。如果你的数据量很大这一步优化能显著降低渲染卡顿。6.2 纹理叠加与直接shade的效果差异对比从最终视觉效果看用二维光照图再贴回三维曲面得到的图片色彩层次更丰富过渡更平滑。因为你在二维环境里可以用更高精度的算法生成RGB图这个精度是三维渲染时无法直接达到的。举个例子你在二维里可以先用高斯滤波把光照强度矩阵模糊一下消除高频噪声然后再映射成颜色这样得到的阴影就很柔和山体表面有一种雾面感。直接对三维曲面做shade的时候没有这个中间处理步骤高频起伏会被渲染成尖锐的锯齿。两种方法没有绝对的好与坏。直接shade简单快捷适合内部快速验证二维纹理叠加则适合最终出图能雕琢的细节多得多。如果你准备把图放进论文或者项目报告里我强烈建议你用二维纹理叠加这条路。6.3 隐藏曲面网格线的三个细节让地形图拥有DEM质感一张高质量的地形图网格线是几乎不可见的。很多初学者的图一看就是matplotlib默认风格就是因为曲面上的黑色网格线没有处理干净。关键在三个参数上linewidth要设成0antialiased要设成True此外还需要把rstride和cstride设置得足够小确保面片足够细密线条间隔自然就消失了。surf ax.plot_surface(X, Y, Z, facecolorsrgb2d, linewidth0, antialiasedTrue, rstride1, cstride1)如果你用了rstride1, cstride1把每个网格节点都画进去而数据分辨率又很高那facecolors数组的维度也要匹配上。这里有一个常见的报错点rgb2d的shape是(M, N, 3)而plot_surface期望facecolors的shape是(M, N, 3)或(M, N, 4)如果维度不匹配画出来要么报错要么颜色错位。遇到这个问题的时候先用print查看一下数据的shape再往下排查不要一上来就怀疑是data的问题。7. 精通颜色映射让地形高程更清晰可读颜色映射在数据可视化里绝对算得上语义学级别的重要环节。地形图的高程色带选得好不好直接决定读者第一眼能不能读对这个区域的地势走向。matplotlib内置的colormap有很多但真正适合地形展示的其实是少数派。7.1 为什么terrain和viridis是地形绘制的两大镇山之宝matplotlib自带的terrain这个colormap设计初衷就是为了地形数据。它的色带从低海拔的绿色过渡到中海拔的黄色、棕色再到高海拔的白色甚至灰色其色彩逻辑跟传统纸质地形图的用色习惯一脉相承。确切说绿色代表低地植被棕色代表山地裸岩白色代表雪线以上的区域。这种色带的好处是观众不需要额外的图例说明凭直觉就能读懂。另一个选择是viridis它是matplotlib的默认colormap最大的特点是色觉障碍友好同时灰度打印出来依然保持较好的区分度。它的逻辑不是模仿自然地理环境而是纯粹为了数据可视化设计的——颜色从暗紫到亮黄渐变适用于需要严格区分数值高低的场景。跟terrain相比viridis少了那种地理直觉但也能提供很好的细节分辨率。7.2 自定义地形色带三段式配色方法很多情况下内置色带没法精确表达特定区域的特征。比如你要显示一个海拔范围特别小的区域蓝色到绿色占主导的地形你想把海洋和陆地的分界点突出出来。这时自定义色带就是必要的。最简单的自定义办法是用LinearSegmentedColormap从一组颜色断点构建色带。拿三段式举例低海拔用蓝色、中海拔用黄绿色、高海拔用褐色。关键是要设计好断点位置让颜色在数值分布不均匀的数据里依然能合理分配。# 自定义一个三段式地形色带 from matplotlib.colors import LinearSegmentedColormap colors [#2c7bb6, #ffffbf, #d7191c] cmap_custom LinearSegmentedColormap.from_list(terrain_custom, colors, N256)不过要提醒一句如果你在高程上叠加了山体阴影自定义色带的饱和度不要太高否则阴影叠加后容易出现脏色。我一般会先把颜色的饱和度降低个20%再考虑要不要加阴影这样出来的颜色层次更干净更有高级感。7.3 颜色条(colorbar)参数的艺术不只是加个刻度很多人在colorbar上栽过跟头。如果你只是简单调用plt.colorbar(surf)颜色条上那些刻度往往是自动生成的跟你的数据分布不一定匹配字太小、位置不合适都是常见问题。在生产环境或者是论文里这种做法会被编辑直接打回来。更可控的做法是先用mappable对象单独生成colorbar然后再精细调整它。mappable在plot_surface返回的Surf对象里可以直接拿到我们只需要把这个对象传给colorbar函数。surf ax.plot_surface(X, Y, Z, facecolorsrgb2d, linewidth0, antialiasedTrue) cbar fig.colorbar(surf, axax, shrink0.6, aspect15, pad0.1) cbar.set_label(Elevation (m), fontsize12) cbar.ax.tick_params(labelsize10)需要注意的坑是当你用facecolors参数时plot_surface返回的surf对象本身并不携带用于colorbar的颜色映射信息因为颜色已经被固定在facecolors里了。这种时候直接调fig.colorbar(surf)可能报错或者根本不显示。解决办法是单独创建一个ScalarMappable对象把你想让colorbar展示的数值范围和colormap对应起来。from matplotlib.cm import ScalarMappable sm ScalarMappable(cmapcmap, normnorm) sm.set_array(Z) cbar fig.colorbar(sm, axax, shrink0.6, aspect15, pad0.1)这个细节知道的人不多但遇到的时候非常卡人。提前了解等于给自己省了半小时排查时间。8. 实战踩坑记录光照方向错误、Z轴比例失调与颜色过曝任何一篇实操文章如果不把自己踩过的坑说出来那这篇内容的信任度要打一个很大的折扣。下面这些坑都是我实际画地形图的时候一个个趟过来的每条都对应一个具体的症状和排查方法。8.1 光照方向跟视角死锁阴影角度不变怎么回事有一次我在调图改了很多次azdeg但图上阴影的角度怎么都不变。当时我怀疑是LightSource的代码出了问题还专门去翻源码后来才发现问题出在作图流程上——我虽然改了azdeg但忘了重新调用shade函数只是重新执行了plot_surface而plot_surface用的还是之前shade已经生成好的facecolors数组。这个错误很低级但很真实它提醒我们shade函数生成的RGB数组在传入plot_surface之后就已经固化了后续任何光照参数的改动都必须重新执行shade才能生效。另一个可能的原因是你用的matplotlib版本太老旧版本LightSource对azdeg的解析有bug升级一下就能解决。8.2 地形被拍扁Z轴刻度跟实际高程不匹配这个问题几乎每个画过真实DEM数据的人都会遇到一次。原因很简单水平坐标通常是经纬度或者平面投影坐标单位是米或者公里范围动辄几十上百公里而高程Z轴的范围只有几百米到几千米。当matplotlib默认使用相等的三轴刻度时Z轴的变化在整体画面中微乎其微看起来就像一块被拍扁的饼。解决办法有两个思路。一是手动设置Z轴的显示范围让它的上下界比真实数据更宽一些比如把Z轴的范围扩大3倍这样能人为放大垂直方向的变化效果二是调整箱体比例用set_box_aspect来控制三个轴的长度比例。ax.set_box_aspect((1, 1, 0.3))我习惯用第二种方式因为它是纯粹控制显示比例不会对数据本身产生困惑。设置box aspect的时候数值的含义是x轴长度:y轴长度:z轴长度的比值。如果你希望山体更突出就把第三个值调大比如(1, 1, 0.8)如果你只想看一个大致的趋势第三个值可以设置得小一点比如(1, 1, 0.2)。这个思路配合前面的vert_exag是双重保险但两者角度不同——vert_exag是让数据在shade计算时就被放大而set_box_aspect是让显示比例改变不会影响光照计算的形态。8.3 画面过曝或者死黑光照强度与颜色饱和度的平衡这种问题很常见尤其是你从网上找了一段代码改参数时。画出来的图要么亮得跟白纸一样要么暗得完全看不清高程颜色。核心原因通常是altdeg设置得太高或太低或者blend_mode选择得不对。altdeg如果超过60度光照会非常强山体表面几乎失去阴影层次颜色变得一片白茫茫altdeg低于15度光照基本是贴着地面走的大面积的阴影会让地形细节完全隐藏在黑暗里。另一个原因是cmap本身的亮度范围跟光照叠加之后出现了溢出。解决方法是先降低vert_exag再调整altdeg到45度附近。如果还是过曝可以手动设置norm的范围让vmax略低于实际最大值这样色彩的亮度会被压缩一点图像看起来更平衡。9. 把地形图玩出进阶花样添加河流、路径与标注当你掌握了基础的地形渲染之后自然会想在上面添加更多信息。matplotlib的三维坐标轴天然支持叠加其他元素比如代表河流的曲线、代表登山路径的折线、关键山峰的标注点等。这些元素放在地形图上能让你的图从展示地形升级为讲述区域故事。9.1 用PolyCollection和plot在三维地形上叠加矢量要素如果要在三维地形上画一条河流你需要把河流的坐标从二维延伸到三维每个点的x、y坐标直接就是这条河在地图上的位置z坐标则应该取自这个位置的对应高程值。关键是要从Z数组中插值出沿着河流路径的高程。from scipy.interpolate import interp2d # 假设river_x, river_y是河流的平面坐标 f interp2d(x, y, Z, kindlinear) river_z np.array([f(xi, yi)[0] for xi, yi in zip(river_x, river_y)]) ax.plot(river_x, river_y, river_z, color#1e88e5, linewidth2.5)河流要素的视觉风格要注意跟地形色带区分开。我一般用深蓝色或者青色线宽稍微粗一点看起来像卡通地图里那种带描边的河流清晰且不失美观。与之类似的还有路径、行政边界等基本逻辑相同只是颜色和线型不同。9.2 典型应用山地徒步路线的高程剖面图叠加这种图在户外运动类的项目报告里非常好用。展示一条徒步路线的三维走向同时在山体地形图上标出起终点和关键转折点受众一眼就能看出攀爬的难度分布和沿途经过的地貌类型。再加上一个二维的高程剖面图就是从路线起点到终点的海拔变化曲线作为子图放在旁边整张图的信息量就非常完整了。实现的方法是在同一张figure里创建两个子图左边是三维地形叠加路线右边是二维的高程剖面折线图。三维地形的视角最好设置在路线的侧上方让路线跟视角保持一个合适的夹角这样路的起伏才能被清晰地看出。二维剖面图则应该把x轴设为路程累积距离y轴是高程颜色跟三维图保持一致。这样两张图可以互相印证读者看起来不会有割裂感。9.3 添加文字标注时如何防止文字被地形遮挡在三维图上做文字标注是个老难题因为文字默认只贴在xy平面上而地形是有起伏的文字很容易被前面的山体遮住。matplotlib提供了text函数支持x、y、z三个自由度的坐标定位你可以把文字的z坐标抬高到指定等高面的上方避免被遮挡。ax.text(peak_x, peak_y, peak_z 0.3, 主峰, colorblack, fontsize12, hacenter)这里把z坐标加了0.3个高度单位实际值取决于你的数据量级总体原则是让文字底部贴在对应地形的上方同时留出足够的安全距离。如果扁平的文字仍然被遮挡可以用ax.text再配合一个浅色的背景框或者加粗的白色轮廓来增强可读性。10. 从静态图到动态交互视角旋转与保存高分辨率的大图跟静态图相比旋转视角可以让你从不同方向观察同一块地形更全面地理解其结构。matplotlib的交互功能虽然在专业3D引擎面前不算强大但对检测地形数据来说已经足够。10.1 用view_init逐帧渲染生成GIF或视频如果你想把图做成简单的旋转动画原理并不复杂创建一个角度序列在每一帧调用ax.view_init设置不同的视角然后保存每一帧的图片最后用PIL或者imageio合成GIF。import imageio import os angles np.arange(0, 360, 5) frames [] for azim in angles: ax.view_init(elev25, azimazim) fig.canvas.draw() frame np.frombuffer(fig.canvas.tostring_rgb(), dtypenp.uint8) frames.append(frame.reshape(fig.canvas.get_width_height()[::-1] (3,))) imageio.mimsave(terrain_rotation.gif, frames, fps15)需要注意的一点是生成动画之前应该把光照的方向一起设好因为在动画里只旋转观察角度光照方向是保持不变的。这样视觉效果是地形在转动而光线的方向恒定看起来非常真实。反过来如果你希望光照也动起来就需要在每一帧重新生成一次shade那样计算量会更大但效果也更炫。10.2 高清大图输出设置dpi、figsize与bbox_inches的配合论文插图对分辨率要求通常很高因此保存图片时要设置合理的参数。默认保存出来的图片清晰度可能不够我们需要显式地控制dpi和bbox_inches。fig.set_size_inches(12, 8) plt.savefig(terrain_map.png, dpi300, bbox_inchestight, facecolorwhite)dpi设置成300是一个比较稳妥的印刷级分辨率如果计算机性能允许450也没问题。bbox_inchestight的作用是自动裁剪掉多余的空白区域让图中的地形主体最大化地充满画面。但剪裁的这个过程会改变图片的边距比例如果你后续还要加外部图例或者其他元素需要考虑这一点。10.3 性能优化建议大数据量地形图的降采样技巧处理真实的超大DEM数据时比如12000x12000的全球尺度数据直接丢给matplotlib画会卡到怀疑人生。这种情况下时间优化的主要手段是降采样。你可以先对Z数组做切片每隔n个点取一个或者是用scipy.ndimage.zoom来做更平滑的缩放。from scipy.ndimage import zoom factor 0.2 Z_small zoom(Z, factor) X_small zoom(X, factor) Y_small zoom(Y, factor)缩放的时候要注意X和Y也要同步缩放并且factor保持一致否则横纵比例会失衡。另外如果你的数据里有非常细碎的纹理需要保留zoom的order参数可以调整order1是线性插值order3是三次样条插值后者更平滑但计算量更大。一般来说用order1就足够了因为后续还要叠加山体阴影高频细节会影响光照计算的稳定性反而让画面显得噪杂。11. 基于matplotlib地形绘制的扩展思考换个库会不会更好写到这里我觉得有必要摊开来聊几句关于工具选择的思考。matplotlib在地形可视化这条路上并不是唯一的选择甚至在某些指标上并不是最优选择。但正是它的不完美让它在很多务实场景中成为最合适的选择。比如plotly的3D地形图支持WebGL渲染拖动起来极其丝滑而且可以在浏览器里交互查看PyVista则更接近专业可视化工具对大规模网格数据处理能力强渲染效果也更接近电影级。但这两者的问题在于它们的依赖更重跟numpy和pandas的集成不如matplotlib直接在论文排版和样式统一性上也需要更多额外工作。我个人的原则是如果目标是把地形图当作数据分析流程中的一环快速出图看图选matplotlib如果目标是要做一个完整的地理信息系统原型或者是一个面向用户的交互式应用那应该毫不犹豫去选更专业的库。工具没有绝对的高下只有是否匹配你的场景。另外值得提一句的是matplotlib社区这几年一直在改进三维渲染能力包括简化API、优化性能、增加更多光照模型选项等。这是好事说明这个库并没有放弃三维可视化这个领域。虽然跟专业3D库仍有很大差距但对我们这种主要做二维数据可视化的用户来说它已经多了一层应急与轻量的价值。也许再过几个大版本这些差距会继续缩小。到那个时候用matplotlib画地形图就更没有心理负担了。12. 几个能直接提升成图质量的小习惯最后分享一些我在长期使用中养成的习惯。这些点单个看都很小但合起来能让整张图的质感明显提升一个台阶而且几乎不需要额外的时间成本。第一始终用set_box_aspect去调整三轴比例不要让matplotlib自动决定。自动决定的结果通常是三轴比例相等这在展示地形时非常吃亏。手动调整为(1, 1, 0.4)左右地形的立体感马上就不一样了。第二用坐标轴标签和刻度字体来暗示图面风格。字体大小统一设成10到12磅之间刻度标签的朝向也统一。如果字体太杂会给人一种拼贴感比较出戏。第三光线方位跟主视角方位错开45度以上这是保证立体感的基本规律。如果光线跟视线完全一致就类似逆光自拍地形细节会被强烈压缩什么纹理都看不清。第四保存图片之前不要怕多试几种cmap。terrain、gist_earth、cividis都是可以用候选每种色带呈现地形的感觉差异很大。多切换几次挑出最能体现你数据特征的那个然后再调光照和视角不要反过来先调完光照再换cmap那样会让你之前调好的参数前功尽弃。第五如果你要把地形图嵌入到更大的图版里面记得在早期就设定好统一的dpi和字体风格。等所有子图画完以后再统一处理往往需要返工浪费的时间比自己预想的多得多。matplotlib画地形图这件事说难也确实有一点门槛主要的难度在于光照模型的理解和参数之间的协同配合说简单也简单一旦把流程跑通后面换数据、换视角、换颜色都只是改几个参数的事。希望这篇内容能帮你把这条路的第一个门槛跨过去。你拿到的数据值得一张足够清晰、足够专业的地形图。