目录
一、基础环境配置
（一）研究区与地图设置
（二）时间范围定义
二、核心函数定义
（一）云去除函数（针对不同 Landsat 传感器）
1. Landsat 4/5/7 云去除（rmL457Cloud）
2. Landsat 8/9 云去除（rmL8Cloud）
（二）遥感指数计算函数（AddIndices）
三、卫星影像数据加载与预处理
（一）波段映射与标准化
（二）Landsat 8/9 数据加载与预处理
（三）Landsat 7 数据加载与预处理
（四）影像集合合并与指数添加
四、结果计算与可视化
（一）NDVI 均值图计算与可视化
（二）频率图计算（绿度、冠层、淹没）
1. 绿度频率（Greenness_frequency）
2. 冠层频率（Canopy_frequency）
3. 淹没频率（Inundation_frequency）
4. 频率计算与可视化
五、关键技术要点总结
六、运行结果
若觉得代码对您的研究 / 项目有帮助，欢迎点击打赏支持！需要完整代码的朋友，打赏后可在后台私信（复制文章标题发给我），我会尽快发您完整可运行代码，感谢支持！
详情请参考：基于 GEE 多源 Landsat 影像识别红树林分布区域-CSDN博客
本代码基于 Google Earth Engine（GEE）平台，针对 2018-2020 年指定研究区域，整合 Landsat 7 和 Landsat 8/9 卫星的 C02 级地表反射率数据，通过云去除、波段标准化、植被与水体指数计算，最终生成 NDVI 均值图及绿度、冠层、淹没频率图，用于生态环境（如植被覆盖、湿地 inundation 等）时空特征分析。
一、基础环境配置
（一）研究区与地图设置
var roi = AOI; // 定义研究区域（AOI为预先设定的矢量边界） Map.addLayer(roi, {'color': 'grey'}, 'StudyArea'); // 在地图上添加研究区图层，灰色显示 Map.centerObject(roi,11); // 地图中心定位到研究区，缩放级别11（1-20，数值越大越精细）- 核心作用：确定分析范围，可视化研究区边界，便于后续数据筛选与结果查看。
- 关键说明：
AOI需提前在 GEE 中导入或绘制（如 shp 矢量文件转换的 GEE 资产），缩放级别可根据研究区大小调整。
（二）时间范围定义
var year_start = 2018; // 起始年份 var year_end = 2020; // 结束年份核心作用：限定卫星影像的时间筛选范围，后续数据加载仅获取该时段内的影像。
二、核心函数定义
（一）云去除函数（针对不同 Landsat 传感器）
Landsat 卫星影像受云、云影、积雪、卷云干扰，需通过质量控制波段（QA_PIXEL）进行掩码处理。
1. Landsat 4/5/7 云去除（rmL457Cloud）
function rmL457Cloud(image) { var qa = image.select('QA_PIXEL'); // 选择质量控制波段 // 基于比特位筛选干扰项（GEE中QA_PIXEL波段比特位定义固定） var cloudShadow = qa.bitwiseAnd(1 << 3); // 3号比特位：云影 var snow = qa.bitwiseAnd(1 << 4); // 4号比特位：积雪 var cloud = qa.bitwiseAnd(1 << 5); // 5号比特位：云 // 生成掩码：排除云影、积雪、云的像素（值为0的保留，1的剔除） var mask = cloudShadow.eq(0).and(snow.eq(0)).and(cloud.eq(0)); // 二次掩码：排除波段不完整的边缘像素 var mask2 = image.mask().reduce(ee.Reducer.min()); // 应用掩码并保留原始影像属性（时间、元数据等） return image.updateMask(mask).updateMask(mask2) .copyProperties(image) .copyProperties(image, ["system:time_start", 'system:time_end']); }2. Landsat 8/9 云去除（rmL8Cloud）
function rmL8Cloud(image) { var qa = image.select('QA_PIXEL'); // 新增6号比特位：卷云（Landsat 8/9新增卷云检测，需额外排除） var cloudShadow = qa.bitwiseAnd(1 << 3); var snow = qa.bitwiseAnd(1 << 4); var cloud = qa.bitwiseAnd(1 << 5); var cirrus = qa.bitwiseAnd(1 << 6); // 掩码：排除云影、积雪、云、卷云 var mask = cloudShadow.eq(0).and(snow.eq(0)).and(cloud.eq(0)).and(cirrus.eq(0)); return image.updateMask(mask) .copyProperties(image, ["system:time_start", 'system:time_end']); }- 核心差异：Landsat 8/9 新增卷云检测比特位，需额外排除卷云干扰；Landsat 7 无卷云波段，需补充边缘像素掩码。
- 关键原理：
bitwiseAnd(1 << n)表示提取 QA_PIXEL 波段的第 n 位比特值，eq(0)表示保留该比特位为 0（无对应干扰）的像素。
（二）遥感指数计算函数（AddIndices）
通过卫星影像的光学波段，计算常用生态环境指数，用于表征植被覆盖、水体分布等。
function AddIndices(image) { // 1. NDVI（归一化植被指数）：(近红外-红)/ (近红外+红)，值域[-1,1]，高值代表高植被覆盖 var ndvi = image.normalizedDifference(["nir", "red"]).rename("NDVI"); // 2. MNDWI（改进型归一化水体指数）：(绿-短波红外1)/ (绿+短波红外1)，高值代表水体 var mNDWI = image.normalizedDifference(["green", "swir1"]).rename("MNDWI"); // 3. LSWI（陆地表面水分指数）：(近红外-短波红外1)/ (近红外+短波红外1)，高值代表高水分 var lswi = image.normalizedDifference(["nir", "swir1"]).rename("LSWI"); // 4. EVI（增强型植被指数）：减少大气和土壤背景干扰，更适合高植被覆盖区 var evi = image.expression( '2.5 * ((B5 - B4) / (B5 + 6 * B4 - 7.5 * B1 + 1))', { "B1": image.select("blue"), // 蓝波段 "B4": image.select("red"), // 红波段 "B5": image.select("nir") // 近红外波段 } ).rename("EVI"); // 5. 衍生指数：水分指数与植被指数的差值，表征植被水分胁迫或湿地特征 var lswi2ndvi = lswi.subtract(ndvi).rename("LSWI2NDVI"); var lswi2evi = lswi.subtract(evi).rename("LSWI2EVI"); // 将计算的指数作为新波段添加到原始影像 return image.addBands(ndvi) .addBands(evi) .addBands(lswi) .addBands(lswi2ndvi) .addBands(lswi2evi); }- NDVI/EVI：植被覆盖度指标，EVI 抗干扰能力更强；
- MNDWI/LSWI：水体 / 水分含量指标；
- LSWI2NDVI/LSWI2EVI：衍生指标，可用于区分湿地（高水分 + 高植被）与普通植被。
三、卫星影像数据加载与预处理
（一）波段映射与标准化
不同 Landsat 传感器的波段编号不同，需统一命名为通用波段名（蓝、绿、红、近红外等），便于后续统一处理。
// Landsat 8/9 C02级波段（SR_B2=蓝，SR_B3=绿，SR_B4=红，SR_B5=近红外，SR_B6=短波红外1，SR_B7=短波红外2） var LC8_BANDS = ['SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B6', 'SR_B7']; // Landsat 7 C02级波段（SR_B1=蓝，SR_B2=绿，SR_B3=红，SR_B4=近红外，SR_B5=短波红外1，SR_B7=短波红外2） var LC7_BANDS = ['SR_B1', 'SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B7']; var STD_NAMES = ['blue', 'green', 'red', 'nir', 'swir1', 'swir2']; // 通用波段名（二）Landsat 8/9 数据加载与预处理
var l8 = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2') // GEE内置Landsat 8 C02级地表反射率数据集 .filterBounds(roi) // 筛选研究区内的影像 .filter(ee.Filter.calendarRange(year_start, year_end, 'year')) // 筛选时间范围 .map(function(image) { // 光谱波段缩放：C02级数据原始值为整数，需转换为真实反射率（GEE官方缩放系数） var opticalBands = image.select('SR_B.').multiply(0.0000275).add(-0.2); // 热红外波段缩放：转换为真实温度值 var thermalBands = image.select('ST_B.*').multiply(0.00341802).add(149.0); // 替换原始波段为缩放后的波段 return image.addBands(opticalBands, null, true) .addBands(thermalBands, null, true); }) .map(rmL8Cloud) // 应用云去除函数 .select(LC8_BANDS, STD_NAMES); // 波段重命名为通用名（三）Landsat 7 数据加载与预处理
var l7 = ee.ImageCollection('LANDSAT/LE07/C02/T1_L2') // Landsat 7 C02级数据集 .filterBounds(roi) .filter(ee.Filter.calendarRange(year_start, year_end, 'year')) .map(function(image) { var opticalBands = image.select('SR_B.').multiply(0.0000275).add(-0.2); // Landsat 7仅1个热红外波段（ST_B6），缩放逻辑与Landsat 8一致 var thermalBands = image.select('ST_B6').multiply(0.00341802).add(149.0); return image.addBands(opticalBands, null, true) .addBands(thermalBands, null, true); }) .map(rmL457Cloud) // 应用Landsat 4/5/7云去除函数 .select(LC7_BANDS, STD_NAMES);（四）影像集合合并与指数添加
// 合并Landsat 7和8的影像集合，按时间排序 var L78Col = ee.ImageCollection(l7.merge(l8)) .sort("system:time_start") .map(AddIndices); // 对每幅影像计算遥感指数 // 释放内存（GEE中大型集合可手动释放，避免内存溢出） l7 = null; l8 = null;四、结果计算与可视化
（一）NDVI 均值图计算与可视化
// 可视化配色方案（从紫色到红色，中间为绿色系，适配NDVI值域） var paletteVis = 'D8BFD8, CE7E45, DF923D, F1B555, FCD163, 99B718, 74A901, 66A000, 529400,' + '3E8601, 207401, 056201, 004C00, 023B01, 012E01, 011D01, FF0000'; // 计算所有影像的NDVI均值，裁剪到研究区 var NDVIMean = L78Col.select("NDVI").mean().clip(roi); // 添加到地图，默认不显示（false），可在GEE界面手动勾选 Map.addLayer(NDVIMean, {'min': 0, 'max': 1, 'palette': paletteVis}, 'NDVIMean', false);意义：NDVI 均值图反映研究区 2018-2020 年平均植被覆盖水平，值越高植被越茂密。
（二）频率图计算（绿度、冠层、淹没）
频率图反映某一特征（如 “高水分 + 高植被”）在研究时段内出现的比例（单位：%），用于分析时空稳定性。
1. 绿度频率（Greenness_frequency）
var Greenness = L78Col.map(function(image) { // 筛选“高水分（LSWI>0）且高植被（EVI>0.2）”的像素（1表示满足，0不满足） var pixels_water = image.select("LSWI").gt(0); var pixel_veg = image.select("EVI").gt(0.2); return pixels_water.and(pixel_veg); // 逻辑与：同时满足两个条件 }).sum(); // 对所有影像求和：得到满足条件的影像次数2. 冠层频率（Canopy_frequency）
var Canopy = L78Col.map(function(image) { // 筛选“高水分（LSWI>0.3）且高植被覆盖（NDVI>0.3）”的像素 var pixels_water = image.select("LSWI").gt(0.3); var pixel_veg = image.select("NDVI").gt(0.3); return pixels_water.and(pixel_veg); }).sum();3. 淹没频率（Inundation_frequency）
var Inundation = L78Col.map(function(image) { // 筛选“水分-植被指数≥0”的像素（表征湿地或淹没植被） var pixels_water = image.select("LSWI2NDVI").gte(0); var pixel_veg = image.select("LSWI2EVI").gte(0); return pixels_water.or(pixel_veg); // 逻辑或：满足任一条件 }).sum();4. 频率计算与可视化
// 计算研究时段内的有效影像总数（按红波段计数，排除云掩膜后的无效影像） var ImgTotal = L78Col.select("red").count().clip(roi); // 频率 = 满足条件的次数 / 总影像数 × 100（转换为百分比） var Greenness_frequency = Greenness.divide(ImgTotal).multiply(100).clip(roi); var Canopy_frequency = Canopy.divide(ImgTotal).multiply(100).clip(roi); var Inundation_frequency = Inundation.divide(ImgTotal).multiply(100).clip(roi); // 添加到地图，默认不显示 Map.addLayer(Greenness_frequency, {'min': 0, 'max': 100, 'palette': paletteVis}, 'Greenness_frequency', false); Map.addLayer(Canopy_frequency, {'min': 0, 'max': 100, 'palette': paletteVis}, 'Canopy_frequency', false); Map.addLayer(Inundation_frequency, {'min': 0, 'max': 100, 'palette': paletteVis}, 'Inundation_frequency', false);- 绿度频率：高值区域为长期保持 “高水分 + 高植被” 的区域（如湿地、森林）；
- 冠层频率：高值区域为长期高植被覆盖且水分充足的区域（如成熟森林、灌溉农田）；
- 淹没频率：高值区域为长期处于淹没或高水分状态的区域（如河流、湖泊、沼泽）。
五、关键技术要点总结
- 数据标准化：Landsat 7/8 的波段编号与缩放系数不同，需统一命名和转换，确保后续指数计算一致性；
- 云去除逻辑：基于 QA_PIXEL 波段的比特位掩码，针对不同传感器的干扰类型（卷云、边缘像素）差异化处理；
- 指数选择：结合 NDVI/EVI（植被）、LSWI/MNDWI（水分）及衍生指数，全面表征生态环境特征；
- 频率分析：通过 “条件筛选 - 求和 - 归一化” 流程，将时序影像转换为频率图，反映特征的时间稳定性。
六、运行结果
若觉得代码对您的研究 / 项目有帮助，欢迎点击打赏支持！需要完整代码的朋友，打赏后可在后台私信（复制文章标题发给我），我会尽快发您完整可运行代码，感谢支持！
