拓十年匠心定制 · 商业建站与技术教学双线并行 咨询热线:400-886-1026 service@lmnt.cn
ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

中国地面气候日值数据集V3.0处理避坑指南

中国地面气候日值数据集V3.0处理避坑指南

处理中国地面气候日值数据集(V3.0),是我这几年做过最锻炼耐心的事之一。这个数据集覆盖全国2400多个国家级地面站、时间跨度从1951年至今,做气候变化、水文模拟、农业气象、极端天气分析的人基本都绕不开它。但问题恰恰出在“好用”这两个字上:文件格式老、编码不统一、缺测值五花八门、单位还习惯用0.1的尺度,第一时间拿到的原始数据基本没法直接塞进模型。这篇文章把我反复踩过的坑和解决方案整理成一份避坑指南,希望能让后来的人少走点弯路。

文章适合的人很明确:正在对V3.0做数据清洗的气象数据处理工程师、气候相关方向的硕博研究生、以及任何需要把这份日值数据接入自己研究流程的从业者。如果你已经卡在某一步好几天了,直接跳到第六节的速查表先对号入座;如果你是第一次接触这个数据集,建议从头到尾过一遍,很多坑其实在第一步读文件的时候就埋下了。

1. 拿到数据先别急着跑:V3.0数据结构和最常见的低级错误

1.1 数据集包里到底有什么

中国地面气候日值数据集V3.0是由国家气象信息中心发布的逐日地面气象观测数据,数据源可以追溯到各个国家级地面气象站的自动站和人工站观测记录。整个数据集解压后,核心内容通常是这几样:大量以台站号命名的文本文件、一个台站信息表、以及一份说明文档(有时候是readme.txt,有时候是PDF或word文档)。

站点文件长这样:54511.txt、58362.txt,文件名就是区站号。每个文件里是该站从建站到更新日期之间的逐日观测记录,一行一天。这个组织方式看起来很清爽,但问题在于它没有数据库导出的那种规范表头,而是纯文本的定宽格式加空格分隔,所以读取方式稍微不对,整张表就废了。

台站信息表非常重要,里面包含每个站的区站号、站名、省份、纬度、经度、海拔高度、资料起止时间等元信息。很多人的处理过程里,数据文件都读对了,结果在合并站点信息的时候翻车,后面我会专门讲这块。

还有一个很多人会忽略的点:数据集会不定期更新增量,比如你年初下载了一份V3.0,年中再下一次,目录里的文件数量和时间范围可能都不一样。建议每次下载后先核对一下说明文档里的版本日期和文件数量,别用旧版本的台站信息表去匹配新版本的数据文件。

注意:说明文档是整份数据集最重要的文件,没有之一。里面的编码规则、缺测标记、质控码含义,不同子版本之间可能有细微差异。网上很多教程写的是V2.0甚至更早的习惯,照搬过来会踩坑。

1.2 编码、分隔符与列名:三个让你第一行代码就崩的细节

先说编码。V3.0的文本文件大多使用GBK或者ANSI编码,不是UTF-8。在Linux服务器上直接pd.read_csv('54511.txt', encoding='utf-8'),大概率会收到一个UnicodeDecodeError。我第一次碰到这个报错还以为是文件损坏了,后来才发现是编码不对。解决方案很简单:

import pandas as pd # 用gbk或者gb18030读取,gb18030兼容性更好 df = pd.read_csv('54511.txt', encoding='gb18030')

如果连gb18030都有部分文件报错,可以用encoding='ISO-8859-1'先兜底。不过这种做法的风险是中文列名可能变成乱码,只适合用来救急。正常的数据文件主体全是数字和空格,不需要依赖中文解析,所以用gb18030就够了。

再说分隔符。V3.0的站点文件不是逗号分隔,也不是严格的Tab分隔,而是空格混合Tab,字段之间有时候还连续好几个空格。如果你直接用默认的sep=','去读,pandas会把整行当成一列,惨不忍睹。正确做法是使用正则表达式\s+把连续的空白字符都当作分隔符:

df = pd.read_csv('54511.txt', sep=r'\s+', header=None, encoding='gb18030')

再说列名。绝大多数站点文件没有表头,header=None这个参数不能少。读进来之后,需要根据说明文档里的字段说明,给每一列赋予有意义的名称。V3.0不同时期的列数不完全一样,有的版本是30列左右,有的扩展到34列,包含平均气温、最高最低气温、气压、水汽压、相对湿度、风速、降水量、日照时数、最大风速极大风速及对应的风向和时间等字段。处理时一定要以你手上那份说明文档的字段顺序为准,先数清楚实际列数再去命名。

2. 数字本身会“撒谎”:解码缺测值、特殊值与单位换算

2.1 32766不是气温,是缺测

V3.0数据里最著名的“陷阱”就是缺测值。气温要素的缺测值通常是32766,气压、相对湿度、风速等要素也有类似的编码。如果你不做任何处理,直接把原始数值除以10还原单位,那一个缺测的日平均气温就会变成3276.6℃,极端值检查直接爆炸。

正确的处理顺序是先识别缺测值、替换成缺失标记,然后再做单位换算。顺序反了,数值就再也回不去了。举个例子:

import numpy as np # 方法一:全局替换(谨慎使用,确保这些值在降水等特殊字段没有特殊含义) df.replace(32766, np.nan, inplace=True) # 方法二:按列精确替换(推荐) temp_columns = ['tem_avg', 'tem_max', 'tem_min'] for col in temp_columns: df.loc[df[col] == 32766, col] = np.nan

在我处理过的数据里,气温缺测值基本是32766,但降水字段的缺测可能出现在32000至32766这个区间,日照时数有时会用到32700附近的编码。稳妥的做法是不要一刀切地把所有大于30000的数值全变成缺失,而是先分类、查文档、再处理。我自己的习惯是写一个函数,把每个字段的取值分布先跑一遍,看到异常集中的数值特征,再去对照说明文档确认它的含义。

注意:有些站点在个别年份会有一整段“数据缺失”,但原始文件里并不是空行,而是一整排缺测值。统计有效天数的时候,一定要先完成缺测替换,否则有效天数会虚高。

2.2 降水的“31”和“32700”家族:别把特殊天气当真实雨量

这是全网讨论最少、但实际影响最大的一个坑。在V3.0的降水要素里,很多数值并不是真实的降水量,而是天气现象的编码。最常见的两个是31和32700附近的数值。

先说“31”。在部分版本的V3.0里,20-20时降水量字段中的31表示“微量降水”,也就是雨量计能捕捉到降水痕迹,但累积量不足0.1毫米。如果你看到某天降水量是31,别急着记成31毫米(或者还原后是3.1毫米),它其实应该按0毫米参与统计,或者单独标记为“痕迹降水”。

再看“32700”这个区间。在降水字段中,32700、32701、32702这一类的数值,可能对应着雾、露、霜、雪暴、雾凇等特殊天气现象编码。这些也不是真实的降水量。如果把这类值当成降水参与总量计算,累计偏差会非常离谱。我遇到过一位做极端降水研究的同事,某年某站的年降水量比周边站高出几十毫米,排查到最后就是这一批特殊值没清理干净。

推荐处理逻辑如下:

# 假设pre_20_20是20-20时降水量字段 # 31及相似编码表示痕迹降水,映射为0 df.loc[df['pre_20_20'] == 31, 'pre_20_20'] = 0.0 # 大于等于32000的值,一般是特殊现象或缺测,统一转为NaN df.loc[df['pre_20_20'] >= 32000, 'pre_20_20'] = np.nan

这套逻辑不是标准答案,因为不同子版本的编码含义有差异,但思路是对的:把“特殊值”和“真实数值”分开处理,绝不能让编码值混进统计口径。

2.3 除以10还是除以100?各单位换算一表理清

V3.0里大量要素为了压缩存储空间,使用了0.1作为计量单位。如果你只对气温做了除以10的操作,而降水量、风速、日照时数都忘了,后续所有的统计结果都是歪的。我把常见要素的换算整理成一个表,处理前对着看一眼:

要素原始单位需要转换成的单位操作
平均气温0.1℃℃除以10
最高/最低气温0.1℃℃除以10
平均气压0.1hPahPa除以10
平均水汽压0.1hPahPa除以10
平均相对湿度1%%不换算
平均风速0.1m/sm/s除以10
20-20时降水量0.1mmmm除以10
日最大降水量0.1mmmm除以10
日照时数0.1hh除以10
日最大风速0.1m/sm/s除以10
日极大风速0.1m/sm/s除以10
风向1(度)度不换算

具体的字段名在不同版本里可能叫法不一样,但原则是一致的:说明文档里写了“0.1”开头的单位,就要除以10。我自己每次新建一个处理脚本,都会把这张表贴在头部当注释,提醒自己别犯低级错误。至于“除以100”的传言,我目前还没有在V3.0里遇到过需要除以100的常规要素,气温、降水、风速这些都遵循除以10的规则。

3. 站点信息里藏着坑:经纬度、站号与台站数据的拼接

3.1 dddmm格式的经纬度:不是十进制度

台站信息表里给出的经纬度,不是常见的十进制度,而是“度分”的拼接写法。比如一个站的纬度写着3954,这里的正确读法是北纬39度54分,而不是十进制的39.54度。经度11623对应的是东经116度23分,也不是116.23度。

如果有人直接用原始数值去做空间插值或者距离计算,站点位置会偏移几十公里甚至更多,这种误差对于区域气候分析来说是致命的。转换公式很简单:

def dms_string_to_decimal(value): """ 将度分拼接格式转换为十进制度。 例如 3954 -> 39 + 54/60 = 39.9 11623 -> 116 + 23/60 = 116.3833 """ value = int(value) degrees = value // 100 minutes = value % 100 return degrees + minutes / 60

注意,有些版本里经纬度是整数,有些版本里可能带了小数或者负号。西南地区一些站点的经度可能超过10000(比如经度100度以上),但严格按照公式处理就没问题。另外,海拔高度在台站信息表里通常是“米”为单位,不需要额外换算,但也不要忽略它,因为很多插值算法需要用到高程修正。

3.2 台站号请当字符串处理

区站号是5位数字,比如北京站是54511,广州站是59287。很多人在读取的时候为了方便直接指定dtype=int,结果导致文件名字符串和数据里的站号匹配不上。看起来是小事,实际非常容易阴沟翻船:当你把站点数据和站点信息表按station_id合并时,一边是整型54511,一边是字符串'54511',pandas会直接给你一个全是NaN的合并结果。

我的建议是从头到尾都把区站号当成字符串处理:

df = pd.read_csv('54511.txt', sep=r'\s+', header=None, dtype={0: str}, encoding='gb18030')

如果已经在某个环节变成了整型,补救方法是临时astype(str).str.zfill(5),把可能丢失的前导零补回来。虽然中国大部分区站号不是0开头,但养成这个习惯能避免在其他数据集上踩同样的坑。

3.3 台站迁移与数据不连续:合并多源数据时的误判

中国的地面气象站不是所有站点都从1951年观测到现在。很多站点经历过迁站、停测、恢复观测,台站信息表里的经纬度和海拔,通常代表的是该站最新的位置。如果你直接用最新经纬度去代表这个站所有历史时段的位置,那早期数据的位置就全错了。

处理方法是把台站信息表和日值数据合并之后,额外做一步检查:对每个站点的经纬度时间序列求变化量,如果某年突然跳变了几十公里,很可能是迁站了。严谨的研究会在论文里特别说明迁站时间和新老位置,或者干脆把迁站前后的数据视为两个不同的站点。

另外,早期站点数量远少于现在,比如1951年可能只有零星几十个站有观测,到20世纪60年代才逐渐形成全国覆盖网络。所以分析长时间序列的时候,不要笼统地说“全国站点平均值”,先看清楚每个时段的站点数量,否则统计结果会被站点分布变化严重影响。

4. 时间序列处理:日界、质控码和边界条件的实战细节

4.1 20-20时降水日界带来的时间归属问题

V3.0里的降水量字段经常标注着“20-20时”,意思是这个日降水量统计的是前一日20时到当日20时之间的累积量。比如“2021年7月20日”这一行的20-20时降水量,实际包含了7月19日晚上的降水。这是气象观测的一个老传统,跟日常生活自然日(0点到24点)的认知不一样。

对于以日为单位的研究,比如做月总量、年总量,这种半天偏移影响不大。但如果你要拿日值降水去匹配小时尺度的天气过程,或者做逐日气象和水文变量的联合分析,这个日界就必须考虑进去。稳妥的做法是保留原始字段名里的“20-20”标记,不要简单地改名为precipitation,免得三个月后回看代码时误以为它是自然日降水。

4.2 质量控制码:0/1/2/3/8分别意味着什么

V3.0给每个核心要素都配了一个质量控制码字段,直接放在要素值的旁边。通用规则大致是:0表示未进行质量控制(或者质量正确,视版本而定),1表示数据正确,2表示数据可疑,3表示数据错误,8表示缺测。不同版本细节会有出入,还是要看说明文档。

处理策略上,我一般用质控码来辅助过滤:

# 假设质量码字段名为 tem_avg_qc # 1为正确,0为未质控但可参考,2以上先谨慎处理 valid_mask = df['tem_avg_qc'].isin([0, 1]) df.loc[~valid_mask, 'tem_avg'] = np.nan

有一点要提醒:不要看到可疑值直接删掉整行,因为一行里可能有多个要素,一个要素的可疑不代表所有要素都不可用。逐字段过滤比整行删除更能保留有效信息。另外,自动站普及之后(2003年前后),数据质量整体有所提高,但在早期人工观测阶段,质控码的含义和使用要更保守一些。

4.3 闰年、缺测月份和59年前的老数据

时间索引是数据分析的第一步,但V3.0的老数据里偶尔会出现一些边缘日期问题。直接对“年-月-日”三列执行pd.to_datetime,遇到2月30日这种非法日期就会报错。用errors='coerce'可以把解析失败的行变成NaT,不至于让整个脚本中断:

dates = pd.to_datetime({ 'year': df['year'], 'month': df['month'], 'day': df['day'] }, errors='coerce') df['date'] = dates df = df.dropna(subset=['date'])

跑完日期解析之后,建议顺手检查一下每个站点每年的记录条数,正常应该接近365或366。如果某年只有100多天或者400多天,那大概率是文件解析出了问题,或者该站当年确实只观测了一部分时段。1951年建站初期的数据尤其要注意,很多站是从年中才开始有记录的,不要简单地在图上画一条断头线就完事了。

5. 一套能直接跑的清洗流程(Python + pandas)

5.1 读取单个站点的标准写法

前面说了很多原则,这里给一套完整可用的示例代码。假设你的数据目录里有一个54511.txt,且说明文档确认了字段顺序和我的示例一致,那么可以这样读取:

import numpy as np import pandas as pd # 按实际版本调整列名和顺序,这里是常见V3.0字段示意 columns = [ 'station_id', 'year', 'month', 'day', 'tem_avg', 'tem_avg_qc', 'tem_max', 'tem_max_qc', 'tem_min', 'tem_min_qc', 'pre_20_20', 'pre_20_20_qc', 'pre_max_day', 'pre_max_day_qc', 'rh_avg', 'rh_avg_qc', 'win_avg', 'win_avg_qc', 'sun_dur', 'sun_dur_qc', 'win_max', 'win_max_qc', 'win_max_dir', 'win_max_dir_qc', 'win_ins', 'win_ins_qc', 'win_ins_dir', 'win_ins_dir_qc', ] df = pd.read_csv( '54511.txt', sep=r'\s+', header=None, names=columns, dtype={'station_id': str}, encoding='gb18030', )

如果你的实际文件列数比这个多或少,以你手上的说明文档为准,把columns列表调整到和文件列数一致再运行,否则pandas会因为列数不匹配直接报错。

5.2 缺测替换、单位还原和质量过滤

读取之后,接下来就是最重要的清洗环节。我的建议是在这一步生成两个版本:一个是处理后的分析版DataFrame,一个是单独保存了原始值的备份版,方便随时回头核查。

# 备份 raw_df = df.copy() # 1. 缺测值替换:气温、风速、湿度等,32766为标准缺测 for col in ['tem_avg', 'tem_max', 'tem_min', 'win_avg', 'rh_avg', 'sun_dur']: df.loc[df[col] == 32766, col] = np.nan # 2. 降水特殊编码处理 # 31 或类似值代表痕迹降水,记0;>=32000的编码视为缺测 df.loc[df['pre_20_20'] == 31, 'pre_20_20'] = 0.0 df.loc[df['pre_20_20'] >= 32000, 'pre_20_20'] = np.nan # 3. 单位还原,一次性处理所有0.1单位字段 for col in ['tem_avg', 'tem_max', 'tem_min', 'pre_20_20', 'pre_max_day', 'win_avg', 'sun_dur', 'win_max', 'win_ins']: df[col] = df[col] / 10.0 # 4. 质量控制码过滤 for col in ['tem_avg', 'tem_max', 'tem_min', 'pre_20_20', 'rh_avg', 'win_avg']: qc_col = col + '_qc' if qc_col in df.columns: df.loc[~df[qc_col].isin([0, 1]), col] = np.nan # 5. 时间索引 df['date'] = pd.to_datetime({ 'year': df['year'], 'month': df['month'], 'day': df['day'] }, errors='coerce') df = df.dropna(subset=['date']) # 6. 排序并设为索引 df = df.sort_values('date').reset_index(drop=True) df = df.set_index('date')

清洗结束后,做一个极值合理性校验。比如气温应该在-60到50℃之间,降水应该大于等于0,相对湿度应该在0到100之间。如果出现离谱的值,说明缺测替换或单位换算有遗漏。

print(df[['tem_avg', 'pre_20_20', 'rh_avg']].describe())

5.3 批量处理全国站点:效率与内存经验

全国两千多个站点,每个站点一个文件,如果全部读入内存再合并,普通笔记本很容易内存爆掉。更现实的做法是循环处理每个文件,先把单个站清洗后的结果持久化到磁盘,最后再统一合并或直接分站使用。

import pathlib import glob data_dir = pathlib.Path('path/to/v3.0/') out_dir = pathlib.Path('path/to/cleaned/') out_dir.mkdir(exist_ok=True) for file_path in glob.glob(str(data_dir / '*.txt')): station_id = pathlib.Path(file_path).stem try: df = read_and_clean(file_path) # 封装上面的流程 except Exception as e: print(f'Failed: {station_id}, {e}') continue if df.empty: continue df.to_parquet(out_dir / f'{station_id}.parquet')

Parquet格式压缩率高、读取快,非常适合这种按站存储的清洗结果。后续不管是用Dask还是直接按需读取,都比反复解析原始文本高效得多。如果你不熟悉Parquet,退一步用to_csv也不丢人,但要做好文件体积膨胀的心理准备。

整个批量过程里,我会额外记录一个日志文件,写明每个站读取了多少行、清洗后剩多少行、有没有异常报错。这样一来,全国站点跑完一遍,我不仅拿到了干净数据,还拿到了一份数据质量清单,后面写文章、答辩都能用上。

6. 常见问题速查表与排查思路

6.1 6个高频报错速查

问题现象直接原因解决方案
读取报UnicodeDecodeError文件是GBK/ANSI编码使用encoding='gb18030'
读出来只有一列,数字挤在一起分隔符不是逗号使用sep=r'\s+'
气温出现327.66℃这类离谱值缺测值32766没替换就做了单位换算先replace(32766, np.nan)再除以10
降水量出现31、32700等数值特殊天气现象编码被当成真实降水按说明文档转化为0或缺失
按站号合并数据全是NaN站号一边是字符串一边是整型统一转为字符串并zfill(5)
日期解析报OutOfBoundsDatetime或非法日期老数据存在2月30日等异常日期使用pd.to_datetime(..., errors='coerce')并过滤

这些都是我真实踩过或者帮别人排查过的问题,频率最高的还是前面三条,基本上每接触一个同事的项目就能看到一次。

6.2 排查方法论:先写探针脚本

与其等到数据跑进模型才发现异常,不如在清洗阶段就建立一个“探针脚本”。所谓探针,就是快速地读取几个大小不同的站文件,统计每列的最小值、最大值、非空数量、唯一值数量,再打印出来给你看。通过这组数值,你能很快定位到是缺测没清干净、单位没换算、还是列顺序本身就对错了。

举个例子,如果某列探针输出的最小值是0,最大值是3276,中间几乎全是几十到几百,那基本可以断定这列是一个0.1单位的物理量,且3276这个最大值对应的原始值就是32766。如果某列降水探针里出现了几千个31,那就说明痕迹降水编码在数据里占比很高,后续统计降水量时必须特殊处理。

处理V3.0数据集这段时间,我最大的感受是:这个数据集本身质量并不差,差的是它的“表达习惯”和现在大家习惯用的CSV、Excel差太远了。只要摸清了编码规则、缺测标记、单位体系和质控码的含义,整个清洗流程是非常机械的。反过来说,很多人在这份数据上翻车,往往不是输在算法上,而是输在对说明文档的敬畏不够。我后来每次新建项目,都会把readme里的字段表单独抽出来放到项目仓库里,处理完数据之后还要留一份处理日志,这样哪怕过半年再回来看,也能很快回忆起当时每一步改了什么、为什么这么改。

返回列表