Python实战：利用天地图API批量下载与拼接遥感瓦片数据

1. 为什么你需要掌握天地图瓦片下载与拼接

如果你正在处理地理信息相关的项目，比如城市规划、环境监测，或者只是想为自己研究的区域制作一张高清的卫星底图，那你很可能遇到过这样的烦恼：网上的在线地图服务虽然方便，但无法直接获取到高清、无偏移、可以自由分析的大范围影像数据。手动截图？那精度和效率简直是一场噩梦。这时候，天地图作为国内权威的地理信息服务，其提供的遥感瓦片数据就成了一个绝佳的选择。

天地图提供了多种类型的地图瓦片，包括电子地图、卫星影像和地形图。我们今天要重点“拿下”的，就是它的卫星影像瓦片。这些瓦片就像一块块标准尺寸的拼图，覆盖了全球的每一个角落。我们的目标，就是用Python写个脚本，像一只勤劳的小蜜蜂，自动飞到天地图的服务器上，把我们指定区域的所有“拼图块”都下载回来，然后无缝地拼接成一张完整的大图。

这个过程听起来复杂，但拆解开来，核心就是三步：第一，把你关心的地理范围（比如一个城市的边界）转换成瓦片地图的行列号坐标；第二，根据这些坐标，批量构造下载链接，把图片一块块抓下来；第三，像玩拼图一样，把这些小图片按照正确的顺序拼接起来。我做过不少类似的项目，实测下来，一旦脚本跑通，获取数据的效率和质量远超手动操作，而且完全可重复、可定制。接下来，我就带你一步步实现这个自动化流程，避开我踩过的那些坑。

2. 核心原理：经纬度如何变成一张张小瓦片

在开始写代码之前，我们得先搞明白，我们熟悉的经纬度（比如北京天安门广场的坐标：116.397, 39.907），是怎么变成地图服务器上一张张编号为(x, y, z)的图片的。这背后是地图学中经典的“瓦片金字塔”模型。

你可以把整个世界想象成一个巨大的正方形。在最顶层（缩放级别z=0），整个世界只用一张256x256像素的图片来表示，信息非常粗略。当我们把地图放大一级（z=1），这张大图就会被切成2行×2列，一共4张小图，每张还是256x256像素，但能显示更多细节。以此类推，每放大一级，瓦片数量就变成上一级的4倍。到了z=18级（非常高的缩放级别），全球的瓦片数量将是惊人的2^18 * 2^18张。天地图提供的卫星影像，最高可以到z=18级，细节足够看清街道和大型建筑物。

那么，一个具体的经纬度点，到底落在哪一张瓦片上呢？这里涉及到一个关键的数学转换。对于采用Web墨卡托投影（EPSG:3857）的瓦片（如Google地图、高德、天地图的vec_w电子地图），其行号（y）的计算公式会稍微复杂一些，因为纬度到y轴的投影不是线性的。公式如下：

import math

def wgs84_to_tile_mercator(lon, lat, zoom):
    """
    将WGS84经纬度转换为Web墨卡托瓦片坐标。
    :param lon: 经度
    :param lat: 纬度
    :param zoom: 缩放级别
    :return: 瓦片的列号x和行号y
    """
    n = 2.0 ** zoom
    # 经度转列号x (这个计算对于经纬度投影和墨卡托投影是通用的)
    tile_x = int((lon + 180.0) / 360.0 * n)
    # 纬度转行号y (墨卡托投影公式)
    lat_rad = math.radians(lat)
    tile_y = int((1.0 - math.log(math.tan(lat_rad) + 1.0 / math.cos(lat_rad)) / math.pi) / 2.0 * n)
    return tile_x, tile_y

但是，请注意！天地图的卫星影像瓦片（img_w服务）采用的坐标系是WGS84经纬度投影（EPSG:4326）。这种投影方式下，经纬度到瓦片行列号的转换是线性的，公式更简单直观。这也是我们原始文章里那个get_tianditu_info函数所用的方法。它的核心是一个分辨率字典，记录了每个层级下，一个像素代表多少度的经纬度。计算行列号的逻辑是：

1. 将经度值加上180度（因为经度范围是-180到180），使其变为0到360的正值。
2. 将这个值除以（当前层级一个瓦片的宽度所代表的度数）。一个瓦片宽256像素，所以“瓦片宽度度数” = resolution[zoom] * 256。
3. 对结果取整，就得到了该点所在的列号（x）。行号（y）的计算类似，用90度减去纬度值（纬度范围是-90到90），再做除法取整。

所以，当我们拿到一个矩形区域（比如一个城市的边界框：左上角经度lon1、纬度lat1，右下角经度lon2、纬度lat2）和一个缩放级别z时，我们就能分别计算出左上角点所在的瓦片行列号(x1, y1)和右下角点所在的瓦片行列号(x2, y2)。那么，我们需要下载的瓦片范围，就是从(x1, y1)到(x2, y2)这个矩形区域内的所有瓦片。理解了这个转换原理，我们才能准确地告诉程序：“嘿，去把这一片‘拼图’都给我拿回来。”

3. 实战第一步：准备你的Python环境和工具箱

工欲善其事，必先利其器。在开始写核心代码之前，我们需要搭建好Python环境并安装必要的库。我强烈建议你使用conda或venv创建一个独立的虚拟环境，这样不会搞乱系统本身的Python包。这里我以使用pip为例。

打开你的终端或命令提示符，依次安装以下库：

pip install requests
pip install pillow
pip install opencv-python
pip install numpy

我来简单解释一下每个库的用途：

• requests: 这是Python中最好用的HTTP库，没有之一。我们将用它来向天地图服务器发送请求，下载瓦片图片。它比标准库的urllib更简洁、更强大。
• Pillow (PIL): Python的图像处理库。我们主要用它来打开从网络下载的图片数据流。
• opencv-python (cv2): 强大的计算机视觉库。在这里，我们主要借用它高效的图像拼接功能（np.vstack和np.hstack）和图像保存功能。
• numpy (np): Pyth