当前位置:首页>python>空间数据炼金术 | 3-[Python] ETMonitor实际蒸散发数据批量处理(一)【日数据转换为月数据,全球】

空间数据炼金术 | 3-[Python] ETMonitor实际蒸散发数据批量处理(一)【日数据转换为月数据,全球】

  • 2026-10-21 13:27:07
空间数据炼金术 | 3-[Python] ETMonitor实际蒸散发数据批量处理(一)【日数据转换为月数据,全球】

01

数据处理介绍

1)原始数据存储方式与格式

① 存储方式

实际蒸散发日数据集存储在以年份命名的文件夹内,每个文件夹包含365或366个数据。同一年份数据,按照日序进行区别,如ETMonitor.GlobalET.2000001.v001.tif。原始数据详细信息和批量下载见:空间数据炼金术 | 2-[FTP]2000—2021年ETMonitor全球1 km陆地表面实际蒸散发数据批量下载。

② 数据格式

原始数据为栅格数据,数据格式为GeoTIFF,文件名后缀为.tif格式。

2)数据处理原则与流程

① 数据处理原则:时间跨度长、电脑内存占用小与数据处理效率高。

② 数据处理流程:

③ 结果数据详情:

a. 全球实际蒸散发月数据,存储在以年份命名的文件夹内。

b. 数据比例系数:0.01→0.1

注:与直接下载自官网的月数据集保持同样的比例系数,便于后期长时序实际蒸散发年和季节数据的获取

3)结果数据

① 2000年1月全球实际蒸散发分布

② 2020年1月全球实际蒸散发分布

③ 2021年1月全球实际蒸散发分布

02

数据批量处理

[Python]代码1

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103

import osimport globimport numpy as npimport rasteriofrom rasterio.windows import Windowfrom datetime import datetime, timedeltafrom tqdm import tqdmimport gc# ================= 配置参数 =================INPUT_ROOT = r'原始数据文件夹绝对路径,存放2000、2020和2021数据文件夹的上层文件夹路径'   # 输入路径OUTPUT_ROOT = r'结果数据文件夹绝对路径'SCALE_FACTOR = 0.1NODATA_VALUE = -1# 512 是最安全的内存块高度,兼顾速度与稳定性BLOCK_HEIGHT = 2048def get_date_from_doy(year, doy):    return datetime(int(year), 1, 1) + timedelta(days=int(doy) - 1)def process_one_month(year, month, file_list, out_path):    """单月处理逻辑:极其严苛的内存管理"""    if os.path.exists(out_path):        return    # 1. 仅打开第一个文件获取参考元数据    with rasterio.open(file_list[0]) as ref:        meta = ref.meta.copy()        width, height = ref.width, ref.height    meta.update(        dtype=rasterio.float32,         nodata=NODATA_VALUE,        compress='lzw', # 建议开启压缩,防止单个文件过大(3.7GB)导致磁盘写入压力        tiled=True,        blockxsize=512,        blockysize=512    )    with rasterio.open(out_path, 'w', **meta) as dst:        # 纵向分块        for y in range(0, height, BLOCK_HEIGHT):            h = min(BLOCK_HEIGHT, height - y)            window = Window(0, y, width, h)            # 预分配缓冲区            month_sum = np.zeros((h, width), dtype=np.float32)            any_valid = np.zeros((h, width), dtype=bool)            # 2. 依次处理每一天:打开 -> 读取 -> 累加 -> 关闭            for f_path in file_list:                with rasterio.open(f_path) as src:                    day_data = src.read(1, window=window)                    valid = (day_data != NODATA_VALUE)                    # 仅在有效区域累加                    np.add(month_sum, day_data, out=month_sum, where=valid)                    any_valid |= valid                    # 强制清理当前循环的日数据内存                    del day_data, valid            # 计算月均值并恢复空值            month_sum *= SCALE_FACTOR            month_sum[~any_valid] = NODATA_VALUE            dst.write(month_sum.astype(np.float32), 1, window=window)            # 手动清理块内存            del month_sum, any_valid            gc.collect()def main():    # 扫描文件    tasks = {}    year_folders = [d for d in os.listdir(INPUT_ROOT) if os.path.isdir(os.path.join(INPUT_ROOT, d))]    for y_str in year_folders:        tif_files = sorted(glob.glob(os.path.join(INPUT_ROOT, y_str, "*.tif")))        os.makedirs(os.path.join(OUTPUT_ROOT, y_str), exist_ok=True)        for f in tif_files:            try:                date_part = os.path.basename(f).split('.')[2]                dt = get_date_from_doy(date_part[:4], date_part[4:])                tasks.setdefault((y_str, dt.month), []).append(f)            except: continue    # 顺序处理每个月    sorted_months = sorted(tasks.keys())    for year, month in tqdm(sorted_months, desc="Overall Progress"):        out_name = f"ETMonitor.GlobalET.{year}{str(month).zfill(2)}.1km.v001.tif"        out_path = os.path.join(OUTPUT_ROOT, year, out_name)        print(f"\nProcessing {year}-{month}...")        process_one_month(year, month, tasks[(year, month)], out_path)        # 每月处理完彻底清理        gc.collect()if __name__ == "__main__":    # 强制限制 GDAL 的缓存大小,防止其私自占用过多内存    with rasterio.Env(GDAL_CACHEMAX=512):         main()

核心代码详解:

1)日序与日期的对应,实现数据按照年-月进行合并

① 提取年份和日序,如: 2020032

date_part = os.path.basename(f).split('.')[2] 

假设文件名:ETMonitor.GlobalET.2020032.v001.tif

split('.') 后得到列表 ['ETMonitor', 'GlobalET', '2020032', 'v001', 'tif']

取索引 [2] 得到 "2020032"

② 将年份和日序转换为年-月-日,调用函数:datetime, timedelta

def get_date_from_doy(year, doy):

用该年的1月1日作为基准,加上 (日序 - 1) 天

return datetime(int(year), 1, 1) + timedelta(days=int(doy) - 1)

具体2020001(2020年第1天): 2020-01-01+0天=1月1日

闰年处理:由于使用了标准库,它能自动识别 2020 年是闰年(2月有29天),因此日序 61 会正确指向 3月1日,而平年(如2021年)日序 60 就会指向 3月1日。

③ 年-月-日返回月份并“归桶”,调用函数:dt.month

dt 是一个 datetime 对象

dt.month 直接返回该日期所属的月份

tasks.setdefault((y_str, dt.month), []).append(f)

扫描完文件夹后,tasks[(2020, 1)] 里面就装满了 1 月份所有的 31 个日数据文件路径。

2)代码运行时显示数据处理进程,实时监控数据处理过程

温馨提示:由于需要通过一个代码实现全球范围数据处理操作,代码长度有些长;代码经作者调试可运行,直接使用即可。此外,该代码核心是完成按年-月求和、空间分块、实时清除中间文件,适用于其他变量栅格数据的操作,数据格式为GeoTIFF,期待大家能够实现代码的多样化用途!如有问题,欢迎各位在留言区留言,互相交流讨论,共同进步!

点个“看一看”吧

最新文章

随机文章