/**** Start of imports. If edited, may not auto-convert in the playground. ****/
var table = ee.FeatureCollection("projects/ee-13727983375/assets/ygao");
/***** End of imports. If edited, may not auto-convert in the playground. *****/
function maskL8sr(image) {
// Bit 0 - Fill
// Bit 1 - Dilated Cloud
// Bit 2 - Cirrus
// Bit 3 - Cloud
// Bit 4 - Cloud Shadow
var qaMask = image.select('QA_PIXEL').bitwiseAnd(parseInt('11111', 2)).eq(0);
var saturationMask = image.select('QA_RADSAT').eq(0);
// Apply the scaling factors to the appropriate bands.
var opticalBands = image.select('SR_B.').multiply(0.0000275).add(-0.2);
var thermalBands = image.select('ST_B.*').multiply(0.00341802).add(149.0);
// Replace the original bands with the scaled ones and apply the masks.
return image.addBands(opticalBands, null, true)
.addBands(thermalBands, null, true)
.updateMask(qaMask)
.updateMask(saturationMask);
}
function SI_cal(img) {
var blue = img.select("SR_B2");
var red = img.select("SR_B4");
var nir = img.select("SR_B5");
var swir1 = img.select("SR_B6");
var swir2 = img.select("SR_B7");
var SI_temp =((swir1.add(red)).subtract(blue.add(nir)))
.divide((swir1.add(red)).add(blue.add(nir)))
//print(SI_temp);
return SI_temp;
}
function IBI_cal(img) {
var green = img.select("SR_B3");
var red = img.select("SR_B4");
var nir = img.select("SR_B5");
var swir1 = img.select("SR_B6");
var IBI_temp =(((swir1.multiply(2.0)).divide(swir1.add(nir)))
.subtract((nir.divide(nir.add(red))).add(green.divide(green.add(swir1)))))
.divide((((swir1.multiply(2.0)).divide(swir1.add(nir)))
.add((nir.divide(nir.add(red))).add(green.divide(green.add(swir1))))))
//print(IBI_temp);
return IBI_temp;
}
function NDVI_cal(img) {
var ndvi_temp = img.normalizedDifference(["SR_B5","SR_B4"]);
//print(ndvi_temp);
return ndvi_temp;
}
function Wet_cal(img) {
var blue = img.select("SR_B2");
var green = img.select("SR_B3");
var red = img.select("SR_B4");
var nir = img.select("SR_B5");
var swir1 = img.select("SR_B6");
var swir2 = img.select("SR_B7");
var wet_temp =blue.multiply(0.1511)
.add(green.multiply(0.1973))
.add(red.multiply(0.3283))
.add(nir.multiply(0.3407))
.add(swir1.multiply(-0.7117))
.add(swir2.multiply(-0.4559))
//print(wet_temp);
return wet_temp;
}
function NDWI_cal(img) {
var nir = img.select("SR_B5");
var green = img.select("SR_B3");
var ndwi_temp = green.subtract(nir).divide(green.add(nir));
return ndwi_temp;
}
function img_normalize(img){
var minMax = img.reduceRegion({
reducer:ee.Reducer.minMax(),
geometry: table,
scale: 1000,
maxPixels: 10e13,
})
var year = img.get('year')
var normalize = ee.ImageCollection.fromImages(
img.bandNames().map(function(name){
name = ee.String(name);
var band = img.select(name);
return band.unitScale(ee.Number(minMax.get(name.cat('_min'))), ee.Number(minMax.get(name.cat('_max'))));
})
).toBands().rename(img.bandNames());
return normalize;
} ;
// Map the function over one year of data.
var dataset = ee.ImageCollection("LANDSAT/LC08/C02/T1_L2")
.filterDate('2021-05-01', '2021-8-31')
.map(maskL8sr)
.median()
var dataset1 = ee.ImageCollection("MODIS/006/MOD11A2")
.filterDate('2021-05-01', '2021-8-31')
.median()
var ndwi=NDWI_cal(dataset);
var NDWI_mask = ndwi.lt(0.1);
var dataset_no_water = dataset.updateMask(NDWI_mask).clip(table);
var dataset1_no_water = dataset1.updateMask(NDWI_mask).clip(table);
var lst = dataset1_no_water.expression(
'a*0.02-273.15',
{
a:dataset1_no_water.select('LST_Day_1km'),
});
var SI = SI_cal(dataset_no_water);
var IBI = IBI_cal(dataset_no_water);
var NDBSI =(SI.add(IBI)).divide(2.0);
var ndvi = NDVI_cal(dataset_no_water);
var wet = Wet_cal(dataset_no_water);
var visParams1 = {
palette: '313695,4575b4,74add1,abd9e9,e0f3f8,ffffbf,fee090,fdae61,f46d43,d73027,a50026'
};
var visParams = {
min: -1,
max: 1,
palette: ['0000FF','ff0000','00ff00']
};
var visParams2 = {
palette: 'a50026,d73027,f46d43,fdae61,fee090,ffffbf,e0f3f8,abd9e9,74add1,4575b4,313695'
};
var unit_ndvi = img_normalize(ndvi);
dataset_no_water=dataset_no_water.addBands(unit_ndvi.rename('NDVI').toFloat())
var unit_NDBSI = img_normalize(NDBSI);
dataset_no_water=dataset_no_water.addBands(unit_NDBSI.rename('NDBSI').toFloat())
var unit_Wet = img_normalize(wet);
dataset_no_water=dataset_no_water.addBands(unit_Wet.rename('Wet').toFloat())
var unit_lst = img_normalize(lst);
dataset_no_water=dataset_no_water.addBands(unit_lst.rename('LST').toFloat())
//print(dataset_no_water);
var bands = ["Wet","NDVI","NDBSI","LST"]
var sentImage =dataset_no_water.select(bands)
var region = table;
var image = sentImage.select(bands);
var scale = 1000;
var bandNames = image.bandNames();
// 图像波段重命名函数
var getNewBandNames = function(prefix) {
var seq = ee.List.sequence(1, bandNames.length());
return seq.map(function(b) {
return ee.String(prefix).cat(ee.Number(b).int());
});
};
//数据平均
var meanDict = image.reduceRegion({
reducer: ee.Reducer.mean(),
geometry: region,
scale: scale,
maxPixels: 1e9
});
var means = ee.Image.constant(meanDict.values(bandNames));
var centered = image.subtract(means);
//主成分分析函数
var getPrincipalComponents = function(centered, scale, region) {
// 图像转为一维数组
var arrays = centered.toArray();
// 计算协方差矩阵
var covar = arrays.reduceRegion({
reducer: ee.Reducer.centeredCovariance(),
geometry: region,
scale: scale,
maxPixels: 1e9
});
// 获取“数组”协方差结果并转换为数组。
// 波段与波段之间的协方差
var covarArray = ee.Array(covar.get('array'));
// 执行特征分析,并分割值和向量。
var eigens = covarArray.eigen();
print(eigens)
// 特征值的P向量长度
var eigenValues = eigens.slice(1, 0, 1);
//计算主成分载荷
var eigenValuesList = eigenValues.toList().flatten()
var total = eigenValuesList.reduce(ee.Reducer.sum())
var percentageVariance = eigenValuesList.map(function(item) {
return (ee.Number(item).divide(total)).multiply(100).format('%.2f')
})
print('特征值',eigenValues )
print("贡献率", percentageVariance)
// PxP矩阵,其特征向量为行。
var eigenVectors = eigens.slice(1, 1);
print(eigenVectors)
// 将图像转换为二维阵列
var arrayImage = arrays.toArray(1);
//使用特征向量矩阵左乘图像阵列
var principalComponents = ee.Image(eigenVectors).matrixMultiply(arrayImage);
// 将特征值的平方根转换为P波段图像。
var sdImage = ee.Image(eigenValues.sqrt())
.arrayProject([0]).arrayFlatten([getNewBandNames('sd')]);
//将PC转换为P波段图像,通过SD标准化。
principalComponents=principalComponents
// 抛出一个不需要的维度,[[]]->[]。
.arrayProject([0])
// 使单波段阵列映像成为多波段映像,[]->image。
.arrayFlatten([getNewBandNames('pc')])
// 通过SDs使PC正常化。
.divide(sdImage);
return principalComponents
};
//进行主成分分析,获得分析结果
var pcImage = getPrincipalComponents(centered, scale, region);
var rsei_un_unit = pcImage.expression(
'constant - pc1' ,
{
constant: 1,
pc1: pcImage.select('pc1')
});
var rsei = img_normalize(rsei_un_unit);
//统
小风飞子
- 粉丝: 374
- 资源: 1961
最新资源
- 搜广推推荐系统中传统推荐系统方法思维导图整理-完整版
- 微藻检测19-YOLO(v5至v11)、COCO、CreateML、Paligemma、TFRecord、VOC数据集合集.rar
- 使用AS的自定义功能块与OS之间WINCC自定义功能块图标,自定义功能块面板教程 1.不是采用西门子APL面板实现 2.AS可以采用LAD或者SCL语言生成功能块 3.实现弹窗功能 4.事件可
- 等发达地区的无穷大无穷大无穷大请问
- Python实现常见排序算法详解
- JWaaaaaaaaaaaaaaaaaaaa
- Python复制重复数据工具.exe
- 2024圣诞节海外消费市场趋势及营销策略分析报告
- 基于Java的网上教务评教管理系统的设计与实现.doc
- EventHandlerError解决办法.md
- NotImplementedError.md
- SecurityException(解决方案).md
- IllegalAccessException(解决方案).md
- NameError.md
- NSRunLoopError如何解决.md
- OSError.md
资源上传下载、课程学习等过程中有任何疑问或建议,欢迎提出宝贵意见哦~我们会及时处理!
点击此处反馈
- 1
- 2
前往页