GEE&Python-demo1：利用Sentinel-2监测北京奥林匹克森林公园2024年NDVI变化(附Python版)
<h1 id="01-说明">01 说明</h1><h2 id="11-逻辑和流程">1.1 逻辑和流程</h2>
<p>简要流程：</p>
<ol>
<li>获取2024年覆盖北京奥林匹克森林公园的所有Sentinel-2影像</li>
<li>对所有不同时间段的影像分别计算NDVI</li>
<li>对于同一时间段的影像，取公园内所有像元NDVI值的中位数作为该时间点的NDVI</li>
<li>将所有时间点的NDVI综合绘制折线图</li>
<li>地图上展示公园的真彩色Sentinel-2图层和NDVI图层</li>
</ol>
<p><code>ps</code>： 提供JS版本和Python版本(在Pycharm中可正常运行, colab未尝试)参考学习，二者逻辑基本保持一致.</p>
<h2 id="12-数据集说明">1.2 数据集说明</h2>
<p>使用的数据集为: <code>Harmonized Sentinel-2 MSl:MultiSpectral Instrument, Level-2A (SR)</code>(SR表示表面反射率即地表反射率, TOA版本为大气层顶反射率<包含大气层的影响>).</p>
<p>使用到的相关波段信息如下：</p>
<table>
<thead>
<tr>
<th>波段名称</th>
<th>描述</th>
<th>分辨率</th>
<th>比例系数</th>
</tr>
</thead>
<tbody>
<tr>
<td>B2</td>
<td>Blue</td>
<td>10 meters</td>
<td>0.0001</td>
</tr>
<tr>
<td>B3</td>
<td>Green</td>
<td>10 meters</td>
<td>0.0001</td>
</tr>
<tr>
<td>B4</td>
<td>Red</td>
<td>10 meters</td>
<td>0.0001</td>
</tr>
<tr>
<td>B8</td>
<td>NIR</td>
<td>10meters</td>
<td>0.0001</td>
</tr>
</tbody>
</table>
<p>其中，比例系数0.0001与像元DN值相乘即可得到真正的Sentinel-2表面反射率(原始值通过整数存储节省存储,通过0.0001缩放回来)</p>
<p>Sentinel-2中存在属性<code>CLOUDY_PIXEL_PERCENTAGE</code>表示影像的云覆盖率(单位为%)。</p>
<h1 id="02-js代码">02 JS代码</h1>
<pre><code class="language-javascript">/*
demo1 利用Sentinel-2监测北京奥林匹克森林公园2024年NDVI变化
北京奥林匹克森林公园经度: 116.388768°E, 纬度: 39.988588°N
*/
// 准备
var start_date = '2024-01-01';
var end_date = '2024-12-31'
var pt = ee.Geometry.Point();// 定义矢量点
var roi = pt.buffer(2000);// 2000m缓冲区
var vis_ture_color = {// 真彩色显示参数
bands: ['B4', 'B3', 'B2'],
min: 0,
max: 1,
gamma: 1.6
}
var vis_ndvi = {// NDVI显示参数
bands: 'NDVI',
min: 0,
max: 1,
palette: ['white', 'green', 'yellow']
}
// 获取sentinel-2的影像集合
var s2_collection = ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED")
.filterBounds(roi)// 按空间范围筛选
.filterDate(start_date, end_date)// 按时间范围筛选
.filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20));// 筛选小于20%云覆盖的影像
// 计算NDVI
var add_ndvi = function(img){
var ndvi = img.normalizedDifference(['B8', 'B4']).rename('NDVI');// 计算NDVI
return img.addBands();// 返回添加了NDVI波段的img
}
var s2_ndvi = s2_collection.map(add_ndvi).select('NDVI');
// 绘制ndvi
var ndvi_chart = ui.Chart.image.series({
imageCollection: s2_ndvi,// 必填
region: roi,// 必填
reducer: ee.Reducer.mean(),// 默认就是该参数, 可不填
scale: 30,// 在30m的分辨率进行reducer(此处即在30m分辨率下进行均值计算), 可不填
xProperty: 'system:time_start'// 默认就是该参数
})
// 为绘制的ndvi添加样式
ndvi_chart.setOptions({
title: '2024年北京奥林匹克森林公园NDVI时间序列',
vAxis: {title: 'NDVI(平均值)', viewWindow: {min: 0, max: null}},
hAxis: {title: '日期', format: 'MM-dd', gridlines: {count: 10}},
series:{
0: {
color: 'green',
lineWidth: 1,
pointSize: 1.5,
pointShape: 'cicrle'
}
}
})
// 在控制台显示该图表
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, '真彩色影像 (median)');
Map.addLayer(ndvi_median, vis_ndvi, 'NDVI (median)');
</code></pre>
<p>结果展示：</p>
<p><img src="https://img2024.cnblogs.com/other/3638406/202508/3638406-20250816222055180-2109493596.png"><br>
<img src="https://img2024.cnblogs.com/other/3638406/202508/3638406-20250816222056323-16442940.png"><br>
<img src="https://img2024.cnblogs.com/other/3638406/202508/3638406-20250816222057471-1718230639.png"></p>
<h1 id="03-python代码">03 Python代码</h1>
<pre><code class="language-python">#%% 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 = '2024-01-01'
end_date = '2024-12-31'
pt = ee.Geometry.Point()# 定义矢量点
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('CLOUDY_PIXEL_PERCENTAGE', 20))# 保留云覆盖率小于20%的影像
# .map(lambda img: img.normalizedDifference(['B8', 'B4']).rename('NDVI'))
# .select('NDVI')
)
#%%
# 计算NDVI
def add_ndvi(img):
ndvi = img.normalizedDifference(['B8', 'B4']).rename('NDVI')# 计算NDVI并修改为波段名称为NDVI
return img.addBands(ndvi)# 添加NDVI波段
s2_ndvi = s2_collection.map(add_ndvi).select('NDVI')
#%%
# 获取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('date', img.date().format('YYYY-MM-dd'))
s2_ndvi_f = s2_ndvi.map(extract_roi_ndvi)
ndvi_list = for p in s2_ndvi_f.getInfo()['features']]
ndvi_df = pd.DataFrame(ndvi_list).sort_values('date').set_index('date')
ndvi_df
#%%
# 绘制NDVI折线图
plt.rcParams['font.family'] = ['Times New Roman', 'SimSun']# 可显示中文(中文宋体-SimSun, 英文新罗马-Times New Roman)
plt.rcParams['axes.unicode_minus'] = False# 可显示负号
fig, ax = plt.subplots(figsize=(12, 6))
ndvi_df.plot(
ax=ax,
style='-o',
color='green',
title='2024年北京奥林匹克森林公园NDVI时间序列'
)
# 优化XY轴信息
ax.set_xlabel('日期', size=16)
ax.set_ylabel('NDVI(中位数)', size=16)
ax.grid(True)
ax.set_ylim(0, None)
# 优化X轴日期显示
ax.xaxis.set_major_formatter(mdates.DateFormatter('%m-%d'))
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, {'min': 0, 'max': 1, 'palette': ['white', 'green', 'yellow']}, 'NDVI (中位数)')
Map.addLayer(s2_true_color, {'bands': ['B4', 'B3', 'B2'], 'min': 0, 'max': 1, 'gamma': 1.6}, 'Sentinel-2 (真彩色)')
#%%
Map
</code></pre>
<p>结果展示：</p>
<p><img src="https://img2024.cnblogs.com/other/3638406/202508/3638406-20250816222058097-1952369917.png"></p>
<p><img src="https://img2024.cnblogs.com/other/3638406/202508/3638406-20250816222100660-761491824.png"></p>
<p><img src="https://img2024.cnblogs.com/other/3638406/202508/3638406-20250816222103256-233164599.png"></p>
<blockquote>
<p>本文由博客一文多发平台 OpenWrite 发布！</p>
</blockquote><br><br>
来源：https://www.cnblogs.com/ChaoQiezi/p/19042735
頁:
[1]