第一次拿到这份“中国1970–2022年碳排放栅格数据集”时我的第一反应就是终于不用再去东拼西凑各个年份的排放清单了。一个数据集把52年、逐年、10km分辨率的CO₂排放全部打包成GeoTIFF还是EDGAR官方源这对做气候变化、城市减排、行业研判的科研党来说基本属于“开箱即用”的宝贝。但话说回来数据好拿不代表上手就容易。我身边有不少人下载完以后解压出来面对几十个tif文件直接懵不清楚单位、不知道坐标系、更不知道怎么裁剪成中国区域、怎么算年际变化。这篇文章我就用一份可复现的实操思路把这个数据集从文件结构到可视化、从裁剪到时序统计一步步拆开讲清楚。不管你是学生、研究员还是相关行业的从业者照着往下走至少能少走一大半弯路。1. 先把这个数据集拆开看它到底给了你什么1.1 纵向52年的连续观测尺度这个数据集覆盖1970年到2022年整整52个年份。很多研究者做长周期分析的时候最头疼的就是“数据断层”——有的数据库只更新到2014年有的干脆只有个别年份。这组数据最值钱的地方是它把改革开放前后、经济高速增长期、双碳目标提出后的几个关键阶段全部串起来了。你想看1970年代的低排放基线行想看2000年以后的快速爬坡行想看2015年之后有没有出现平台期甚至拐点也行。时间维度拉满以后做趋势分析、突变检测、政策效果评估底气和单一年份的数据完全不是一回事。1.2 10km栅格的空间结构所谓10km分辨率在这个数据集里实际对应的是一个约0.1度经纬度网格。整个中国区域大概会划分成几十万个有效格网。别小看这个尺度它比省市级统计公报精细得多又不像百米级数据那样动辄几个GB、处理起来直接把内存吃满。10km精度能支撑什么城市群级别的排放差异能看得很清楚一个省内不同地市的高排放聚簇也能分辨出来。同时它又足够“轻便”用一台普通笔记本就能跑动全流程。比如山东省和河北省之间、长三角和珠三角之间到底谁的排放强度更高、空间上呈现什么形态用这个数据直接能画出来。1.3 EDGAR官方源背后的可靠逻辑EDGAR是国际权威机构长期维护的全球温室气体排放数据库被大量顶刊文章引用方法论相对稳定且透明。它的逻辑是通过部门活动数据、燃烧技术参数和排放因子来估算排放所以世界各国的数据口径基本一致这一点非常重要——做全球对比或者跨国分析时最怕的就是各国上报数据的标准都不统一。直接用EDGAR源意味着你不需要自己再去做一堆数据匹配、单位换算和缺失年份插值。它已经把“逐年、逐网格、逐部门”的排放数据整理到位你拿到手之后更多是聚焦自己的研究问题而不是把时间浪费在数据清洗上。当然它不是没有缺点后面我会讲到边界和局地偏差的问题。2. 拿到数据之后文件结构、坐标系和单位才是第一道门槛2.1 一探GeoTIFF文件内部GeoTIFF是遥感领域最常见的一种栅格格式本质是在普通TIFF图像里嵌入地理位置信息。你下载完解压后大概率看到一堆以年份命名的tif文件名字可能类似CO2_1970.tif、CO2_1971.tif这种形式。每个文件只有一个波段波段里的每个像元值就是那个格网当年的CO₂排放量。注意这里有个常见误区很多新手以为打开文件直接就是一张“好看的图”实际上默认拉伸出来往往一片黑或者一片白。因为原始数据是浮点型值的范围极大必须做符号化设置或者数值裁剪才能正常显示。也就是说拿到数据的第一步不是“看图”而是“读数值”。2.2 EDGAR坐标系与CRS处理这个数据集用的坐标系通常是WGS84地理坐标系也就是我们常说的经纬度坐标系EPSG代码是4326。这意味着数据的行列位置是跟经纬度直接对应的而不是像投影坐标系那样用米作单位。正因为是经纬度网格0.1度在赤道附近大概是11km左右但到了中国北方高纬度地区经度方向的实际距离会明显收缩。所以严格说10km只是一个近似称呼。你在做面积统计或者距离计算的时候最好还是投影到适合中国区域的Albers等积投影或者UTM分区投影里不然算出来的面积会有偏差。这个问题我遇到不止一次后面会细说。2.3 单位到底怎么理解单位是整个数据集里最容易翻车的环节。很多栅格产品的单位是“吨/格网/年”或者“千克/格网/年”甚至可能是“吨/年/平方公里”。如果不先确认单位直接做后续计算最后总和跟官方统计一对比可能会差好几个数量级。我的习惯是拿到文件先读一下元数据再随机取几个像元值做个数量级判断。如果像元值是几百上千大概率是吨如果是零点几、几十可能是千克或者经过面积换算后的强度值。最稳妥的办法是读取整个栅格求和再跟EDGAR官方发布的该年中国排放总量做对比误差在个位数百分比以内就说明单位判断正确。3. 实操用Python从零读取一张年度CO₂栅格3.1 环境搭建与依赖工欲善其事必先利其器。读取GeoTIFF最推荐的工具是rasterio它是在GDAL之上封装的Python库API对新手特别友好。再配合geopandas做矢量裁剪matplotlib做可视化numpy做数值计算一套组合拳基本够用了。pip install rasterio geopandas matplotlib numpy装完以后建议先用rasterio打开一个文件感受一下数据集的基本属性。下面这段代码可以读取2020年的栅格打印出CRS、分辨率、数据范围、NoData值并做一个最简单的统计。import rasterio import numpy as np path CO2_2020.tif with rasterio.open(path) as src: print(CRS:, src.crs) print(宽度:, src.width, 高度:, src.height) print(空间范围:, src.bounds) print(NoData值:, src.nodata) data src.read(1).astype(float64) valid data[data 0] print(有效像元个数:, valid.size) print(排放总和(吨):, valid.sum()) print(平均值:, valid.mean())如果你运行以后发现有效像元个数少得可怜别急着怀疑人生——很可能是NoData值参与了统计被当成0处理了这时就要用掩膜把无效值剔除。3.2 处理NoData和异常值在栅格数据里NoData一般用-9999、-3.402823e38等特殊值表示部分文件内部还会直接把它标为nodata。但当你用普通方式读取时这些值不会自动变成NaN而是原样保留直接计算结果就有问题了。处理办法很简单读取的时候设一个掩码。with rasterio.open(path) as src: data src.read(1) nodata src.nodata if nodata is not None: data np.where(data nodata, np.nan, data)还有一种情况是数据里出现了很大的负值这可能是海洋掩膜或者无数据区的残留。严格来说CO₂排放量不可能为负所以遇到负值可以直接统一替换成NaN。完成这步之后再算总和就有意义了。3.3 快速可视化第一眼看到空间格局数据读进来以后直接画图是不可能的因为值域跨度过大。最常用的处理方式是取对数或者用分位数把异常高值压下去。我个人倾向用对数拉伸因为排放分布本身就是偏态分布胡焕庸线两侧的差异非常明显。import matplotlib.pyplot as plt from matplotlib.colors import LogNorm with rasterio.open(CO2_2020.tif) as src: data src.read(1).astype(float64) data[data 0] np.nan fig, ax plt.subplots(figsize(10, 8)) im ax.imshow(data, cmapRdYlBu_r, normLogNorm(vmin1, vmax1e4)) plt.colorbar(im, labelCO2 emissions (t/grid/year)) ax.set_title(China CO2 Emissions 2020 (EDGAR)) plt.show()画完以后你基本能一眼看到东部沿海排放强度明显高于西部京津冀、长三角、珠三角几个热点区域非常醒目。这种直观展示是统计数据做不到的。4. 进阶中国区域裁剪、多年度时序分析与统计4.1 用GeoPandas裁剪出中国边界EDGAR是全球数据文件里肯定包含了周边国家和地区。做中国研究时第一步就是把中国区域裁剪出来排除境外排放对统计的影响。裁剪的方法无外乎两种一是直接用经纬度范围做粗裁二是用中国矢量边界做精确裁剪。粗裁适合快速预览但会把印度、尼泊尔、俄罗斯、蒙古等国家的排放一并框进来。精确裁剪需要一份中国国界矢量文件用geopandas读取后配合rasterio的mask功能实现。import geopandas as gpd from rasterio.mask import mask china gpd.read_file(china_boundary.shp) with rasterio.open(CO2_2020.tif) as src: out_image, out_transform mask(src, china.geometry, cropTrue) out_meta src.meta.copy() out_meta.update({ driver: GTiff, height: out_image.shape[1], width: out_image.shape[2], transform: out_transform }) with rasterio.open(CO2_2020_china.tif, w, **out_meta) as dst: dst.write(out_image)需要注意这份矢量边界最好与你研究时段匹配因为国界线在几十年间有调整。如果只是做粗略趋势研究用最新边界问题也不大但如果是做国界敏感分析务必查清楚边界数据的时效性。另外裁剪以后务必检查一下格网数量和总和的变化以防裁掉不该裁的区域。4.2 52年逐年年际动态循环读取与折线图有了单年份的处理流程做多年度时序就简单了。写一个循环逐年读取裁剪后的栅格计算总排放量然后输出成一条折线。这样你就可以看到中国从1970到2022年碳排放的大致变化曲线。years range(1970, 2023) totals [] for year in years: with rasterio.open(fCO2_{year}_china.tif) as src: data src.read(1).astype(float64) data[data 0] np.nan totals.append(np.nansum(data)) plt.figure(figsize(12, 5)) plt.plot(years, totals, linewidth2) plt.xlabel(Year) plt.ylabel(Total CO2 emissions) plt.title(China annual CO2 emissions from EDGAR) plt.grid(alpha0.3) plt.show()画出来以后你大概率会看到1970年代低位平稳1990年代中期开始明显爬升2000年后进入陡增阶段近几年增速放缓和波动。这种长时序曲线放在报告里比任何文字描述都有说服力。另外你还可以算同比增长率、五年均值、以及不同阶段年均增速都是论文里的加分项。4.3 省域汇总让栅格数据与行政区域结合如果仅仅看全国总量10km分辨率就有点大材小用了。把栅格数据分区统计到省级甚至地市级才能真正释放空间数据的价值。做法是先把行政区划矢量文件读进来通过空间连接给每一个格网匹配省份代码再按省份聚合。一个省覆盖多个格网部分格网会跨越省界所以严格起见要做面积加权或者重心归属。简单学术研究可以直接用“栅格中心点落在哪个省就归哪个省”的规则误差可以接受。import pandas as pd from shapely.geometry import Point # 将栅格转为点或直接使用经纬度提取省份属性 x, y np.meshgrid(np.linspace(src.bounds.left, src.bounds.right, src.width), np.linspace(src.bounds.top, src.bounds.bottom, src.height)) points [Point(xi, yi) for xi, yi in zip(x.ravel(), y.ravel())]不过说实话这种点归类方法对10km网格来说够用但对边境省市可能会有边界效应。真要发布严谨结果还是用gis软件里的Zonal Statistics as Table工具更靠谱。Python里也有rasterstats库函数名叫zonal_stats一行就能搞定。5. 常见问题与排查技巧实录5.1 内存爆炸全图读取直接死机这个数据集虽然只有0.1度分辨率但全球范围的tif文件尤其某些年份还自带多波段或高压缩比一次性read(1)可能占用几百MB甚至上GB内存。如果是老笔记本很容易直接卡死。解决方法是窗口读取。rasterio支持按地理范围读入比如我只是想看华东区域就先设定范围用src.read(1, window...)。或者用rasterio.windows.from_bounds把经纬度范围转换成窗口索引再读取局部数据。这样内存占用马上降下来处理速度反而更快。from rasterio.windows import from_bounds west, south, east, north 110, 20, 130, 40 with rasterio.open(path) as src: w from_bounds(west, south, east, north, src.transform) data_local src.read(1, windoww)5.2 数值对不上跟官方总量差太多我踩过最大的坑就在这。用EDGAR源做全国汇总算出来总和跟某些年度统计公报的数值怎么都对不上。原因基本有三类一是单位理解错了二是裁剪范围包含了周边国家三是没有排除NoData或者负值。排查步骤很简单先看单格网的数值量级再看有效像元总数和覆盖面积最后和iea或者ceads等参考数据库做对比。只要偏差在合理范围内一般是5%以内就不用太担心。如果差了一倍以上优先检查NoData处理。5.3 可视化时一片空白或颜色诡异拿到tif直接imshow经常是一片空白这不代表数据有问题。浮点栅格的值域跨度太大默认线性拉伸会把你关注的低值区域全部压到同一色阶里。用LogNorm或者Percentile stretch就能解决。如果用了对数变换还是有大片空白那就要看数据是否在某个值以下全是NaN。另一个小技巧是设置合适的颜色条范围比如在主要人口密集区排放量通常在几十到几千吨/格网/年你只要卡在1到100之间细节就会很丰富。5.4 边界锯齿明显、海岸线粗糙10km分辨率数据画出来海岸线会有明显的锯齿感这是栅格数据的天性在放大输出的时候尤其扎眼。如果你做地图出版级别的图件建议叠加矢量海岸线作为底图或者把栅格图层的透明度调低配合自然地理底图一起展示科学感立刻提升。常见问题速查表问题描述常见原因处理办法总量比官方统计大一截未排除周边国家或NoData精确裁剪、掩膜删除NoData总量小好几个数量级单位判断错误与官方总量对比校准打开图像全黑/全白值域跨度大或比例尺不当使用对数拉伸或分位数拉伸内存不足一次性读取全图窗口读取或分块处理边界有锯齿分辨率偏低叠加矢量底图或重采样经纬度范围对不上坐标系没搞清检查CRS统一到4326写在最后的小技巧我个人在实际操作中发现这份数据最适合做的第一张图不是全国平面分布而是“1970–2022年各省排放热力趋势图”把每年各省排放总量标准化以后按年份排成一条横向热力带颜色的冷暖变化能把整个减排进程讲得非常清楚。另外因为这组数据时间跨度长做逐年格网排放变化率的时候尽量别直接算相邻两年的比值年际波动太大会掩盖真实趋势建议做五年滑动平均或者用Theil-Sen趋势估计结果稳定很多。还有一个小建议如果你要发表论文记得把数据处理流程和精度评估写清楚EDGAR数据集本身也是有版本差异的年份不同可能算法微调过这就需要前后保持版本一致性。这份数据集的潜力远不止做个折线图。配合人口格网、GDP格网、夜间灯光数据你可以做排放的空间驱动因子分析识别高排放低产值区域配合气象数据可以做污染扩散的排放清单输入配合未来情景数据还能做减排路径模拟。它所提供的这种长序列、高分辨率、标准化口径的“底座”才是真正拉开研究差距的地方。
看内容编号体系：设计逻辑与实操方法)
