CSIF（Contiguous Solar-induced Fluorescence， 日光诱导叶绿素荧光）数据时间分辨率转换，在月尺度上聚合求均值（4天转为monthly）

问题：CSIF数据，月尺度上只有0.5°空间分辨率的

目标：将0.05°空间分辨率、4d时间分辨率的数据转成0.05° monthly

思路：

• 首先，读取路径下的文件
• 其次，将每个月的tif文件先归拢在一起（涉及到将DOY转为datetime
• 然后，基于每个月的tif文件求均值（注意使用np.nanmean，剔除nan值
• 最后，将均值输出为tif，每月一幅

代码如下：

# -*- coding: utf-8 -*-
"""
- 将每个月的tif归拢到一起，求平均值
- 每个月生成一幅tif输出
20240316, 15:45, YMJ
"""

import os
from datetime import datetime, timedelta
from osgeo import gdal
import numpy as np

# 文件路径
path = 'D:\\CSIF\\SIF_daily'
# 获取文件列表
files = os.listdir(path)

# 创建一个字典来存储每个月份的数据
monthly_data = {}

# 读取每个文件
for file in files:
    if file.endswith('.tif'):
        # 解析文件名以获取年份和DOY
        year = int(file[0:4])
        doy = int(file[4:7])

        # 使用datetime模块计算具体日期
        date = datetime(year=year, month=1, day=1) + timedelta(days=doy - 1)
        # 获取月份
        month = date.month

        # 使用GDAL读取.tif文件
        dataset = gdal.Open(os.path.join(path, file))
        array = dataset.ReadAsArray()

        # 将数据添加到月份字典中
        key = (year, month)
        if key not in monthly_data:
            monthly_data[key] = [array]
        else:
            monthly_data[key].append(array)

# 计算每个月份的平均值
for key, data in monthly_data.items():
    year, month = key
    monthly_mean = np.nanmean(data, axis=0)
    # print(monthly_mean[200,450])
    monthly_mean[monthly_mean < -10000] = np.nan

    # 创建一个新的.tif文件来保存平均值
    driver = gdal.GetDriverByName('GTiff')
    if month < 10:
        month = str("0") + str(month)
    else:
        month = str(month)
    output_path = f'D:\\CSIF\\SIF_monthly\\{year}{month}.tif'
    print(output_path)

    output_dataset = driver.Create(output_path, dataset.RasterXSize, dataset.RasterYSize, 1, gdal.GDT_Float32)
    output_dataset.GetRasterBand(1).WriteArray(monthly_mean)
    output_dataset.SetProjection(dataset.GetProjection())
    output_dataset.SetGeoTransform(dataset.GetGeoTransform())
    output_dataset.FlushCache()
    output_dataset = None