GEE下载的数据有空值,把空值设为固定值,比如-9999,然后再进行镶嵌,空值就不会压着分块数据了。
也适用于其他有空值范围的tif数据。
直接镶嵌
在这里插入图片描述

from osgeo import gdal
import numpy as np
import glob
import os

def modify_nodata_and_mosaic_gdal(input_files, output_file):
    # 修改每个文件的NoData值
    modified_files = []
    for file in input_files:
        # 读取数据
        ds = gdal.Open(file)
        band = ds.GetRasterBand(1)
        data = band.ReadAsArray()
        
        # 获取地理信息
        geotransform = ds.GetGeoTransform()
        projection = ds.GetProjection()
        
        # 将NaN替换为-9999
        data = np.nan_to_num(data, nan=-9999)
        
        # 创建临时文件
        temp_file = file.replace('.tif', '_temp.tif')
        modified_files.append(temp_file)
        
        # 创建新文件
        driver = gdal.GetDriverByName('GTiff')
        out_ds = driver.Create(temp_file, 
                             ds.RasterXSize, 
                             ds.RasterYSize, 
                             1, 
                             gdal.GDT_Float32,
                             options=['COMPRESS=LZW', 'TILED=YES'])
        
        # 设置地理信息
        out_ds.SetGeoTransform(geotransform)
        out_ds.SetProjection(projection)
        
        # 写入数据
        out_band = out_ds.GetRasterBand(1)
        out_band.SetNoDataValue(-9999)
        out_band.WriteArray(data)
        
        # 清理
        out_ds = None
        ds = None

    # 创建VRT进行镶嵌
    vrt_options = gdal.BuildVRTOptions(srcNodata='-9999', VRTNodata='-9999')
    gdal.BuildVRT('temp.vrt', modified_files, options=vrt_options)
    
    # 转换为最终的TIFF文件,添加压缩和分块选项
    translate_options = gdal.TranslateOptions(format='GTiff', 
                                            creationOptions=['COMPRESS=LZW', 
                                                           'TILED=YES',
                                                           'BLOCKXSIZE=256',
                                                           'BLOCKYSIZE=256'])
    gdal.Translate(output_file, 'temp.vrt', options=translate_options)
    
    # 添加金字塔
    ds = gdal.Open(output_file, 1)  # 1表示可写入模式
    ds.BuildOverviews("NEAREST", [2,4,8,16,32,64])
    ds = None
    
    # 清理临时文件
    os.remove('temp.vrt')
    for file in modified_files:
        os.remove(file)

# 设置输入输出路径
input_folder = r'E:\云量数据\其他数据\GD分析数据\MODIS_LAI_FPAR_2021\MODIS_MeanLAI_2021'
input_files = glob.glob(os.path.join(input_folder, '*.tif'))
output_file = os.path.join(input_folder, 'mosaic_result.tif')

# 执行镶嵌
modify_nodata_and_mosaic_gdal(input_files, output_file)
Logo

DAMO开发者矩阵,由阿里巴巴达摩院和中国互联网协会联合发起,致力于探讨最前沿的技术趋势与应用成果,搭建高质量的交流与分享平台,推动技术创新与产业应用链接,围绕“人工智能与新型计算”构建开放共享的开发者生态。

更多推荐