Matplotlib/Cartopy/Python 无法检索高缩放级别（即缩小）的地图图块

Question

我正在使用 Python 和 Cartopy 创建一个地图应用程序，并尝试使用开源地图图块作为背景，以获得比默认 Cartopy 地图更多的选项。

它非常适合放大到相当近的地图，但是当我尝试从更高的高度获取视图时，出现了问题。如果我将缩放设置为 11，它就可以工作。如果我将它设置为 12，它会无限期地挂起并且不会给出回溯。

OSM 和 Stamen 地图服务器的结果相同。

这是一个简短的自包含示例（请注意，一两行可能是我尝试过的各种方法的产物）

import matplotlib as mpl
mpl.use('Agg')
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.io.img_tiles as cimgt


def custom_background(source_point):

    source_point = source_point.split(" ")
    source_point = (float(source_point[0]), float(source_point[1]))
    dx = 1.5
    dy = 1.5
    lon_min, lon_max = source_point[0]-dx, source_point[0]+dx
    lat_min, lat_max = source_point[1]-dy, source_point[1]+dy
    zoom = 7
    map_url = "https://www.openstreetmap.org/#map={}/{}/{}".format(zoom,source_point[0],source_point[1])
    tile = cimgt.OSM(url=map_url)
    tile = cimgt.StamenTerrain()
    ax = plt.axes(projection=ccrs.PlateCarree())
    ax.set_extent([lat_min, lat_max, lon_min, lon_max])
    ax.add_image(tile, zoom)
    #~ ax.add_image(tile)
    return ax

custom_background("45.068466 -66.45477")
plt.savefig("tile.png")

结果，缩放 = 7:

但如果我将缩放比例更改为 14，则无论我允许多长时间，程序都不会完成运行。

传递给 cimgt.OSM() 的 url 参数是可选的。无论有没有它，我都会得到相同的结果。（参见：https://scitools.org.uk/cartopy/docs/v0.16/cartopy/io/img_tiles.html#cartopy.io.img_tiles.OSM）

我是不是漏掉了什么？任何帮助将不胜感激，谢谢。

Answer 1

"zoom" levels are based on a Quadtree。从本质上讲，增加 "zoom" 分辨率会使图块数量增加四次方。

所以：

zoom level 0: 4^0 = 1 tile(s) to cover the globe
zoom level 1: 4^1 = 4 tile(s) to cover the globe
...
zoom level 7: 4^7 = 16,384 tile(s) to cover the globe
...
zoom level 14: 4^14 = 268,435,456 tile(s) to cover the globe

因此，如果您请求大面积的高缩放级别的图块，您最终可能会请求很多个图块。

在名为 find_images 的 Tiler 对象上有一个有用但未记录的方法。它的实现不太复杂：
https://github.com/SciTools/cartopy/blob/v0.16.0/lib/cartopy/io/img_tiles.py#L103-L122.

通过这种方法，我们实际上可以看到将用于给定范围的图块。重要的是，需要在平铺机的坐标系中提供范围（几乎完全是 Web Mercator）。

In [1]: import cartopy.io.img_tiles as cimgt

In [2]: import shapely.geometry as sgeom

In [3]: import cartopy.crs as ccrs

In [4]: tiler = cimgt.OSM()

In [5]: pt = 45.068466, -66.45477

In [6]: target = sgeom.box(pt[0] - 1.5, pt[1] - 1.5, pt[0] + 1.5, pt[1] + 1.5)

In [7]: target_mercator = tiler.crs.project_geometry(target, ccrs.Geodetic()).geoms[0]

所以所有的部分都就位后，我们可以开始找出我们需要在特定缩放级别绘制哪些图块：

对于您指定的目标，在缩放级别 0，我们需要以下图块 (x, y, z)：

In [8]: list(tiler.find_images(target_mercator, 0))
Out[8]: [(0, 0, 0)]

对于z=1，它仍然只有一个瓦片：

In [9]: list(tiler.find_images(target_mercator, 1))
Out[9]: [(1, 1, 1)]

但是对于 z=2，我们显然跨越了一个 tile 边界，因为我们现在需要两个 tile 来覆盖目标域：

In [10]: list(tiler.find_images(target_mercator, 2))
Out[10]: [(2, 2, 2), (2, 3, 2)]

Naturally, the list grows as we increase the zoom level:

In [11]: list(tiler.find_images(target_mercator, 3))
Out[11]: [(4, 5, 3), (5, 5, 3), (4, 6, 3), (5, 6, 3)]

In [12]: list(tiler.find_images(target_mercator, 6))
Out[12]: [(39, 47, 6), (40, 47, 6), (39, 48, 6), (40, 48, 6)]

当我们点击 z=7 时，我们发现我们需要 8 个图块来表示相关区域：

In [13]: list(tiler.find_images(target_mercator, 7))
Out[13]: 
[(79, 94, 7),
 (79, 95, 7),
 (80, 94, 7),
 (80, 95, 7),
 (79, 96, 7),
 (79, 97, 7),
 (80, 96, 7),
 (80, 97, 7)]

我相信您可以看到这是怎么回事，但让我们尝试将缩放级别设置为 14 以了解我们需要多少图块。为了节省电力和阅读时间，我们只打印该列表的长度...

In [14]: len(list(tiler.find_images(target_mercator, 14)))
Out[14]: 47334

是的，因此在您想要的范围内请求缩放级别 14 需要您以 256x256（8 位彩色映射 PNG）下载约 47,000 个图块。 z=0 磁贴 (https://a.tile.openstreetmap.org/0/0/0.png) 压缩后约为 8882 字节，因此假设这是典型情况，您最终会下载 ~420420588 字节 (400MB)。为了将这些数据保存在内存中，您还需要大约 2.9GB 的 RAM。最后，要将此数据重新投影到 PlateCarree，您需要至少将此 RAM 量加倍，并且假设一个高效的重新投影实现（cartopy 不是）。

希望这会告诉您为什么 您的代码似乎需要很长时间 - 您要求它做很多工作。关于当请求的图块数量过多时 cartopy 是否应该发出警告已经有一些讨论，但它总是归结为一个不合理的数字（你实际上可能想要获取那个许多瓷砖！）。我们也谈到了自动缩放级别选择——如果有足够的需求，这是really feasible

的东西。

HTH

Matplotlib/Cartopy/Python 无法检索高缩放级别（即缩小）的地图图块

Unable to retrieve map tiles at high zoom level (i.e. zoomed out) with Matplotlib/Cartopy/Python

python

gis

matplotlib

openstreetmap

cartopy