Autor
Guía paso a paso sobre cómo recopilar y preprocesar imágenes de satélite data para algoritmos de aprendizaje automático
Este artículo es la primera parte de una serie de dos que explora los usos de los algoritmos de aprendizaje automático en imágenes de satélite. Aquí nos centraremos en los pasos de recopilación y preprocesamiento de data , mientras que la segunda parte mostrará el uso de varios algoritmos de aprendizaje automático.
RESUMEN
Este artículo:
Utilizaremos como ejemplo la detección y clasificación de campos agrícolas. Se utilizan Python y jupyter notebook para los pasos de preprocesamiento. El cuaderno jupyter para reproducir los pasos está disponible en github.
Este artículo asume fundamentos básicos en data science y python.
La visión por ordenador es una fascinante rama del aprendizaje automático que se ha desarrollado rápidamente en los últimos años. Con el acceso público a muchas bibliotecas de algoritmos(Tensorflow, OpenCv, Pytorch, Fastai...), así como a conjuntos de datos de código abierto(CIFAR -10, Imagenet, IMDB...), es muy fácil para un científico de data adquirir experiencia práctica en este tema.
Pero ¿sabía que también puede acceder a imágenes de satélite de código abierto? Estos tipos de imágenes pueden utilizarse en multitud de casos de uso (agricultura, logística, energía...). Sin embargo, suponen retos adicionales para un científico de data debido a su tamaño y complejidad. He escrito este artículo para compartir algunos de mis principales aprendizajes trabajando con ellas.
Paso 1 - Seleccione el satélite abierto data fuente
En primer lugar, debe elegir el satélite que desea utilizar. Las principales características a tener en cuenta son:
Estos son los satélites cuya data está disponible gratuitamente :
Sistemas de satélites de referencia
Dado que nuestro objetivo es detectar e identificar cultivos agrícolas, el espectro visible (que se utiliza para proporcionar imágenes en color estándar) debería ser bastante útil. Además, existen correlaciones entre el nivel de nitrógeno de los cultivos y las bandas espectrales NIR (infrarrojo cercano). En cuanto a las opciones libres, Sentinel 2 es actualmente el mejor en términos de resolución con 10 m para esas bandas espectrales, por lo que vamos a utilizar este satélite.
Por supuesto, la resolución de 10 m es muy inferior a la de las imágenes de satélite más avanzadas, con satélites comerciales que alcanzan una resolución de 30 cm. El acceso a los satélites comerciales data puede desbloquear muchos casos de uso, pero a pesar de la creciente competencia con la entrada de nuevos actores en este mercado, siguen teniendo un precio elevado y, por tanto, son difíciles de integrar en un caso de uso a escala.
Paso 2 - Recoger las imágenes de satélite data de Sentinel 2
Ahora que hemos seleccionado la fuente de nuestro satélite, el nuevo paso es descargar las imágenes. Sentinel 2 forma parte del programa Copérnico de la Agencia Espacial Europea. Podemos acceder a su data utilizando la API python de SentinelSat.
Primero tenemos que crear una cuenta en el Copernicus Open Access Hub. Alternativamente, tenga en cuenta que este sitio web también ofrece una forma de consultar imágenes con una GUI.
Una forma de especificar las zonas que desea explorar es proporcionar a la API un archivo geojson, que es un json que contiene las coordenadas GPS de posición (latitud y longitud) de una zona. Geojson puede crearse rápidamente en este sitio web.
De este modo, podrá consultar todas las imágenes de satélite que se crucen con la zona indicada.
Aquí se enumeran todas las imágenes disponibles en una zona y franja horaria determinadas:
api = SentinelAPI(
credenciales["nombre de usuario"],
credenciales["contraseña"],
"https://scihub.copernicus.eu/dhus"
)
shape = geojson_to_wkt(read_geojson(geojson_path))
imágenes = api.query(
forma,
date=(date(2020, 5, 1), date(2020, 5, 10)),
platformname="Centinela-2″,
processinglevel="Level-2A",
porcentaje de nubosidad=(0, 30)
)
imágenes_df = api.to_dataframe(imágenes)
Obsérvese el uso de argumentos:
imágenes_df
Imágenes_df es un marco de datos pandas que contiene todas las imágenes que coinciden con nuestra consulta. Seleccionamos una con poca cloud y la descargamos usando la api :
api.download(uuid)
Paso 3 - Saber qué bandas espectrales utilizar y gestionar la cuantización
Una vez descargado el archivo zip podemos ver que en realidad contiene varias imágenes.
Esto se debe a que como dijimos en el paso 1, Sentinel 2 captura muchas bandas espectrales y no sólo Rojo Verde y Azul y cada banda es capturada en una sola imagen.
Además sus instrumentos no tienen todos la misma resolución espacial, por lo que hay una subcarpeta para cada resolución diferente 10, 20 y 60 metros.
If we check the GRANULE/<image_id>/IMG_DATA folder we see that there is a subfolder for each different resolution: 10, 20 and 60 meters as the spectral bands. We are going to use those with max resolution 10m as they provide Band 2, 3, 4 and 8 which are Blue, Green, Red and NIR. Band 7, 8, 9 which are vegetation red edge (spectral bands between red and NIR that mark a transition in reflectance value for vegetation) could also be interesting for detecting crops but since they are lower resolution we are going to leave them aside now
Aquí cargamos cada banda como una matriz numpy 2D:
def get_band(image_folder, band, resolution=10):
subfolder = [f for f in os.listdir(image_folder + “/GRANULE”) if f[0] == “L”][0]
image_folder_path = f”{image_folder}/GRANULE/{subfolder}/IMG_DATA/R{resolution}m”
image_files = [im for im in os.listdir(image_folder_path) if im[-4:] == “.jp2”]
selected_file = [im for im in image_files if im.split(“_”)[2] == band][0]
with rasterio.open(f”{image_folder_path}/{selected_file}”) as infile:
img = infile.read(1)
devolver img
band_dict = {}
para banda en ["B02", "B03", "B04", "B08"]:
diccionario_bandas[banda] = get_bandas(carpeta_imagen, banda, 10)
Una prueba rápida para asegurarse de que todo está funcionando es tratar de recrear la imagen en color verdadero (es decir, RGB) y mostrarla usando opencv y matplotlib:
img = cv2.merge((band_dict["B04"], band_dict["B03"], band_dict["B02"]))
plt.imshow(img)
Imagen de satélite - Problema de intensidad de color (Copernicus Sentinel data 2020)
Sin embargo, hay un problema con la intensidad del color. Esto se debe a que la cuantización (número de valores posibles) difiere del estándar para matplotlib que asume un valor int 0-255 o un valor float 0-1. Las imágenes estándar suelen utilizar una cuantización de 8 bits (256 valores posibles), pero las imágenes de Sentinel 2 utilizan una cuantización de 12 bits y un reprocesado convierte los valores a un entero de 16 bits (65536 valores).
Si queremos escalar un entero de 16 bits a un entero de 8 bits, en teoría deberíamos dividirlo por 256, pero en realidad esto da como resultado una imagen muy oscura, ya que no se utiliza toda la gama de valores posibles. El valor máximo de nuestra imagen es en realidad 16752 de 65536 y pocos píxeles alcanzan valores superiores a 4000, por lo que dividir por 8 da una imagen con un contraste decente.
img_processed = img / 8
img_processed = img_processed.astype(int)
plt.imshow(img_processed)
Imagen de satélite - Intensidad de color derecha (Copernicus Sentinel data 2020)
Ahora tenemos un tensor de 3 dimensiones de forma (10980, 10980, 3) de enteros de 8 bits en forma de array numpy
Tenga en cuenta que :
Paso 4 - Dividir las imágenes en tamaños aptos para el aprendizaje automático
Podríamos intentar aplicar un algoritmo a nuestra imagen actual, pero en la práctica eso no funcionaría con la mayoría de las técnicas debido a su tamaño. Incluso con una gran potencia de cálculo, la mayoría de los algoritmos (y especialmente el aprendizaje profundo) quedarían fuera de juego.
Así que debemos dividir la imagen en fragmentos. Una pregunta es cómo establecer el tamaño del fragmento.
Decidimos dividir cada imagen en una cuadrícula de 45 * 45 (siendo 45 convenientemente un divisor de 10980)
Pasada esta fijación de parámetros, este paso es sencillo utilizando la manipulación de matrices numpy . Almacenamos nuestros valores en un dict usando tuplas (x, y) como claves.
frag_count = 45
frag_size = int(img_processed.shape[0] / frag_count)
frag_dict = {}
for y, x in itertools.product(range(frag_count), range(frag_count)):
frag_dict[(x, y)] = img_processed[y*tamaño_frag: (y+1)*tamaño_frag,
x*tamaño_frag: (x+1)*tamaño_frag, :]
plt.imshow(frag_dict[(10, 10)])
Fragmento de la imagen del satélite (Copernicus Sentinel data 2020)
Paso 5 - Vincule su data a las imágenes de satélite mediante la conversión de GPS a UTM
Por último, tenemos que vincular nuestras imágenes con el sitio data que tenemos sobre los campos agrícolas para permitir un enfoque de aprendizaje automático supervisado. En nuestro caso, tenemos data sobre los campos agrícolas en forma de lista de coordenadas GPS (latitud, longitud), similar al geojson que utilizamos en el paso 2, y queremos ser capaces de encontrar la posición de sus píxeles. Sin embargo, en aras de la claridad, voy a empezar por mostrar una conversión de píxeles a coordenadas GPS y luego mostrar la inversa
Los metadatos recogidos al solicitar la API de SentinelSat proporcionan las coordenadas GPS de las esquinas de la imagen (columna de la huella) y nosotros queremos las coordenadas de un píxel concreto. La latitud y la longitud, que son valores angulares, no evolucionan linealmente dentro de una imagen, por lo que no podemos hacer una simple relación utilizando la posición del píxel. Hay maneras de resolver esto utilizando ecuaciones matemáticas, pero las soluciones ya existen en python.
Una solución rápida es convertir la posición del píxel a UTM (Universal Transverse Mercator) en el que una posición está definida por una zona, así como una posición (x, y) (en unidad de metros) que evoluciona linealmente con la imagen. Podemos entonces utilizar una conversión de UTM a Latitud / Longitud proporcionada en la biblioteca utm. Para ello, primero tenemos que obtener la zona UTM y la posición de la esquina superior izquierda de la imagen que se puede encontrar en los metadatos de la imagen en color verdadero.
#Obtención de meatadatos
transform = rasterio.open(tci_ruta_archivo, driver='JP2OpenJPEG').transform
número_zona = int(tci_ruta_archivo.split("/")[-1][1:3])
letra_zona = tci_file_path.split("/")[-1][0].
utm_x, utm_y = transform[2], transform[5]
# Conversión de la posición del píxel a utm
este = utm_x + pixel_columna * 10
norte = utm_y + pixel_row * - 10
# Convertir UTM a latitud y longitud
latitud, longitud = utm.to_latlon(este, norte, número_zona, letra_zona)
La conversión de posición GPS a píxel se realiza utilizando las fórmulas inversas. En ese caso tenemos que asegurarnos de que la zona obtenida coincide con nuestra imagen. Esto puede hacerse especificando la zona. Si nuestro objeto no está presente en la imagen, se obtendrá un valor de píxel fuera de límites.
# Convertir latitud y longitud a UTM
este, norte, número_zona, letra_zona = utm.from_latlon(
latitud, longitud, forzar_número_zona=número_zona
)
# Convertir UTM a columna y fila
pixel_columna = round((este - utm_x) / 10)
pixel_row = round((norte - utm_y) / -10)
ver en bruto
Conclusión
Ya estamos listos para utilizar las imágenes de satélite para el aprendizaje automático.
En el próximo artículo de esta serie de 2 partes veremos cómo podemos utilizar esas imágenes para detectar y clasificar superficies de campos utilizando aprendizaje automático supervisado y no supervisado.
Gracias por leer y no dude en seguir el blog Artefact tech si desea que le avisemos cuando se publique el próximo artículo.
