Geotiffから非ピクセル情報を抽出するには？

Question

これはGeotiff画像のNDVI計算スクリプトの一部です。各ピクセルのNDVI値は、ラスターバンドに対するwirttenとして計算されます。 gdalやその他のモジュールでgisの処理に関する知識を持っている人は、出力されたイメージからNDVIを検索する方法を手伝ってください。

for i in range(len(red_tuple)): 
       print nir_tuple[i] 
       # Calculate the NDVI for the current pixel. 
       ndvi_lower = (nir_tuple[i] + red_tuple[i]) 
       ndvi_upper = (nir_tuple[i] - red_tuple[i]) 
       ndvi = 0 
       # Becareful of zero divide 
       if ndvi_lower == 0: 
        ndvi = 0 
       else: 
        ndvi = ndvi_upper/ndvi_lower 

       # Add the current pixel to the output line 
       outputLine = outputLine + struct.pack('f', ndvi) 
       #print(ndvi) 
      # Write the completed line to the output image    
outDataset.GetRasterBand(1).WriteRaster(0, line,red_band.XSize, 1,outputLine, buf_xsize=red_band.XSize, 
              buf_ysize=1, buf_type=gdal.GDT_Float32) 

Answer 1

スクリプトにはいくつか問題があります。まず、画像全体の計算を一度に行うことができます（十分なメモリがあれば、数GBの大きさの画像を処理できるはずです）。次に、書き込み後にファイルを閉じることはありません。部）。

from osgeo import gdal 

red_ds = gdal.Open("pathroredband.tif") 
nir_ds = gdal.Open("pathtonirband.tif") 

# Assuming that these are reflectance, float type. 
red = red_ds.ReadAsArray() 
nir = nir_ds.ReadAsArray() 

# To save memory, we will use in-place operations 
num = red.copy() # Numerator 
num -= nir 
nir += red # Denominator 
num /= nir # num now contains NDVI 

# Write the output 
driver = gdal.GetDriverByName('GTiff') 
ds_out = driver.Create(path, xsize, ysize, 1, gdal.GDT_Float32) 
ds_out.GetRasterBand(1).WriteArray(num) 
ds_out = None # Close the dataset; important!