python 拼接 遥感影像_如何用Python| 制作遥感影像拼接
本文的文字及图片来源于网络,仅供学习、交流使用,不具有任何商业用途,版权归原作者所有,如有问题请及时联系我们以作处理
以下文章来源于腾讯云,作者:bugsuse。
0.前言
因为没有喝上“秋天的第一份奶茶”,准备来更新一篇推送。
在上一篇推文中,我展示了如何使用Python结合Landsat制作遥感影像图(Python干货 | 制作遥感影像图)。
对于Landsat数据来说,对某个区域的重访周期为16天,每个位置使用全球参考系(WRS)进行索引,即每一个位置都会对应一个Path和Row,相邻的影像之间会有部分区域是重叠的。
Fig.1 World Reference System
在某些遥感影像的应用场景中,如果我们关注的区域正好处于两景影像的交界处,如下图中的象山港,那我们就需要将影像拼接起来才可以使用。
单张影像是这样。
本文合并后是这样。
1.准备工作
相较于上一篇推送,我们这次为了实现遥感影像的镶嵌拼接,我们使用到了两个库, rasterio和gdal。
import rasterio as rio
import gdal
先介绍一下我们实现两组遥感影像拼接的思路,首先选取两景相邻的影像,分别得到他们的空间范围,再得到两景组合到一起之后的空间范围,使用gdal新建一个tif文件(数据中转用),分别得到原来两景影像在新建的tif文件中的起始位置,将对应的数据写入新的tif文件中,即实现镶嵌拼接。
上面说的是两景影像的拼接,如果是更多影像拼接同样适用,但是现阶段的方法如果拼接多的影像的话,需要的内存空间很大,容易导致内存溢出,感兴趣的朋友可以思考一下如何高效实现多景影像的拼接。
其中还有两处关键处理,一是如何去除重叠区域的无效信息,二是重叠区域的数据如何选择。希望各位看官能从代码里面找到答案。
2.动起手来
得到输入影像的四个角点。
def tiffileList2filename(tiffileList):
filename = []
prefix = []
for ifile in tiffileList:
file0 = ifile.split("\\")[-1]
prefix.append(os.path.join(ifile, file0))
filename.append(os.path.join(ifile, file0) + "_B1.TIF")
return filename, prefix
def get_extent(tiffileList):
filename, prefix = tiffileList2filename(tiffileList)
rioData = rio.open(filename[0])
left = rioData.bounds[0]
bottom = rioData.bounds[1]
right = rioData.bounds[2]
top = rioData.bounds[3]
for ifile in filename[1:]:
rioData = rio.open(ifile)
left = min(left, rioData.bounds[0])
bottom = min(bottom, rioData.bounds[1])
right = max(right, rioData.bounds[2])
top = max(top, rioData.bounds[3])
return left, bottom, right, top, filename, prefix
得到新建tif文件的size,这里已知Landsat空间分辨率为30m,如果是其他遥感数据,需对应进行修改。
if __name__ == '__main__':
tiffileList = [r'PathofLandsat8\LC08_L1TP_118039_20160126_20170330_01_T1',
r'PathofLandsat8\LC08_L1TP_118040_20160126_20170330_01_T1']
left, bottom, right, top, filename, prefix = get_extent(tiffileList)
cols, rows= getRowCol(left, bottom, right, top)
bands = ['B7', 'B5', 'B3']
n_bands = len(bands)
arr = np.zeros((n_bands, rows, cols), dtype=np.float)
# 打开一个tif文件
in_ds = gdal.Open(filename[0])
for i in range(len(bands)):
ibands = bands[i]
# 新建一个tif文件
driver = gdal.GetDriverByName('gtiff')
out_ds = driver.Create(ibands + 'mosaic.tif', cols, rows)
# 设置tif文件的投影
out_ds.SetProjection(in_ds.GetProjection())
out_band = out_ds.GetRasterBand(1)
# 设置新tif文件的地理变换
gt = list(in_ds.GetGeoTransform())
gt[0], gt[3] = left, top
out_ds.SetGeoTransform(gt)
# 对要拼接的影像进行循环读取
for ifile in prefix:
in_ds = gdal.Open(ifile + '_' + ibands + '.tif')
# 计算新建的tif文件及本次打开的tif文件之间的坐标漂移
trans = gdal.Transformer(in_ds, out_ds, [])
# 得到偏移起始点
success, xyz = trans.TransformPoint(False, 0, 0)
x, y, z = map(int, xyz)
# 读取波段信息
fnBand = in_ds.GetRasterBand(1)
data = fnBand.ReadAsArray()
# 写入tif文件之前,最大值设置为255,这一步很关键
data = data / 65535 * 255
data[np.where(data == 255)] = 0
# 影像重合部分处理,重合部分取最大值
xSize = fnBand.XSize
ySize = fnBand.YSize
outData = out_band.ReadAsArray(x, y, xSize, ySize)
data = np.maximum(data, outData)
out_band.WriteArray(data, x, y)
del out_band, out_ds
file2read = ibands + 'mosaic.tif'
arr[i, :, :] = tiff.imread(file2read)
os.remove(file2read)
plot_rgb(arr, rgb=(0, 1, 2))
3.小结 大功告成!
python 拼接 遥感影像_如何用Python| 制作遥感影像拼接相关推荐
- python 矩阵运算 for循环_如何用 Python 科学计算中的矩阵替代循环
展开全部 因为在Mathematica中使用循环确实是低效的.32313133353236313431303231363533e78988e69d8331333361313961..... 深层次的原 ...
- python搭建自动化测试平台_如何用python语言搭建自动化测试环境
原标题:如何用python语言搭建自动化测试环境 技术分享:基于Python语言的Web自动化测试环境搭建 近期发现很多初学者在学习自动化的过程当中,在环境安装环节总是出现问题,所以详细的出一篇环境搭 ...
- 如何制作python检查小软件_如何用Python制作整蛊小程序
原标题:如何用Python制作整蛊小程序 下面的整蛊程序,千万不要发代码,否则就实现不了你整蛊的目的了.完成后一定要打包成一个exe程序,再发给朋友使用 . 1. 使用 pip install pyi ...
- python rest api 测试_如何用Python编写REST API的单元测试
在过去的几个月中,正在从事一个名为B的项目.它是带有简单Web UI的徽章生成器,用于添加数据并生成PDF可打印徽章.B后端现在已转移到REST-API并测试REST-API中使用的功能,我们需要一些 ...
- python的out模式_如何用python中的DataFrame列的模式替换NA值?
我对Python(和本网站)完全陌生,目前正试图用它们的模式替换特定数据帧列中的NA值.我试过了各种不起作用的方法.请帮我看看我做错了什么:如何用python中的DataFrame列的模式替换NA值? ...
- 用python做一张图片_如何用python下载一张图片
如何用python下载一张图片 这里要用到的主要工具是requests这个工具,需要先安装这个库才能使用,该库衍生自urllib这个库,但是要比它更好用.多数人在做爬虫的时候选择它,是个不错的选择. ...
- 用python处理excel表格_如何用python处理excel数据 | 用python处理excel表格数据类型
python 读取EXCEL文件中的数据格式 扩展库 xlrd 读excle xlwt 写excle 直上搜就能下载 下载后使用 import xlrd 就可以读excle了 打开文件: xls = ...
- python turtle画动物_如何用python画简单的动物
首先来看一下实现效果,如下图:程序猿的生活:Python入门到精通资料大汇总,不啰嗦,全是珍藏资料!zhuanlan.zhihu.com 具体实现代码请看: # -*- coding:utf-8 - ...
- 用python画机器猫代码_如何用Python画一只机器猫?| 原力计划
原标题:如何用Python画一只机器猫?| 原力计划 作者 | 人邮异步社区 责编 | 胡巍巍 出品 | CSDN博客 自信心是成功的源泉,对刚入门编程行业的初级程序员来说,多敲代码多做项目就是构建自 ...
- python 读取excel图片_如何用Python读取Excel中图片?
公众号: 早起Python 作者:刘早起 大家好,在使用Python进行办公自动化操作时,一定少不了与Excel表格的交互,我们通常是用pandas处理表格数据,但大多数情况下,都是读取表格中的数值进 ...
最新文章
- ggplot01:R语言坐标轴离散、连续与图例离散连续的区分
- KeilKill.bat删除keil编译生成的过程文件
- 计算机视觉与深度学习 | TensorMask: A Foundation for Dense Object Segmentation(何凯明团队新作)近5年目标检测综述
- Springboot整合Hikari数据库连接池,密码加密
- Jerry带您了解Restful ABAP Programming模型系列之二:Action和Validation的实现
- 利用bootstrap框架做了一个采摘节节日活动网页
- 网站安全之XSS漏洞攻击以及防范措施
- idea 使用 Gradle 构建过程中控制台中文显示乱码解决
- Python菜鸟入门:day05列表
- 蓝桥杯2020年第十一届Python省赛第一题-门牌制作
- c#textBox控件限制只允许输入数字及小数点
- CRMEB 知识付费模版消息修改教程
- 【图像配准】基于sift算法实现图像配准matlab源码
- 阿里云企业邮箱标准版/集团版/尊享版区别对比
- win10系统中如何查看wifi密码
- linux 即时通讯,Linux即时通讯Pidgin简洁漂亮的插件Screenlets[图文]
- 横河变送器EJA430E-JCS4G-917DB
- 这是一个关于五台山关于点化顿悟的真实事件,源起偶然,其为必然,一趟旅程获得一份机缘。
- vue实现考勤排班日历(备忘)
- 台式计算机多少g的显卡怎么看,怎样看电脑配置|怎样看电脑显卡配置?
热门文章
- 重邮python实验课之华氏温度转摄氏温度速查表
- 吴恩达深度学习笔记(40)-指数加权平均数优化算法
- File /py-faster-rcnn/tools/../lib/datasets/imdb.py, line 108, in append_flipped_images assert
- 杰理之低延时无线麦功能支持以下三种组合配置【篇】
- Eighth Week's ARST
- 按键精灵的5级开发认证,笔试题参考
- 分享一款超棒的jQuery Google地图插件:Gmaps
- c#关于GMap离线地图加载的问题
- 全国计算机四级薪资,全国计算机四级通过率有多少
- 笔记:全网最详细jQuery教程