# 在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程序员能更自然地去触碰那些带着经纬度信息的“像素世界”。它不试图包办一切，而是在自己的领域里把事情做得干净利落。下次当你再看到卫星图时，或许可以想想，它背后的那些数字，正安静地等着被这样的代码读取和解读呢。
