CLCD、tiff、30米分辨率、41年连续——这几个关键词同时出现基本说明你正在为长时间序列土地利用/覆被变化分析找数据。2025年武大CLCD数据集又更新了时间跨度从1985年拉到2025年累计41期每期都是30米分辨率的GeoTIFF既提供全国整幅也提供分省版本。做城市扩张、耕地流失、生态修复、国土空间规划或者只想给自己论文加一条可信的地表覆被时间序列这套数据都是值得优先试的方案。这篇文章不打算复制官网说明书而是把数据落地后的读取、裁剪、统计、变化分析、报错排查全部过一遍让你拿到手就能用。1. 项目概述41年连续图斑背后的CLCD数据集1.1 数据集来源与技术路线CLCD全称China Land Cover Dataset官方中文名通常叫“中国土地覆盖数据集”由武汉大学相关遥感团队开发并持续更新。早期版本基于Landsat系列卫星的30米多光谱影像结合训练样本和随机森林分类算法经过时间滤波、质量控制等步骤产出逐年覆盖中国全境的土地覆盖分类结果。2021年前后团队将1990到2019年的数据正式对外发布之后每年滚动更新到2025年这次发布已经累积到1985-2025年共41期。如果你写论文需要引用通常会引用发表在Earth System Science Data期刊上的相关文章。和动不动就几百米分辨率的全球土地覆盖产品相比CLCD的30米空间分辨率能清楚区分街区尺度的建设用地、小面积坑塘以及零散耕地。时间连续性则是它最大的卖点之一41年逐年更新没有断层这对年际变化检测、趋势分析和政策评估来说特别友好。毕竟土地利用变化不是一年两年能看出来的有了这个长序列就不用来回拼多源数据了。1.2 九类分类体系与编码CLCD的分类体系一共9类所有的tif文件里每个像元的值都直接对应一个地表类型编码类别典型地物1农田水田、旱地、种植园2森林乔木林、混交林、竹林3灌木灌木林、灌草丛4草地草原、草甸、稀疏草地5水体河流、湖泊、水库、坑塘6冰雪永久冰川、积雪7裸地戈壁、沙地、裸土8不透水面城市建筑、道路、工矿用地9湿地沼泽、滩涂、红树林这个编码规则在CLCD内部基本是统一的但不同渠道发布的版本偶尔会把类别顺序调整所以分析前先看下随数据发布的说明文件别直接把值当类别名称用。后面5.2里我还会专门讲怎么快速验证编码。理解这套编码是后续所有统计的基础。例如要提取城市扩张范围只需要把值为8的像元筛出来要算耕地面积就把值为1的像元全部统计一遍。对于不熟悉遥感的人来说这比自己去训练一个分类模型门槛低很多因为它已经是一份成品。1.3 为什么选CLCD而不是其他产品市面上的土地覆盖数据不少GlobeLand30精度高但主要发布2000、2010、2020等几个年份做逐年分析不够用ESA WorldCover有10米分辨率但只是单年份或多期产品MODIS MCD12Q1有逐年序列但分辨率500米做县域尺度分析太粗。CLCD的优势正好卡在“30米分辨率逐年连续”这个位置它既有足够细的空间颗粒度也具备完整的年际序列。当然它不是万能的。CLCD是自动分类结果在林地与灌木、旱地与草地交错区域个别年份可能有错分受Landsat影像质量和季节影响某些像元会出现年际抖动。所以我在实际项目中一般把它当作基础底图和趋势分析主数据在关键区域再叠加高分辨率影像抽样验证。下面的实操内容全部基于CLCD的GeoTIFF文件展开不管你是研究城市、流域还是农田这些步骤都能复用到你自己的数据上。2. 数据文件形态与选型思路2.1 全国版与分省版的区别和选择标题里注明“全国/分省”说明这次发布提供两种粒度一个文件覆盖整个中国一类是每个省级行政区单独裁剪一份。全国版每个年份是一整幅大tif通常几百MB到几个GB不等好处是范围完整、坐标系一致、从头到尾不会出现跨省边界拼接问题坏处是打开和处理时比较吃内存如果只要某个省仍然得裁剪一次。分省版则是提前按行政区边界裁剪好的每个省一个文件体积小、打开快适合只研究一个省的情况。但这类文件有几个潜在问题一是裁剪边界可能直接用省界矢量海岸线、湖岸、小岛在裁剪时容易被截掉二是省界附近如果有飞地或者河流跨越省界边缘会出现不自然的割裂三是不同省版本之间如果制作批次不同坐标系和投影也可能不完全一致。我给一个很具体的选型建议如果你只做单个省的大尺度统计可以用分省版节省时间如果你做全国对比、跨省流动分析或者涉及省界两侧的连续地块最好还是下全国版自己按研究区裁剪。这样能把“边界漂移”这一类问题控制在自己手里后面4.1会给出裁剪脚本。2.2 文件命名与内部结构CLCD的tif文件命名虽然不同渠道略有区别但常见规律是包含数据集名、年份、范围三个信息。例如全国版可能是CLCD_1985.tif、CLCD_2025.tif这样的形式分省版可能类似CLCD_1985_Guangdong.tif这样的结构具体以你拿到的实际文件名为准。文件内部是GeoTIFF标准结构影像元数据里记录了坐标系、仿射变换参数、像元尺寸、四至边界和无效值信息。打开文件后你看到的“图像”其实是一个单波段栅格每个像元存一个整数分类码而不是一张RGB彩色图。很多朋友第一次打开CLCD发现是黑白的以为是文件坏了实际上这是正常现象——它本身就是单波段分类数据需要通过“分类渲染/唯一值渲染”给它赋予颜色才能直观查看。2.3 坐标系与像元属性CLCD官方常见版本一般使用WGS84地理坐标系也就是经纬度坐标像元尺寸约0.00027度大约对应30米。有些处理好的版本也可能被重投影到Albers等积投影或UTM投影。这里最关键的一点是拿到文件后第一件事就是用工具查看它的CRS和像元大小不要想当然。因为后续裁减、面积计算、变化检测都依赖坐标系的一致性。我常用的查看方法是在QGIS里双击图层看“图层属性-信息”或者在Python里用rasterio打印crs、transform、res这些字段。如果发现是WGS84做面积统计时最好先投影到等积投影否则按像素直接乘900平方米会引入变形误差。这个坑太常见了后面5.1专门再讲。3. 环境准备与读取操作3.1 安装Python依赖处理CLCD这种GeoTIFF最顺手的工具还是Python。核心依赖就那么几个rasterio负责读写栅格geopandas负责矢量边界numpy负责数据计算。如果你还没有环境直接创建虚拟环境并安装conda create -n clcd python3.10 conda activate clcd conda install -c conda-forge rasterio geopandas numpy matplotlibrasterio自带GDAL所以不需要单独装GDAL命令行。QGIS可以作为快速可视化工具后面3.3里也会提到配置方法。这里解释下为什么选rasterio而不是纯GDAL命令行在脚本里做掩膜裁剪、随机采样、属性读取、批量循环rasterio的API更贴近Python习惯代码量小错误信息也友好。3.2 用rasterio读取tiff基础信息先写一段最小代码看看这一帧tif的“体检报告”import rasterio tif_path CLCD_2025.tif with rasterio.open(tif_path) as src: print(文件名:, tif_path) print(宽x高:, src.width, x, src.height) print(波段数:, src.count) print(坐标系:, src.crs) print(地理范围:, src.bounds) print(仿射变换:, src.transform) print(有效值/NoData:, src.nodata) print(数据类型:, src.dtypes[0])输出里重点看三个字段坐标系决定了后面投影转换地理范围能帮你确认数据是否覆盖目标区域NoData值则直接影响统计结果如果文件里把背景值设成了255之类的特殊值统计类别面积前必须把它过滤掉。接着读取整幅影像import numpy as np with rasterio.open(tif_path) as src: data src.read(1) print(唯一值分布:, np.unique(data, return_countsTrue))这一步就能快速看到分类编码和每个类别的像元数量。如果文件很大一次性read整个数组可能内存爆炸后面5.3会说处理大文件的办法。3.3 QGIS快速可视化不会写代码也可以直接用QGIS处理。把tif拖进QGIS窗口后默认显示可能是灰色渐变右键图层选择“属性-符号系统”把渲染类型改为“单波段伪彩色-唯一值”按1到9分别设置颜色。比如不透水面用红色、农田用黄色、森林用绿色、水体用蓝色这样一张直观的土地覆盖图就出来了。QGIS里还可以用“栅格计算器”快速提取某一类。比如提取不透水面表达式写成CLCD_202518输出结果会变成布尔栅格1表示不透水面0表示其他之后就能做面积统计或转矢量。对于只想临时看一个区域的情况QGIS已经够用但如果要做批量处理和可复现分析还是建议用Python脚本。4. 实操裁剪、统计与变化分析4.1 按行政区掩膜裁剪研究区往往不是全国第一步通常是把目标区域裁剪出来。我这里以任意研究区shp为例给出标准掩膜裁剪流程import geopandas as gpd import rasterio from rasterio.mask import mask # 1. 读取研究区边界并统一坐标系 shp_path research_area.shp gdf gpd.read_file(shp_path) with rasterio.open(CLCD_2025.tif) as src: # 2. 与栅格坐标系统一 if gdf.crs ! src.crs: gdf gdf.to_crs(src.crs) # 3. 按边界裁剪cropTrue表示只保留范围内 out_image, out_transform mask(src, shapesgdf.geometry, cropTrue) out_meta src.meta.copy() out_meta.update({ height: out_image.shape[1], width: out_image.shape[2], transform: out_transform, nodata: src.nodata }) with rasterio.open(CLCD_2025_clip.tif, w, **out_meta) as dst: dst.write(out_image)这段代码里最关键的是坐标统一如果研究区shp是CGCS2000或UTM投影而tif是WGS84不转换直接裁剪很可能得到空结果或者错位结果。另外裁剪时mask函数默认会把研究区外的像元设为NoData如果你只想保留研究区内完整像元还需要检查研究区边界与栅格边缘的对齐情况。裁剪之后记得检查输出文件范围和原数据是否重叠with rasterio.open(CLCD_2025_clip.tif) as clip: print(clip.bounds)4.2 各类别面积统计面积统计是最常见的需求。先提醒一句如果tif的坐标系是WGS84经纬度直接按像元计数乘900得到的面积只能当估算不能用于正式报告和论文因为不同纬度每个像元的实际地表距离不同。正确做法是先投影到适合中国的等积投影。以中国常用的Albers等积投影为例用rasterio的reproject把原数据投影过去再统计面积import rasterio from rasterio.warp import calculate_default_transform, reproject, Resampling from pyproj import CRS import numpy as np dst_crs CRS.from_proj4(projaea lat_00 lon_0105 lat_125 lat_247 datumWGS84 unitsm no_defs) with rasterio.open(CLCD_2025.tif) as src: transform, width, height calculate_default_transform( src.crs, dst_crs, src.width, src.height, *src.bounds ) kw src.profile.copy() kw.update({ crs: dst_crs, transform: transform, width: width, height: height, compress: lzw }) with rasterio.open(CLCD_2025_albers.tif, w, **kw) as dst: for i in range(1, src.count 1): reproject( sourcerasterio.band(src, i), destinationrasterio.band(dst, i), src_transformsrc.transform, src_crssrc.crs, dst_transformtransform, dst_crsdst_crs, resamplingResampling.nearest ) # 统计各类别面积 with rasterio.open(CLCD_2025_albers.tif) as src: data src.read(1) classes, counts np.unique(data, return_countsTrue) pixel_area 30 * 30 # 等积投影下每个像元约30米 for cls, cnt in zip(classes, counts): area_km2 cnt * pixel_area / 1_000_000 print(f分类{cls}: {cnt}像元, {area_km2:.2f}平方公里)注意投影后像元大小不一定是精确的30米但等积投影下每个像元代表的地表面积基本稳定所以按900平方米估算足够。如果做更严谨的面积核算可以用transform计算每个像元的实际面积不过CLCD做常规统计时差别不大。4.3 多年份转移矩阵计算转移矩阵是土地利用变化分析里最常用的工具。比如想知道1985-2025年有多少农田变成了建设用地需要把两个年份的分类栅格叠在一起统计。首先需要确保两个tif的坐标系、范围、像元网格完全一致。如果都是从同一个CLCD系列下载的全国版通常是一样的可以直接计算import numpy as np import rasterio def read_cls(path): with rasterio.open(path) as src: return src.read(1) cls_1985 read_cls(CLCD_1985.tif) cls_2025 read_cls(CLCD_2025.tif) # 保证维度一致 assert cls_1985.shape cls_2025.shape, 两个文件网格不一致需要先配准 n_class 10 transition np.zeros((n_class, n_class), dtypenp.int64) for y in range(cls_1985.shape[0]): for x in range(cls_1985.shape[1]): a cls_1985[y, x] b cls_2025[y, x] if a 0 and b 0 and a n_class and b n_class: transition[a, b] 1这段双重循环可读性好但太慢了。实际上可以用NumPy的广播直接统计flat_1985 cls_1985.ravel() flat_2025 cls_2025.ravel() valid (flat_1985 1) (flat_1985 9) (flat_2025 1) (flat_2025 9) idx flat_1985[valid] * 10 flat_2025[valid] counts np.bincount(idx, minlength100).reshape(10, 10) # counts[a][b]: 1985年为a类、2025年为b类的像元数量转移矩阵出来后可以进一步计算各类的净变化量行和减列和各类的稳定性对角线占比两两转移的主要去向每行非对角的最大值再把它转成CSV方便画桑基图import pandas as pd df pd.DataFrame(counts[:9, :9], index[f1985_class{i1} for i in range(9)], columns[f2025_class{j1} for j in range(9)]) df.to_csv(transition_1985_2025.csv)4.4 从tiff中提取点位的类别或高程信息有时候我们不只是要整幅图而是想知道一系列样点的土地利用类型。比如野外调查点了100个坐标想快速匹配每个点所在位置的CLCD类别。rasterio提供了很方便的采样方法import rasterio # 点位坐标注意顺序是(x经度, y纬度) points [(113.25, 23.15), (114.08, 29.60), (120.20, 30.30)] with rasterio.open(CLCD_2025.tif) as src: samples list(src.sample(points)) for pt, val in zip(points, samples): print(f点{pt} - 类别{val[0]})这里有个容易错的地方src.sample要求输入点的坐标与栅格坐标系统一。如果tif是WGS84而你的点是UTM坐标就需要先转坐标。另外顺带说一句最近有朋友问我“能不能用CLCD的tif提取区域内高程点”。CLCD不是高程数据它存的是地表分类代码没有海拔信息。要从高程数据里提取点的高程值应该用SRTM、ALOS DEM或者NASA DEM这类产品但读取方法完全一样换成dem.tif路径后用src.sample()就能拿到每个点的高程值。所以如果你拿到了DEM的tif上面这段代码同样适用只是返回值从分类编码变成了海拔米数。4.5 批量处理41期文件的脚本思路做长时序研究时最烦的是一个个年份手动处理。CLCD总共41期建议写一个循环批量提取不透水面面积或者批量裁剪然后拼接成一张时间序列表。下面这个例子统计全国版文件里1985到2025年每年不透水面的像元数import glob import numpy as np import rasterio years list(range(1985, 2026)) results [] for year in years: path fCLCD_{year}.tif try: with rasterio.open(path) as src: data src.read(1) impervious_count np.sum(data 8) results.append((year, int(impervious_count))) print(year, 完成) except FileNotFoundError: results.append((year, None)) print(year, 缺失)如果你是用分省版把路径换成省名即可。跑完后保存import pandas as pd df pd.DataFrame(results, columns[year, impervious_pixels]) df[area_km2] df[impervious_pixels] * 900 / 1_000_000 df.to_csv(impervious_trend.csv, indexFalse)这种批量脚本的实际意义在于把41年的数据变成一条可计算的趋势线。你在分析报告里画出的城市扩张曲线、耕地减少曲线本质上就是这样一个循环统计出来的。5. 常见问题与避坑清单5.1 坐标系混用导致面积严重失真我见过最多的问题就是拿WGS84经纬度的tif直接去算面积还乘900平方米。这样算出来的结果在高纬度地区会明显偏大。原因很简单经纬度坐标下经度方向的真实距离随纬度变化一个0.00027度的像元在低纬度可能接近30米在高纬度会缩水不少乘出来的面积自然不是真实面积。解决方式就两条一是做面积敏感分析前先把数据投到适合中国的Albers等积投影再统计二是用geopandas先把研究区矢量投到等积投影再配合栅格统计。CLCD本身提供了坐标信息重投影不复杂4.2里的代码可以直接复制用。5.2 分类编码容易弄混CLCD大部分版本采用1到9的编码但不同平台上重新发布的剪辑版可能把背景值设成了0或255甚至把分类顺序调整过。统计之前建议先打印np.unique看看最大值、最小值再用4.2的方法算一下每类像元数量是否合理。如果发现第8类不透水面面积占了全国一半那基本都是编码理解错了赶紧去查随文件的分级说明。5.3 大tiff内存不足全国版tif动辄几GBread(1)直接把整个数组载入内存8GB内存的电脑很容易卡死或溢出。推荐几种办法用窗口读取window参数分批处理先用rasterio.warp.reproject做重投影时顺便降低分辨率比如聚合到100米做探索性分析或者先用QGIS/GDAL把文件裁剪到研究区再读入内存。import rasterio from rasterio.windows import from_bounds # 只读取指定地理范围的窗口 with rasterio.open(CLCD_2025.tif) as src: win from_bounds(west110, south20, east115, north25, transformsrc.transform) data src.read(1, windowwin)用窗口读取对内存非常友好适合只关注某一小块区域的情况。5.4 省级版数据边界不干净分省版数据在省界处经常有几个问题一是沿海岛屿可能因为边界分割被单独删掉二是跨省水体在省界两侧被硬生生切开导致同一水面分成两半三是如果行政区边界更新过旧版本裁剪的省文件可能和新版行政边界不完全匹配。我的经验是分省版适合快速看个大概正式成果里如果涉及边界区域一定再用标准行政边界重新裁剪全国版。如果只能拿到分省版至少要检查目标区域的边界和参考影像是否一致。5.5 年份与文件对应关系CLCD每年发布时文件名里的年份代表数据年份而不是发布日期。比如2025年发布的最新一期可能是CLCD_2025.tif也可能因为生产周期问题实际只更新到2024年。拿到数据后先检查文件数量和时间范围确认41期是否齐了。另外如果发现某一年缺失不要直接跳过先查一下官方说明确认是“该年份没有产品”还是“下载漏了”。5.6 下载后快速验证数据完整性文件下载完不一定会立即打开。我习惯用几行代码做快速体检确认文件没损坏、范围正确import os, rasterio for f in [CLCD_1985.tif, CLCD_2025.tif]: filesize os.path.getsize(f) / 1024 / 1024 with rasterio.open(f) as src: print(f, f{filesize:.1f}MB, src.crs, src.bounds)如果文件大小异常小比如只有几KB基本就是下载失败或不完整。如果地理范围显示为0到0说明坐标信息缺失这类文件最好重新获取。除了上面这些技术问题还有一个容易被忽略的点CLCD的客观局限。它是基于Landsat反射率自动分类的在山区阴影、云覆盖残留和地物混匀区域偶尔会把森林错分成灌木或者把水田错分成湿地。所以项目方法论里最好写明“以CLCD为基础数据并抽样验证”而不是把CLCD当成绝对真值。做小区域高精度研究时我通常还会叠加高分影像做参考。6. 写在最后41个tif文件带来的教训6.1 先定好坐标系再动手我在实际项目中最早踩过一回拿到全国版CLCD直接按WGS84统计某省会城市的建设用地面积出来的结果比统计公报多了两成排查半天才发现是经纬度坐标系直接乘平方米导致的。后来我给自己定了一条规矩凡是涉及面积指标一律先投影到等积投影再做统计。这个习惯看着不起眼却省掉了大量返工。6.2 不要迷信自动分类结果CLCD虽然好用但它毕竟是被动遥感自动分类的产物。我在做某区域耕地变化分析时候对照过高分影像发现个别年份的坡耕地和草地确实存在混淆尤其是靠近山地丘陵的地带。因此我现在拿到数据后一定会先抽三五个验证点检查分类是否符合当地常识再往下做转移矩阵和趋势分析。如果你是学生写论文这个“抽样验证”步骤还能直接写进研究方法部分提升可信度。最后再分享一个我自己的小技巧无论你用的是全国版还是分省版都保留一个“原始未处理”的备份目录不要直接在源文件上修改。CLCD年份多、文件大万一处理过程把坐标系写坏或者把值减掉了重新下载41个tif可太折腾了。愿这41期数据能帮你在城市蔓延、生态变化、耕地保护这些老问题里找到一个更清晰的答案。