【GDAL基础教程】多张二维tif数据转三维tif数据

【GDAL基础教程】多张二维tif数据转三维tif数据

今天分享一下多张二维单波段tif数据合并为一张三维多波段tif数据的脚本,话不多说,详见代码。

原数据

请添加图片描述

# -*- encoding: utf-8 -*-
'''
@File    :   tif_3d.py
@Time    :   2022/04/26 21:10:58
@Author  :   HMX 
@Version :   1.0
@Contact :   kzdhb8023@163.com
'''

# here put the import lib
from osgeo import gdal
import os
import glob
import numpy as np
import time

t1 = time.time()
tiflist = []
filepath = r'E:\Project\XINAN\WWLLN_NEW\CG\YEAR'
for file in glob.glob(os.path.join(filepath,'*.tif')):
    print(file)
    tiflist.append(file)
inds = gdal.Open(tiflist[0])
inband = inds.GetRasterBand(1)
indriver = gdal.GetDriverByName('GTiff')
outds = indriver.Create(r'D:\公众号\NO.10\out.tif',inband.XSize,inband.YSize,len(tiflist),inband.DataType)
outds.SetProjection(inds.GetProjection())
outds.SetGeoTransform(inds.GetGeoTransform())
ds = inds.ReadAsArray()
m,n = ds.shape[0],ds.shape[1]
ds_new = np.zeros((len(tiflist),m,n))
for r in range(len(tiflist)):
    rds = gdal.Open(tiflist[r]).ReadAsArray()
    nodatavalue = rds[0,0]
    rds[rds==nodatavalue] = np.nan
    ds_new[r,:,:] = rds[:,:]
    outband = outds.GetRasterBand(r+1)
    outband.WriteArray(ds_new[r,:,:])
    outband.SetDescription(os.path.basename(tiflist[r]))
inds = None
rds = None
t2 = time.time()
print('共计用时:{:.2f}s'.format(t2-t1))

结果

利用gis打开可以发现已经成功合并。
请添加图片描述

如果对你有帮助的话,请‘点赞’、‘收藏’,‘关注’,你们的支持是我更新的动力。欢迎关注【森气笔记】。

  • 2
    点赞
  • 6
    收藏
    觉得还不错? 一键收藏
  • 2
    评论
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值