GEE随记（一）：利用NDVI和NDWI阈值提取区域裸土用地类型

1. 首先根据数据源、时间、掩膜、云量等选取适合的影像

// 裁剪和加载矢量
var roi = ee.FeatureCollection("projects/ee-****/assets/zhuhai_boundary");
Map.centerObject(roi,13);

var empty =ee.Image().toByte();
var outline = empty.paint({
  featureCollection:roi,
  color:0,
  width:3
});
Map.addLayer(outline, {palette: "ff0000"}, "outline");

//设置去云影像函数(TOA)
//包括云周边、卷云、云、云阴影和雪
function maskL8toa(image){
  var cloudShadowBitMask = 1 << 3;
  var cloudBitMask = 1<<4;
  //获取pixel QA band，这里需要注意，有些产品是‘piexl_qa’
  var qa = image.select('QA_PIXEL');
  // 明确条件，设置两个值都为0
  var mask = qa.bitwiseAnd(cloudShadowBitMask).eq(0)
               .and(qa.bitwiseAnd(cloudBitMask).eq(0));
  return image.updateMask(mask);
}
//加载TOA影像
var l8 = ee.ImageCollection('LANDSAT/LC08/C02/T1_TOA')
           .filterBounds(roi)
           .filterDate('2018-01-01','2018-12-31')
           .map(maskL8toa)
           //.map(applyScaleFactors)
           .median().clip(roi);
print(l8);

2. 构建NDVI和NDWI指数

//构建遥感指数
//NDVI
function NDVI(img){
  var ndvi = img.normalizedDifference(["B5","B4"]);
  return ndvi}
var visndvi = {
 min: 0,
 max: 1,
 palette: 'FFFFFF, CE7E45, DF923D, F1B555, FCD163, 99B718, 74A901, 66A000, 529400,' +
   '3E8601, 207401, 056201, 004C00, 023B01, 012E01, 011D01, 011301'
};
var ndvi =NDVI(l8);
//Map.addLayer(ndvi, visndvi, "NDVI");

//NDWI
function NDWI(img){
  var ndwi = img.normalizedDifference(['B3','B5']);
  return ndwi}
var visndwi = {
  min:-1,
  max:1,
  palette:['00FFFF','0000FF']
};
var ndwi = NDWI(l8);
//Map.addLayer(ndwi,visndwi,"NDWI");

3. 阈值提取裸地

NDVI和NDWI的计算结果范围均在-1到1之间，其中在NDVI，<0.1一般为岩石、沙等，0.2~0.3为灌木和草甸，0.6~0.8为森林等；而在NDWI中>0.2为水体，0~0.2为一般洪水，-0.3~0为中度干旱区、非水表面，-1~0.3为干旱、非水性表面。利用二者之间关系构建阈值提取。

//阈值提取裸地
// 通过两个指数的阈值设定来监测裸地范围
// 该阈值旨在检测裸地。
var bare1 = ndvi.lt(0.1).and(ndwi.lt(0));
// 定义光谱索引的可视化参数。
var ndviViz = {min: -1, max: 1, palette: ['FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', 
                       '99B718','74A901', '66A000', '529400', '3E8601', '207401', '056201',
                       '004C00', '023B01', '012E01', '011D01', '011301']};

var ndwiViz = {min: 0.2, max: 1, palette: ['00FFFF', '0000FF']};
// 掩膜和镶嵌可视化图像
var mosaic = ee.ImageCollection([
  // NDWI > 0.2 设置成蓝色为水域
  ndwi.updateMask(ndwi.gte(0.2)).visualize(ndwiViz),
  // NDVI > 0.2 设置绿色为植被
  ndvi.updateMask(ndvi.gte(0.2)).visualize(ndviViz),
  //将带有阴影（bare2 但不是bare1）的裸露区域可视化为灰色。
  //bare2.updateMask(bare2.and(bare1.not())).visualize({palette: ['AAAAAA']}),
  // 将其他裸露区域可视化为白色。
  bare1.updateMask(bare1).visualize({palette: ['FFFFFF']}),
]).mosaic();
Map.addLayer(mosaic, {}, 'Visualization mosaic');

4. 结果呈现

由于较为笼统的划分，可以看到城镇建设用地等也是白色显示，可以根据类似方法结合不同指数进行更进一步的裸土或者不同水体的阈值提取。