Загружайте, экспериментируйте и скачивайте геотифы, оптимизированные для облака (COG), используя Python с Google Colab.

Для доступа к спутниковым изображениям требуется много места для хранения, когда вы загружаете их локально. Этот процесс может занять все ваше хранилище, если вам нужно несколько графиков в течение длительного периода. Хуже того, вам может понадобиться только часть всего изображения.

Чтобы справиться с такими проблемами, последние достижения в области геопространственных данных используют то, что мы Cloud Optimized Geotiffs (COG).

Оптимизированный для облака GeoTIFF (COG) - это обычный файл GeoTIFF, предназначенный для размещения на файловом сервере HTTP с внутренней организацией, которая обеспечивает более эффективные рабочие процессы в облаке. Это достигается за счет использования возможности клиентов, отправляющих запросы диапазона HTTP GET, запрашивать только те части файла, которые им нужны. - https://www.cogeo.org/

А теперь представьте, что вы можете перематывать, перематывать и останавливать большое видео, чтобы посмотреть нужную вам часть видео. COG позволяет делать это с помощью спутниковых снимков. Вы можете получить плитку с отверстиями, визуализировать ее с низким разрешением, подмножество и замаскировать область и даже выполнить любую обработку, даже не загружая ее.

Из этого руководства вы узнаете, как получить доступ к изображениям Landsat, хранящимся в хранилище AWS s3, прямо в Google Colab с помощью Python. В первой части рассказывается, как найти подходящее изображение для интересующей вас области, а во второй части показано, как получить доступ, визуализировать и обработать спутниковое изображение в Python.

Поиск спутниковых изображений для интересующей вас области

Я получаю много запросов о том, как получить подходящие изображения для интересующей области. В этом случае мы используем изображения Landsat, поэтому я расскажу вам, как получить путь и строку спутникового изображения. Перейдите на сайт USGS explorer и найдите интересующее вас место, указав адрес или место (выделено красным на рисунке ниже). Как только вы получите результаты, нажмите на адрес. Вы можете отфильтровать по дате, если хотите.

Следующим шагом является щелчок по набору данных и выбор ваших наборов данных. В этом уроке мы используем Landsat 8, поэтому выбираем его (выделено красным ниже)

Щелкните Дополнительные критерии, чтобы отфильтровать облачное покрытие. Как показано на рисунке ниже, мы отфильтровываем облачное покрытие до менее 10%.

Наконец, нажмите на результаты, чтобы узнать все изображения, доступные в вашем районе. Как только вы узнаете изображение, которое хотите использовать для анализа, вы можете просто скопировать идентификатор или где-нибудь записать его (как показано ниже).

Вам также необходимо записать путь и ряд на спутниковом снимке. Вы можете нажать кнопку метаданных (красный прямоугольник) или просто записать третью часть идентификатора: 0030065.

В следующем разделе мы увидим, как мы можем получить доступ к данным непосредственно в Google Colab с помощью Python. Landsat хранится как в AWS, так и в Google Cloud Platform, но в этом руководстве мы получаем данные с помощью AWS.

Доступ к Landsat с AWS и Google Colab

Давайте сначала импортируем библиотеки, которые мы собираемся использовать. Нашим основным инструментом является Rasterio, который предоставляет простой в использовании API для обработки спутниковых изображений.

import numpy as np
import matplotlib.pyplot as plt
import rasterio as rio
from rasterio.plot import show

Вы можете построить свой URL-путь следующим образом.

fpath = ‘http://landsat-pds.s3.amazonaws.com/c1/L8/003/065/LC08_L1TP_003065_20190925_20191017_01_T1/LC08_L1TP_003065_20190925_20191017_01_T1_B4.TIF'

Мы можем разделить указанный выше URL на разные части:

http://landsat-pds.s3.amazonaws.com/c1/

Первая часть URL-адреса всегда одинакова и показывает URL-адрес хранилища AWS.

L8/003/065/

L8 означает Landsat 8. Путь: 003. Строка: 065. Эта часть URL-адреса будет меняться в зависимости от набора данных, к которому вы обращаетесь, и интересующей области.

LC08_L1TP_003065_20190925_20191017_01_T1

Эта часть является идентификатором изображения, которое вы можете получить, как показано на изображении выше.

LC08_L1TP_003065_20190925_20191017_01_T1_B4.TIF

Наконец, вы повторяете идентификатор изображения, но на этот раз с Band, к которому вы хотите получить доступ (здесь B4), а также с расширением изображения (.TIFF)

Теперь, когда мы настроили URL-адрес доступа, мы можем написать простую функцию для открытия изображения с помощью Rasterio.

def rasterio_open(f):
return rio.open(f)
src_image = rasterio_open(fpath)

Как только мы открываем растровое изображение; вы можете построить его, используя любую библиотеку визуализации или только методы визуализации данных Rasterio.

fig, ax = plt.subplots(1, figsize=(12, 10))
show(src_image, ax=ax)
plt.show()

Здесь отображается все изображение, к которому мы получили доступ. Мы получили только Band 4. Позже мы увидим, как мы можем комбинировать разные полосы для создания изображения RGB.

Как видите, у нас есть значения NaN, которые отображаются черным цветом. Мы можем избавиться от них, используя функциональность Numpy. Сначала мы конвертируем изображение в массивы Numpy, просто читая его с помощью Rasterio.

src_image_array = src_image.read(1)
src_image_array = src_image_array.astype(“f4”)
src_image_array

Изображение преобразуется в массивы, как показано ниже.

array([[0., 0., 0., …, 0., 0., 0.], [0., 0., 0., …, 0., 0., 0.], [0., 0., 0., …, 0., 0., 0.], …, [0., 0., 0., …, 0., 0., 0.], [0., 0., 0., …, 0., 0., 0.], [0., 0., 0., …, 0., 0., 0.]], dtype=float32)

Теперь вы можете легко удалить эти нулевые массивы, назначив им np.nan.

src_image_array[src_image_array==0] = np.nan
fig, ax = plt.subplots(1, figsize=(12, 10))
show(src_image_array, ax=ax)
plt.show()

У нас больше нет нулевых значений по краям.

Как это часто бывает, вас может интересовать только часть изображения. В следующем разделе мы расскажем, как разбить изображения на части с помощью функции окна Rasterio.

Подмножество изображений

Чтобы получить доступ только к определенной части изображения, вы можете отфильтровать по строкам, столбцам, ширине и высоте изображения. Скажем, нам нужно не все изображение, а изображение размером 750 X 850 (ширина и высота), обрезанное на 1200 столбцов и 1200 строк.

# Window(col_off, row_off, width, height)
window = rio.windows.Window(1200, 1200, 750, 850)
subset = src_image.read(1, window=window)
fig, ax = plt.subplots(1, figsize=(12, 10))
show(subset, ax=ax)
ax.set_axis_off()
plt.show()

Теперь изображение уменьшено до подмножества около города Эйрунепе, как показано ниже.

Наконец, в следующем разделе показано, как создать изображение RGB и загрузить его локально.

Создать RGB и скачать

Сначала мы обращаемся к каждому бэнду, предоставляя каждому отдельный URL. Обратите внимание, что все URL-адреса одинаковы, за исключением последней части перед расширением .TIF. Это полосы, и в этом случае, поскольку мы хотим создать изображение RGB, мы получаем доступ к полосе 4 (красный), полосе 3 (зеленый) и полосе 2 (синий).

rpath = ‘http://landsat-pds.s3.amazonaws.com/c1/L8/003/065/LC08_L1TP_003065_20190925_20191017_01_T1/LC08_L1TP_003065_20190925_20191017_01_T1_B4.TIF'
gpath = ‘http://landsat-pds.s3.amazonaws.com/c1/L8/003/065/LC08_L1TP_003065_20190925_20191017_01_T1/LC08_L1TP_003065_20190925_20191017_01_T1_B3.TIF'
bpath = ‘http://landsat-pds.s3.amazonaws.com/c1/L8/003/065/LC08_L1TP_003065_20190925_20191017_01_T1/LC08_L1TP_003065_20190925_20191017_01_T1_B2.TIF'
red = rio.open(rpath)
green = rio.open(gpath)
blue = rio.open(bpath)

Мы создаем изображение RGB, открывая пустое изображение и заполняя его тремя полосами, которые мы открыли выше.

# Create an RGB image
with rio.open(‘RGB.tiff’,’w’,driver=’Gtiff’, width=red.width, height=red.height,count=3,crs=red.crs,transform=red.transform, dtype=red.dtypes[0]) as rgb:
rgb.write(blue.read(1),1)
rgb.write(green.read(1),2)
rgb.write(red.read(1),3)
rgb.close()

Изображение теперь сохранено. И вы можете открыть его с помощью Rasterio или загрузить локально.

Заключение

В этом руководстве мы рассмотрели, как вы можете искать и находить нужные данные для интересующей вас области. Мы также увидели, как получить доступ к изображениям Landsat с помощью хранилища Python и AWS.

Код для этого руководства доступен в этом блокноте Colab.