John Cunha Universidade Federal de Campina Grande Paraíba, Brasil

Utilização da detecção remota para estimar a distribuição espacial do

balanço de energia e a evapotranspiração em regiões semi-

áridasJohn CunhaUniversidade Federal de Campina GrandeParaíba, Brasil

Lisboa, 10 de novembro de 2010.

Projeto de Cooperação Internacional do Semi-Árido (CISA)

• O projeto CISA - Cooperação Internacional do Semi-árido é financiado pelo Ministério de Ciência e Tecnologia, com recursos do CTHIDRO - Fundo Setorial de Recursos Hídricos, repassados através da FINEP e CNPq.

• O projeto tem duração de 3 anos e objetiva promover cooperação entre instituições de pesquisa voltadas para a solução dos problemas gerados pela escassez de água nas regiões semi-áridas

Países parceiros no projeto

Brasil, Argentina, Cabo Verde, França, Portugal, Estados Unidos

Mais informaçõeswww.cisa.ufpe.br/twiki/bin/view/Cisa/WebHomewww.hidro.ufcg.edu.br/twiki/bin/view/Cisa/

Universidade Federal de Campina Grande (UFCG)

Figura 1 – Localização da cidade de Campina Grande e UFCG

Mais Informações:www.ufcg.edu.brwww.hidro.ufcg.edu.br/

• Geoprocessamento de dados termais e suas aplicações em estudos do regime pluviométrico no estado da Paraíba.

• O Projeto tem como objetivo identificar ilhas de calor nos principais centros urbanos do estado da Paraíba a partir de imagens do infravermelho termal (temperatura de superfície), em diferentes épocas e estabelecer correlações com possíveis alterações do regime pluviométrico destas áreas (e próximo a estas), utilizando técnicas de sensoriamento remoto e geoprocessamento

Projeto GeoTermal

• Escala local• Lisímetros;• Método das correlações turbulentas;• Método do balanço de energia baseado na Razão de Bowen;• Método aerodinâmico;• Balanço hídrico no solo.

• Escala regional• Utilização de técnicas de sensoriamento remoto baseado em

imagens de satélites.

Métodos de estimativa da evapotranspiração (ET)

• Sensores remotos

• Sensores aerotransportados

• Técnicas de Processamento Digital de Imagens (PDI)

Princípios básicos

• Satélites que possuam sensores que contenham as bandas de refletividade e infravermelho termal

• Landsat 5 sensor TM• Landsat 7 sensor ETM+• Satélites NOAA sensor AVHRR• Satélite TERRA sensores ASTER e MODIS• Satélite AQUA sensor MODIS• Satélite GOES

Satélites

• Brasileira - Instituto Nacional de Pesquisas Espaciais-INPE• Dados Gratuitos

• Imagens dos satélites • Produtos de sensoriamento remoto• Produtos de Modelos Climáticos• Softwares de Processamento Digital de Imagens (PDI) e de Sistema

de Informações Geográfica (SIG)

Política de dados

Mais informações:www.inpe.br

• Europa – EUMETSAT• Projetos EUMETCAST e GEONETCAST

• Dados Gratuitos• Imagens de satélites • Produtos de sensoriamento remoto• Produtos de Modelos Climáticos• Softwares de Processamento Digital de Imagens (PDI) e de Sistema

de Informações Geográfica (SIG)

Política de dados

Mais informações:http://www.eumetsat.int

• SPRING – INPE- Brasileiro

• ILWIS – ITC- Holandês

• ERDAS – Software comercial

Softwares

Algoritmo SEBAL

No SEBAL a evapotranspiração é obtida como resíduo da equação clássica do balanço de energia à superfície, qual seja:

onde Rn é o saldo de radiação, LE a densidade de fluxo de calor latente, H a densidade de fluxo de calor sensível e G a densidade de fluxo de calor no solo, todos em W.m-2.

LE = Rn - H - G

Saldo de Radiação

Radiação de onda curta

Radiação de onda Longa

onda curtaIncidente

onda curtarefletida

onda longaincidente

onda longarefletida

onda longaemitida

Superfície vegetada

Saldo de radiação = ganhos - perdas

Figura 2 – Balanço de Radiação na superfície (adaptado de Allen et al., 2002)

Saldo de Radiação

Fluxograma 1 – Etapas para obtenção do Saldo de Radiação (Rn)

Imagens de satélite (DN)

Banda 1 Banda 2 Banda 3 Banda 4 Banda 5 Banda 6 Banda 7

Saldo de Radiação

Figura 3 – Diferentes bandas do satélite LANDSAT

.

.

.

Imagens de satélite (DN)

Saldo de Radiação

Figura 4 - Composição colorida R3G2B1 e extração dos pixels para uma planilha


Saldo de Radiação

ND255

abaL iiiλi

Bandas Comprimento de Onda(μm)

Coeficientes de Calibração

a b1 (azul) 0,45 – 0,52 -1.52 193.0

2 (verde) 0,52 – 0,60 -2.84 365.0

3 (vermelho) 0,63 – 0,69 -1.17 264.0

4 (IV-próximo) 0,76 – 0,79 -1.51 221.0

5 (IV-médio) 1,55 – 1,75 -0.37 30.2

6 (IV-termal) 10,4 – 12,5 1.2378 15.303

7 (IV-médio) 2,08 – 2,35 -0.15 16.5

)μmsr(Wm 112

Saldo de Radiação

Radiância espectral monocromática

Tabela 1 - Descrição das bandas do Mapeador Temático (TM) do Landsat 5, com os correspondentes intervalos de comprimento de onda, coeficientes de calibração (radiância mínima – a e máxima – b), após 5 de maio de 2003. (Chander & Markham, 2003)

Saldo de Radiação

Radiância espectral monocromática

Figura 5 – Model Maker do Software ERDAS 9.2

$n4_Custom_Float + ($n5_Custom_Float - $n4_Custom_Float) * $n20_memory / 255

Saldo de Radiação


rλi

λiλi d.cos.k

L.πρ

Z

)365/2.cos(033,01 DSAdr

)90cos(cos Ez

Bandas Irradiância Espectral no Topo da Atmosfera

1 (azul) 1957

2 (verde) 1826

3 (vermelho) 1554

4 (IV-próximo) 1036

5 (IV-médio) 215,0

6 (IV-termal) -

7 (IV-médio) 80,67

)μm(Wm 12

Saldo de Radiação

Reflectância planetária

Tabela 2 - Descrição das bandas do Mapeador Temático (TM) do Landsat 5, com as correspondentes irradiâncias espectrais no topo da atmosfera terrestre (TOA) .(Chander & Markham, 2003)

Saldo de Radiação


754321toa 0,011ρρ0,033ρ0,157ρ0,233ρ0,274ρ0,293α

2sw

toa

τ

ααα atm

z2.100,75τ 5sw

Saldo de Radiação

Albedo Planetário

Albedo da Superfície

é a da radiação solar refletida pela atmosfera, que varia entre 0,025 e 0,04, mas para o modelo SEBAL é recomendado o uso do valor de 0,03, com base em Bastiaanssen (2000)

atmα

swτ é a transmissividade atmosférica que para condições de céu claro, pode ser obtida por (Allen et al., 2002)

Saldo de Radiação


VIV

VIV

ρρ

ρρNDVI

)ρρ(L

)ρL)(ρ(1SAVI

VIV

VIV

0,91

0,59SAVI0,69

ln

IAF

Saldo de Radiação

Índices de Vegetação

Índice de Vegetação da Diferença Normalizada (NDVI)

Índice de Vegetação Ajustado para os Efeitos do Solo (SAVI)

Índice de Área Foliar (IAF)

O seu valor mais frequente de L = 0,5 (Huete &Warrick, 1990; Accioly et al., 2002; Boegh et al., 2002)

Saldo de Radiação


IAF0,003310,97εNB IAF0,010,95ε0

Como cada pixel não emite radiação eletromagnética como um corpo negro, há a necessidade de introduzir a emissividade de cada pixel no domínio espectral da banda termal , qual seja: 10,4 – 12,5 μm. Por sua vez, quando do cômputo da radiação de onda longa emitida por cada pixel, há de ser considerada a emissividade no domínio da banda larga (5 – 100 μm). Segundo Allen et al. (2002), as emissividades e podem ser obtidas, para NDVI > 0 e IAF < 3, segundo:

Saldo de Radiação

Emissividades

Para pixels com , . Para corpos de água (NDVI < 0), no caso do lago de Sobradinho e do leite do Rio São Francisco, Silva & Cândido (2004) utilizaram os valores de 0,99 e 0,985, conforme Allen et al. (2002).

3IAF 0,98εε 0NB

Saldo de Radiação


1L

Kεln

KT

λ,6

1NB

2s

Onde e são constantes de calibração da banda termal do Landsat 5 –TM. (Allen et al., 2002; Silva et al., 2005)

Saldo de Radiação

Temperatura da Superfície

1121 μmsrWm607,76 k K1260,562 k

Saldo de Radiação


4s0L

T.σ.εR

é a emissividade de cada pixel; é a constante de Stefan-Boltzman; é a temperatura da superfície (K).

0εσ )KWm5,67.10(σ 428

sT

Saldo de Radiação

Radiação de Onda Longa Emitida

Saldo de Radiação


swrsτ.d.cos.SR Z

S é a constante solar (1367 w/m²)

Saldo de Radiação

Radiação de Onda curta incidente

Saldo de Radiação


é a emissividade atmosférica, obtida por:

4aaL T.σ.εR

aε

0,09swa )lnτ0,85.(ε (Allen et al., 2002)

é a temperatura do ar (K)aT

σ é a constante de Stefan-Boltzman

Saldo de Radiação

Radiação de onda longa incidente

Saldo de Radiação


LoLLss RRRRRRn )1(

Algoritmo SEBAL



LE = Rn - H - G

Fluxo de calor no solo

Rn)NDVI0,98α)(10.0074(0,0038G 4Ts

O fluxo de calor no solo (G) pode ser obtido segundo equação empírica desenvolvida por Bastiaanssen (2000), que representa valores próximos ao meio-dia:

Onde é a temperatura da superfície (°C), o albedo da superfície e NDVI o índice de vegetação da diferença normalizada, todos computados pixel a pixel. Para efeito de correção dos valores do fluxo de calor no solo para corpos de água (NDVI<0), pode ser utilizada a seguinte expressão: G = 0,3Rn, usada por Silva & Cândido (2004) ou G = 0,5Rn, segundo Allen et al.(2002).

Algoritmo SEBAL



LE = Rn - H - G

Fluxo de calor sensível

Equação do Fluxo de calor sensível

ah

p rρcH

dT

ρ é a massa específica do ar (1,15 kg m-3)Cp é o calor específico do ar a pressão constante (1004 J kg-1 K-1)dT é a diferença da temperatura entre dois níveis Z1 e Z2

rah é a resistência aerodinâmica ao transporte de calor (s/m)



Dados de uma estação agrometeorológica

h0,12z m0

A velocidade de fricção é computada usando o perfil logaritmo do vento para a condição de estabilidade neutra, qual seja:

m0

x

x*

zz

ln

kuu

Dessa forma, obtém-se o coeficiente de rugosidade local em função da altura média da vegetação segundo equação de Brutsaert (1982):



Velocidade do vento a 100 metros

Considerando-se, ainda, a atmosfera em equilíbrio neutro, é estimada a velocidade do vento ao nível de z = 100 m [ u100 (m/s)], chamada de blending height, onde se assume que os efeitos da rugosidade da superfície são desprezíveis, e que é dada pela equação:

k

z100

ln

uu m0*100


Fluxo de calor sensível Velocidade de fricção (u

*)

m0

100*

z100

ln

kuu

Com a hipótese de que é constante em toda a cena estudada, pode ser obtida a velocidade de fricção para cada pixel da imagem, através da equação:

onde pode ser obtido em função do SAVI segundo equação desenvolvida por Bastiaanssen (2000):

SAVI)5,625,809exp(z m0


k.u

zz

ln

r*

1

2

ah


Resistência aerodinâmica (rah)

A resistência aerodinâmica é computada admitindo-se a atmosfera em condição de estabilidade neutra, pela seguinte expressão:

Onde Z1 e Z2 são as alturas (em metros) acima da superfície.



Pixels “Frio” e “Quente”O cômputo da diferença de temperatura próxima à superfície dT para cada pixel é computada pelo SEBAL através de uma relação linear entre dT e (temperatura da superfície):

sbTadT

onde os coeficientes a e b são obtidos através dos pixels âncoras (quente e frio). O pixel “frio” da imagem é escolhido admitindo-se que este se encontra numa área bem irrigada, onde se assume que o fluxo de calor sensível é nulo e o fluxo de calor latente LE é dado por:

GRnfrio LE


Pixels “Frio” e “Quente”Por sua vez, o pixel “quente” é escolhido numa área com grande fração de solo exposto, onde se assume que o fluxo de calor latente é nulo (LE = 0) e o fluxo de calor sensível H quente (w/m²) será, portanto, dado por:

ah

spquente r

)bT(aρcGRnH


pahs /ρG)(RnrbTa cdT

0bTa s dT



ah

p rρcH

dT

Pixel “quente”

Pixel “frio”




ah

p rρcH sbTa

No entanto, os valores obtidos não representam adequadamente o H de cada pixel e servem, tão somente, como valores iniciais de um processo iterativo, e que nas etapas seguintes são consideradas, efetivamente, a condição de estabilidade de cada pixel. Em virtude dos efeitos turbulentos afetar as condições atmosféricas e a resistência aerodinâmica, aplica-se a teoria da similaridade de Monin-Obukhov.



Comprimento de Monin-Obukhov

O comprimento de Monin-Obukhov L é utilizado para definir as condições de estabilidade da atmosfera e é computado em função dos fluxos de calor sensível e de momentum pela seguinte expressão:

kgH

TuρcL s

3*p

onde ρ é a densidade do ar, Cp é o calor especifico do ar, u* é a da velocidade de fricção de cada pixel , Ts é a temperatura da superfície (K), g é o módulo do campo gravitacional terrestre (9,81 m/s) e H é o fluxo de calor sensível (w/m²) , obtido inicialmente considerando a condição de neutralidade.Os valores de L definem as condições de estabilidade da seguinte forma: se L < 0, a atmosfera é considerada instável; se L > 0, a atmosfera é considerada estável e se L = 0 a atmosfera é considerada neutra.



1) Se L<0 (condição de instabilidade):

π0,5)arctg(x22

x1ln

2

x1ln2ψ m)001(

2m)001(m)001(

m)001(m

2

x1ln2ψ

2m)(2

m)(2h

2

x1ln2ψ

2m)(0,1

m)(0,1h

0,25

m)001( L

001161x

0,25

m)(2 L

2161x

0,25

m)(0,1 L

0,1161x

Onde:



2) Se L>0 (condição de estabilidade):

L

1005ψ m)001m(

L

25ψ m)m(2

L

0,15ψ m)m(0,1

3) Se L=0 (condição de neutralidade):

0ψm 0ψh

.

e



Velocidade de fricção corrigida

m)001(mm0

001*

ψz

001ln

k.uu

Onde u100 é a velocidade do vento a 100 m , k é a constante de von Karman (0,41), Zom é o coeficiente de rugosidade de cada pixel (m) e é a correção da estabilidade para o transporte de momentum a 100m.

m)001(mψ


Fluxo de calor sensível Resistência aerodinâmica ao transporte de calor corrigida

k.u

ψψz

zln

r*

)h(z)h(z1

2

ah

12

onde Z1= 2,0 m, Z1 = 0,1 m, e e são as correções de estabilidade para o transporte de calor a 2,0 m e 0,1 m, respectivamente.

Após obtidos os valores desses parâmetros, retorna-se ao cômputo da função da diferença de temperatura, repetindo-se os cálculos mencionados anteriormente até que se obtenha a convergência dos valores da resistência aerodinâmica (rah ). Para isto, são em geral necessárias de cinco a sete iterações para se atingir a convergência do processo.

)h(z2ψ )h(z1

ψ

Fluxo de calor Latente

Fluxo de calor Latente

O fluxo de calor latente LE (w/m²) corresponde ao fluxo de massa (água sob a forma de vapor) que deixa a superfície, através dos processos evaporação mais transpiração vegetal, e representa uma grande quantidade de energia invisível, daí o termo latente. Sua determinação é obtida como resíduo da equação do Balanço de Energia, qual seja:

LE = Rn – G – H

Onde os fluxos correspondem ao instante da passagem do satélite.

Evapotranspiração Diária

Fração da evapotranspiração (EF instantânea )

GRn

LEFET

ainstantâne

hFETFET 24ainstantâne

hhh RnFETLE 242424

Fluxo de calor Latente (LE24h)

Saldo de radiação (Rn24h), da acordo com De Bruin, 1987.

hwhSh RRn 242424 110)1(

(Brutsaert and Sugita, 1992; Crago, 1996),


L

LEET horar .3600,

Onde L = 2,45 x 106 J kg-1 3600, é um fator de conversão do valor instantâneo para valor horário, conforme equação (Allen et al., 2002a, Trezza, 2002):

Evapotranspiração real horária (ETr, hora)

Fração da evapotranspiração de referência horária (FET0,hora )

horao

horarhora ET

ETFET

,

,_0

Onde ET0,hora é a evapotranspiração de referência horária.


De acordo com Trezza (2002) é relativamente constante em todo o período diurno.

ho

hr

horao

horarhhorao ET

ET

ET

ETFETFET

24,

24,

,

,24,0,

Dessa forma, a evapotranspiração real diária é calculada pela seguinte equação:

24,024,0 .ETFETETr

• Área de estudo

Aplicações

Está área caracteriza-se por apresentar vegetação mais susceptível às mudanças climáticas, respondendo rapidamente à presença ou falta de chuvas. O clima da região de estudo segundo a classificação de Köppen, é do tipo Aw’- Tropical quente e úmido.

• Duas datas de passagem do satélite LANDSAT 5 sensor TM, final do período chuvoso e período seco, capturadas nos dias 29 de agosto de 2008 e 01 de novembro de 2008, e a essas imagens foi aplicado o algoritmo SEBAL.

• Os produtos MODIS, MOD11 A2 e MOD13 Q1, para todo o ano de 2008.

Matérias

Matérias• Elevação - SRTM

Elevação (em metros)

PontosP1 P2 P3 P4 P5 P6 P7 P8 P9 P10 P11 P12

709 350 477 284 319 327 420 310 350 288 358 228

Matérias• Pluviometria para o ano de 2008.

Posto Pluviométrico

Jan Fev Mar Abr Maio Jun Jul Ago Set Out Nov Dez Total

Engenheiro Ávido

16,8 176,9 462,0 167,5 96,8 13,0 0,0 3,2 0,0 0,0 0,0 0,0 936,2

Lagoa de Arroz 275,3 242,4 481,9 295,8 157,7 19,2 28,6 1,1 0,0 0,0 0,0 89,6 1591,6Antenor Navarro

124,2 71,6 677,1 383,6 137,3 12,0 29,8 6,8 0,0 0,0 0,0 71,4 1513,8

São Gonçalo 249,8 131,5 714,6 326,8 327,3 24,8 19,8 5,5 0,0 5,6 0,0 31,7 1837,4Aparecida 214,0 258,8 529,4 191,7 133,2 37,5 54,4 6,6 0,0 0,0 0,0 0,0 1425,6Bernadino

Batista87,2 149,4 261,6 323,6 95,8 51,0 61,4 11,6 0,0 0,0 0,0 0,0 1041,6

Bom Jesus 241,4 84,1 590,2 348,0 132,1 37,7 19,6 0,0 0,0 0,0 0,0 43,4 1496,5Cachoeira dos

Índios289,0 245,2 445,0 313,6 187,3 3,4 20,0 0,2 0,0 0,0 0,0 0,0 1503,7

Cajazeiras 210,9 205,1 590,4 411,4 252,4 38,9 42,0 351,9 0,0 0,0 0,0 53,4 2156,4Poço José de

Moura86,7 49,7 378,7 227,6 67,2 0,0 48,0 9,7 0,0 0,0 0,0 63,5 931,1

Santa Cruz 52,1 62,3 482,9 315,8 181,4 45,2 59,6 29,4 0,0 0,0 0,0 4,3 1233,0Santa Helena 259,0 200,3 572,5 288,8 159,8 5,8 22,9 2,6 0,0 0,0 0,0 50,8 1562,5

Santarém 161,2 49,3 301,8 267,4 112,5 27,0 49,8 12,7 0,0 0,0 0,0 35,8 1017,5São Francisco 79,2 69,6 384,7 261,6 94,9 15,1 55,8 18,4 1,6 0,0 0,0 5,1 986,0

Souza 132,2 162,7 516,6 229,2 282,6 28,2 37,5 3,0 0,0 0,0 0,0 7,0 1399,0Triunfo 166,6 111,5 559,4 293,3 96,0 0,0 36,3 8,2 7,2 0,0 0,0 20,6 1299,1

Figura 5 – Mapa de NDVI (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008

Resultados

Figura 6 – Mapa de temperatura da superfície (A) em 29 de agosto de 2008 e (B) em 01 de Novembro de 2008

Resultados

Resultados

Figura 8 – Mapa de albedo da superfície (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008.

Figura 9 – Mapa de saldo d radiação (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008.

Resultados

Figura 10 – Mapa de fluxo de calor no solo (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008

Resultados

Figura 11 – Mapa de fluxo de calor sensível (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008.

Resultados

Figura 12 – Mapa de fluxo de calor Latente (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008.

Resultados

Figura 12 – Mapa de evapotranspiração real (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008.

Resultados

Figura 5 – Mapa de NDVI (A) em 29 de agosto de 2008 e (B) em 01 de novembro de 2008

Resultados

Figura 4 – Gráfico do NDVI obtido a partir do MODIS e precipitação no tempo, em (A) ponto 1, em (B) ponto 2, em (C) ponto 3 e em (D) ponto 4

Resultados

Figura 6 – Gráfico da Temperatura da superfície obtido a partir do MODIS e precipitação no tempo, em (A) ponto 1, em (B) ponto 2, em (C) ponto 3 e em (D) ponto 4.

Resultados

• Imagens de satélite com sensores de média resolução.

• Imagens de satélite com alta resolução temporal.

• Investigação de outros algoritmos que possibilite a obtenção da evapotranspiração.

• Implementação dos algoritmos em Software Livres.

• Validações do resultados obtidos pela técnicas de processamento digital de imagens.

Comentários finais

Obrigado a todos!

Perguntas?

Documents

John Cunha Universidade Federal de Campina Grande Paraíba, Brasil