2

I am trying to calculate Multispectral Vegetation Indices from multi spectral photos, using Python.

So far, I've managed to calculate the NDVI, which is easy to implement as there are many examples online.

#Read the image file
image_file = "my_tif_image.tif"


with rasterio.open(image_file) as src:
    band_red = src.read(3)

with rasterio.open(image_file) as src:
    band_nir = src.read(4)
   

# Allow division by zero
numpy.seterr(divide='ignore', invalid='ignore')

band_nir = numpy.array(band_nir)
band_red = numpy.array(band_red)

NDVI = (band_nir.astype(float) - band_red.astype(float)) / (band_nir + band_red)

fig1, ax = plt.subplots(figsize=(12, 12))
fig1 = plt.gcf()
ep.plot_bands(NDVI,
              cmap='PiYG',
              scale=False,
              cbar = False,  #Change this to TRUE to display the color bar
              ax=ax,
              vmin=-1, vmax=1,
              title="NDVI")
plt.show()
fig1.savefig('NDVI.png',bbox_inches='tight')

However, apart from the above index, I couldn't manage to find anywhere examples on how to calculate other indexes (at least the ones I am interested in).

For example, instead of the NDVI index, I tried the formula of Burned Area Index:

BIR = 1.0 / (((0.1 + band_red) ** 2) + ((0.06 + band_nir) ** 2))

but I get an empty plot and the calculated array numbers are very small (in the numbers of 10^-9).

My question is, why is it so hard to find examples that calculate and plot indices other than the NDVI? And lastly, if I can't do it with python, is there any other programmatic way to do it?

I've stumbled across the Arcpy, but I wasn't able to import that module. Do I have to have the ArcGIS software in order to use that library?

3

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.