У меня есть два растровых изображения, одно из группы 4 с B4 в конце и другое из группы 5 с B5 в конце.Я хочу настроить растр B5 на 800x600, затем отобразить его и сохранить как GeoTiff.Затем я хочу вычислить NDVI (я предполагаю, что для этого мне понадобятся и B4, и B5, но я не уверен).Затем я хочу отобразить подмножество NDVI растра B5.Отобразите его и сохраните как GeoTiff.
Как бы я создал что-то вроде подмножества растрового изображения TIFF размером 800 x 600 пикселей?Я также хочу взять этот TIFF и сгенерировать NDVI-образ для этого подмножества.
ПРИМЕЧАНИЕ. Я работаю с изображением Landsat.На изображении в конце заголовка файла стоит буква B5.
Что я сделал до сих пор:
import rasterio
from rasterio.windows import Window
import matplotlib.pyplot plt # for later use
with rasterio.open('MyRasterImage.tif') as src:
w = src.read(1, window=Window(0, 0, 800, 600))
Я хочу отобразить это, используя записные книжки Spyder или Jupyter.Поэтому я подумал использовать matplotlib и выполнил текущий код:
# Plot
plt.imshow(w)
plt.show()
При этом создается окно matplotlib размером 800x600, но оно все фиолетовое, не знаю, зачем оно это создает.
Теперь я хочучтобы иметь возможность отображать это изображение 800x600.Затем, после этого, я хочу предварительно сформировать NDVI для этого изображения подмножества 800x600.Затем отобразите подмножество 800x600 с отображением NDVI.
Я знаю, что формула: NDVI = (NIR - красный) / (NIR + красный)
Но как мне извлечь NIR и красныйвот из этого единственного изображения Landsat?
Моя попытка сделать это:
band1 = dataset.read(1)
band2 = dataset.read(2)
band3 = dataset.read(3)
print(band[2])
Когда я запускаю этот код для полос, я получаю ошибку:
rasterio indexerror: band index 2 out of range (not in (1,))
КогдаЯ запускаю этот код:
print(w.count)
Возвращает '1'.
Значит, это означает, что изображение Landsat имеет только одну полосу?Но для того, чтобы сделать NDVI, мне не нужны 3 полосы?
Я думаю написать некоторый код, подобный этому, чтобы получить NDVI из этого растра.Но я не уверен, что делать с извлечением полос:
# We handle the connections with "with"
with rasterio.open(bands[0]) as src:
b3 = src.read(1)
with rasterio.open(bands[1]) as src:
b4 = src.read(1)
# Allow division by zero
numpy.seterr(divide='ignore', invalid='ignore')
# Calculate NDVI
ndvi = (b4.astype(float) - b3.astype(float)) / (b4 + b3)
Этот код не работает, потому что полосы не определены как что-либо, поэтому я не знаю, как определить полосы для получения NDVI.
После этого я не уверен, что делать с отображением и сохранением отрендеренного изображения.