GEE | 22秒带你看完密云水库2015年-2021年水体变化

22秒带你看完密云水库2015年-2021年水体变化

1.背景介绍

在北京的朋友可能感觉到了,北京貌似整个夏天一直在下雨。偶然间我看到这样的一个新闻:

这就激起了我的兴趣了,我想用遥感看看这几年密云水库是如何变化的。

2.工具与数据

工具:Google Earth Engine

数据:Sentinel-2 MSI Level 1C

时间:2015年9月—2021年9月

3.研究方法

4.方法实现

数据筛选:

//时间
var startDate = ee.Date('2015-9-1'); 
var endDate = ee.Date('2021-9-30'); 
//筛选
var sentinelcollection = SENTINEL
                  .filterDate(startDate, endDate)//时间过滤
                  .filterBounds(roi)//位置过滤
                  .filterMetadata("MGRS_TILE", "equals", "50TMK")//条带号过滤
                  .filter(ee.Filter.lte('CLOUDY_PIXEL_PERCENTAGE',5))//云量过滤

格式转换

//待筛选的波段
var bands = ['B3','B4','B8']
//数据格式转换
sentinelcollection = sentinelcollection.map(function(image){
      //影像像素深度转为Unit8
      image=image.add(32768).divide(65536).multiply(255).uint8()
      //波段重命名
      image=image.select(['B3', 'B4', 'B8'],['B8', 'B4', 'B3'])
      return image.select(bands)
    })

影像拉伸

//影像拉伸
var sentinelcollection2 = sentinelcollection.map(function(image){
      image=(image.subtract(120)).divide(40).multiply(255)
      return image.uint8()
    })

导出视频

// 导出视频
Export.video.toDrive({
  collection: sentinelcollection2,
  description: 'water',
  dimensions: 720,
  framesPerSecond: 3,//帧率
  region: roi,
})

完整代码:

//批量处理函数
Map.centerObject(roi, 11)
//时间
var startDate = ee.Date('2015-9-1'); 
var endDate = ee.Date('2021-9-30'); 
//筛选
var sentinelcollection = SENTINEL
                  .filterDate(startDate, endDate)//时间过滤
                  .filterBounds(roi)//位置过滤
                  .filterMetadata("MGRS_TILE", "equals", "50TMK")//条带号过滤
                  .filter(ee.Filter.lte('CLOUDY_PIXEL_PERCENTAGE',5))//云量过滤
//待筛选的波段
var bands = ['B3','B4','B8']
//数据格式转换
sentinelcollection = sentinelcollection.map(function(image){
      //影像像素深度转为Unit8
      image=image.add(32768).divide(65536).multiply(255).uint8()
      //波段重命名
      image=image.select(['B3', 'B4', 'B8'],['B8', 'B4', 'B3'])
      return image.select(bands)
    })
//影像拉伸
var sentinelcollection2 = sentinelcollection.map(function(image){
      image=(image.subtract(120)).divide(40).multiply(255)
      return image.uint8()
    })
​
Map.addLayer(sentinelcollection.first())
Map.addLayer(sentinelcollection2.first())
// 导出视频
Export.video.toDrive({
  collection: sentinelcollection2,
  description: 'water',
  dimensions: 720,
  framesPerSecond: 3,//帧率
  region: roi,
})

链接:
完整代码链接:
https://code.earthengine.google.com/5993bea8a031d03ad353d38a36c7dc07

打开链接,直接在Google earth中运行。如果你需要换成其他区域,该代码只需要在筛选中去掉“条带号”的指定语句。

​5.图像对比图


2020年9月9日密云水库的水体分布图


2021年9月9日密云水库的水体分布图

可以前往“地信遥感数据汇”(https://www.gisrsdata.com/)获取更多数据。

要基于Google Earth Engine (GEE) 平台实现密云水库2024消落范围的提取,可以参考文献中的方法,具体步骤如下: ### 1. **准备工作** - **注册和配置GEE账户**:确保你已经注册了GEE账户并安装了GEE Python API。 - **导入必要的库**:在Python环境中导入GEE库和其他必要的库。 ```python import ee import geemap ``` - **初始化GEE**: ```python ee.Initialize() ``` ### 2. **数据准备** - **选择研究区域**:定义密云水库的地理范围。 ```python # 密云水库的经纬度范围 region = ee.Geometry.Polygon( [[[116.4, 40.4], [116.8, 40.4], [116.8, 40.8], [116.4, 40.8]]] ) ``` - **选择Landsat 8 OLI数据**:获取2024的Landsat 8 OLI影像数据。 ```python start_date = '2024-01-01' end_date = '2024-12-31' l8 = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2').filterDate(start_date, end_date).filterBounds(region) ``` ### 3. **图像预处理** - **云掩膜处理**:去除云和云阴影。 ```python def cloud_mask(image): qa = image.select('QA_PIXEL') cloud_shadow_bit_mask = 1 << 3 clouds_bit_mask = 1 << 5 mask = qa.bitwiseAnd(cloud_shadow_bit_mask).eq(0) \ .And(qa.bitwiseAnd(clouds_bit_mask).eq(0)) return image.updateMask(mask) l8_clean = l8.map(cloud_mask) ``` - **大气校正**:Landsat 8 OLI数据已经进行了大气校正,但可以进一步处理以提高精度。 ```python def apply_scale_factors(image): optical_bands = image.select('SR_B.').multiply(0.0000275).add(-0.2) thermal_bands = image.select('ST_B.*').multiply(0.00341802).add(149.0) return image.addBands(optical_bands, None, True).addBands(thermal_bands, None, True) l8_preprocessed = l8_clean.map(apply_scale_factors) ``` ### 4. **水湿指数计算** - **计算多种水湿指数**:包括NDWI、MNDWI、NDMI、NMDI、TCWI等。 ```python def calculate_indices(image): ndwi = image.normalizedDifference(['SR_B3', 'SR_B5']).rename('NDWI') mndwi = image.normalizedDifference(['SR_B3', 'SR_B6']).rename('MNDWI') ndmi = image.normalizedDifference(['SR_B5', 'SR_B6']).rename('NDMI') nmdi = image.normalizedDifference(['SR_B5', ('SR_B6', 'SR_B7')]).rename('NMDI') tcwi = image.expression( '0.2651 * SR_B2 + 0.2367 * SR_B3 + 0.1296 * SR_B4 + 0.059 * SR_B5 - 0.7506 * SR_B6 - 0.5386 * SR_B7', {'SR_B2': image.select('SR_B2'), 'SR_B3': image.select('SR_B3'), 'SR_B4': image.select('SR_B4'), 'SR_B5': image.select('SR_B5'), 'SR_B6': image.select('SR_B6'), 'SR_B7': image.select('SR_B7')} ).rename('TCWI') return image.addBands(ndwi).addBands(mndwi).addBands(ndmi).addBands(nmdi).addBands(tcwi) l8_indices = l8_preprocessed.map(calculate_indices) ``` ### 5. **图像合成** - **确定高水位和低水位影像**:使用75百分位数和25百分位数合成影像。 ```python high_water = l8_indices.select('MNDWI').reduce(ee.Reducer.percentile([75])) low_water = l8_indices.select('MNDWI').reduce(ee.Reducer.percentile([25])) ``` ### 6. **背景均一化** - **使用MFCM算法**:增强湿地与非湿地的背景反射率差异。 ```python def mfcn(image): # MFCM算法的具体实现 # 这里假设你已经有了MFCM算法的实现 return image high_water_mfcn = mfcn(high_water) low_water_mfcn = mfcn(low_water) ``` ### 7. **动态阈值分割** - **使用OTSU算法**:确定分割阈值。 ```python def otsu(image): # OTSU算法的具体实现 # 这里假设你已经有了OTSU算法的实现 return image high_water_otsu = otsu(high_water_mfcn) low_water_otsu = otsu(low_water_mfcn) ``` ### 8. **消落范围提取** - **计算消落范围**:通过高水位和低水位影像的差异来确定消落范围。 ```python drawdown_zone = high_water_otsu.subtract(low_water_otsu).gt(0) ``` ### 9. **结果可视化** - **将结果导出或可视化**: ```python # 导出结果 task = ee.batch.Export.image.toDrive( image=drawdown_zone, description='drawdown_zone_2024', folder='GEE_Results', region=region.getInfo()['coordinates'], scale=30, maxPixels=1e13 ) task.start() # 可视化结果 Map = geemap.Map() Map.addLayer(drawdown_zone, {'palette': ['blue']}, 'Drawdown Zone') Map.centerObject(region, 10) Map ``` ### 10. **精度验证** - **选择验证样本点**:使用高分辨率影像(如高分一号、WorldView、Google Earth)进行精度验证。 - **计算精度指标**:包括总体精度、生产者精度、用户精度和Kappa系数。 通过以上步骤,你可以基于GEE平台实现密云水库2024消落范围的提取。每个步骤的具体实现细节可以根据实际需求和数据情况进行调整。
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

叫我锐多宝

请我喝杯啤酒吧

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值