回复
#冲刺创作新星#PIE-Engine APP：计算水体面积提取 原创
此星光明
发布于 2022-9-26 10:26
浏览
0收藏
本次我们来查看进行水体处理的分布，我们这里首先对数据进行预处理，先进行NDWI，AWEI、MNDWI等计算函数和去云函数，第二部分市机器学习部分这里有三个机器学习模型，分别是随机森林、贝叶斯和支持向量机
NDWI=（Band2-Band4）/（Band2+Band4），式中Band2表示绿光波段的反射率，Band4表示近红外波段的反射率。该方法尽管已经较为古老，但其是最为常用的水体提取方法（部分高分辨率数据仅有4个波段），并且目前很多的水体指数法都是在该方法地基础上进行地变化。该方法对于大部分的常规水体均可有效提取，但是同样受到其他因素的影响较大。
MNDWI=（Band2-Band5）/（Band2+Band5）式中Band2表示绿光波段的反射率，Band5表示中红外波段的反射率。该方法是NDWI的变种，但是总体提取精度要优于NDWI。但是对于部分高分辨率数据不适用（不存在中红外波段），因此在进行中低分辨率的水体自动提取时常用。且由于该方法对于部分湖泊湿地中的富营养化的水体提取也较为适用，因此也常用来提取此类水体。
randomColumn(columnName,seed,distribution)
向FeatureCollection中添加一列确定性伪随机数。
方法参数：
- featureCollection(FeatureCollection)
FeatureCollection实例
- columnName(String)
新增列的名称，默认为 ‘random’
- seed(Long)
随机种子，默认为0
- distribution(String)
生成随机数的分布类型。赋值为’uniform’ 、'normal’之一
返回值：FeatureCollection
confusionMatrix()
计算监督分类分类器结果的混淆矩阵
方法参数：
- Classifier(Classifier)
监督分类分类器实例
返回值：ConfusionMatrix
代码：
* @Name : 基于PIE-Engine的水体频率变化长时序遥感监测自动计算平台
* @Time : 2021/06/30
* @Author : 中国地质大学（武汉）水体频率小组
* @Desc : 基于水体指数或监督分类方法的水体频率计算
* @Source : 航天宏图第四届 “航天宏图杯”PIE软件二次开发大赛云开发组二等奖获奖作品
*/
//设定预定的的变量备用
var layerKey = null;
var roiKey = null;
var selectDates = ["2016-1-1", "2016-12-31"];
var selectCode = "";
var selectTag = "NDWI";
var selectcity = "wuhan";
var selectyear = "2018";
var selectThreshold = "";
var k = 0,
y = 0;
var selectYZ = "随机序列作为验证样本"
//这里是将图层
var layerTF = [false, false, false, false, false, false];
////////////////////////////////////////////////对数据的预处理部分/////////////////////////////////////////////////
//获得感兴趣的研究区域
function getROI(cityname) {
var ChinaCity = pie.FeatureCollection('user/pieadmin/ChinaCity');
var roi = ChinaCity.filter(pie.Filter.eq('市', cityname)).first().geometry();
Map.centerObject(roi, 8);
roiKey = Map.addLayer(roi, { color: "#ff0000", fillColor: "#00000000" }, "roi");
return roi;
}
//计算水体指数
function NDWI(image) {
var ndwi = image.normalizedDifference(['B3', 'B5']).rename('NDWI');
return ndwi.gt(pie.Number(selectThreshold));
};
function AWEI(image) {
var awei = image.select(["B2", "B3", "B5", "B6", "B7"]).expression(
'B2+2.5*B3-1.5*(B5+B6)-0.25*B7', {
B2: image.select("B2"),
B3: image.select("B3"),
B5: image.select("B5"),
B6: image.select("B6"),
B7: image.select("B7"),
}).rename('AWEI');
return awei.gt(pie.Number(selectThreshold));
};
function MNDWI(image) {
var mndwi = image.normalizedDifference(['B3', 'B6']).rename('mNDWI');
return mndwi.gt(pie.Number(selectThreshold));
}
//训练样本波段范围0-5000，LC08/02/SR数据集范围0-50000，除以10处理
function divide10(image) {
return imgd10 = image.divide(10);
}
//SR去云处理
function maskL8sr(image) {
var qa = image.select('QA_PIXEL');
var mask = qa.bitwiseAnd(1 << 3).eq(0)
.and(qa.bitwiseAnd(1 << 4).eq(0))
.and(qa.bitwiseAnd(1 << 5).eq(0));
return image.updateMask(mask);
}
///////////////////////////////////////////////机器学习分类水体/////////////////////////////////////////
function Machinelearning(images, samcity, k) {
// 添加训练样本
var TrainingPoints = pie.FeatureCollection('user/pieadmin/' + samcity + "ALL");
// 预测使用的波段
var bands = ['B2', 'B3', 'B4', 'B5', 'B6', 'B7', 'NDVI', 'mNDWI', 'AWEI'];
// 分类标签
var label = 'waterclass';
//选择训练训练器
switch (k) {
case 0:
var classifer = pie.Classifier.rTrees().train(TrainingPoints, label, bands);
break;
case 1:
var classifer = pie.Classifier.normalBayes().train(TrainingPoints, label, bands);
break;
case 2:
var classifer = pie.Classifier.svm().train(TrainingPoints, label, bands);
break;
}
function water_index(img) {
var image = img.select(["B2", "B3", "B4", "B5", "B6", "B7"]);
var ndvi = image.normalizedDifference(['B5', 'B4']).rename('NDVI');
var mndwi = image.normalizedDifference(['B3', 'B6']).rename('mNDWI');
var awei = image.select(["B2", "B3", "B5", "B6", "B7"]).expression(
'B2+2.5*B3-1.5*(B5+B6)-0.25*B7', {
B2: image.select("B2"),
B3: image.select("B3"),
B5: image.select("B5"),
B6: image.select("B6"),
B7: image.select("B7"),
}).rename('AWEI');
return img.addBands(ndvi).addBands(mndwi).addBands(awei);
}
var image = images.map(water_index);
var resultImage = image.map(function(image) {
var Rfiamge = image.select(bands).classify(classifer);
return Rfiamge;
});
//随机序列的生成
var sampleFeatureCollection = TrainingPoints.randomColumn('random');
if (y == 0) var sampleTestingFeatures = sampleFeatureCollection.filter(pie.Filter.gt("random", 0.7));
if (y == 1) var sampleTestingFeatures = pie.FeatureCollection('user/pieadmin/' + samcity + "ALL");
var checkM = classifer.confusionMatrix();
print("训练矩阵-ACC系数：", checkM.acc(), "训练矩阵-Kappa系数：", checkM.kappa());
var predictResult = sampleTestingFeatures.classify(classifer, "classification");
var errorM = predictResult.errorMatrix("waterclass", "classification")
print("验证矩阵-ACC系数：", errorM.acc(), "验证矩阵-Kappa系数：", errorM.kappa());
return resultImage;
}
//计算有效像元
function validPixel(image) {
return image.select('B2').gte(0);
};
//计算水体频率并分类
function FrequencyC(roi, startDate, endDate, way, samcity) {
l8_images = pie.ImageCollection("LC08/02/SR")
.filterDate(startDate, endDate)
.filterBounds(roi)
.filter(pie.Filter.lt('cloud_cover', 30))
.select(['B2', 'B3', 'B4', 'B5', 'B6', 'B7', 'QA_PIXEL'])
.map(maskL8sr);
print("影像数量为：");
print(l8_images.size());
//计算有效像元个数
var pixel_validNumber = l8_images.map(validPixel).sum().clip(roi);
//分方法计算水体个数和水体频率
var water_validNumber;
if (way == "NDWI") {
water_validNumber = l8_images.map(NDWI).sum().clip(roi);
var waterFrequency = water_validNumber.divide(pixel_validNumber).rename('frequency');
}
if (way == "AWEI") {
water_validNumber = l8_images.map(AWEI).sum().clip(roi);
var waterFrequency = water_validNumber.divide(pixel_validNumber).rename('frequency');
}
if (way == "MNDWI") {
water_validNumber = l8_images.map(MNDWI).sum().clip(roi);
var waterFrequency = water_validNumber.divide(pixel_validNumber).rename('frequency');
}
if (way == "随机森林") {
water_validNumber = Machinelearning(l8_images.map(divide10), samcity, 0).sum().clip(roi);
var waterFrequency = water_validNumber.divide(pixel_validNumber).rename('frequency');
}
if (way == "正态贝叶斯") {
water_validNumber = Machinelearning(l8_images.map(divide10), samcity, 1).sum().clip(roi);
var waterFrequency = water_validNumber.divide(pixel_validNumber).rename('frequency');
}
if (way == "支持向量机") {
water_validNumber = Machinelearning(l8_images.map(divide10), samcity, 2).sum().clip(roi);
var waterFrequency = water_validNumber.divide(pixel_validNumber).rename('frequency');
}
var vis2 = {
min: 0,
max: 1,
palette: ['ffffff', 'f8f8ff', '87cefa', '1e90ff', '4169e1', '000080']
};
var data = {
title: "水体频率",
colors: ["#ffffff", "#f0f8ff", "#add8e6", "#87ceeb", "#6495ed", "#4169e1", "#0000ff", "#0000cd", "#000080"],
labels: ["0", "0.25", "0.5", "0.75"],
step: 2
};
var style = {
right: "150px",
bottom: "10px",
height: "70px",
width: "350px"
};
var legend = ui.Legend(data, style);
Map.addUI(legend);
var PermanentWater = waterFrequency.gte(0.75).rename("waterclass");
var mask1 = waterFrequency.gte(0.25);
var mask2 = waterFrequency.lt(0.75);
var SensonalWater = pie.ImageCollection.fromImages([mask1, mask2]).sum().eq(2).rename('waterclass');
var Land = waterFrequency.lt(0.25).rename("waterclass");
var MaskPS = waterFrequency.gte(0.25);
Map.addLayer(waterFrequency.updateMask(MaskPS), vis2, '水体频率图层');
//显示各水体类别的图层
Map.addLayer(PermanentWater.setMaskValue(0), { uniqueValue: '1', palette: ['000080'] }, 'permanent waterbody', layerTF[3]);
Map.addLayer(SensonalWater.setMaskValue(0), { uniqueValue: '1', palette: ['FF0000'] }, 'seasonal waterbody', layerTF[4]);
Map.addLayer(Land.setMaskValue(0), { uniqueValue: '1', palette: ['FFFF00'] }, 'land', layerTF[5]);
Map.addLayer(l8_images.mosaic().clip(roi), { min: 0, max: 30000, bands: ["B4", "B3", "B2"] }, 'B4/B3/B2', layerTF[0]);
Map.addLayer(l8_images.mosaic().clip(roi), { min: 0, max: 30000, bands: ["B5", "B4", "B3"] }, 'B5/B4/B3', layerTF[1]);
Map.addLayer(l8_images.mosaic().clip(roi), { min: 0, max: 30000, bands: ["B5", "B6", "B4"] }, 'B5/B6/B4', layerTF[2]); //432,543,564
return waterFrequency;
}
///////////////////////////////////////////////UI设计部分/////////////////////////////////////////
var label1 = ui.Label("基于PIE-engine的水体频率变化长时序遥感监测自动计算平台", { "font-size": "20px" })
var label2 = ui.Label("一、区域地表水水体频率计算：", { "font-size": "18px" })
var label3 = ui.Label("请选择研究的时间范围，选择起始时间和终止时间(2013.07——现在数据可用)：", { "font-size": "15px" });
var label4 = ui.Label("请选择水体提取方法，目前包括指数法和监督分类法：", { "font-size": "15px" });
//选择研究区模块
var textBox1 = ui.TextBox({
placeholder: "请输入行政市名字(XX市)",
value: selectCode,
onChange: function(value) {
selectCode = value;
},
disabled: false
})
var textboxName1 = ui.Label("请设置研究区范围，区域类别为市级行政区划范围，如武汉市：", { "font-size": "15px" });
var textboxPanel1 = ui.Panel({
widgets: [textboxName1, textBox1],
layout: ui.Layout.flow("vertical")
});
//选择时间范围模块
var selectdateName = ui.Label("日期范围：", { "font-size": "15px" });
var dateSelect = ui.DateSelect({
type: "daterange",
placeholder: "请输入数值",
value: selectDates,
onChange: function(value) { selectDates[0] = value[0];
selectDates[1] = value[1]; },
disabled: false,
});
var selectdatePanel = ui.Panel({
widgets: [selectdateName, dateSelect],
layout: ui.Layout.flow("horizontal")
});
//选择提取水体方法模块
var select1 = ui.Select({
items: ['AWEI', 'NDWI', "MNDWI", "随机森林", "正态贝叶斯", "支持向量机"],
placeholder: "请选择",
value: selectTag,
multiple: false,
onChange: function(value) {
selectTag = value;
}
})
var selectName = ui.Label("水体提取方法：", { "font-size": "15px" });
var selectPanel = ui.Panel({
widgets: [selectName, select1],
layout: ui.Layout.flow("horizontal")
});
//下一步按钮
var btn1 = ui.Button({
label: "下一步",
type: "success",
onClick: clickBtn1,
style: { left: "250px" }
});
//第一层界面
var panel = ui.Panel({
widgets: [
label1, label2,
textboxPanel1,
label3,
selectdatePanel,
label4,
selectPanel,
btn1
],
style: {
width: "350px",
backgroundColor: "#fff"
}
});
var panel1;
function clickBtn1() {
print("选择的参数是：", "市名：" + selectCode, "时间范围：" + selectDates[0] + "_" + selectDates[1], "计算方法:" + selectTag);
var label3 = ui.Label("选择的提取方法为指数法");
var label4 = ui.Label("选择的提取方法为监督分类法");
//选择阈值模块
var textBox2 = ui.TextBox({
placeholder: "阈值(推荐值:0-0.2)",
value: selectThreshold,
onChange: function(value) {
selectThreshold = value;
},
disabled: false
})
var textboxName2 = ui.Label("请输入指数法计算阈值：", { "font-size": "13px" });
var textboxPanel2 = ui.Panel({
widgets: [textboxName2, textBox2],
layout: ui.Layout.flow("horizontal")
});
//计算水体频率按钮
var btn2 = ui.Button({
label: "计算水体频率",
type: "success",
onClick: clickBtn2,
style: { left: "50px" }
});
var btn2N = ui.Button({
label: "上一步",
onClick: function() { ui.root.add(panel);
ui.root.remove(panel1); },
style: { left: "25px" }
});
var btnPanel = ui.Panel({
widgets: [btn2N, btn2],
layout: ui.Layout.flow("horizontal")
});
//选择训练样本模块
var selectF = ui.Select({
items: ['wuhan', 'suzhou', "tianjin", "ALL"],
placeholder: "请选择",
value: selectcity,
multiple: false,
onChange: function(value) {
selectcity = value;
}
})
var selectCName = ui.Label("选择训练样本城市名：", { "font-size": "13px" });
var selectFPanel = ui.Panel({
widgets: [selectCName, selectF],
layout: ui.Layout.flow("horizontal")
});
//选择允许显示的图册模块
var ckboxName = ui.Label("请选择需要加载显示的图层：");
var checkbox = ui.Checkbox({
label: ["B4/B3/B2", "B5/B4/B3", "B5/B6/B4", "永久性水体", "季节性水体", "陆地"],
value: [],
disabled: [],
onChange: function(value) {
layerTF[0] = false;
layerTF[1] = false;
layerTF[2] = false;
for (var i = 0; i < value.length; i++) {
if (value[i] == "B4/B3/B2") { layerTF[0] = true; }
if (value[i] == "B5/B4/B3") { layerTF[1] = true; }
if (value[i] == "B5/B6/B4") { layerTF[2] = true; }
if (value[i] == "永久性水体") { layerTF[3] = true; }
if (value[i] == "季节性水体") { layerTF[4] = true; }
if (value[i] == "陆地") { layerTF[5] = true; }
}
}
});
//选择为指数法的界面
if (selectTag == "AWEI" || selectTag == "NDWI" || selectTag == "MNDWI") {
panel1 = ui.Panel({
widgets: [
label3,
textboxPanel2,
ckboxName,
checkbox,
btnPanel
],
style: {
width: "350px",
backgroundColor: "#fff"
}
});
k = 0;
}
//选择为机器学习法的界面
if (selectTag == "随机森林" || selectTag == "正态贝叶斯" || selectTag == "支持向量机") {
panel1 = ui.Panel({
widgets: [
label4,
selectFPanel,
ckboxName,
checkbox,
btnPanel
],
style: {
width: "350px",
backgroundColor: "#fff"
}
});
k = 1;
}
ui.root.add(panel1);
ui.root.remove(panel);
}
function clickBtn2() {
var str = "";
if (k == 0) { print("计算的阈值为：" + selectThreshold);
str = "计算的阈值为：" + selectThreshold; }
if (k == 1) { print("选择的分类器为:" + selectcity);
str = "选择的分类器为:" + selectcity; }
if (selectYZ == "随机序列作为验证样本") y = 0;
if (selectYZ == "训练样本作为验证") y = 1;
var roi = getROI(selectCode);
var waterFrequency = FrequencyC(roi, selectDates[0], selectDates[1], selectTag, selectcity);
var printlabel = ui.Label("用户已选的参数为：", { "font-size": "20px" });
var printlabel1 = ui.Label("市名：" + selectCode);
var printlabel2 = ui.Label("时间范围：" + selectDates[0] + "_" + selectDates[1]);
var printlabel3 = ui.Label("计算方法:" + selectTag);
var printlabel4 = ui.Label(str);
var selectTag1 = "B4/B3/B2";
var label5 = ui.Label("请用户选择需要下载的影像，包括原始影像和水体频率：", { "font-size": "18px" })
function clickBtn3() {
print("导出影像:" + selectTag1);
function exportImageToDrive(image, region) { //导出影像
Export.image({
image: image,
description: "ExportImage",
region: region,
scale: 30
});
}
if (selectTag1 == "B4/B3/B2") {
timages = pie.ImageCollection("LC08/02/SR")
.filterDate(selectDates[0], selectDates[1])
.filterBounds(roi)
.filter(pie.Filter.lt('cloud_cover', 20))
.select(['B4', 'B3', 'B2'])
.map(maskL8sr).mosaic().clip(roi);
exportImageToDrive(timages, roi);
}
if (selectTag1 == "水体频率影像") { exportImageToDrive(waterFrequency, roi); }
}
var select2 = ui.Select({
items: ['B4/B3/B2', '水体频率影像'],
placeholder: "请选择",
value: selectTag1,
multiple: false,
onChange: function(value) { selectTag1 = value; }
});
var btn12 = ui.Button({
label: "确定",
type: "success",
onClick: clickBtn3,
style: { left: "0px" }
});
var selectName2 = ui.Label("输出影像：");
var selectPanel1 = ui.Panel({
widgets: [selectName2, select2, btn12],
layout: ui.Layout.flow("horizontal")
});
var panel2 = ui.Panel({
widgets: [printlabel, printlabel1, printlabel2, printlabel3, printlabel4, label5, selectPanel1],
style: {
width: "350px",
backgroundColor: "#fff"
}
});
ui.root.add(panel2);
ui.root.remove(panel1);
}
ui.root.add(panel);
这个APP的训练集中只有武汉，天津等几个城市的训练集，而且计算时间会非常长，并且结果依旧没有出现。
©著作权归作者所有，如需转载，请注明出处，否则将追究法律责任
已于2022-10-27 15:56:12修改
赞
1
收藏
回复
相关推荐
