GEE | 单波段、双波段、多波段水深反演尝试

文摘   2024-06-02 15:56   山西  

LXX

读完需要

2
分钟

速读仅需 1 分钟

前言:

水深是水文研究和水生环境变化中的重要参数,常用的测量手段是声纳。近年来,随着遥感技术的发展,一些研究者尝试利用遥感技术进行水深测量,即基于遥感的水深反演,是为了弥补声纳测量的缺陷而发展起来的。遥感水深反演是指利用遥感数据测量水深,从而估算水深的方法。本文通过 GEE 平台,对常用的三种水深反演方法(单波段、双波段、多波段)进行尝试。


希望各位同学点个关注,点个小赞,这将是更新的动力,不胜感激❥(^_-)


1

   

原理

水深反演的基本原理基于水体对光的吸收和散射特性。当太阳辐射进入水体后,会经历反射、吸收和散射,传感器接收到的能量包括以下几部分:

在大气传输过程中散射的太阳辐射。在水-空气界面反射的辐射。进入水体后再次散射返回大气的辐射。到达水底后反射的辐射。

进入水体的光能随着深度的增加而逐渐衰减,最终只有一小部分光能到达水底并被反射回传感器。这部分光携带了水下地形的信息,是遥感水深反演的主要来源。

光在水体中的衰减可以表示为:

在湖泊清澈且基底均匀的情况下,可以直接使用上述公式估算水深。对于基底分布不均且悬浮物浓度高的湖泊,如洞庭湖,使用多光谱反射率模拟水深更为适宜。


2

   

方法

2.1. 单波段模型 (Single-band model)

单波段模型基于水深反演的基本原理,使用单一波段进行水深反演:

    其中,z水深,RE 是水体的反射率,a0 和a1是通过现场测量数据和相对水深      关系求解的待定系数。

2.2. 多波段模型 (Multi-band model)

多波段模型集成了多个波段的信息,相比单波段模型更为准确,适用于底质分布不均且悬浮物浓度高的湖泊:

其中, z是水深,a是待定系数,R是第n个波段的反射率。

2.3. 双波段比值模型 (Dual-band ratio model)

双波段比值模型通过两个波段的底反射率比值消除基底类型和水体类型的影响,用于计算水深:

其中,Ri和Rj分别是第i和j个波段的反射率。


3

   

研究区

洞庭湖(北纬 28°30′–30°20′,东经 111°40′–113°10′)是中国第二大淡水湖(面积 2691 平方公里),位于湖南省北部,长江以南盆地内,通过支流河道与长江相连(图 1)。湖的东、南、西三大部分已被《拉姆萨尔公约》指定为国际重要湿地。湖面平均水深 6.7 米,最深处 30.8 米。年平均径流量 3126 亿立方米。年水位落差约为 12–14 米,最大出现在 8 月,最小出现在 1 月或 2 月。这构成了维持大面积滨海湿地的基本水文状况,而滨海湿地是洞庭湖地区最具活力的景观。东洞庭湖是洞庭湖最大的区域,直接与长江相连,在调蓄洪涝灾害中发挥着重要作用。由于东洞庭湖具有重要的生态价值,尤其是多样性,1994 年起被设立为重要自然保护区——东洞庭湖国家级自然保护区。

4

   

GEE代码

其中需要通过NDWI来确定水体范围。

// 加载 Landsat 8 图像var landsat8 = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR');
// 定义感兴趣区域 (ROI)var roi = ee.FeatureCollection('users/geestudy2/DTH'); // table为自己上传的矢量边界
// 第一步:使用 QA 波段去云function maskClouds(image) { var qa = image.select('pixel_qa'); // 获取云和云影的位信息 var cloud = qa.bitwiseAnd(1 << 5).eq(0); var cloudShadow = qa.bitwiseAnd(1 << 3).eq(0); return image.updateMask(cloud.and(cloudShadow)) .select(['B1', 'B2', 'B3', 'B4', 'B5', 'B6', 'B7']) .divide(10000);}
// 第二步:图像合成var composite = landsat8 .filterDate('2016-11-21', '2016-12-21') .filterBounds(roi) .map(maskClouds) .median() .clip(roi);
// 显示合成图像Map.addLayer(composite, {bands: ['B4', 'B3', 'B2'], max: 0.3}, 'Composite Image');
// 检查合成波段的范围var compositeStats = composite.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, bestEffort: true});print('Composite band ranges:', compositeStats);
// 第三步:计算 NDWIvar ndwi = composite.normalizedDifference(['B3', 'B5']).rename('NDWI');
// 显示 NDWI 图像Map.addLayer(ndwi, {min: -1, max: 1, palette: ['00FFFF', '0000FF']}, 'NDWI');
// 检查 NDWI 的范围var ndwiStats = ndwi.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, bestEffort: true});print('NDWI range:', ndwiStats);
// 第四步:使用 NDWI 进行水体遮罩var waterMask = ndwi.gt(0);var waterOnlyComposite = composite.updateMask(waterMask);
// 显示仅包含水体的合成图像Map.addLayer(waterOnlyComposite, {bands: ['B4', 'B3', 'B2'], max: 0.3}, 'Water-only Composite');
// 检查仅包含水体的合成波段的范围var waterOnlyCompositeStats = waterOnlyComposite.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, bestEffort: true});print('Water-only composite band ranges:', waterOnlyCompositeStats);
// 第五步:阳光反射校正
// 第六步:水深反演计算// 单波段方法var depthSingleBand = waterOnlyComposite.expression( '16.968 * log(GREEN) + 50.975', { 'GREEN': waterOnlyComposite.select('B3').add(1).log() });
// 多波段方法var depthMultiBand = waterOnlyComposite.expression( '(-5.469) * log(BLUE) + 24.356 * log(GREEN) + 3.063 * log(RED)+60 ', { 'BLUE': waterOnlyComposite.select('B2').add(1).log(), 'GREEN': waterOnlyComposite.select('B3').add(1).log(), 'RED': waterOnlyComposite.select('B4').add(1).log() });
// 双波段比值方法var depthDualBand = waterOnlyComposite.expression( '32.176 * log(GREEN / RED) + 22.725', { 'GREEN': waterOnlyComposite.select('B3').add(1).log(), 'RED': waterOnlyComposite.select('B4').add(1).log() });
// 检查计算的水深范围var singleBandStats = depthSingleBand.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, bestEffort: true});
var multiBandStats = depthMultiBand.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, bestEffort: true});
var dualBandStats = depthDualBand.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, bestEffort: true});
print('Single-band depth range:', singleBandStats);print('Multi-band depth range:', multiBandStats);print('Dual-band ratio depth range:', dualBandStats);
// 第七步:可视化和图例创建var depthVisParams = { min: 0, max: 30, palette: ['blue', 'cyan', 'green', 'yellow', 'red']};
Map.centerObject(roi, 10);Map.addLayer(depthSingleBand, depthVisParams, 'Depth (Single-band)');Map.addLayer(depthMultiBand, depthVisParams, 'Depth (Multi-band)');Map.addLayer(depthDualBand, depthVisParams, 'Depth (Dual-band ratio)');
// 图例创建var legend = ui.Panel({style: {position: 'bottom-left'}});var legendTitle = ui.Label({ value: 'Water Depth (m)', style: {fontWeight: 'bold'}});legend.add(legendTitle);
var makeColorBar = function(palette) { var colorBar = ui.Thumbnail({ image: ee.Image.pixelLonLat().select(0), params: { bbox: [0, 0, 1, 0.1], dimensions: '100x10', format: 'png', min: 0, max: 1, palette: palette, }, style: {stretch: 'horizontal', margin: '0px 8px'} }); return colorBar;};
var colorBar = makeColorBar(depthVisParams.palette);legend.add(colorBar);
var legendLabels = ui.Panel({ widgets: [ ui.Label('0', {margin: '4px 8px'}), ui.Label('15', {margin: '4px 8px', textAlign: 'center', stretch: 'horizontal'}), ui.Label('30', {margin: '4px 8px'}) ], layout: ui.Panel.Layout.flow('horizontal')});legend.add(legendLabels);
Map.add(legend);

结果:

参考文献水深反演结果:


【参考】

采用的模型系数参考以下文献:Yang Nan, Li Jianhui, Mo Wenbo, Luo Wangjun, Wu Di, Gao Wanchao, Sun Changhao,Water depth retrieval models of East Dongting Lake, China, using GF-1 multi-spectral remote sensing images,Global Ecology and Conservation,Volume 22,2020,e01004,ISSN 2351-9894,https://doi.org/10.1016/j.gecco.2020.e01004.


遥感小屋
分享遥感相关文章、代码,大家一起交流,互帮互助
 最新文章