本節主要利用python 和 GDAL來計算一些常用的遙感指數，如，NDVI（歸一化植被指數）、RVI（比值植被指數）、NDWI （歸一化水體指數）、SAVI（土壤調節植被指數）。
1.讀取波段數據
import numpy as np from osgeo import gdal in_ds = gdal.Open(r"D:\data\test.tif") # 打開樣本文件 xsize = in_ds.RasterXSize # 獲取行列數 ysize = in_ds.RasterYSize bands = in_ds.RasterCount # 獲取波段數 geotransform = in_ds.GetGeoTransform() # 獲取投影信息 projection = in_ds.GetProjectionRef() block_data = in_ds.ReadAsArray(0,0,xsize,ysize).astype(np.float32)# 獲取影像信息 B = block_data[0, :, :] G = block_data[1,:, :] R = block_data[2,:, :] NIR = block_data[3,:, :]
2.表達式計算
#歸一化植被指數 NDVI=(NIR-R)/(NIR+R) #歸一化水體指數 NDWI=(G-NIR)/(G+NIR) #比值植被指數 RVI=NIR/R #土壤調節植被指數 L=0.5 #L是隨着植被密度變化的參數，取值范圍從0-1，當植被覆蓋度很高時為0，很低時為1。一般取0.5 SAVI=(NIR-R)*(1+L)/(NIR+R+L)
3.保存成果
#以NDVI為例 driver = gdal.GetDriverByName('GTiff') out_dataset=driver.Create("NDVI.tif",xsize,ysize,1,gdal.GDT_Float32) out_band1 = out_dataset.GetRasterBand(1) out_band1.WriteArray(NDVI) out_dataset.SetGeoTransform(geotransform) # 寫入仿射變換 out_dataset.SetProjection(projection)