python如何使用gdal扩充影像

最新推荐文章于 2024-06-14 17:09:53 发布

rs_gis

最新推荐文章于 2024-06-14 17:09:53 发布

阅读量363

点赞数

分类专栏： python GIS 文章标签： python

本文链接：https://blog.csdn.net/rs_gis/article/details/120591196

版权

python 同时被 2 个专栏收录

26 篇文章

订阅专栏

GIS

8 篇文章

订阅专栏

本文介绍了使用Python的GDAL库进行遥感影像读取、处理和扩充的方法，包括打开tif文件、获取元数据、创建和写入新的tif文件，以及计算并应用偏移量。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

import os
import gdal
import numpy as np
import pandas as pd
import datetime
import glob

def readTif(fileName):
    dataset = gdal.Open(fileName)
    if dataset == None:
        print(fileName + "文件无法打开")
        
    width = dataset.RasterXSize 
    height = dataset.RasterYSize 
    bands = dataset.RasterCount 

    data = dataset.ReadAsArray()
    geotrans = dataset.GetGeoTransform()
    proj = dataset.GetProjection()
    del dataset
    return width, height, bands, data, geotrans, proj

def writeTiff(im_data,input_data,im_geotrans,xOffset,yOffset,im_proj, path):
    if 'int8' in im_data.dtype.name:
        datatype = gdal.GDT_Byte
    elif 'int16' in im_data.dtype.name:
        datatype = gdal.GDT_UInt16
    else:
        datatype = gdal.GDT_Float32
        
    im_bands, im_height, im_width = im_data.shape
        
    if len(input_data.shape) == 3:
        input_bands, input_height, input_width = input_data.shape
    elif len(input_data.shape) == 2:
        input_data = np.array([input_data])
        input_bands, input_height, input_width = input_data.shape
    else:
        print("image is error")
            
    #创建文件
    driver = gdal.GetDriverByName("GTiff")
    dataset = driver.Create(path, int(im_width), int(im_height), int(input_bands), datatype)
    if(dataset!= None):
        dataset.SetGeoTransform(im_geotrans) #写入仿射变换参数
        dataset.SetProjection(im_proj) #写入投影
        
    for i in range(input_bands):   
        
        dataset.GetRasterBand(i + 1).WriteArray(input_data[i],xOffset,yOffset)

    del dataset

if __name__ == '__main__':
    start = datetime.datetime.now()

    ref_tif = r"G:\xdn\class_result\sy_index_result_pro.tif"
    samplefile = r'G:\xdn\sample\sy_sample_pro211.tif'
    path = r'G:\xdn\sample\extecd_sample5.tif'

    width_ref, height_ref, bands_ref, data_ref, geotrans_ref, proj_ref = readTif(ref_tif)
    width_s, height_s, bands_s, data_s, geotrans_s, proj_s = readTif(samplefile)

    print("开始扩充影像。。。。。。。。")
    #计算待扩充影像在扩充影像左上角的偏移量
    xOffset=int((geotrans_s[0]-geotrans_ref[0])/geotrans_ref[1])
    yOffset=int((geotrans_s[3]-geotrans_ref[3])/geotrans_ref[5])

    im_data = data_ref
    input_data = data_s
    im_geotrans = geotrans_ref
    im_proj = proj_ref
    print(xOffset,yOffset)
    print(width_ref,height_ref)
    print(width_s,height_s)

    writeTiff(im_data,input_data,im_geotrans,xOffset,yOffset,im_proj, path)
    end = datetime.datetime.now()
    print(end-start)