GEE第一段程序：沙洲面积提取
打开GEE，先选定一块区域（使用ROI在地图上选定）可以看到此时代码框中出现了Imports (1 entry) ，我们可以对这个变量var geometry改名，例如这里我们将之改为roi,执行下边的代码：
Map.centerObject(roi, 13);
var Area_roi=roi.area().divide(1000000);
print(Area_roi);
var bounds= ee.Image().toByte()
.paint({
featureCollection:
ee.FeatureCollection([ee.Feature(roi)]),
color: null,
width: 1
});
Map.addLayer(bounds, {palette: "red"}, "bounds");
// 选择遥感影像的时间段（startDate 包含输入时间，endDate 不包含输入时
//间。下图中的遥感影像选择时段实际上是 2018 年 12 月 28 日-2019 年 5 月 26
//日）
var startDate = ee.Date("2018-12-28");
var endDate =ee.Date("2019-5-27");
// 引入所有遥感影像源（这些遥感影像已经经过了预处理，可以直接利用）
var l5col=ee.ImageCollection("LANDSAT/LT05/C01/T1_SR");
var l7col= ee.ImageCollection("LANDSAT/LE07/C01/T1_SR");
var l4col=ee.ImageCollection("LANDSAT/LT04/C01/T1_SR");
var l8col =ee.ImageCollection("LANDSAT/LC08/C01/T1_SR");
//选择某一遥感影像源，进行影像的区域，云量筛选，并进行去云处理,计算
// MNDWI(即影像二值化).
var Col=ee.ImageCollection(l8col)
.filterBounds(roi)
.filterDate(startDate, endDate)
.filter(ee.Filter.lte("CLOUD_COVER", 20))
.map(
function(image) {
var cloudShadowBitMask = 1 << 3;
var cloudsBitMask = 1 << 5;
var snowBitMask = 1 << 4;
var qa = image.select('pixel_qa');
var mask = qa.bitwiseAnd(cloudShadowBitMask).eq(0)
.and(qa.bitwiseAnd(cloudsBitMask).eq(0))
.and(qa.bitwiseAnd(snowBitMask).eq(0));
return image.updateMask(mask);
}
)
.map(function(image) {
var time_start = image.get("system:time_start");
image = image.multiply(0.0001);
image = image.set("system:time_start", time_start);
return image;
})
.map(function(image) {
return image.addBands(image.normalizedDifference(["B3", "B6"])
.rename("MNDWI"));
})
.select("MNDWI");
print("Col", Col);
//大津算法（otsu）选择最优阈值方法，并进行影像的分割，水体和陆地分
// 离，提取陆地面积。
function otsu(histogram) {
var counts = ee.Array(ee.Dictionary(histogram).get('histogram'));
var means = ee.Array(ee.Dictionary(histogram).get('bucketMeans'));
var size = means.length().get([0]);
var total = counts.reduce(ee.Reducer.sum(), [0]).get([0]);
var sum = means.multiply(counts).reduce(ee.Reducer.sum(),
[0]).get([0]);
var mean = sum.divide(total);
var indices = ee.List.sequence(1, size);
var bss = indices.map(function(i) {
var aCounts = counts.slice(0, 0, i);
var aCount = aCounts.reduce(ee.Reducer.sum(), [0]).get([0]);
var aMeans = means.slice(0, 0, i);
var aMean = aMeans.multiply(aCounts)
.reduce(ee.Reducer.sum(), [0]).get([0])
.divide(aCount);
var bCount = total.subtract(aCount);
var bMean = sum.subtract(aCount.multiply(aMean)).divide(bCount);
return aCount.multiply(aMean.subtract(mean).pow(2)).add(
bCount.multiply(bMean.subtract(mean).pow(2)));
});
return means.sort(bss).get([-1]);
}
var sCol1 = Col.map(function(image) {
var histogram = image.reduceRegion({
reducer: ee.Reducer.histogram(),
geometry: roi,
scale: 2,
maxPixels: 1e13,
tileScale: 10
});
var threshold = otsu(histogram.get("MNDWI"));
var mask = image.lte(threshold);
var water = mask.updateMask(mask).rename("water");
var newWater = water.addBands(image);
//转矢量
var vectors = newWater.reduceToVectors({
geometry: roi,
scale: 2,
geometryType: 'polygon',
reducer: ee.Reducer.mean(),
maxPixels: 1e13,
tileScale: 10
});
vectors = vectors.filterBounds(roi);
water = water.clip(vectors);
//计算陆地面积
var areaImage = water.multiply(ee.Image.pixelArea());
//统计指定区域中所有陆地像素和，也就是沙体的面积
var dict = areaImage.reduceRegion({
reducer: ee.Reducer.sum(),
geometry: vectors.geometry(),
scale: 2,
maxPixels: 1e13,
tileScale: 10
});
water = water.set("area",
ee.Number(dict.get("water")).divide(1000000));
water = water.set("threshold", threshold);
return water.toByte();
});
print("sCol1", sCol1);
Map.addLayer(sCol1.first(), {palette: "red"}, "BAR");
// 整合所有影像的时间和计算出的沙洲面积并输出
Export.table.toDrive({collection: sCol1,
description: "guanzhou",
fileNamePrefix: "guanzhou",
selectors:["system:index","area"]
});