Como extrair valores de pixels de uma coleção de imagens no Google Earth Engine com Python?

Confira nosso curso online de Word

Aprenda como obter os valores dos pixels de uma coordenada dentro de uma coleção de imagens no Google Earth Engine usando Python e crie gráficos deles.

O Google Earth Engine é uma ferramenta com um extenso banco de imagens de satélite. Além da possibilidade de trabalhar com imagens do mundo inteiro, também temos acesso à imagens históricas.

Diante disso, pode ser necessário para alguns estudos ambientais conhecer a evolução de determinados indicadores (e.g. NDVI ou Temperatura de Superfície) para uma determinada área.

Por isso, nessa postagem vamos mostrar como filtrar uma coleção de imagens usando datas específicas (e.g. Datas sem um intervalo contínuo); como criar uma função em Python; como criar data frames com a biblioteca pandas; e como criar gráficos com os dados obtidos.

No final, teremos um gráfico mostrando a variação da temperatura de superfície de uma área em recuperação ambiental no município de Treviso/SC.

Ajustando as bibliotecas e dados de entrada

Inicialmente, vamos importar as bibliotecas necessárias para rodar nosso código, isto é, as bibliotecas ee (Google Earth Engine), matplotlib.pyplot, datetime e pandas.

import ee
# ee.Authenticate() # Descomentar se não foi feito a autenticação ainda
ee.Initialize()

import matplotlib.pyplot as plt
import datetime as dt
import pandas as pd

Além dessas bibiotecas, também vamos usar uma função, já apresentado em outra postagem, para obter as datas de uma coleção de imagens.

def ymdList(imgcol):
    def iter_func(image, newlist):
        date = ee.Number.parse(image.date().format("YYYYMMdd"));
        newlist = ee.List(newlist);
        return ee.List(newlist.add(date).sort())
    ymd = imgcol.iterate(iter_func, ee.List([]))
    return list(ee.List(ymd).reduce(ee.Reducer.frequencyHistogram()).getInfo().keys())

Agora vamos definir alguns dados de entrada, tais como coordenadas geográfica do nosso ponto de interesse e o intervalo das datas.

ee_coord = ee.Geometry.Point([-49.4548067,-28.5375093])
startDate = '2013-06-01'
endDate = '2020-06-01'

As coordenadas geográficas se referem a uma área em recuperação ambiental, sendo as datas arbitrárias.

É importante lembrar que, quando maior for o intervalo da data ou quanto maior for o número de imagens, pode ser que seu código fique rodando por muito tempo (em função da grande quantidade de dados).

Filtrar coleção de imagens com datas específicas

Iremos buscar os dados de temperatura de superfície do sensor MODIS, o qual obtem dados diários da superfície terrestre com uma resolução de 1.000 m. Para não termos uma grande quantidade de dados, vamos filtrar as datas usando o intervalo de tempo de passagem do satélite LANDSAT, que é de 15 dias.

Vamos, primeiro, criar uma variável para armazenar a coleção LANDSAT e filtrá-la com a nossa localização (ref_coll). Em seguida, vamos criar outra variável para filtrar essa nossa coleção com o nosso intervalo de dadas (ref_date).

ref_coll = ee.ImageCollection("LANDSAT/LC08/C01/T1_SR")\
    .filterBounds(ee_coord)

ref_date = ref_coll.filterDate(startDate, endDate)

Com isso em mãos, vamos utilizar a função ymdList() para ter as datas das imagens armazenadas na nossa coleção. Depois, iremos passar uma rotina para adicionar um traço entre o ano, mês e dia e por fim, vamos aplicar uma função para converter essas datas em milisegundos (de forma que o Google Earth Engine possa usar elas para filtragem).

Note que usamos duas formas diferentes na programação deste bloco. Primeiro, usamos um loop for em apenas uma linha e segundo, usamos o lambda do Python para passar uma função com apenas uma linha também.

ref_coll = ee.ImageCollection("LANDSAT/LC08/C01/T1_SR")\
    .filterBounds(ee_coord)

ref_date = ref_coll.filterDate(startDate, endDate)

dt_0 = ymdList(ref_date)
dates = [dt_0[i][0:4]+'-'+dt_0[i][4:6]+'-'+dt_0[i][6:8] for i in range(len(dt_0))]
date_list = ee.List(dates).map(lambda d: ee.Date(d).millis())

Com tudo isso em mãos, agora vamos passar nossa lista de datas em milisegundos (date_list) para nossa coleção de imagens MODIS.

MODIS_coll = ee.ImageCollection("MODIS/006/MOD11A1")\
    .select("LST_Day_1km")\

coll = MODIS_coll.filter(ee.Filter.inList("system:time_start", date_list))

Agora já temos nossa coleção de imagens MODIS filtrada com as datas das imagens LANDSAT.

Como criar uma função em Python?

A definição de funções em Python é realizada da seguinte forma: Usamos def para indicar ao python que estamos iniciando uma função, damos um nome à ela em seguida e colocamos entre parênteses as variáveis de entrada.

Em seguida, finalizamos essa primeira linha com dois pontos e começamos a função.

No nosso caso, vamos criar uma função chamada get_point_coll com os seguintes dados de entrada: coleção de imagens, coordenada e escala.

Na nossa função, iremos obter a quantidade de imagens disponíveis na coleção (num_img) e vamos converter a coleção de imagens para uma lista (coll_list). Com esses dados salvos, vamos criar um loop for para pegar uma imagem da lista (img_temp), obter a data dela (e converter ela de milisegundos para ano, mês e dia) (tmstpA, tmstpB e tmstpC), extrair o valor do pixel na nossa coordenada de interesse (eeDict_value e value_pixel) e salvar esses valores (data e valor do pixel) numa variável para retornar ela como resultado da função (value_output).

def get_point_coll(img_coll, coord, escala):
    num_img = img_coll.size().getInfo()
    coll_list = img_coll.toList(num_img)
    values_output = []
    
    for i in range(0, num_img, 1):
        img_temp = ee.Image(coll_list.get(i))
        
        tmstpA = ee.Date(img_temp.get('system:time_start'))
        tmstpB = tmstpA.getInfo()['value']/1000
        tmstpC = dt.datetime.utcfromtimestamp(tmstpB).strftime('%Y-%m-%d')

        eeDict_value = img_temp.reduceRegion(reducer = ee.Reducer.first(), 
                                             geometry = coord, scale = escala)
        value_pixel = list(eeDict_value.getInfo().values())[0]
        
        values_output.append((tmstpC, value_pixel))
    
    return values_output

Com a nossa função pronta, vamos fornecer nossos dados para ela e salvar os valores dos pixels da nossa coleção de imagens.

tabela = get_point_coll(coll, ee_coord, 1000)

Criando um data frame com Pandas no Python

Agora que já temos os dados salvos na variável tabela, vamos criar um data frame com duas colunas, uma para as datas e outra para as temperaturas. Também iremos remover os valores em branco com a função dropna().

É importante comentar que os valores dos pixels, muitas vezes não são as unidades físicas propriamente ditas, sendo necessário realizar algumas operações para convertê-las. Para os dados do MODIS, na banda LST_Day_1km, é necessário multiplicar o valor por 0,02.

df_cols = ['Data', 'LST_K']
df = pd.DataFrame(tabela, columns = df_cols).dropna()
df['LST_K'] = 0.02*df['LST_K'] # Converte valores MODIS para K

Como criar gráficos no Python?

E com os nossos dados organizados, iremos criar um gráfico de linha. Você terá que fornecer alguns dados básicos, como quem é o eixo X, eixo Y, cor, mas também iremos configurar para que o eixo Y sempre fique próximo dos valores obtidos (parâmetro ylim).

No final do código, também vamos salvar os dados obtidos em um arquivo CSV.

grafico = df.plot(kind = 'line', x = 'Data', y = 'LST_K', color = 'darksalmon', rot = 90, 
                  ylim = (round((min(df['LST_K'])*0.95),0), round((max(df['LST_K'])*1.05),0)),
                  grid = True, marker = '.', label='LST (K)')
grafico.set_ylabel('Temperatura (K)')
grafico.grid(linestyle='dotted')
grafico.set_axisbelow(True)
plt.show()

df.to_csv(r'C:\Users\ferna\Desktop\dados_point.csv')

Ao rodar todo o nosso código, você irá ter um gráfico conforme a figura abaixo.

Variação de temperatura com imagem MODIS (LST_Day_1km).

E com isso terminamos nosso tutorial. Você pode conferir o código completo no meu GitHub. Caso tenha ficado com alguma dúvida, é só perguntar nos comentários.



Clique na figura abaixo e assine nossa lista de emails para receber nosso ebook "Como criar mapas de localização com ArcGIS 10.x".

Apostila Mapa de Localização Banner

Author: Fernando BS

Engenheiro Ambiental e de Segurança do Trabalho. Atua nas áreas de geoprocessamento, mineração e hidrologia. Busca soluções utilizando softwares como QGIS, R e Python.

Deixe uma resposta

O seu endereço de e-mail não será publicado. Campos obrigatórios são marcados com *