Python GDAL缩放栅格文件各波段数值

简介: 本文介绍基于Python中的gdal模块,批量读取大量多波段遥感影像文件,分别对各波段数据加以数值处理,并将所得处理后数据保存为新的遥感影像文件的方法。首先,看一下本文的具体需求。我们现有一个文件夹,其中含有大量.tif格式的遥感影像文件;其中,这些遥感影像文件均含有4个波段,每1个波段都表示其各自的反射率数值。而对于这些遥感影像文件,有的文件其各波段数值已经处于0至1的区间内(也就是反射率数据的正常数值区间),而有的文件其各波段数值则是还没有乘上缩放系数的(在本文中,缩放系数是0.0001)。

本文介绍基于Python中的gdal模块,批量读取大量多波段遥感影像文件,分别对各波段数据加以数值处理,并将所得处理后数据保存为新的遥感影像文件的方法。


首先,看一下本文的具体需求。我们现有一个文件夹,其中含有大量.tif格式的遥感影像文件;其中,这些遥感影像文件均含有4个波段,每1个波段都表示其各自的反射率数值。而对于这些遥感影像文件,有的文件其各波段数值已经处于0至1的区间内(也就是反射率数据的正常数值区间),而有的文件其各波段数值则是还没有乘上缩放系数的(在本文中,缩放系数是0.0001)。


例如,如下图所示,即为文件夹中某一景遥感影像。可以看到其各波段数值都是大于1的,这是因为其数值都是还没有乘上缩放系数的,即是真实的反射率数值的10000倍。

1718710166145.jpg

我们希望实现的是,对于这些遥感影像中,还没有乘上缩放系数0.0001的遥感影像,将其像元值乘上这个缩放系数;而对于已经缩放过(也就是像元数值已经落在0至1区间内)的遥感影像,则不加以任何处理。最后,将经过上述操作后的所有图像(无论是否执行缩放)均保存至指定的输出结果文件夹中。


本文所需代码如下:

# -*- coding: utf-8 -*-
"""
Created on Thu Apr 18 12:37:22 2024
@author: fkxxgis
"""
import os
from osgeo import gdal
original_folder = r"E:\04_Reconstruction\99_MODIS\new_data\GF_Original"
output_folder = r"E:\04_Reconstruction\99_MODIS\new_data\GF_Rec"
for filename in os.listdir(original_folder):
    if filename.endswith('.tif'):
        dataset = gdal.Open(os.path.join(original_folder, filename), gdal.GA_ReadOnly)
        width = dataset.RasterXSize
        height = dataset.RasterYSize
        
        band_count = dataset.RasterCount
        driver = gdal.GetDriverByName('GTiff')
        output_dataset = driver.Create(os.path.join(output_folder, "New_" + filename), width, height, band_count, gdal.GDT_Float32)
        
        for band_index in range(1, band_count + 1):
            band = dataset.GetRasterBand(band_index)
            data = band.ReadAsArray()
            
            if band_index == 1:
                data = data.astype(float)
                data[data > 1] /= 10000
            elif band_index == 2:
                data = data.astype(float)
                data[data > 1] /= 10000
            elif band_index == 3:
                data = data.astype(float)
                data[data > 1] /= 10000
            elif band_index == 4:
                data = data.astype(float)
                data[data > 1] /= 10000
            output_band = output_dataset.GetRasterBand(band_index)
            output_band.WriteArray(data)
            output_band.FlushCache()
        output_dataset.SetGeoTransform(dataset.GetGeoTransform())
        output_dataset.SetProjection(dataset.GetProjection())
        dataset = None
        output_dataset = None

首先,我们使用os.listdir()函数遍历原始数据文件夹中的所有文件,并使用if语句筛选出以.tif结尾的文件;随后,使用gdal.Open()函数打开原始影像数据集,并指定只读模式;接下来,使用dataset.RasterXSize和dataset.RasterYSize获取影像数据集的宽度和高度。


随后,使用dataset.RasterCount获取波段数量,并使用gdal.GetDriverByName()创建输出数据集的驱动程序对象;紧接着,通过Create()方法创建输出数据集,并指定输出文件的路径、宽度、高度、波段数量和数据类型(gdal.GDT_Float32表示浮点型)。


接下来,就可以开始使用循环,对每个文件每个波段进行处理。首先,使用dataset.GetRasterBand()方法获取当前波段对象,然后使用band.ReadAsArray()将波段数据读取为数组;根据波段索引的不同,对波段数据进行处理。在本文中,对4个波段进行的其实是相同的处理,即将大于1的像素值除以10000。


其次,使用output_dataset.GetRasterBand()方法获取输出数据集中的当前波段对象,并使用output_band.WriteArray()方法将处理后的数据写入输出数据集。


再次,使用dataset.GetGeoTransform()和dataset.GetProjection()分别获取原始数据集的地理转换和投影信息,并使用output_dataset.SetGeoTransform()和output_dataset.SetProjection()设置输出数据集的地理转换和投影信息。


最后一步,关闭数据集对象。至此,代码就完成了对每个.tif文件的处理,并将处理后的数据保存到输出文件夹中。


此时,打开本文开头展示的那1景遥感影像,可以看到其像素数值已经是乘上缩放系数之后的了,也就是落在了0至1的区间内;如下图所示。

1718710264362.jpg

至此,大功告成。

作者:疯狂学习GIS

链接:https://juejin.cn/post/7381423581018751016

相关文章
|
6天前
|
数据格式 Python
【Python】已解决:Excel无法打开文件test.xIsx“,因为文件格式或文件扩展名无效。请确定文件未损坏,并且文件扩展名与文件的格式匹配。
【Python】已解决:Excel无法打开文件test.xIsx“,因为文件格式或文件扩展名无效。请确定文件未损坏,并且文件扩展名与文件的格式匹配。
34 0
|
2天前
|
安全 网络安全 数据安全/隐私保护
`certifi`是一个Python包,它提供了一个包含Mozilla证书颁发机构(CA)Bundle的PEM文件。
`certifi`是一个Python包,它提供了一个包含Mozilla证书颁发机构(CA)Bundle的PEM文件。
19 7
|
2天前
|
存储 对象存储 Python
`openpyxl`是一个用于读写Excel 2010 xlsx/xlsm/xltx/xltm文件的Python库。它不需要Microsoft Excel,也不需要.NET或COM组件。
`openpyxl`是一个用于读写Excel 2010 xlsx/xlsm/xltx/xltm文件的Python库。它不需要Microsoft Excel,也不需要.NET或COM组件。
6 1
|
5天前
|
Python
python存取yaml文件
python存取yaml文件
|
2天前
|
存储 Python
`tempfile`模块在Python中用于创建临时文件和目录。
`tempfile`模块在Python中用于创建临时文件和目录。
4 0
|
5天前
|
Unix Linux C++
相对路径vs绝对路径 python文件的添加与删除
相对路径vs绝对路径 python文件的添加与删除
|
6天前
|
开发者 iOS开发 MacOS
【Python】已解决:Pycharm设置Python每个文件开头自定义模板(带上作者名、时间等)
【Python】已解决:Pycharm设置Python每个文件开头自定义模板(带上作者名、时间等)
9 0
|
1月前
|
存储 Python 内存技术
python WAV音频文件处理—— (1)读写WAV文件
python WAV音频文件处理—— (1)读写WAV文件
55 14
|
2月前
|
数据挖掘 索引 Python
Python 读写 Excel 文件
Python 读写 Excel 文件
124 0
|
2月前
|
存储 JSON JavaScript
Python中读写(解析)JSON文件的深入探究
Python中读写(解析)JSON文件的深入探究
92 0