AI Earth ——开发者模式案例8：利用Landsat-8数据进行地表温度反演

利用 Landsat-8 数据进行地表温度反演¶

初始化环境¶

import aie

aie.Authenticate()
aie.Initialize()

Landsat-8 数据检索¶

指定区域、时间、云量检索 Landsat-8 ，并对数据进行去云处理。

region = aie.FeatureCollection('China_Province') \
            .filter(aie.Filter.eq('province', '上海市')) \
            .geometry()

dataset = aie.ImageCollection('LANDSAT_LC08_C02_T1_L2') \
             .filterBounds(region) \
             .filterDate('2019-06-01', '2019-08-31') \
             .filter(aie.Filter.lte('eo:cloud_cover', 10.0))

print(dataset.size().getInfo())

image = dataset.median()   
map = aie.Map(
    center=image.getCenter(),
    height=800,
    zoom=7
)
rgb_params = {
    'bands': ['SR_B4', 'SR_B3', 'SR_B2'],
    'min': 8000,
    'max': 13000
}
map.addLayer(
    image,
    rgb_params,
    'raw_img',
    bounds = image.getBounds()
)
map

计算 NDVI¶

ndvi = image.normalizedDifference(['SR_B5', 'SR_B4']).rename(['NDVI'])
ndvi_params = {
    'min': -1.0,
    'max': 1.0,
    'palette': [
        '#FFFFFF', '#CE7E45', '#DF923D', '#F1B555', '#FCD163', '#99B718', '#74A901',
        '#66A000', '#529400', '#3E8601', '#207401', '#056201', '#004C00', '#023B01',
        '#012E01', '#011D01', '#011301'
    ]
}

map.addLayer(
    ndvi,
    ndvi_params,
    'NDVI',
    bounds = image.getBounds()
)
map

计算 Fractional Vegetation¶

min = aie.Image(-0.38)
max = aie.Image(0.68)
fv = (ndvi.subtract(min).divide(max.subtract(min))).pow(aie.Image(2)).rename(['FV'])

fv_params = {
    'min': 0,
    'max': 1
}
map.addLayer(
    fv,
    fv_params,
    'fv',
    bounds = image.getBounds()
)
map

计算 Thermal¶

th = image.select('ST_B10').multiply(aie.Image(0.00341802)).add(aie.Image(149)).rename(['TH'])

th_params = {
    'min': 260,
    'max': 320,
    'palette': [
        '#0000FF', '#FFFFFF', '#008000'
    ]
}
map.addLayer(
    th,
    th_params,
    'thermal',
    bounds = image.getBounds()
)
map

计算比辐射率 Emissivity¶

a = aie.Image(0.004)
b = aie.Image(0.986)
em = fv.multiply(a).add(b).rename(['EM'])

em_params = {
    'min': 0.9865619146722164,
    'max': 0.989699971371314
}
map.addLayer(
    em,
    em_params,
    'emissivity',
    bounds = image.getBounds()
)
map

计算地表温度 ( LST )¶

tb = th.select(['TH'])
eb = em.select(['EM'])
lst = tb.divide(tb.multiply(aie.Image(0.00115)).divide(aie.Image(1.4388)).multiply(eb.log()).add(aie.Image(1))).rename(['LST'])

lst_params = {
    'min': 291,
    'max': 330,
    'palette': ['#040274', '#040281', '#0502a3', '#0502b8', '#0502ce', '#0502e6',
                '#0602ff', '#235cb1', '#307ef3', '#269db1', '#30c8e2', '#32d3ef',
                '#3be285', '#3ff38f', '#86e26f', '#3ae237', '#b5e22e', '#d6e21f',
                '#fff705', '#ffd611', '#ffb613', '#ff8b13', '#ff6e08', '#ff500d',
                '#ff0000', '#de0101', '#c21301', '#a71001', '#911003']
}

map.addLayer(
    lst,
    lst_params,
    'lst',
    bounds = image.getBounds()
)
map