我有一个文件夹，内有1961-2020年逐年平均气温栅格数据，格式为dat，dat文件不是以年份命名的，而是一长串名字的最后4位是年份，现在我要对所有年份进行波段合成，合成一个dat，dat每个波段为不同年份的数据，且波段顺序按照年份排列，生成的dat文件能够在ENVI中打开。

import os
from osgeo import gdal# 设置工作目录和输出文件路径
input_folder = r"D:\itm\平均气温"
output_file = r"D:\itm\平均气温.dat"# 创建一个列表用于存储每年的栅格数据路径
input_files = []# 遍历文件夹内的文件
for file in os.listdir(input_folder):if file.endswith(".dat"):file_path = os.path.join(input_folder, file)input_files.append(file_path)# 对文件进行按年份排序
input_files.sort(key=lambda x: int(x[-8:-4])) # 根据你文件命名的年份位置进行调整# 创建虚拟栅格
vrt_options = gdal.BuildVRTOptions(separate=True)  # 设置为separate以确保每个文件都在单独的波段中
vrt_file = output_file + ".vrt"
gdal.BuildVRT(vrt_file, input_files, options=vrt_options)# 打开虚拟栅格
vrt_dataset = gdal.Open(vrt_file)if vrt_dataset is not None:# 创建输出栅格driver = gdal.GetDriverByName("ENVI")output_dataset = driver.Create(output_file, vrt_dataset.RasterXSize, vrt_dataset.RasterYSize, len(input_files), gdal.GDT_Float32)# 将每年的数据按顺序写入对应的波段for i, file_path in enumerate(input_files):year = file_path[-8:-4]  # 从文件路径中提取年份band = vrt_dataset.GetRasterBand(i + 1)output_band = output_dataset.GetRasterBand(i + 1)output_band.WriteArray(band.ReadAsArray())output_band.SetDescription(year)  # 设置波段描述为年份output_band.SetNoDataValue(9999)  # 设置无效值output_band.FlushCache()# 设置地理参考信息output_dataset.SetGeoTransform(vrt_dataset.GetGeoTransform())output_dataset.SetProjection(vrt_dataset.GetProjection())# 关闭数据集vrt_dataset = Noneoutput_dataset = Noneos.remove(vrt_file)  # 删除虚拟栅格文件
else:print(f"Failed to open the virtual raster file: {vrt_file}")