简介这份全球地形数据以 GeoTIFF 栅格形式提供基于 ETOPO1_Ice_g 数据集分辨率约 2km采用 WGS1984 坐标覆盖全球陆地与海底地形适合 GIS 制图、遥感影像处理、三维地形可视化、地学分析与教学演示等场景面向地信、测绘、海洋与地质相关专业的学习者和研究人员。压缩包共 3 个文件核心为 tif 地形栅格另含 ovr 金字塔索引与 xml 元数据整体约 266.01MBovr 可提升大范围缩放浏览效率xml 便于识别数据组织与投影信息。已有 3036 人学习下载。借助它可快速搭建全球尺度 DEM 底图进行高程查询、坡度起伏观察、海陆过渡带分析与可视化渲染也能为论文插图、课程实验和算法验证提供统一坐标与分辨率的基础数据。1. 全球地形数据 tif 格式一份文件里到底装了什么有人把一份号称覆盖全球的地形文件拖进看图工具渲染出来的青藏高原惨白一片太平洋漆黑一块中间还嵌着几块反相的高亮斑——那些不是山峰是-32768的 nodata 被当成真实高程参与了色带映射。全球地形数据 tif 格式最容易被误解的地方就在这里它是 GeoTIFF但 GeoTIFF 只是容器真正决定这份数据能不能用的是里面的像素值、数值类型、nodata 约定、投影和垂直基准这五样东西。搞清楚这五样你才能把它用在淹没分析、坡度计算、三维地形底图或者切片服务上搞不清楚后面每一个结果都会带着几十米甚至上百米的系统性偏差。这篇写给已经拿到文件、准备真正动手的人。2. 选源分辨率、垂直基准和体积的三笔账在写第一行代码之前先花十分钟把“我要哪一份”定下来比后面调三天参数都值。全球地形数据的选型从来不是“分辨率越高越好”而是三笔账互相牵制格网越细体积按平方级上涨覆盖范围越广越可能混进不同来源的拼缝垂直基准不统一再高的分辨率也白搭。下面按这三笔账逐一拆。2.1 分辨率、覆盖范围与体积的换算表先建立一个数量级直觉。全球 1 弧秒格网约 30 米在赤道附近的像元总数是 360×3600 列、180×3600 行合计约 8.4×10¹¹ 个像元。用 16 位整型存每像元 2 字节未压缩体积约 1.68 TB换成 32 位浮点直接翻倍到 3.4 TB。这个数字决定了后面所有的工程选择——你不可能在笔记本上一次性读进内存。产品代号类型标称格网覆盖范围典型数值类型全球未压缩体积主要短板雷达测绘 1 弧秒1″约 30 m南北纬 60°/56° 之间Int16约 1.7 TB高纬缺失雷达阴影区有空洞雷达测绘 3 弧秒3″约 90 m同上Int16约 190 GB同上精度更粗光学立体像对 1 弧秒1″南北纬 83° 之间Int16约 1.7 TB噪声与条带明显需编辑商业雷达 DSM 30 m1″全球陆地Int16 / Float32 视分发1.7–3.4 TB是 DSM 不是 DTM含建筑植被商业雷达 DSM 90 m3″全球陆地Int16约 190 GB同上陆海一体 15 弧秒15″约 450 m全球含海底Int16 / Float32约 7.5 GB只适合宏观尺度选型的实用判据是这样做全国或大洲尺度的水文、气候、可视域分析3 弧秒完全够用190 GB 的量级一台普通工作站能扛做城市级坡度、光伏选址、无人机航线预规划必须上 1 弧秒做海底地形、海啸传播就得换陆海一体的 15 弧秒产品因为陆域 DEM 在海洋部分全是 nodata 或者 0。至于 1 弧秒全球全量落地除非你要做在线切片服务否则没必要——按需取窗口才是常态。还有一点常被忽略不同产品的覆盖边界不一样。雷达测绘类产品在高纬地区没有覆盖如果你拼一张全球图却只用了它北极圈附近会整片空白看起来像代码 bug其实是数据本来就没有。选源时先把覆盖范围画出来对一遍。2.2 垂直基准同一座山为什么能差出上百米这是全球地形数据里最隐蔽的一类错误。高程值本身是个数字但它“相对于什么面”量的决定了这个数字的含义。三种基准最常见大地水准面正高、参考椭球面椭球高、以及区域平均海平面。多数全球 DEM 给的是正高也就是相对于某个大地水准面模型的高程而 GPS 接收机、很多三维引擎、部分仿真软件默认用的是椭球高。两者之间的差值叫高程异常全球范围内大致在 -105 米到 85 米之间波动。也就是说如果你把一份正高的地形直接当成椭球高去和 GPS 点位比对在某些区域会稳定地差出几十米而且这个偏差不是随机的是有空间规律的——它会让你的误差图看起来像一张光滑的大尺度波纹非常有迷惑性。判断方法很直接查元数据。合规的分发一般会在 GeoTIFF 的元数据段或者随附说明里写明垂直基准。用命令行看一眼# 只看元数据段和高程统计不加载全图 gdalinfo -nomd -mm SRTMGL1_N39E116.tif | head -20 gdalinfo SRTMGL1_N39E116.tif | grep -i -E datum|geoid|vertical|EGM第一行输出会告诉你尺寸、像元大小、数值类型、nodata第二条专门捞基准相关的键值。如果什么都没捞到那就按“正高、大地水准面未知”处理在跨源叠加前先做一次控制点比对找几个已知水准点或者 GNSS 静态观测点反算差值再决定要不要整体加一个偏移量。这个动作花不了半小时但能省掉后面反复怀疑代码的时间。提示不同大地水准面模型之间本身也有一米左右的差异跨产品拼接时如果精度要求到米级这一步不能省。2.3 三步验明一份 tif 的真实身份拿到任何一份全球地形 tif我习惯先跑三个检查一分钟内把它的底细摸清。第一步看结构第二步看数值分布第三步看物理像元大小。三步都过才敢往流水线里放。import math import numpy as np import rasterio path SRTMGL1_N39E116.tif with rasterio.open(path) as ds: print(CRS:, ds.crs) print(dtype:, ds.dtypes, nodata:, ds.nodata) print(size:, ds.width, x, ds.height) print(transform:, ds.transform) # 该景中心纬度处的实际像元尺寸经向要乘 cos(纬度) lat (ds.bounds.top ds.bounds.bottom) / 2.0 dy abs(ds.transform.e) * 111320.0 dx abs(ds.transform.a) * 111320.0 * math.cos(math.radians(lat)) print(中心纬度 %.3f 处像元: 南北 %.1f m, 东西 %.1f m % (lat, dy, dx)) data ds.read(1) valid data[data ! ds.nodata] if ds.nodata is not None else data print(nodata 占比: %.2f%% % (100.0 * (1 - valid.size / data.size))) print(有效高程: min%d max%d % (valid.min(), valid.max()))这段代码里三个地方值得说。ds.transform是仿射变换a是经向像元步长度e是纬向步长度负值表示自北向南。乘以 111320 是把度换算成米经向再乘cos(纬度)才是真实地面距离——在纬度 40° 处1 弧秒的经向地面距离只有约 23.6 米而不是 30.9 米。这一步是后面所有面积、坡度、缓冲区计算的根基做错了整条链路都是错的。nodata判断要小心有些产品的 nodata 是-32768有些是-9999还有些干脆没有声明 nodata海洋部分直接写 0。最后一种最坑因为 0 是一个合法高程值脚本无法自动区分“海平面”和“无数据”。所以“有效高程”和“nodata 占比”这两行统计一定要看如果 nodata 占比是 0% 但最小值是 -32768说明声明缺失得手动补。如果最大值出现 32767 这类边界值那基本是填充或溢出同样要处理。第三步的物理像元大小除了用于计算还能验证数据有没有被误重投影。如果一份号称 30 米的数据算出来南北像元是 0.000833 度3 弧秒那说明你手上是 90 米版本后面所有精度预期都要跟着下调。3. 拼接从上千个 1° 瓦片到一块可用的区域 DEM全球地形数据几乎都是以 1°×1° 分幅分发的一格 3601×3601 个采样点多一个点是为了相邻瓦片共享边线。要用就得先拼。拼接这件事看着简单实际上内存、接缝、坐标范围三个问题会同时找上门。3.1 为什么不该一上来就 gdalwarp 一把梭很多人的第一反应是gdalwarp把所有瓦片一次性投到一个目标坐标系里。这个做法在小范围、单投影的场景下没问题用在全球数据上会撞三堵墙。第一堵是内存。gdalwarp默认会尝试把整个输出范围纳入处理全球 1 弧秒输出就是 8.4 亿×1000 个像元中间缓冲区直接把内存撑爆。第二堵是不必要的重采样。源数据本来就是 EPSG:4326目标也是 4326重投影这一步是恒等变换但gdalwarp依然会按目标分辨率做一次重采样白白损失精度还拖慢速度。第三堵是边界外插值。当目标范围超出源数据覆盖时gdalwarp会用 nodata 填充但如果 nodata 没设对就会填出 0 值接缝处出现莫名其妙的“海平面台阶”。正确顺序是反过来的先做零拷贝的虚拟拼接把上千个文件在逻辑上合成一个数据集再按需取窗口。虚拟拼接不产生新文件、不复制数据、不重采样几乎瞬时完成。3.2 用 gdalbuildvrt 做零拷贝虚拟镶嵌# 把所有瓦片路径写进清单一行一个 find ./tiles -name *.tif | sort tiles.txt wc -l tiles.txt # 虚拟镶嵌优先最高分辨率统一 nodata隐藏 nodata 使其透明 gdalbuildvrt \ -resolution highest \ -srcnodata -32768 \ -vrtnodata -32768 \ -hidenodata \ -input_file_list tiles.txt \ mosaic.vrt # 验证结果应能看到正确的范围、分辨率与波段数 gdalinfo -stats mosaic.vrt | head -30-resolution highest是关键。不同来源的瓦片分辨率可能不一致有的 1 弧秒有的 3 弧秒默认行为是取平均或者取最低都会让最终结果变糊。highest保证按最细的那一档建格网粗分辨率的瓦片在读取时由 GDAL 自动上采样。-srcnodata和-vrtnodata成对出现作用是告诉 GDAL“源里这个值是空的输出里也标成这个值”。-hidenodata则让这个值在渲染和统计时被忽略避免出现第 1 章说的那种反相高亮块。-input_file_list比在命令行里堆几百个路径稳得多不受 shell 参数长度限制。生成出来的.vrt是个几百 KB 的 XML记录的是瓦片路径和地理变换读它等于读全部瓦片但内存占用只有你实际读取的那一块。如果只关心某个区域可以顺手加-te 73 3 136 54限定输出范围VRT 就只覆盖这个框后面的窗口读取也更快。-tap可以让输出对齐到整像元避免出现 0.4999 个像元这种偏移。3.3 用窗口读取只碰你需要的那一块有了 VRT取任意区域就是一次窗口读。下面这段是可直接抄的脚本输入 VRT输出一块带压缩和分块的目标 tif。import rasterio from affine import Affine from rasterio.enums import Resampling from rasterio.windows import from_bounds SRC mosaic.vrt DST region_3arcsec.tif BOUNDS (73.0, 3.0, 136.0, 54.0) # 左, 下, 右, 上单位度 DECIM 3 # 3 表示从 1 弧秒降采样到 3 弧秒 with rasterio.open(SRC) as src: nodata src.nodata if src.nodata is not None else -32768 # 按地理范围切出窗口并取整到整像元避免半个像元的偏移 win from_bounds(*BOUNDS, transformsrc.transform).round_offsets().round_lengths() h max(1, int(win.height // DECIM)) w max(1, int(win.width // DECIM)) # 关键带 out_shape 读取GDAL 只解码需要的块并按平均值聚合 data src.read( 1, windowwin, out_shape(h, w), resamplingResampling.average, boundlessTrue, fill_valuenodata, ) # 输出变换 窗口变换 × 缩放 out_transform src.window_transform(win) * Affine.scale(win.width / w, win.height / h) prof src.profile.copy() prof.update( driverGTiff, heighth, widthw, count1, dtypedata.dtype, transformout_transform, nodatanodata, tiledTrue, blockxsize512, blockysize512, compressdeflate, predictor2, BIGTIFFIF_SAFER, ) with rasterio.open(DST, w, **prof) as dst: dst.write(data, 1) print(done:, DST, w, x, h)逻辑上有三处需要展开。第一处是from_bounds加round_offsets().round_lengths()。地理范围换算成行列号一般不是整数不取整会让输出栅格与源栅格错开半个像元后续和多源数据叠加时会出现系统性偏移——这种偏移在地形上表现为整体轻微“滑动”肉眼很难发现但做坡度分析会放大成明显的条纹。第二处是out_shape配合Resampling.average。这是整个脚本性能的核心GDAL 不会把窗口内的原始像元全部读进内存再降采样而是在块解码阶段就完成聚合内存占用只和输出尺寸相关。对 1 弧秒到 3 弧秒的降采样average比nearest平滑得多但会削平尖峰做洪水淹没这类关心极值的分析用min或max更合适具体看你要保哪一端。第三处是boundlessTrue, fill_valuenodata。当请求范围超出 VRT 的覆盖边界时不加这两个参数会直接抛异常。加了之后越界部分填 nodata边界处理干净。3.4 参数怎么调块大小、重采样和 nodata这三个参数直接决定输出文件的质量和后续读取性能值得单独说。块大小blockxsize/blockysize必须是 16 的整数倍常用 256 或 512。块太小文件里块索引的开销占比升高块太大单次点查询要读的数据变多随机访问变慢。经验值是 512×512一个 Int16 单波段块正好 512 KB既能一次读满一条网络请求的量又不至于让点查询付出太多代价。重采样方法的选择取决于用途没有通用最优场景推荐方法原因降采样做底图渲染average平滑无摩尔纹计算流域最低点、淹没范围min保住谷底不被平均抬高计算山顶点、可视域max保住山脊线上采样做精细坡度cubic / bilinear连续曲面避免阶梯分类后地形如土地利用叠加nearest不产生新类别nodata 的处理要贯穿全流程。读的时候要maskedTrue或者手动比较写的时候要显式声明nodata转格式的时候要带-srcnodata/-dstnodata。只要有一环漏了这个值就会在下游变成真实高程参与运算。我见过最典型的事故是统计某区域平均海拔结果 nodata 的 -32768 拉低了平均值近千米排查了半天才发现是中间某一步gdal_translate没带 nodata 参数。4. 转成 COG压缩、分块与金字塔的参数取舍拼接好的 tif 如果只是本地放着看普通 GeoTIFF 就够。但只要涉及多进程并发读、网络按需取块、或者给地图服务做数据源就该转成 COGCloud Optimized GeoTIFF。它不是新格式是 GeoTIFF 的一种布局约定——让读取端能用一次 HTTP Range 请求只拿到需要的那几百 KB。4.1 COG 必须满足的四个硬条件第一内部必须分块tiled不能是条带strip布局。条带布局意味着要读一行就得扫过大半文件完全失去随机访问的意义。第二必须有金字塔overview。缩放到全球视图时如果不走金字塔GDAL 只能全分辨率读取再降采样全球 1 弧秒数据一次请求就是 TB 级。金字塔把多个降采样层级预存进同一个文件读哪层取哪层。第三主图像数据与 IFD影像文件目录的排布必须让头部先到。这样客户端读前 16 KB 就能拿到全部元数据和块索引再决定去取哪些块。第四压缩。不压缩的 COG 在网络传输上毫无优势。压缩算法与 predictor 的搭配直接影响体积和解码速度。4.2 一条 gdal_translate 命令完成转换gdal_translate mosaic.vrt global_dem_cog.tif \ -of COG \ -co COMPRESSDEFLATE \ -co PREDICTOR2 \ -co BLOCKSIZE512 \ -co OVERVIEW_RESAMPLINGAVERAGE \ -co OVERVIEW_COUNT6 \ -co BIGTIFFIF_SAFER \ -co NUM_THREADSALL_CPUS \ --config GDAL_NUM_THREADS ALL_CPUS-of COG是驱动名GDAL 3 之后内置不需要额外插件。COMPRESSDEFLATE是无损压缩里兼容性最好的一个几乎所有客户端都支持。PREDICTOR2是针对整型数据的差分预测对连续变化的地形高程特别有效通常能把压缩率再提升 20%–40%如果源数据是浮点用PREDICTOR3。BLOCKSIZE512前面解释过。OVERVIEW_RESAMPLINGAVERAGE决定金字塔怎么生成做底图选它没问题如果这份 COG 主要用途是极值查询改成NEAREST更保险避免金字塔层把山峰削平。OVERVIEW_COUNT6表示生成 6 级金字塔从 2 倍到 64 倍降采样覆盖从区域到全球的缩放范围一般够用。BIGTIFFIF_SAFER是必须的。全球尺度文件很容易超过 4 GB 的经典 TIFF 上限不加这个参数会在写到一半时直接报错而且报错信息往往只说“写入失败”不告诉你是尺寸问题。NUM_THREADSALL_CPUS让压缩并行实测在 16 核机器上能把转换时间压到单线程的三分之一左右。4.3 用 rio-cogeo 做校验与批量重写转完之后不能只看文件生成了就完事要验证它真的是 COG。rio-cogeo这个工具做校验很省事pip install rio-cogeo # 校验会逐项检查 tiled、overview、IFD 位置、压缩等 rio cogeo validate global_dem_cog.tif # 批量转换把目录下所有 tif 重写为 COG rio cogeo create input.tif output_cog.tif \ --cog-compress DEFLATE \ --cog-blocksize 512 \ --overview-resampling average \ --overview-level 6validate的输出会逐条列出通过与否。常见的失败项是 “Image is not tiled” 和 “Overviews are not internal”前者说明源文件是条带布局后者说明金字塔被写成了外部的.ovr附属文件——附属文件在单机没问题但传到对象存储上就读不到了必须内部化。create支持大批量重写配合 shell 循环或者 Python 的concurrent.futures就能并行处理整个瓦片目录。注意并行度不要开太高压缩本身是 CPU 密集的开满核反而会因为 IO 争抢变慢一般设成核数的 70% 比较稳。4.4 压缩算法与块大小的实测取舍不同压缩算法在全球地形数据上的表现差异不小下面是大致量级的经验值具体数字随地形复杂度波动山区和平原能差一倍以上压缩配置相对未压缩体积解码速度适用场景无压缩100%最快仅本地临时文件LZW约 45%–55%快兼容性优先的老客户端DEFLATE PREDICTOR2约 25%–35%中通用首选ZSTD level 9约 22%–30%快新版 GDAL 环境ZSTD level 15约 18%–25%中归档存储读少写多JPEG有损约 8%–12%最快纯渲染底图不做分析选的时候记住一条分析用无损渲染可有损。做坡度、体积、淹没计算绝不能用 JPEG 压缩因为它会在高频地形区域引入块状伪影坡度一阶差分直接把伪影放大成假陡坡。只做视觉底图的话JPEG 能把全球 1 弧秒压到一百多 GB代价可以接受。块大小方面512 是通用答案但如果你的典型查询是“给一个经纬度返回高程”256 反而更好因为单次读取的字节数减半。如果是“给一个省的范围做统计”1024 能减少块索引跳转次数。这个参数没有绝对最优取决于你的访问模式上线前拿真实查询压一遍最靠谱。5. 避坑全球地形数据 tif 处理里最容易翻车的五件事下面这五条都是我在实际项目里真实撞过的按“现象 → 原因 → 解决”记录每条都给出可操作的判定方法。5.1 高程整体差几十米垂直基准没对齐现象把地形数据和一批 GNSS 实测点比对发现所有点的残差都是同一个方向的几十米而且残差在空间上呈大尺度平滑变化不是随机散布。原因两边用的垂直基准不同。全球 DEM 多为正高相对于某个大地水准面而 GNSS 直接给的是椭球高。两者之间的高程异常在全球范围内从 -105 米到 85 米不等在你关注的那个区域内可能恰好是个稳定值所以看起来像“整体偏移”。解决第一步查元数据确认基准第二步找 3–5 个同时具备水准和 GNSS 观测的控制点反算差值。如果差值在区域内变化不大比如都在 40–45 米之间整体加一个常数就能满足米级需求如果变化超过你的精度要求就得引入高程异常格网做逐像元改正。5.2 出现大片空洞或极端值nodata 没被识别现象统计出来的最低高程是 -32768或者渲染时出现规则的黑块、白块又或者某个区域的均值莫名其妙地偏低。原因nodata 在某一步丢失了。常见于三种情况源文件本来就没声明 nodata中间用了不支持 nodata 传递的工具做了格式转换把数据转成 NumPy 数组后做了运算再写回写回时忘了重新声明。解决读取时统一用maskedTrue拿到掩膜数组或者手动data ds.nodata建掩膜任何一次写出都显式指定nodata跨工具链时用gdalinfo抽查每一步的中间产物。养成“每步都看一眼 nodata 声明”的习惯比事后排查省事得多。5.3 拼接处出现台阶分辨率与重采样不一致现象在瓦片边界处能看到一条笔直的亮线或暗线放大后是一两米的台阶沿着经纬度整数度分布。原因相邻瓦片来自不同批次、不同版本甚至是不同产品分辨率或高程基准有细微差异虚拟镶嵌时没有统一重采样策略GDAL 用了默认的最近邻边界处直接取了一侧的值。解决建 VRT 时明确-resolution highest并在读取时指定Resampling.average或bilinear让边界处平滑过渡如果差异较大可以在拼接前对每块做一次直方图匹配把系统性偏移拉平。另外要注意 1° 瓦片本身有 1 个像元的重叠边处理不当会在重叠处产生双倍权重。5.4 面积和距离算错拿度当米用现象在 EPSG:4326 的栅格上直接算面积得到的数字和实际相差几倍坡度计算出来在低纬正常、高纬离谱。原因经纬度坐标下的像元不是正方形。1 弧秒在南北方向约 30.9 米在东西方向是 30.9×cos(纬度) 米纬度 60° 处只有 15.4 米。直接用度做单位算面积等于默认所有像元都是正方形高纬地区会被严重高估。解决要么重投影到等积坐标系再算要么在计算时逐像元乘上cos(纬度)修正系数。坡度计算更要注意x方向和y方向的地面距离必须分别用不同分母否则在纬度 45° 以上误差会超过 40%。5.5 DSM 当 DTM 用树和楼被算成地形现象在城市或者森林区域算出来的地形起伏明显偏大建筑物轮廓清晰可见做洪水淹没模拟时水流被“挡”在树冠和屋顶外面。原因部分全球产品本质上是 DSM数字表面模型包含植被冠层和人工建筑的高度不是裸地 DTM。产品说明里通常会写明但很多人只看分辨率不看类型。解决选源阶段就确认是 DSM 还是 DTM。如果只有 DSM 而需要 DTM可以用公开的建成区、森林覆盖掩膜做粗略剔除或者用滤波算法形态学开运算、渐进式三角网加密做点云级的表面滤波——后者对栅格数据精度有限只能算权宜之计。真正做城市级水文还是得用专门加工的 DTM 产品。6. 把全球 DEM 压到毫秒级金字塔层级 按块读取的组合数据落成 COG 之后真正决定查询体验的只剩一件事你让它读了多少字节。全球 1 弧秒 COG 即使压缩后也在几百 GB 量级但只要金字塔层级选对、读取窗口收窄到一个块单次点查询完全可以做到几十毫秒。这个技巧我自己用了两年是整套方案里性价比最高的一步。核心思路是把“缩放级别”映射到金字塔层级。假设屏幕上一千像素宽要显示全球 360 度经度那每像素对应 0.36 度约 1296 弧秒最接近的金字塔层级是 2¹⁰1024 倍降采样。此时 GDAL 只读那一层的对应块一次 512×512 的 Int16 块是 512 KB网络请求一次搞定如果误用第 0 层同样一个屏幕需要读的原始数据是它的百万倍。import rasterio from rasterio.windows import Window from rasterio.enums import Resampling with rasterio.open(global_dem_cog.tif) as ds: # 打印所有金字塔层级及其对应分辨率选层级时直接看这个 print(overviews:, ds.overviews(1)) lon, lat 116.39, 39.90 row, col ds.index(lon, lat) # 点查询只读 1 个像元 z ds.read(1, windowWindow(col, row, 1, 1), maskedTrue) print(高程:, z[0, 0]) # 缩略图查询走第 5 层金字塔一次只解码一小块 small ds.read(1, out_shape(256, 256), resamplingResampling.average)ds.overviews(1)返回的是各层的降采样倍数比如[2, 4, 8, 16, 32, 64]。做 Web 地图时把缩放级别和这个列表对齐就能保证每一级只读一层金字塔。index()把经纬度换成行列号注意它按左上角原点计算纬度是自北向南递增的行号别把顺序写反——这个错误我在早期代码里犯过结果是查出来的高程总是隔壁纬度的。还有一个容易被忽略的细节点查询用Window(col, row, 1, 1)时GDAL 实际读的是包含这个像元的整个块512×512。如果服务 QPS 很高块大小从 512 降到 256 能把单次读取量压到四分之一代价是文件块数变四倍、索引变大。这个权衡要靠压测决定我给的经验值是点查询为主选 256范围统计为主选 512两者都重要就复制一份不同块大小的副本分别服务不同接口。最后说个我踩过的血泪教训金字塔层级不是越多越好。早期我给一份全球数据生成了 12 级金字塔结果文件头部索引膨胀到几十 MB客户端光读索引就花了半秒。后来砍到 6 级覆盖从 2 倍到 64 倍降采样实际够用头部也回到了 KB 级。金字塔的层级数应该由你的最小显示比例决定而不是顺手写个大数字。希望帮到你。本文还有配套的精品资源点击获取