python rasterio

# 在Python里摆弄卫星照片：聊聊rasterio这个库

如果你曾经对着一张卫星地图或者遥感影像发呆，好奇这些数据在代码里究竟长什么样，那么rasterio这个库或许能帮你打开一扇窗。它不是什么新潮的框架，但在处理地理空间栅格数据这个特定领域里，算得上是个得心应手的工具。

他是什么

简单来说，rasterio是Python里专门用来读写地理栅格数据文件的库。所谓栅格数据，你可以把它想象成一张特殊的图片，但每个像素点除了颜色信息，还绑定了具体的地理位置坐标。常见的.tif格式的卫星影像、数字高程模型（就是那种表示地形起伏的数据），大多属于这个范畴。

这个库的底层其实是C语言写的GDAL库，但rasterio用Python的方式把它包装了一遍，接口更符合Python开发者的习惯。你不用去记那些冗长的GDAL命令行参数，也不用和复杂的内存管理打交道，用几行Python代码就能把数据读进来，像操作NumPy数组一样去处理它。

他能做什么

rasterio干的事情很专注，主要就是三件：读、写、处理地理栅格数据。

比如你手头有一份某地区的多光谱卫星影像，想看看近红外波段的情况，用rasterio可以轻松提取出这个波段，进行运算。或者你想把几张相邻区域的影像拼接成一张完整的大图，它也能帮上忙。再比如，你需要把处理好的数据，按照原有的地理坐标信息，重新写成一个新的GeoTIFF文件，交给GIS软件去打开，这个过程用rasterio来实现也非常顺畅。

它更像是一个数据搬运工和格式转换器，把复杂的、带有地理外衣的数据，变成程序员熟悉的数组，处理完再原样包装回去。至于更复杂的空间分析、高级渲染，那可能需要借助其他专门的库了。

怎么使用

用rasterio的第一步，永远是打开一个数据源。这通过一个上下文管理器来完成，感觉很Pythonic。

importrasteriowithrasterio.open('path/to/your/image.tif')asdataset:# 数据现在就在dataset对象里了data=dataset.read()# 读出所有波段的数据，是个三维数组profile=dataset.profile# 获取数据的元信息，比如坐标系、数据类型等bounds=dataset.bounds# 获取图像的地理范围

读出来的数据，本质上就是NumPy的ndarray。你可以用任何熟悉的NumPy或者SciPy方法去处理它，比如做归一化、滤波、或者波段运算。处理完之后，如果要保存，就需要参考原始数据的“profile”来创建一个新的文件。

# 假设我们对数据做了一些处理，得到了new_datanew_profile=profile.copy()new_profile.update(dtype=rasterio.float32,count=1)# 更新一些参数，比如数据类型和波段数withrasterio.open('output.tif','w',**new_profile)asdst:dst.write(new_data,1)# 将数据写入新文件的第一个波段

整个过程就像是在操作一个带有地理标签的NumPy数组，读写步骤被清晰地分开了。

最佳实践

用久了会发现，有些小习惯能让代码更稳当。比如，总是使用上下文管理器（with语句）来打开文件。这能确保文件句柄被正确关闭，尤其是在写数据的时候，避免文件损坏或不完整。

处理大型栅格数据是家常便饭，这些文件动辄几个GB。一股脑用dataset.read()读到内存里，很容易就把机器搞卡了。这时候可以看看数据的分块（block）信息，用dataset.block_window等方法进行分块读取和处理，或者直接使用rasterio.windows模块来读写数据的一个子区域（窗口）。这有点像吃一个大蛋糕，一次只切一小块来吃。

数据的地理信息是它的灵魂。做任何坐标相关的操作，比如裁剪、重投影，务必先搞清楚数据的坐标系（CRS）。dataset.crs属性会告诉你答案。如果要把数据从一个坐标系转到另一个（比如从WGS84经纬度转到Web墨卡托），rasterio.warp模块里的函数是专门干这个的，虽然计算有点耗时，但精度有保障。

元数据（profile）的传递也值得留心。当你基于一个源数据创建新数据时，最好先复制一份它的profile，然后只修改需要变动的部分（比如数据类型、压缩方式），而不是全部自己从头设置。这样能最大程度保留原始数据的“上下文”。

和同类技术对比

在Python的地理空间世界里，处理栅格数据的库不止一个。最常被拿来和rasterio比较的，可能是GDAL自己的Python绑定（osgeo.gdal）和xarray。

直接用GDAL的Python绑定，功能无疑是最全、最底层的，因为它就是GDAL本身。但它的API设计更接近C语言风格，用起来不那么“Python”，错误信息有时也比较晦涩。rasterio相当于是它的一个友好外壳，覆盖了80%的常用场景，代码写起来更简洁直观。如果你需要调用GDAL某个非常生僻的功能，那可能还得回去找原生的绑定。

xarray是另一个强大的库，它专注于带标签的多维数组，在处理气候、海洋等科学数据时大放异彩。它也能通过rioxarray这样的扩展来读写地理栅格数据，并且提供了非常优雅的、基于标签的数据选择和运算接口。如果你处理的数据本身维度多、结构复杂，或者你需要做复杂的沿维度聚合运算，xarray的抽象可能更合适。而rasterio则更“质朴”一些，它聚焦于地理栅格数据的核心读写和简单处理，和NumPy数组的对接更直接，在只需要处理二三波段影像、做做裁剪重投影的场景下，反而更轻快。

所以，选哪个有点像选工具。rasterio像一把趁手的螺丝刀，专攻一点，用起来直接；xarray像一套多功能工具箱，能应付更庞杂的数据结构。没有绝对的优劣，就看手头的活儿具体是什么了。

说到底，rasterio解决的是一个很具体但也很基础的问题：让Python程序员能更自然地去触碰那些带着经纬度信息的“像素世界”。它不试图包办一切，而是在自己的领域里把事情做得干净利落。下次当你再看到卫星图时，或许可以想想，它背后的那些数字，正安静地等着被这样的代码读取和解读呢。