这篇笔记是GEE学习笔记的第29篇。前面我写过Sentinel-2的单期NDVI、写过水体指数提取，这次要解决一个特别高频的需求：把一整年的Sentinel-2影像处理成一张高精度NDVI年均值数据，并直接导出下载。所谓“高精度”，在这里指的是使用Level-2A表面反射率产品、逐像元云掩膜和10米空间分辨率——不是拿DN值随手一算，也不是用粗分辨率搞近似。适用场景很明确，比如区域植被年际变化分析、退耕还林评估、农业长势监测、生态红线的本底摸查，以及作为后续机器学习分类的输入特征。无论你是刚接触GEE的新手，还是正在做时间序列分析的老手，这份笔记都可以当作模板直接改。源码我会完整贴出来，并把它拆开讲清楚每一步在做什么、为什么这样做。

1. 项目背景与整体思路

1.1 为什么选Sentinel-2做NDVI年均值

这些年大家聊植被指数，绕不开MODIS、Landsat、Sentinel-2这三个数据源。MODIS的NDVI产品时间序列长、处理方便，但250米分辨率在县域尺度上会有明显的混合像元问题；Landsat 30米分辨率可用，但单星重访16天，再加上云，一年内干净影像数量并不算多。Sentinel-2A/B双星组网后5天重访，10米分辨率，还有Level-2A表面反射率成品，非常适合做中高分辨率的年度植被监测。

我自己的体感是，在南方丘陵地区，一块几平方公里的破碎耕地，MODIS基本是一片混合值，Sentinel-2能把田埂和小水塘都大致分开。所以做区域尺度的NDVI年均值，我首选S2。下面是几个常见数据源的对比，可以直观看到差距。

这里说的“高精度”其实有两层含义。第一层是空间精度，10米像元能把地物边界看得比较清楚。第二层是光谱精度，我们用的是COPERNICUS/S2_SR地表反射率产品，而不是顶层大气反射率COPERNICUS/S2，这样可以减少大气散射对红光和近红外波段的影响，NDVI本身也会更加可信。

1.2 从影像集合到年均值的核心链路

整体链路可以概括为四步：选数据、洗数据、算指数、汇总导出。选数据是确定时间范围、空间范围、云量阈值；洗数据是用QA60波段做像元级云/卷云掩膜；算指数是归一化差分；汇总导出是把一年内所有合格NDVI像元做均值，再按区域裁剪导出。

这个流程听起来简单，真正决定成果质量的其实是第2步和第4步之间的配合。云掩膜不彻底，年均值会被云污染像元拉低；掩膜太激进，有效观测数不足，合成结果出现空洞。所以我在源码里特意加了一个有效观测数图层，肉眼检查比事后找原因快得多。你看到红色区域的Valid observations值很低，就要警觉了。

2. 完整源码与逐段解析

2.1 可直接运行的GEE源码

先整理成一段能直接跑的完整源码。我默认把研究区设在北京中心点10公里缓冲区，年份写的是2023年，你只需要替换geometry和year两个变量就能用在别处。

// GEE学习笔记 29：一键下载 Sentinel-2 高精度 NDVI 年均值 // 平台：Google Earth Engine Code Editor（JavaScript API） // 1. 设置研究区：默认示例为北京中心10km缓冲区 // 如果有自己的矢量，可以直接写：var geometry = table; var geometry = ee.Geometry.Point([116.4074, 39.9042]).buffer(10000); // 2. 设置年份 var year = 2023; var startDate = ee.Date.fromYMD(year, 1, 1); var endDate = ee.Date.fromYMD(year + 1, 1, 1); // 3. Sentinel-2 Level-2A 表面反射率集合 var s2 = ee.ImageCollection('COPERNICUS/S2_SR') .filterBounds(geometry) .filterDate(startDate, endDate) .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20)); // 4. 云/卷云掩膜（基于 QA60 位掩膜，只保留干净像元） function maskS2clouds(image) { var qa = image.select('QA60'); var cloudBitMask = 1 << 10; var cirrusBitMask = 1 << 11; var mask = qa.bitwiseAnd(cloudBitMask).eq(0) .and(qa.bitwiseAnd(cirrusBitMask).eq(0)); return image.select(['B8', 'B4']) .updateMask(mask) .divide(10000) .copyProperties(image, ['system:time_start']); } // 5. 计算 NDVI function addNDVI(image) { var ndvi = image.normalizedDifference(['B8', 'B4']).rename('NDVI'); return image.addBands(ndvi); } // 6. 一年内全部合格像元求平均 var ndviAnnual = s2 .map(maskS2clouds) .map(addNDVI) .select('NDVI') .mean() .clip(geometry) .rename('NDVI_' + year); // 有效像元数统计（排查云掩膜过度或影像缺失） var validCount = s2 .map(maskS2clouds) .map(addNDVI) .select('NDVI') .count(); // 7. 地图可视化 var visParam = { min: -0.2, max: 1.0, palette: ['#0000ff', '#8b4513', '#d2b48c', '#ffff00', '#008000'] }; Map.centerObject(geometry, 11); Map.addLayer(ndviAnnual, visParam, 'NDVI annual mean ' + year); Map.addLayer(validCount, {min: 1, max: 60, palette: ['#ff0000', '#ffff00', '#00ff00']}, 'Valid observations'); // 8. 导出到 Google Drive（大范围、正式结果用这个） Export.image.toDrive({ image: ndviAnnual, description: 'S2_NDVI_Annual_Mean_' + year, folder: 'GEE_Exports', fileNamePrefix: 'S2_NDVI_Annual_Mean_' + year, region: geometry, scale: 10, maxPixels: 1e13 }); // 小区域也可以直接生成下载链接（控制台 print 后点击） print(ndviAnnual.getDownloadURL({ name: 'S2_NDVI_Annual_Mean_' + year, scale: 10, region: geometry }));

运行方式很简单：打开 Google Earth Engine Code Editor，粘贴这段代码，点Run。如果一切正常，地图上会出现两个图层：一个是带颜色的NDVI年均值，一个是有效观测数。右侧Tasks面板里会出现S2_NDVI_Annual_Mean_2023的导出任务，点击Run即可把GeoTIFF存到你的Google Drive里。

2.2 核心函数逐段拆解

源码里有几个关键点值得单独拿出来讲，尤其是maskS2clouds这个函数，很多人直接复制但看不懂，一旦结果不对就不知道从哪里排查。

先说影像集合筛选。filterBounds(geometry)是空间过滤，只保留覆盖研究区的影像；filterDate(startDate, endDate)是时间过滤，取年初到次年年初；filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20))是一个场景级云量过滤，意思是整景影像的云量百分比不能超过20%。这个阈值可以根据地区调整，后面我会展开讲。

再说云掩膜。QA60波段是Sentinel-2的质量波段，它用一个整数编码了很多位信息。1 << 10是1左移10位，也就是二进制第10位；bitwiseAnd(cloudBitMask)就是把这个位提取出来。如果结果是0，说明这个像元没有被云覆盖；如果结果是1，说明有云。卷云同理用第11位。最后把干净像元保留，云和卷云像元变成空值。

divide(10000)是因为Sentinel-2 Level-2A产品的表面反射率是以10000倍存储的整数，除以10000后得到0到1之间的反射率。这个步骤不是必须的，因为normalizedDifference在求比值时上下都同一个尺度，最终NDVI不变。但保留这个步骤更规范，也方便你后续输出反射率做别的分析。

addNDVI里的核心就一行：normalizedDifference(['B8', 'B4'])，公式是(NIR - Red) / (NIR + Red)，对应Sentinel-2就是(B8 - B4) / (B8 + B4)。注意B8是近红外，B4是红波段，顺序一旦写反，NDVI会变成负值。

3. 关键参数精讲：精度从何而来

3.1 场景级云量阈值与像元级QA60掩膜

很多初学者只做filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20))，以为这样就处理干净了，其实这是两个层次的事情。场景级云量是一个全局统计量，代表整景影像里有百分之多少的云；而QA60掩膜是在像元级别，把影像里残留的云和卷云像元单独抠掉。

比如一景影像云量10%，看起来不高，但云可能正好集中在你的研究区上空。如果你不做像元级掩膜，那这10%的云区域就会参与NDVI平均，直接把结果污染。反过来，QA60也不是万能的，它对薄云和云阴影的识别不完全，冬季高海拔地区的雪也可能误判。这里我给一个经验值：干燥地区场景级云量阈值可以放到20到30，湿润地区建议10到15。阈值太小容易把整年影像筛没了，阈值太大残留云多。

如果你所在区域地形复杂、阴影和雪比较麻烦，可以考虑用S2的SCL场景分类波段做二次过滤。SCL中的低值通常代表阴影或云，但这会让代码更复杂，暂不展开。

3.2 均值合成、中位数合成与最大值合成怎么选

“年均值”这三个字听着简单，但合成策略并不唯一。我用同一景区域分别跑过mean、median、max三种方法，结果差异其实很大。均值反映全年平均绿度状态，对农作物的多茬种植和常绿林地比较友好；中位数更稳健，能消除残留云和极端值的影响；最大值反映的是年内绿度峰值，常用于估算植被覆盖潜力。

这个任务既然标题叫年均值，我默认用mean()。如果你只是想初步看一眼区域植被生产力，中位数其实更稳。我个人的习惯是：正式跑趋势之前，至少同时算一份mean和一份median，看看结论是否稳健。如果两者趋势方向一致，说明你不是被某个季节的极端影像带偏的。

3.3 分辨率、坐标系与导出设置

scale参数直接决定导出像元大小。源码里写的10米对应Sentinel-2原始分辨率，是“高精度”的关键。但需要认清一个现实：10米分辨率是以像素数量成平方增长的。你的研究区如果是一个县，10米导出大约会有几万乘几万像元；如果是整个省，文件几个GB甚至更大，导出任务很容易失败。我建议根据范围调整：

坐标系方面，源码里我没有显式指定crs，这样GEE会使用影像默认投影，通常是原始S2瓦片所在的UTM投影。好处是输出像元与原始数据对齐，缺点是如果你一个大范围跨多个UTM分带，默认投影可能不是你想要的。如果你需要统一的WGS84经纬度坐标，可以在导出设置里加crs: 'EPSG:4326'。但注意，4326是地理坐标系，跨带变形会比UTM明显，所以常规做法还是默认投影，本地再用GIS转。

maxPixels是导出任务的兜底限制。10米分辨率大范围跑的时候，一定要把这个值调大，否则会报Error: Number of pixels too high。但也不要盲目填一个天文数字，先算一下研究区像元数量，再乘以1.2倍左右比较稳妥。

4. 实操过程：从零到一张年均NDVI图

4.1 三步快速跑通默认例子

拿到代码，最直接的方式就是先跑通默认参数，不要一上来就换成自己的陌生区域。具体操作分三步。

第一步，打开GEE Code Editor，新建脚本，粘贴源码，点击Run。如果代码没有语法错误，右侧Console会输出NDVI影像的信息和下载链接，地图也会自动加载图层。

第二步，检查地图上的NDVI年均值图层。植被区域应该呈绿色，裸地呈黄色或棕色，水体和建设用地偏蓝。这个直观检查能提前暴露波段组合或云掩膜的问题。再看Valid observations图层，如果一片区域全是红色，说明那一年有效观测数非常少，后面即使导出，结果也没有可信度。

第三步，点开Tasks面板，找到导出任务，点击右侧的Run按钮。第一次使用时GEE会要求你授权Google Drive，跟着提示走即可。等待任务完成，去你的Google Drive里找到GEE_Exports文件夹，下载GeoTIFF到本地。下载后可以在QGIS或ArcGIS里打开，查看波段属性，坐标参考是否正确。

4.2 把自己的研究区和年份替换进去

默认示例始终只是示例，实际项目里肯定要换成自己的区域。最快的方式是用GEE左侧的绘图工具画一个多边形，画完后在代码里把开头两句替换掉。

如果你画完多边形后GEE生成的是一个变量table，代码就写成：

var geometry = table.geometry();

如果你是通过Assets上传了Shapefile或GeoJSON，代码写成：

var geometry = ee.FeatureCollection('projects/你的用户名/assets/你的矢量路径').geometry();

还有一种常见做法是从其他公开数据集里筛出行政边界，像我源码里演示的其实是一行注释掉的写法：

// var geometry = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017') // .filter(ee.Filter.eq('country_na', 'China')).geometry();

注意，导出时region参数需要的是一个Geometry对象，而不是FeatureCollection。如果控制台打印出来是FeatureCollection，记得加.geometry()，否则导出任务会报 invalid region。

年份修改更简单，直接把var year = 2023;改成想要的年份。年份区间是自动计算的，startDate是当年1月1日，endDate是下一年1月1日，所以不需要手动改日期。

4.3 下载链接方式与Drive导出方式取舍

源码里我同时加入了getDownloadURL()，它会在Console输出一个可点击的下载链接。这个方案看似方便，但不是对所有区域都适用。GEE对于直接下载链接有像素数量限制，小范围比如十几平方公里以内没问题，一旦研究区变大，或者你想要10米高精度，链接会直接报错。

我的建议是：正式项目一律用Export.image.toDrive，因为任务式导出没有交互式链接那么严格的像素限制，而且可以异步排队，不容易断掉。直接下载链接更适合快速验证、临时拿一小块数据。如果你一定要大范围直接用链接下载，可以把大范围切成格网，循环生成小块的URL。但那样会让代码复杂不少，不如直接用Export然后等Drive同步。

5. 常见问题与排查技巧实录

5.1 导出任务成功但影像全黑或全空白

这个现象我遇到不止一次，最常见的原因是图层可视化参数和实际数值范围不匹配。比如你把visParam的min设成了0，max设成了1，但如果年均NDVI在某些高海拔区域普遍偏低，低于0就全部显示成蓝色，看起来像“空白”。这种情况下先用Map.addLayer配合一个更宽的范围比如min: -1, max: 1试试，而不是急着怀疑导出数据。

另一个原因是geometry类型不对。Export.image.toDrive的region参数只接受Geometry，如果你传入了FeatureCollection，GEE会在任务提交阶段报Invalid region，但有时候你看到的错误信息不够明确。我在排查时习惯先加一行print(geometry)，看看类型是什么，再做转换。

还有一种是几何边界太复杂。如果矢量边界有几十万个顶点，GEE处理这个边界就需要大量计算，结果任务虽然能导出，但研究区边缘可能会出现奇怪的黑边或空洞。对复杂边界可以先用geometry.simplify()做简化，或者对影像先clip再selfMask()。

5.2 NDVI数值范围异常或出现大量空值

如果你看到NDVI整片都是负值，大概率是波段顺序写反了，把normalizedDifference(['B8', 'B4'])写成了['B4', 'B8']。这在所有GEE新手翻车原因里排前三，我自己也干过。另一个容易踩的坑是把B8A当成了B8，Sentinel-2的窄近红外波段是B8A，但NDVI公式里用的是宽波段B8，两者数值有差异，导出后与别人产品对比时也会不一致。

出现大量空值，常见原因是云掩膜过于激进。比如你在湿润地区把场景级云量阈值设成5，一年下来可能只剩十来景影像，再抠掉云和卷云，很多像元在所有影像里都是被掩膜状态，合成结果自然为空。应对方式是放宽场景级云量阈值，或者改用中位数合成。如果只是背景区域空值，可以在导出前加一个unmask(-9999)，但要注意，这样会把没有观测值的区域变成-9999，而不是真正的NODATA，后续分析前一定要处理。

5.3 任务排队时间太长或内存溢出的应对

GEE的导出任务排队时间取决于当前服务器负载，尤其是下午时段，全球用户都在跑任务，排队一两个小时很正常。如果等了很久还是Queued，不要反复提交同一个任务，先把代码停掉，等一会儿再看。小任务如果也长时间排队，可以换Export.image.toCloudStorage，但需要配置云存储，不是所有用户都有条件。

内存溢出通常是大范围加上10米分辨率导致的。之前有个朋友想导出一个省份的NDVI年均值，scale直接写10，任务提交后很快报User memory limit exceeded。我的建议是先把maxPixels提高，然后看是否还报错；还报错就降分辨率，比如20或30米。大范围高分辨率的需求最好分块处理，不要指望一台浏览器就能合成全国10米数据，这不是代码能力问题，是地球引擎的资源限制。

6. 扩展思路与个人经验

6.1 改造成逐月或逐季节NDVI合成

年均值是基础，实际项目里很多需求要的是物候特征，比如生长季峰值、峰值出现时间、春季NDVI增速。把源码改造成逐月合成其实不难，用一个map遍历月份即可：

var months = ee.List.sequence(1, 12); var monthlyMean = ee.ImageCollection(months.map(function(m) { var start = ee.Date.fromYMD(year, m, 1); var end = start.advance(1, 'month'); var monthly = s2.filterDate(start, end) .map(maskS2clouds) .map(addNDVI) .select('NDVI') .mean() .set('month', m); return monthly; }));

得到的结果是一个以“月”为单位的ImageCollection，你可以直接导出每一张，或者用monthlyMean.toBands()合成一个12波段的影像。后面做物候分析时，这个12波段的堆叠影像会非常方便。

6.2 用NDVI年均值快速估算FVC植被覆盖度

很多研究需要植被覆盖度，而不仅仅是一个NDVI值。最常用的遥感估算方式是二分法模型，公式是：

FVC = (NDVI - NDVI_soil) / (NDVI_veg - NDVI_soil)

其中NDVI_soil是纯裸土NDVI，NDVI_veg是纯植被NDVI。这两个阈值有地区差异，干旱区通常取0.05和0.70，湿润地区可以取0.10和0.85。代码实现非常简短：

var fvc = ndviAnnual.subtract(0.10).divide(0.85 - 0.10).clamp(0, 1).rename('FVC');

注意，不要直接把0和1以外的值扔掉，有些区域NDVI比土壤值还低，比如水体，会被clamp成0；高密度植被NDVI超过植被阈值，会被clamp成1。这样做出来的FVC在0到1之间，后续统计覆盖度百分比会方便很多。

6.3 构建多年NDVI序列，衔接深度学习和地物分类

如果你想做2019到2024年连续6年的植被趋势，可以写一个循环把每年的ndviAnnual放到一个列表里，最后用toBands()叠加成多波段影像。波段名类似NDVI_2019、NDVI_2020，这个堆叠结果可以直接作为随机森林或深度学习的输入特征。

我自己的经验是，单靠年均NDVI做地物分类，区分草地和灌丛的效果往往一般；如果加上年内最大值、峰值时间和夏季均值，分类精度会有明显提升。这也是为什么我建议先把源码跑通，后续扩展的空间真的很大。目前GEE也支持把影像导出成TFRecord格式，配合TensorFlow做深度学习，前提是你对数据尺度、样本量有清醒的规划。刚开始不用追求太复杂的模型，先把手头NDVI时序数据整理干净，比任何算法都重要。

最后讲一点我自己踩出来的经验。年均NDVI并不是年份越长越适合直接比较，不同年份的有效观测数不同，云覆盖差异会造成合成结果系统偏低。所以每次导出一张年均值图，我都会同时导出有效观测数，跑趋势时把观测数小于10的像元剔除。这个习惯帮我避开了很多假信号。源码里的Valid observations图层就是干这个的，别嫌它占存储，建议一并导出。祝大家跑图顺利。