©2020-2024 ioDraw All rights reserved,
Privacy Policy
刚刚开始想着手学习Python，决定从最基础的读写遥感影像开始。学习借鉴了网上很多前辈们的经验，自己出现了一些小问题写在这里，算是记录一下。
这是USGS上下载的一景Landsat8影像，地点在青海。
用ENVI截取了扎陵湖的一部分，真彩色合成图如下：
想利用Python把nir,red,green三个波段合成假彩色，代码如下：
import os import numpy as np from osgeo import gdal class IMAGE: # 读图像文件 def
read_img(self,filename): dataset = gdal.Open(filename) # 打开文件 im_width =
dataset.RasterXSize # 栅格矩阵的列数 im_height = dataset.RasterYSize # 栅格矩阵的行数 #
im_bands = dataset.RasterCount # 波段数 im_geotrans = dataset.GetGeoTransform() #
仿射矩阵，左上角像素的大地坐标和像素分辨率 im_proj = dataset.GetProjection() # 地图投影信息，字符串表示 im_data
= dataset.ReadAsArray(0, 0, im_width, im_height) del dataset #关闭对象dataset，释放内存
# return im_width, im_height, im_proj, im_geotrans, im_data,im_bands return
im_proj, im_geotrans, im_data, im_width,im_height,im_bands # 遥感影像的存储 #
写GeoTiff文件 def write_img(self,filename, im_proj, im_geotrans, im_data): #
判断栅格数据的数据类型 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 # 判读数组维数 if len(im_data.shape) == 3: # 注意数据的存储波段顺序：im_bands,
im_height, im_width im_bands, im_height, im_width = im_data.shape else:
im_bands, (im_height, im_width) = 1, im_data.shape #没看懂 # 创建文件时 driver =
gdal.GetDriverByName("GTiff")，数据类型必须要指定，因为要计算需要多大内存空间。 driver =
gdal.GetDriverByName("GTiff") dataset = driver.Create(filename, im_width,
im_height, im_bands, datatype) dataset.SetGeoTransform(im_geotrans) # 写入仿射变换参数
dataset.SetProjection(im_proj) # 写入投影 if im_bands == 1:
dataset.GetRasterBand(1).WriteArray(im_data) # 写入数组数据 else: for i in
range(im_bands): dataset.GetRasterBand(i + 1).WriteArray(im_data[i]) del
dataset if __name__ == "__main__":
os.chdir(r'F:\ProfessionalProfile\DataRelevant\testZone') #切换路径到待处理图像所在文件夹
run=IMAGE() # 第一步 # proj, geotrans, data1, row1, column1,band =
run.read_img(r'F:\ProfessionalProfile\DataRelevant\testZone\testZone1_Bnad5.tif')
# 读数据 # proj, geotrans, data2, row2, column2,band =
run.read_img(r'F:\ProfessionalProfile\DataRelevant\testZone\testZone1_Bnad4.tif')
# 读数据 # proj, geotrans, data3, row3, column3,band =
run.read_img(r'F:\ProfessionalProfile\DataRelevant\testZone\testZone1_Bnad3.tif')
# 读数据 proj, geotrans, data1, row1, column1, band =
run.read_img('testZone1_Bnad5.tif') # 读数据 proj, geotrans, data2, row2, column2,
band = run.read_img('testZone1_Bnad4.tif') # 读数据 proj, geotrans, data3, row3,
column3, band = run.read_img('testZone1_Bnad3.tif') # 读数据 #
第二步:将上述读取的3个波段放到一个数组中 # 出错（ 'str' object has no attribute 'datatype'） data =
np.array((data1,data2,data3), dtype=data1.dtype) # 按序将3个波段像元值放入 # 第三步
run.write_img('jcs543.tif', proj, geotrans, data) # 写为3波段数据,假彩色，nir,red,green
运行完以后，把结果在ENVI中打开，如下：
出现了以下两个错误：
1.ValueError: too many values to unpack (expected
5)。报错显示在72行。找了一下，发现定义的read_img函数返回6个参数，而后续调用只用到了5个参数。去掉read_img函数中返回的im_bands返回值，或者在后续调用read_img函数的时候加上一个band参数就可以了。
2.解决上述问题之后，又报错AttributeError: 'str' object has no attribute
'dtype'。发现在调用read_img函数的时候，参数顺序和返回值顺序不一致。调整顺序以后再运行就没有问题了。
这也侧面告诉了自己：光拼凑、借鉴别人写的不行，还是有很多小问题需要注意的。