rgee尝鲜Landsat9，以及Landsat8对比
看到吴秋生老师发布的推文，Landsat 9数据已经加入Google Earth Engine，我也迫不及待想要试一试。
试验内容
- 使用
rgee查看Landsat9数据 - 交互式地图划定感兴趣区ROI
- Landsat8和Landsat9数据对比
使用rgee查看Landsat9数据
这部分代码基本上是根据吴秋生老师的推文和我以前写的Landsat的文章中的代码改的：Landsat Collection2 Surface Reflectance产品介绍、数据获取、去云
主要实现以下功能：
- 引用LANDSAT9数据，按云量检索排序
- 对光学波段进行校正（需要校正热红外波段的同学请参考吴老师推文）
- 影像合成
median()，在这我和吴老师的一样，都是取得中值 - 标准假彩色合成可视化
具体代码如下：
library(rgee)
library(mapedit)
library(sf)
ee_Initialize(drive = T)
LT9 = ee$ImageCollection('LANDSAT/LC09/C02/T1_L2')$
sort('CLOUD_COVER')
#应用校正
applyScaleFactor = function(image){
opticalBands = image$select('SR_B.')$multiply(0.0000275)$add(-0.2)
return(opticalBands)
}
LT9mosaic = LT9$select(c('SR_B5', 'SR_B4', 'SR_B3', 'SR_B2'))$
map(applyScaleFactor)$median()
imageVisParam <- list(bands = c('SR_B5', 'SR_B4', 'SR_B3'), min = 0, max = 0.3)
Map$addLayer(LT9mosaic, imageVisParam, 'LT9mosaic')
交互式地图划定感兴趣区ROI
我根据参考文献2的指导使用editMap()进行了交互式制图，划定感兴趣区，这个函数来自mapedit包，需要支持的包比较多，记得安装。
交互式绘图，拖动绘制矩形，绘制完成点击Done：
代码如下：
Map$addLayer(LT9mosaic, imageVisParam, 'LT9mosaic') %>%
editMap() -> my_roi
ee_roi = my_roi$drawn %>% #交互式绘图ROI转gee矢量
sf_as_ee()
LT9_ROI = LT9mosaic$clip(ee_roi)
Map$centerObject(ee_roi)
Map$addLayer(LT9_ROI, imageVisParam, 'LT9_ROI')
Landsat8和Landsat9数据对比
最后我选取了2021年10月至2022年2月的Landsat8影像和Landsat9进行一个对比。
对比一下Landsat9和8的预览图，差别很小，看起来Landsat9的亮度更高一些。
对比代码如下：
LT8 = ee$ImageCollection('LANDSAT/LC08/C02/T1_L2')$
filterDate('2021-10-01', '2022-02-14')$
sort('CLOUD_COVER')
LT8mosaic = LT8$select(c('SR_B5', 'SR_B4', 'SR_B3', 'SR_B2'))$
map(applyScaleFactor)$median()
LT8_ROI= LT8mosaic$clip(ee_roi)
Map$centerObject(ee_roi)
Map$addLayer(LT8_ROI, imageVisParam, 'LT8_ROI')
m1 = Map$addLayer(LT9_ROI, imageVisParam, 'LT9_ROI')
m2 = Map$addLayer(LT8_ROI, imageVisParam, 'LT8_ROI')
m1 |m2 #卷帘式预览
m1 + m2 #上下图层预览
参考文献
- GEE（谷歌地球引擎）已经开始收录Landsat 9数据
- https://csaybar.github.io/blog/2020/06/10/rgee_01_worldmap/
- https://r-earthengine.github.io/image_02/