Python基于GDAL镶嵌拼接遥感影像
没啥好说的,处理高分辨率影像时,数据高达几十G。用ENVI或者ArcGIS进行影像的拼接时,往往会出现未响应的情况。出现未响应的话,运气好等个一晚上可能会动一动,运气不好就等着强制关闭重做吧。
所以搞了一个Python进行拼接操作的代码,虽然速度不算快,但至少不会未响应。同时如果对代码进行一些改进,还可以进行批量拼接的操作,百利而无一害。
一、导入GDAL库
from osgeo import gdal
二、查看影像信息
为了凑字数的,可以查看影像的投影、宽度、高度、波段数等信息。不过需要注意的是在ENVI中没有投影坐标系,只有地理坐标系是做不了镶嵌拼接的。
这个代码我还不太清楚能不能不要投影坐标系进行拼接,你们可以自己试试。但最好还是用包含投影坐标系的影像进行拼接。所以在拼接之前就可以用这段代码先看一看。
def Get_data(filepath): ds = gdal.Open(filepath) # 打开数据集dataset ds_width = ds.RasterXSize # 获取数据宽度 ds_height = ds.RasterYSize # 获取数据高度 ds_bands = ds.RasterCount # 获取波段数 ds_geo = ds.GetGeoTransform() # 获取仿射地理变换参数 ds_prj = ds.GetProjection() # 获取投影信息 print("影像的宽度为:" + str(ds_width)) print("影像的高度为:" + str(ds_height)) print("仿射地理变换参数为:" + str(ds_geo)) print("投影坐标系为:" + str(ds_prj)) # data = ds.ReadAsArray(0, 0, ds_width, ds_height) # 以数组的形式读取整个数据集
三、镶嵌模块
这里用到了GDAL中的Warp函数,是不是有点熟悉。之前裁剪也是用的这个函数,不得不说这个函数是真的强大,之前有文章介绍了其中的函数,大家有兴趣可以去看下,同时给个赞吧!【Python&RS】GDAL批量裁剪遥感影像/栅格数据
代码中的srcSRS,dstSRS分别是输入投影和输出投影,这里用一样的就行了。因为我们做的是镶嵌操作,肯定是不用动原始坐标系的。其他的参数都在代码中表明了,这里就不介绍了,如果大家有什么问题,可以留言交流。
def Mosaic_GDAL(path_image1, path_image2, path_out): """ :param path_image1: 需要镶嵌的影像 :param path_image2: 需要镶嵌的影像 :param path_out: 镶嵌后输出的影像路径 :return: None """ image1 = gdal.Open(path_image1, gdal.GA_ReadOnly) # 第一幅影像 input_proj = image1.GetProjection() image2 = gdal.Open(path_image2, gdal.GA_ReadOnly) # 第二幅影像 options = gdal.WarpOptions(srcSRS=input_proj, dstSRS=input_proj, format='GTiff', resampleAlg=gdal.GRA_NearestNeighbour,callback=Show_Progress) # 输入投影,输出投影,输出格式,重采样方法 gdal.Warp(path_out, [image1, image2], options=options) # 输出路径,需要镶嵌的数据,参数配置 ds = gdal.Open(path_out, gdal.GA_ReadOnly) ds.BuildOverviews(overviewlist=[2, 4, 8, 16]) # 创建金字塔 del image1, image2, ds
四、回调函数(没啥用)
上一篇博文已经介绍过了,就两点:1.水字数,2.记录进度。省的看着代码发呆。
def Show_Progress(percent, msg, tag): """ :param percent: 进度,0~1 :param msg: :param tag: :return: """ if 25 <= percent*100 <= 26: print("进度:" + "%.2f" % (percent*100) + "%") if 50 <= percent*100 <= 51: print("进度:" + "%.2f" % (percent*100) + "%") if 75 <= percent*100 <= 76: print("进度:" + "%.2f" % (percent*100) + "%")
五、完整代码
# -*- coding: utf-8 -*- """ @Time : 2023/6/25 16:28 @Auth : RS迷途小书童 @File :Raster Data Mosaic.py @IDE :PyCharm @Purpose:遥感影像、栅格数据的镶嵌拼接(有明显拼接线) """ from osgeo import gdal def Get_data(filepath): ds = gdal.Open(filepath) # 打开数据集dataset ds_width = ds.RasterXSize # 获取数据宽度 ds_height = ds.RasterYSize # 获取数据高度 ds_bands = ds.RasterCount # 获取波段数 ds_geo = ds.GetGeoTransform() # 获取仿射地理变换参数 ds_prj = ds.GetProjection() # 获取投影信息 print("影像的宽度为:" + str(ds_width)) print("影像的高度为:" + str(ds_height)) print("仿射地理变换参数为:" + str(ds_geo)) print("投影坐标系为:" + str(ds_prj)) # data = ds.ReadAsArray(0, 0, ds_width, ds_height) # 以数组的形式读取整个数据集 def Show_Progress(percent, msg, tag): """ :param percent: 进度,0~1 :param msg: :param tag: :return: """ if 25 <= percent*100 <= 26: print("进度:" + "%.2f" % (percent*100) + "%") if 50 <= percent*100 <= 51: print("进度:" + "%.2f" % (percent*100) + "%") if 75 <= percent*100 <= 76: print("进度:" + "%.2f" % (percent*100) + "%") def Mosaic_GDAL(path_image1, path_image2, path_out): """ :param path_image1: 需要镶嵌的影像 :param path_image2: 需要镶嵌的影像 :param path_out: 镶嵌后输出的影像路径 :return: None """ image1 = gdal.Open(path_image1, gdal.GA_ReadOnly) # 第一幅影像 input_proj = image1.GetProjection() image2 = gdal.Open(path_image2, gdal.GA_ReadOnly) # 第二幅影像 options = gdal.WarpOptions(srcSRS=input_proj, dstSRS=input_proj, format='GTiff', resampleAlg=gdal.GRA_NearestNeighbour, callback=Show_Progress) # 输入投影,输出投影,输出格式,重采样方法 gdal.Warp(path_out, [image1, image2], options=options) # 输出路径,需要镶嵌的数据,参数配置 ds = gdal.Open(path_out, gdal.GA_ReadOnly) ds.BuildOverviews(overviewlist=[2, 4, 8, 16]) # 创建金字塔 del image1, image2, ds if __name__ == "__main__": path1 = "B:/proj_clip" path2 = "B:/2_proj_clip" path_output = "warp1.tif" print("开始拼接......") Mosaic_GDAL(path1, path2, path_output) print("拼接完成......")
这里有一点需要注意的就是,用这个方法进行镶嵌拼接操作时,影像有一条明显的拼接线,不知道是不是我数据的问题。你们可以自己尝试一下。只要修改主函数中的路径即可。
如果大家在学习Python或者RS时有什么问题,可以随时留言交流!
到此这篇关于Python基于GDAL镶嵌拼接遥感影像的文章就介绍到这了,更多相关Python遥感影像镶嵌内容请搜索脚本之家以前的文章或继续浏览下面的相关文章希望大家以后多多支持脚本之家!
相关文章
Python Pexpect库自动化交互式进程控制的expect_list方法解析
Pexpect是一个Python库,为自动化和交互式进程控制提供了丰富的功能,而expect_list方法是其功能强大且灵活的一部分,将详细探讨如何使用这一方法,并提供多个示例来说明其应用场景和功能2024-01-01
最新评论