GDAL 实现影像合并

0 阅读5分钟

^ 关注我,带你一起学GIS ^

前言

GDAL作为地理空间数据处理的核心工具,其影像合并功能为多源栅格数据的集成与分析提供了高效、灵活的解决方案。无论是遥感影像镶嵌、地图瓦片拼接,还是时间序列数据的融合,该功能能够帮助用户将分散的影像片段整合为具有统一地理参考和连贯信息表达的完整数据集,为后续的空间分析、可视化及应用构建可靠的数据基础。

由于本文由一些前置知识,在正式开始之前,需要你掌握一定的Python开发基础和GDAL的基本概念。在之前的文章中讲解了如何使用GDAL或者ogr2ogr工具将txt以及csv文本数据转换为Shp格式,可以作为基础入门学习。本篇教程在之前一系列文章的基础上讲解如何使用GDAL 实现影像合并

如果你还没有看过,建议从以上内容开始。

1. 开发环境

本文使用如下开发环境,以供参考。

时间:2025年

系统:Windows 11

Python:3.11.7

GDAL:3.11.1

2. 数据准备

俗话说巧妇难为无米之炊,数据就是软件的基石,没有数据,再美好的设想都是空中楼阁。因此,第一步需要下载遥感影像数据。

但是,影像数据在哪里下载呢?别着急,本文都给你整理好了。

数据下载可参考文章:GIS 影像数据源介绍

如下,这是我在【地理空间数据云】平台下载的landsat8遥感影像。

3. 导入依赖

GeoTIFF作为一种栅格数据格式,可以使用GDAL直接进行处理,以实现影像数据的合并操作。影像合并涉及到大量数学运算,所以还需要导入numpy模块进行处理。

from osgeo import gdal
import numpy as np

如果你没有安装numpy模块的话,可执行命令pip install numpy进行安装。

4. 影像合并

定义一个方法merge_raster(input_files, output_file, target_resolution=None)用于实现栅格数据的合并。

"""
说明:GDAL 影像合并
参数:
    -input_files:输入合并的影像文件列表
    -output_file:输出影像文件
    -target_resolution:影像合并分辨率
"""
def merge_raster(input_files, output_file, target_resolution=None):

在获取到影像数据后,需要对数据范围进行合并。

"""
需要合并图像尺寸后进行输出
"""

# 获取所有影像的边界和分辨率
min_x, max_x, min_y, max_y = float('inf'), -float('inf'), float('inf'), -float('inf')
res_x, res_y = None, None

for file in input_files:
    # 打开数据集
    ds = gdal.Open(file)
    # 获取地理变换信息
    gt = ds.GetGeoTransform()

    # 获取行、列数
    x_size = ds.RasterXSize
    y_size = ds.RasterYSize

    # 计算实际边界
    x_min = gt[0]
    y_max = gt[3]
    x_max = gt[0] + x_size * gt[1]
    y_min = gt[3] + y_size * gt[5]

    min_x = min(min_x, x_min)
    max_x = max(max_x, x_max)
    min_y = min(min_y, y_min)
    max_y = max(max_y, y_max)

    if res_x is None or abs(gt[1]) < abs(res_x):
        res_x = abs(gt[1])
        res_y = abs(gt[5])

    # 关闭数据源
    ds = None

计算影像合并分辨率。

# 使用目标分辨率或计算出的分辨率
if target_resolution:
    res_x = res_y = target_resolution

获取影像大小,输出合并影像的总行数和总列数。

# 计算输出尺寸
out_cols = int((max_x - min_x) / res_x + 0.5)
out_rows = int((max_y - min_y) / res_y + 0.5)

调用gdal对象方法Warp进行影像合并,该函数第一个参数destNameOrDestDS 为输出数据集名称或者数据源,第二个参数srcDSOrSrcDSTab为源数据,第三个参数options为可选项描述,用于定义合并影像信息。

# 设置Warp参数
warp_options = gdal.WarpOptions(
    format='GTiff',
    outputBounds=[min_x, min_y, max_x, max_y],
    xRes=res_x,
    yRes=res_y,
    resampleAlg='bilinear',  # 根据需要选择:near, bilinear, cubic, lanczos等
    dstNodata=0,
    targetAlignedPixels=True,  # 确保像素对齐
    multithread=True,
    creationOptions=['COMPRESS=LZW''TILED=YES''BIGTIFF=IF_SAFER']
)

# 执行合并
gdal.Warp(output_file, input_files, options=warp_options)

main函数中调用合并方法。

if __name__ == "__main__":

    # 输入影像文件
    input_files = [ 
        "E:\ArcGIS\bandmerge\lc8_band_432.tif",
        "E:\ArcGIS\bandmerge\lc8_band_432_d.tif"
    ]
    # 输出影像文件
    output_file"E:\ArcGIS\bandmerge\merge_result.tif"

    merge_raster(input_files, output_file, target_resolution=10)

图片效果

OpenLayers示例数据下载,请在公众号后台回复:vector

全国信息化工程师-GIS 应用水平考试资料,请在公众号后台回复:GIS考试

GIS之路 公众号已经接入了智能 助手,可以在对话框进行提问,也可以直接搜索历史文章进行查看。

都看到这了,不要忘记点赞、收藏 + 关注

本号不定时更新有关 GIS开发 相关内容,欢迎关注 


    

GeoTools 开发合集(全)

OpenLayers 开发合集

小小声说一下GDAL的官方API接口

《云南省加快构建现代化产业体系推进产业强省建设行动计划》发布

ArcGIS Pro 添加底图的方式

为什么每次打开 ArcGIS Pro 页面加载都如此缓慢?

ArcGIS Pro 实现影像波段合成

自然资源部党组关于苗泽等4名同志职务任免的通知

GDAL 创建矢量图层的两种方式

GDAL 实现矢量数据转换处理(全)

GDAL 实现投影转换

国产版的Google Earth,吉林一号卫星App“共生地球”来了

2026年全国自然资源工作会议召开

日本欲打造“本土版”星链系统