geemap学习笔记022：如何找出一年中最绿的一天？

本文链接：https://blog.csdn.net/qq_41951977/article/details/134775664

前言

这虽然是一个问题，但是解决这个问题之后，就会学习到很多的内容。包括如何计算NDVI、如何进行镶嵌、如何获取影像的时间等等。

1 导入库并显示地图

import ee
import geemap

ee.Initialize()
Map = geemap.Map()
Map

2 定义一个感兴趣区（ROI）

countries = ee.FeatureCollection('users/giswqs/public/countries') #这是一个公共的数据集
Map.addLayer(countries, {}, 'coutries')

roi = countries.filter(ee.Filter.eq('ISO_A3', 'USA')) #过滤USA矢量区域
Map.addLayer(roi, {}, 'roi')

3 筛选图像集合

start_date = '2019-01-01'
end_date = '2019-12-31'

l8 = (
    ee.ImageCollection('LANDSAT/LC08/C02/T1_TOA')
    .filterBounds(roi)
    .filterDate(start_date, end_date)
)

# print(l8.size().getInfo()) #可以看一下一共有多少景数据，数据量越大时间越长

4 创建一个合成的中值图像

median = l8.median() #计算每个像素所有值中的中值来减少图像集合

visParams = {
    'bands': ['B4', 'B3', 'B2'],
    'min': 0,
    'max': 0.4,
}

Map.addLayer(median, visParams, 'Median')

合成的中值图像效果也并不是很好，后续也没用到，只是作为一个对比吧

5 定义函数添加NDVI、时间波段

def addNDVI(image): #定义NDVI函数
    ndvi = image.normalizedDifference(['B5', 'B4']).rename('NDVI')
    return image.addBands(ndvi)

def addDate(image): #定义获取时间的函数，格式为20191222
    img_date = ee.Date(image.date())
    img_date = ee.Number.parse(img_date.format('YYYYMMdd'))
    return image.addBands(ee.Image(img_date).rename('date').toInt())

def addMonth(image): #定义获取月份的函数，格式为12
    img_date = ee.Date(image.date())
    img_doy = ee.Number.parse(img_date.format('M'))
    return image.addBands(ee.Image(img_doy).rename('month').toInt())

def addDOY(image):  #定义获取这是一年中第几天获取的，格式为356
    img_date = ee.Date(image.date())
    img_doy = ee.Number.parse(img_date.format('D'))
    return image.addBands(ee.Image(img_doy).rename('doy').toInt())

6 利用map执行函数计算

withNDVI = l8.map(addNDVI).map(addDate).map(addMonth).map(addDOY)

7 进行快速镶嵌

greenest = withNDVI.qualityMosaic('NDVI') #使用NDVI波段大小作为每个像素的排序函数，合成集合中的所有图像。
greenest.bandNames().getInfo()

8 展示最大值的波段

ndvi = greenest.select('NDVI') #将合成的NDVI添加到地图中
palette = [
    '#d73027',
    '#f46d43',
    '#fdae61',
    '#fee08b',
    '#d9ef8b',
    '#a6d96a',
    '#66bd63',
    '#1a9850',
]
Map.addLayer(ndvi, {'palette': palette}, 'NDVI')

Map.addLayer(greenest, visParams, 'Greenest pixel')
Map

NDVI值最大的合成影像

9 展示时间波段

Map.addLayer(
    greenest.select('month'),
    {'palette': ['red', 'blue'], 'min': 1, 'max': 12},
    'Greenest month',
) #显示合成结果最好的月份是哪个月份

Map.addLayer(
    greenest.select('doy'),
    {'palette': ['brown', 'green'], 'min': 1, 'max': 365},
    'Greenest doy',
) #显示合成结果最好的一天是第几天