找回密码
 立即注册
首页 资源区 代码 利用Sentinel-2监测北京奥林匹克森林公园2024年NDVI变化 ...

利用Sentinel-2监测北京奥林匹克森林公园2024年NDVI变化

跟尴 5 天前
01 说明

1.1 逻辑和流程

简要流程:

  • 获取2024年覆盖北京奥林匹克森林公园的所有Sentinel-2影像
  • 对所有不同时间段的影像分别计算NDVI
  • 对于同一时间段的影像,取公园内所有像元NDVI值的中位数作为该时间点的NDVI
  • 将所有时间点的NDVI综合绘制折线图
  • 地图上展示公园的真彩色Sentinel-2图层和NDVI图层
ps: 提供JS版本和Python版本(在Pycharm中可正常运行, colab未尝试)参考学习,二者逻辑基本保持一致.
1.2 数据集说明

使用的数据集为: Harmonized Sentinel-2 MSl:MultiSpectral Instrument, Level-2A (SR)(SR表示表面反射率即地表反射率, TOA版本为大气层顶反射率<包含大气层的影响>).
使用到的相关波段信息如下:
波段名称 描述 分辨率 比例系数
B2Blue10 meters0.0001
B3Green10 meters0.0001
B4Red10 meters0.0001
B8NIR10meters0.0001
其中,比例系数0.0001与像元DN值相乘即可得到真正的Sentinel-2表面反射率(原始值通过整数存储节省存储,通过0.0001缩放回来)
Sentinel-2中存在属性CLOUDY_PIXEL_PERCENTAGE表示影像的云覆盖率(单位为%)。
02 JS代码
  1. /* demo1 利用Sentinel-2监测北京奥林匹克森林公园2024年NDVI变化 北京奥林匹克森林公园经度: 116.388768°E, 纬度: 39.988588°N */ // 准备 var start_date = &#39;2024-01-01&#39;; var end_date = &#39;2024-12-31&#39; var pt = ee.Geometry.Point([116.388768, 39.988588]); // 定义矢量点 var roi = pt.buffer(2000); // 2000m缓冲区 var vis_ture_color = { // 真彩色显示参数 bands: [&#39;B4&#39;, &#39;B3&#39;, &#39;B2&#39;], min: 0, max: 1, gamma: 1.6 } var vis_ndvi = { // NDVI显示参数 bands: &#39;NDVI&#39;, min: 0, max: 1, palette: [&#39;white&#39;, &#39;green&#39;, &#39;yellow&#39;] } // 获取sentinel-2的影像集合 var s2_collection = ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED") .filterBounds(roi) // 按空间范围筛选 .filterDate(start_date, end_date) // 按时间范围筛选 .filter(ee.Filter.lt(&#39;CLOUDY_PIXEL_PERCENTAGE&#39;, 20)); // 筛选小于20%云覆盖的影像 // 计算NDVI var add_ndvi = function(img){ var ndvi = img.normalizedDifference([&#39;B8&#39;, &#39;B4&#39;]).rename(&#39;NDVI&#39;); // 计算NDVI return img.addBands([ndvi]); // 返回添加了NDVI波段的img } var s2_ndvi = s2_collection.map(add_ndvi).select(&#39;NDVI&#39;); // 绘制ndvi var ndvi_chart = ui.Chart.image.series({ imageCollection: s2_ndvi, // 必填 region: roi, // 必填 reducer: ee.Reducer.mean(), // 默认就是该参数, 可不填 scale: 30, // 在30m的分辨率进行reducer(此处即在30m分辨率下进行均值计算), 可不填 xProperty: &#39;system:time_start&#39; // 默认就是该参数 }) // 为绘制的ndvi添加样式 ndvi_chart.setOptions({ title: &#39;2024年北京奥林匹克森林公园NDVI时间序列&#39;, vAxis: {title: &#39;NDVI(平均值)&#39;, viewWindow: {min: 0, max: null}}, hAxis: {title: &#39;日期&#39;, format: &#39;MM-dd&#39;, gridlines: {count: 10}}, series:{ 0: { color: &#39;green&#39;, lineWidth: 1, pointSize: 1.5, pointShape: &#39;cicrle&#39; } } }) // 在控制台显示该图表 print(ndvi_chart) // 地图显示 var img_median = s2_collection.median().clip(roi).multiply(0.0001); // 中位数合成 var ndvi_median = s2_ndvi.median().clip(roi); // 中位数合成ndvi Map.centerObject(roi, 14); Map.addLayer(img_median, vis_ture_color, &#39;真彩色影像 (median)&#39;); Map.addLayer(ndvi_median, vis_ndvi, &#39;NDVI (median)&#39;);
复制代码
结果展示:



03 Python代码
  1. #%% md # demo1 利用Sentinel-2监测北京奥林匹克森林公园2024年NDVI变化 北京奥林匹克森林公园经度: 116.388768°E, 纬度: 39.988588°N #%% import geemap import ee from matplotlib import pyplot as plt from matplotlib import dates as mdates import pandas as pd #%% # 准备 start_date = &#39;2024-01-01&#39; end_date = &#39;2024-12-31&#39; pt = ee.Geometry.Point([116.388768, 39.988588]) # 定义矢量点 roi = pt.buffer(2000) # 做2000m缓冲区视为北京奥林匹克森林公园区域 #%% # 获取Sentinel-2影像集合 s2_collection = (ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED") .filterDate(start_date, end_date) # 筛选时间范围 .filterBounds(roi) # 筛选空间范围 .filter(ee.Filter.lt(&#39;CLOUDY_PIXEL_PERCENTAGE&#39;, 20)) # 保留云覆盖率小于20%的影像 # .map(lambda img: img.normalizedDifference([&#39;B8&#39;, &#39;B4&#39;]).rename(&#39;NDVI&#39;)) # .select(&#39;NDVI&#39;) ) #%% # 计算NDVI def add_ndvi(img): ndvi = img.normalizedDifference([&#39;B8&#39;, &#39;B4&#39;]).rename(&#39;NDVI&#39;) # 计算NDVI并修改为波段名称为NDVI return img.addBands(ndvi) # 添加NDVI波段 s2_ndvi = s2_collection.map(add_ndvi).select(&#39;NDVI&#39;) #%% # 获取roi的NDVI均值(中位数) def extract_roi_ndvi(img): stats = img.reduceRegion( reducer=ee.Reducer.median(), # 计算的统计量 geometry=roi, # 统计的区域 scale=30 # 允许在30m分辨率进行统计而非此处s2的10m分辨率 ) return ee.Feature(None).set(stats).set(&#39;date&#39;, img.date().format(&#39;YYYY-MM-dd&#39;)) s2_ndvi_f = s2_ndvi.map(extract_roi_ndvi) ndvi_list = [p[&#39;properties&#39;] for p in s2_ndvi_f.getInfo()[&#39;features&#39;]] ndvi_df = pd.DataFrame(ndvi_list).sort_values(&#39;date&#39;).set_index(&#39;date&#39;) ndvi_df #%% # 绘制NDVI折线图 plt.rcParams[&#39;font.family&#39;] = [&#39;Times New Roman&#39;, &#39;SimSun&#39;] # 可显示中文(中文宋体-SimSun, 英文新罗马-Times New Roman) plt.rcParams[&#39;axes.unicode_minus&#39;] = False # 可显示负号 fig, ax = plt.subplots(figsize=(12, 6)) ndvi_df.plot( ax=ax, style=&#39;-o&#39;, color=&#39;green&#39;, title=&#39;2024年北京奥林匹克森林公园NDVI时间序列&#39; ) # 优化XY轴信息 ax.set_xlabel(&#39;日期&#39;, size=16) ax.set_ylabel(&#39;NDVI(中位数)&#39;, size=16) ax.grid(True) ax.set_ylim(0, None) # 优化X轴日期显示 ax.xaxis.set_major_formatter(mdates.DateFormatter(&#39;%m-%d&#39;)) plt.xticks(size=14) plt.tight_layout() # 自动调整布局(紧凑) plt.show() #%% # 地图Map上显示NDVI和真彩色底图 ndvi_median = s2_ndvi.median().clip(roi) s2_true_color = s2_collection.median().clip(roi).multiply(0.0001) Map = geemap.Map() Map.centerObject(roi, 14) # roi居中显示 Map.addLayer(ndvi_median, {&#39;min&#39;: 0, &#39;max&#39;: 1, &#39;palette&#39;: [&#39;white&#39;, &#39;green&#39;, &#39;yellow&#39;]}, &#39;NDVI (中位数)&#39;) Map.addLayer(s2_true_color, {&#39;bands&#39;: [&#39;B4&#39;, &#39;B3&#39;, &#39;B2&#39;], &#39;min&#39;: 0, &#39;max&#39;: 1, &#39;gamma&#39;: 1.6}, &#39;Sentinel-2 (真彩色)&#39;) #%% Map
复制代码
结果展示:



本文由博客一文多发平台 OpenWrite 发布!

来源:程序园用户自行投稿发布,如果侵权,请联系站长删除
免责声明:如果侵犯了您的权益,请联系站长,我们会及时删除侵权内容,谢谢合作!

相关推荐

您需要登录后才可以回帖 登录 | 立即注册