网站哪家做得好,网页游戏平台51,wordpress顶部加广告,国内做视频课程的网站有哪些遥感工作者离不开栅格数据#xff0c;有时候我们可能需要修改栅格数据的值#xff0c;但ENVI和ArcGIS中并没有直接修改DN值的工具#xff0c;只有栅格计算器、Band math这些工具去计算整个波段的值#xff0c;或者Edit Classification Image工具可以修改ENVI分类后的像元值… 遥感工作者离不开栅格数据有时候我们可能需要修改栅格数据的值但ENVI和ArcGIS中并没有直接修改DN值的工具只有栅格计算器、Band math这些工具去计算整个波段的值或者Edit Classification Image工具可以修改ENVI分类后的像元值但这个工具只对分类格式有效博主整不明白怎么把单波段数据转为分类格式所以这些工具都无法满足我们的需求。 今天跟大家分享一下如何使用Python的GDAL库将栅格数据中特定的DN值修改成我们想要的。 一、安装库
import os
import numpy as np
from osgeo import gdal
二、读取栅格数据信息 这一步可有可无只是让你了解一下栅格数据的基本信息如投影信息、长度、宽度等。
def Get_data(filepath):ds gdal.Open(filepath) # 打开数据集datasetds_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) # 以数组的形式读取整个数据集
三、修改DN值 我这里写的函数只适用于修改单波段的栅格影像如果需要修改多波段就自己加一个循环。其实原理都一样就是将波段读成数组然后通过索引读取第几行第几列像元的值然后判断这个值是否为你想要修改的值如果是就将其赋予新的值。
def Modify_value(filepath, out_path, original, target):print(-----正在进行DN值的修改-----)ds gdal.Open(filepath) # 打开数据集datasetds_width ds.RasterXSize # 获取数据宽度ds_height ds.RasterYSize # 获取数据高度ds_geo ds.GetGeoTransform() # 获取仿射地理变换参数ds_prj ds.GetProjection() # 获取投影信息array_band ds.GetRasterBand(1).ReadAsArray(0, 0, ds_width, ds_height).astype(np.float64)# 读取第一个波段全部for row in range(0, ds_height):# 循环当前波段的行for col in range(0, ds_width):# 循环当前波段的列if array_band[row][col] original:# 判断第row行第col列的DN值是否为需要修改的值array_band[row][col] target# 修改该值else:continuedriver gdal.GetDriverByName(GTiff) # 载入数据驱动用于存储内存中的数组ds_result driver.Create(out_path, ds_width, ds_height, bands1, eTypegdal.GDT_Float64)# 创建一个数组宽高为原始尺寸ds_result.SetGeoTransform(ds_geo) # 导入仿射地理变换参数ds_result.SetProjection(ds_prj) # 导入投影信息ds_result.GetRasterBand(1).SetNoDataValue(0) # 将无效值设为0ds_result.GetRasterBand(1).WriteArray(array_band) # 将结果写入数组del ds_result# 删除内存中的结果否则结果不会写入图像中print(计算完成)
四、完整代码
# -*- coding: utf-8 -*-Time 2023/9/8 8:51
Auth RS迷途小书童
File Modifying Raster Data DN Values.py
IDE PyCharm
Purpose修改栅格数据DN值import os
import numpy as np
from osgeo import gdaldef Get_data(filepath):ds gdal.Open(filepath) # 打开数据集datasetds_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 Modify_value(filepath, out_path, original, target):print(-----正在进行DN值的修改-----)ds gdal.Open(filepath) # 打开数据集datasetds_width ds.RasterXSize # 获取数据宽度ds_height ds.RasterYSize # 获取数据高度ds_geo ds.GetGeoTransform() # 获取仿射地理变换参数ds_prj ds.GetProjection() # 获取投影信息array_band ds.GetRasterBand(1).ReadAsArray(0, 0, ds_width, ds_height).astype(np.float64)# 读取第一个波段全部for row in range(0, ds_height):# 循环当前波段的行for col in range(0, ds_width):# 循环当前波段的列if array_band[row][col] original:# 判断第row行第col列的DN值是否为需要修改的值array_band[row][col] target# 修改该值else:continuedriver gdal.GetDriverByName(GTiff) # 载入数据驱动用于存储内存中的数组ds_result driver.Create(out_path, ds_width, ds_height, bands1, eTypegdal.GDT_Float64)# 创建一个数组宽高为原始尺寸ds_result.SetGeoTransform(ds_geo) # 导入仿射地理变换参数ds_result.SetProjection(ds_prj) # 导入投影信息ds_result.GetRasterBand(1).SetNoDataValue(0) # 将无效值设为0ds_result.GetRasterBand(1).WriteArray(array_band) # 将结果写入数组del ds_result# 删除内存中的结果否则结果不会写入图像中print(计算完成)if __name__ __main__:os.environ[PROJ_LIB] G:/Anaconda/envs/pyDL/Lib/site-packages/osgeo/data/projos.environ[GDAL_DATA] G:/Anaconda/envs/pyDL/Lib/site-packages/osgeo/data# 添加PROJ至环境变量消除警告file_path rB:\1m_xiugai.tif # 输入的栅格数据路径out_path1 rB:\1m_xiugai1.tif # 导出的文件路径data1 int(input(请输入需要修改的DN值:))data2 int(input(请输入目标DN值 :))Get_data(file_path) # 执行函数获取影像基本信息print(\n)print(--------------DN值修改--------------)Modify_value(file_path, out_path1, data1, data2) # 执行函数修改DN值今天的分享就到这里了大家需要注意的是我这段代码只适用于单波段数据且想要修改的值只有一种时如你想要将所有DN值等于1的像元全部改成0就可以直接使用我的点吗改数据路径然后再输入1和0就可以了因为我的任务就是将分类数据DN值为01234中分错的部分改成正确的。不要问为什么不使用工具因为我手里的分类数据不是ENVI支持的分类格式泪目。 如果大家在学习Python或者RS时有什么问题可以随时留言交流如果大家对批量处理有兴趣同样可以留言给博主博主会分享相关代码以供学习
