资讯详情

GEE遥感影像预处理实战:影像加载、去云与波段系数转换

📅 2026/10/8 3:29:21 | 华诺云谱 👁 阅读
GEE遥感影像预处理实战:影像加载、去云与波段系数转换
如果你做过一段时间遥感数据处理大概会有这种感觉一张影像从天上下到本地硬盘还没开始做大气校正、去云、波段运算光是下载、裁剪、配准、再上传到服务器跑模型这一套流程就能耗掉大半天。而Google Earth EngineGEE这个平台把影像存储、计算资源和常用算法都放在云端直接在浏览器里写JavaScript或Python就能完成原本要本地折腾很久的活。这篇文章就从一个最常见的需求切入影像加载、去云、波段系数转换把我在GEE里实际操作时用到的思路、代码和踩坑点完整梳理一遍。这篇文章适合刚接触GEE、想快速上手做遥感影像预处理的同学也适合已经会写基础代码但被云掩膜、系数转换这些细节绕晕的人。我尽量把每一步为什么这么做讲清楚而不只是扔一段能跑的代码。1. 核心思路为什么这三件事要在GEE里一起做1.1 从“下载再处理”到“云上直接算”传统遥感处理流程里我们往往是先去USGS或者Copernicus Hub下载影像然后打开ENVI或者QGIS做辐射定标、大气校正、去云、裁剪最后再算植被指数。这个流程在地块少、时相少的时候还能接受一旦遇到“研究区覆盖几百平方公里、需要连续三年每月一景影像”这种需求数据下载量就是几百GB起步本地磁盘和电脑性能立刻变成瓶颈。GEE的做法完全不同。它把所有公开遥感数据Landsat、Sentinel、MODIS等都存成了云端数据集你写代码时只需要指定数据集名称、时间范围、空间范围GEE会在云端把筛选好的影像返回给你而且波段运算、去云掩膜、指数计算这些操作都在云端完成最终只把结果比如一张分类图或一段时间序列统计下载下来。这意味着本地不需要同时存储几百GB原始影像处理效率也可以大幅提升。我个人的使用感受是GEE把遥感的“数据工程”部分压缩到了极致。以前花在下载、解压、转格式、裁剪上的时间现在全部转化为写筛选条件和波段运算代码的时间。对于一个需要在多个地块反复实验的研究项目来说这种模式真的会让人上瘾。1.2 方案选型Landsat还是Sentinel-2做影像加载时首先面临的选择就是数据源。标题里的“影像加载”看似简单实际选哪种影像直接决定了后续去云方案和波段系数转换的逻辑。数据源空间分辨率重访周期常用产品适用场景Landsat 8/930米全色15米16天Collection 2 Level-2大范围、长时间序列分析Sentinel-210米部分波段20米/60米5天COPERNICUS/S2_SR中等尺度、需要较细空间细节MODIS250米~1000米1天MOD09GA等大尺度、植被物候分析做去云处理时如果选Landsat Collection 2的Surface Reflectance产品可以直接使用自带的QA_PIXEL波段来生成云掩膜如果选Sentinel-2 SR数据可以使用SCL波段或者QA60波段。我在实际项目中更倾向于Sentinel-2因为它的10米分辨率在小地块农业分析中优势太明显了而且5天重访周期让云遮挡的影响相对可控。选好数据源之后再去考虑“用什么波段、怎么去云、系数如何转换”整个处理链条就会清晰很多。2. 影像加载从数据集筛选到目标区域裁剪2.1 核心APIImageCollection的筛选与排序GEE中影像加载的本质是从ImageCollection影像集合里筛选出符合条件的一系列Image。最常用的筛选方法是filterBounds空间范围、filterDate时间范围和filterMetadata元数据筛选。下面这段JavaScript代码演示了一个最基础的影像加载流程// 加载研究区矢量 var roi ee.FeatureCollection(users/your_username/your_roi); // 加载Sentinel-2 SR影像集合 var s2 ee.ImageCollection(COPERNICUS/S2_SR) .filterBounds(roi) .filterDate(2023-01-01, 2023-12-31) .filterMetadata(CLOUDY_PIXEL_PERCENTAGE, less_than, 20); print(符合条件的影像数量:, s2.size());这段代码看起来很简单但有几个细节值得展开。filterMetadata这里用了影像自带的云量元数据表示只保留云量低于20%的影像。这是第一次粗筛真正精确的逐像元去云要在后面的环节里做。很多初学者把元数据云量当成最终的去云结果导致合成影像的边缘区域依然有云残留这里先埋个伏笔。另外加载完影像集合之后最好打印一下集合里的影像数量和时间分布确认数据是否充足。如果一年内影像数量太少可能要考虑扩大时间范围或者同时使用Landsat做补充。2.2 时间范围选择策略影像加载的时间范围不是随便写的它决定了最终合成的影像代表哪一段时期的地表状态。比如我要分析夏季植被长势通常选6月到9月如果做全年土地利用分类最好每个季节都选取影像避免只用夏季影像导致落叶林和常绿林难以区分。一个比较实用的技巧是把时间范围稍微放宽留出余量让云掩膜筛选后有足够的剩余像元。比如目标时间段是7月到8月实际加载时可以选6月1日到9月30日之后再通过去云和合成步骤提取中间时段的信息。这个逻辑有点像拍合影时多拍几张备选最后挑最好的用。2.3 研究区裁剪与统一坐标系加载影像后通常要clip到研究区范围。这一步不是为了省存储而是防止后续计算量爆炸var s2_clipped s2.map(function(img) { return img.clip(roi); });需要注意GEE在计算时是按需处理的也就是说你不主动裁剪后面做reduceRegion或者export时它也会根据你设定的区域去计算但提前裁剪可以减少中间结果的数据量让后续统计和可视化都更直观。坐标系方面GEE默认的投影是EPSG:4326WGS84但做面积统计时最好用等面积投影比如EPSG:32650UTM zone 50N否则高纬度地区的面积统计会有偏差。这个细节等到波段系数转换和面积统计阶段就会体现出来。3. 去云实现从逐景掩膜到合成去云3.1 单景影像去云QA波段与位运算去云是整个流程里最容易出问题也最值得深入理解的部分。GEE里不同数据源去云方式不同核心思路都是利用影像自带的云检测波段把这个波段解析成一个二值掩膜云像元标记为0非云像元标记为1然后把掩膜乘到影像各波段上去。以Sentinel-2为例COPERNICUS/S2_SR数据里包含一个SCL波段Scene Classification其中赋值分别为0代表无数据1代表饱和2代表暗影3代表云影4代表云5代表晴空6代表水7代表未分类8代表云边缘9代表薄卷云10代表高云11代表雪。SCL波段看起来很好用但实际项目中我更喜欢用QA60波段因为QA60是位编码波段直接用位运算就可以提取云和卷云掩膜。下面这个函数是我常用的去云方法用的是Sentinel-2官方推荐方式但做了一些调整function maskS2clouds(image) { var qa image.select(QA60); // 第10位和第11位分别代表云和卷云 var cloudBitMask 1 10; var cirrusBitMask 1 11; // 生成云掩膜云和卷云的位置设为0 var mask qa.bitwiseAnd(cloudBitMask).eq(0) .and(qa.bitwiseAnd(cirrusBitMask).eq(0)); return image.updateMask(mask); }这段代码是整个去云步骤的核心拆开解释一下为什么这样写。QA60是16位整型波段每一位代表一种检测结果。1 10是位运算左移即十进制值1024。如果把QA60波段的值和1024做bitwiseAnd运算结果不为0说明第10位是1也就是该像元被识别为云结果等于0说明第10位是0像元不是云。cirrusBitMask同理。最后用.and把两个条件合并得到同时不是云也不是卷云的像元作为有效像元。我第一次用这个函数时犯过一个典型错误把mask条件写反了导致云的位置变成1晴空变成0。Emm这个错误不仔细看很难发现因为最终合成的影像看起来只是“暗”了一点直到和原始影像对比才发现云的位置反而被保留了。这个教训让我养成了一个习惯做完掩膜一定要先单独可视化一下mask图层确认云区域是黑色、晴空是白色再往下走。3.2 时序合成去云中值合成法逐景去云只能处理单张影像但如果一个区域在目标时间段内有多景影像更推荐的做法是先逐景做云掩膜然后把它们合成为一张无云影像。GEE里最常用的合成方法是median()也就是对每个像元取所有影像在该位置的中值。为什么要用中值而不是平均值因为即便做了云掩膜仍然可能残留一些薄云或阴影这些异常值在中值统计中会被弱化。想象一下一景影像的某个像元因为薄云覆盖而反射率异常偏高平均值会被这个异常拉高而中值可以忽略这些离群值保留大多数晴空像元的正常值。这就像几个人统计身高如果一个数据明显是错误值取中位数比取平均数更稳妥。合成代码非常简洁var composite s2_clipped .map(maskS2clouds) .select([B2,B3,B4,B8]) .median();小提示select里面我只选了蓝、绿、红和近红外四个波段因为后续计算NDVI只需要这些波段没必要把所有波段都保留。GEE里养成按需选择波段的习惯能显著降低内存和计算压力。3.3 合成结果太暗或太亮怎么办用中值合成后可能会出现整体偏暗的问题尤其是目标时间段内影像数量较少时。这多半是因为薄云或者气溶胶的影响没有完全消除中值虽然减弱了异常值的影响但无法彻底清除。解决办法包括加入云概率波段作为额外条件进一步筛选比如要求CLOUD_PROBABILITY低于某个阈值改用“先选云量最少的一景然后只用这一景”的策略前提是该景影像在研究区内云覆盖确实少在合成前先用QA波段生成云掩膜然后用一次mosaic或者qualityMosaic替代median。我实际用下来最稳定的组合是filterMetadata粗筛云量小于30%然后map去云最后用qualityMosaic根据NDVI值选择每个像元的最优观测。这个方法在后面实操演示部分我会给出完整代码。3.4 Landsat 8/9去云和Sentinel-2的差异如果你的项目选的是Landsat去云方式和Sentinel-2完全不同。Collection 2 Level-2数据提供了一个QA_PIXEL波段它同样是位编码但每一位的含义不同。最常用的去云方法是判断第3位云和第4位云影function maskLandsat( image ) { var qa image.select(QA_PIXEL); // QA_PIXEL的第3位是云第4位是云影第5位是雪 var clouds qa.bitwiseAnd((1 3).or(1 4)).eq(0); return image.updateMask(clouds); }这里有一个容易忽略的点QA_PIXEL里第1位表示“云置信度高”如果只判断第1位会导致很多边缘薄云没有被过滤掉所以官方建议同时检查第3位和第4位。处理Landsat时还建议把辐射定标后的数据先转换成表面反射率这一步在Collection 2里已经默认做完了所以加载后直接就能用。4. 波段系数转换定标、反射率与指数换算4.1 三种“系数转换”到底指什么标题里的“波段系数转换”在不同语境下语义会不同我在这个章节把常见三种场景都讲一遍你大概率能在其中找到自己需要的部分。第一种是缩放系数转换。Sentinel-2的SR数据每个波段的像元值范围是0到10000而标准反射率范围是0到1。处理时要把像元值乘以0.0001才能得到反射率。Landsat Collection 2的SR数据已经自动应用了缩放系数但部分旧数据集没有需要手动乘。第二种是辐射定标系数转反射率。如果你加载的是原始数字量化值DN需要利用M常量、增益系数等做辐射定标把DN值转换为大气顶部反射率这一步常被称为TOA反射率转换。Landsat数据集的元数据中包含RADIANCE_MULT_BAND_x等系统参数在GEE里用multiply和add方法实现。第三种是波段组合运算时的理化意义系数。计算NDVI、NDWI、EVI等指数时不仅要用不同波段做差值和比值还要加上土壤调节系数、大气校正系数等常量参数。比如EVI计算公式里就包含L土壤调节系数、C1和C2大气修正系数这些系数是通过大量实验验证过的固定值在不同的软件实现中略有差别。这三种场景在GEE里处理逻辑不同下面分开讲。4.2 Sentinel-2波段缩放从整数到反射率Sentinel-2 Level-2A的SR波段在GEE里是以整数存储的取值范围为0到10000对应的物理意义是反射率乘以10000后的结果。如果想得到0到1的反射率需要除以10000也就是乘以0.0001。GEE官方推荐直接在指数计算时统一处理比如var ndvi composite.normalizedDifference([B8, B4]) .multiply(10000).int16();这里normalizedDifference计算出来的结果是-1到1之间的浮点数如果要存储成整型以节省空间可以用multiply(10000).int16()把范围扩大到-10000到10000。很多教程不解释这一步导致新手看到别人的代码里生成NDVI数据之后往往要多乘一个10000不知道是为什么。我自己更推荐反过来做先把原始影像乘以0.0001转成反射率之后再算所有指数都更直观var compositeRef composite.multiply(0.0001); print(反射率可视化参数:, compositeRef.select([B4,B3,B2]));转换成反射率后的优势在于你可以直接设定可视化范围为0到0.3真彩色合成的颜色会比较自然而不是一整个区域都是灰白色。而且后续做阈值分割比如提取水体时阈值也可以按照物理意义来设置比如水体在近红外波段的反射率通常低于0.05。4.3 Landsat定标系数与TOA反射率转换如果需要自己处理Landsat Level-1数据辐射定标这一步就躲不过。Landsat 8 OLI的每个热红外波段有一套定标参数都存储在影像元数据里。GEE中获取并应用这些参数的方法是var l8 ee.Image(LANDSAT/LC08/C01/T1_TOA/...); var multi l8.select(B4).multiply(l8.get(RADIANCE_MULT_BAND_4)); var offset multi.add(l8.get(RADIANCE_ADD_BAND_4));不过在实际应用中绝大多数场景直接用Collection 2的SR产品就好了GEE已经把定标和大气校正都处理完直接进入波段运算阶段。只有在研究算法开发或者需要原始DN值做辐射传输模拟时才需要走手动定标的流程。4.4 指数计算中的固定系数系数转换的另一个常见场景是计算EVI、SAVI这类带参数的指数。看一个EVI的标准公式EVI 2.5 × (NIR - Red) / (NIR 6 × Red - 7.5 × Blue 1)其中2.5、6、7.5、1都是经验系数。在GEE中实现如下var evi compositeRef.expression( 2.5 * ((NIR - RED) / (NIR 6 * RED - 7.5 * BLUE 1)), { NIR: compositeRef.select(B8), RED: compositeRef.select(B4), BLUE: compositeRef.select(B2) } );这里需要特别注意expression里做除法时分母可能为零或者极接近零会导致输出值出现极端异常。实际项目中我通常在做表达式计算前先对影像做一次clamp把反射率限制在合理范围内比如0.0001到0.9避免异常值干扰后续应用。5. 完整实操流程一个可直接复制的示例5.1 准备研究区和哨兵2影像下面这个示例整合了前面所有内容实现一个完整流程加载研究区、获取Sentinel-2影像、逐景去云、合成、转反射率、计算NDVI并导出。你可以把这个代码段当作一个模板替换研究区路径和时间范围后直接使用。// 1. 研究区 var roi ee.FeatureCollection(users/your_username/your_roi); // 2. 影像集合加载与粗筛 var s2 ee.ImageCollection(COPERNICUS/S2_SR) .filterBounds(roi) .filterDate(2023-01-01, 2023-12-31) .filterMetadata(CLOUDY_PIXEL_PERCENTAGE, less_than, 30); // 3. 去云函数 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.updateMask(mask); } // 4. 应用去云并选波段 var s2Masked s2 .map(maskS2clouds) .select([B2, B3, B4, B8]); // 5. 中值合成 var medianComposite s2Masked.median().clip(roi); // 6. 转反射率 var compositeRef medianComposite.multiply(0.0001); // 7. 计算NDVI var ndvi compositeRef.normalizedDifference([B8, B4]) .rename(NDVI); // 8. 可视化 Map.centerObject(roi, 10); Map.addLayer(compositeRef.select([B4, B3, B2]), {min: 0, max: 0.3}, True color); Map.addLayer(ndvi, {min: -0.2, max: 0.8, palette: [blue, white, green]}, NDVI); // 9. 导出 Export.image.toDrive({ image: ndvi, description: NDVI_2023, folder: GEE_exports, region: roi, scale: 10, crs: EPSG:32650, maxPixels: 1e10 });这个代码看起来不长但每一步都包含了前面讲的原理。运行起来后如果研究区不大几秒钟就能在Map面板看到结果。导出任务会在右侧的Tasks选项卡里出现点击Run后GEE云服务器会在后台处理处理完自动存到你的Google Drive。5.2 让合成效果更稳的qualityMosaic方案如果你发现median合成后的影像在部分区域还是不够干净试试qualityMosaic它允许在合成时根据某个波段的值挑选每个像元的最佳观测。比如我们希望优先选择NDVI最高的那次观测对于植被区域来说NDVI高通常意味着更少云干扰和更好的植被信号var withNdvi s2Masked.map(function(img) { var nd img.normalizedDifference([B8, B4]).rename(NDVI); return img.addBands(nd); }); var best withNdvi.qualityMosaic(NDVI).clip(roi);qualityMosaic的工作原理是对每个像元在所有影像中选指定波段值最高的那一景然后用这一景的所有波段覆盖该像元。这种合成方式在植被分析中比median更“锐利”细节更丰富但也可能引入单景影像的噪声需要自己对效果做对比评估。一个值得记住的经验是不要盲目追求复杂方案先跑一次简单median目视检查云遮挡情况。如果没问题就继续用效率最高。如果不行再升级到qualityMosaic。多数项目里filterMetadata粗筛加map去云加median合成已经够用了没必要一上来就上复杂的时空融合算法。5.3 导出时的投影、缩放与像素限制导出代码里有一个细节容易被忽略crs参数。我在示例里指定了EPSG:32650UTM zone 50N这是为了让NDVI影像以米制投影导出避免默认WGS84在赤道地区以外的变形问题。如果你不确定自己研究区的UTM分带可以在GEE里用roi.geometry().projection()查看默认投影或者手动算一下经度对应的带号公式是带号 floor((经度 180) / 6) 1。Export.image.toDrive函数的scale参数决定了导出分辨率。Sentinel-2最好用的波段是10米如果导出其他波段可能需要用20米或60米。如果你不确定当前图像的波段分辨率可以使用图像band的nominalScale方法查询。maxPixels限制导出的最大像元数小地块通常1e7就够大地块记得调大到1e10甚至1e11。6. 常见问题与排查技巧6.1 云边缘残留与掩膜空洞做去云处理时最常遇到的困惑是明明用了QA60掩膜结果影像里还是有一些白花花的地方。原因通常是薄卷云没有被QA60完全识别出来或者云影区域在QA60里被标记为“未分类”。解决办法有两条路一是试用SCL波段把云影类值3也加进掩膜条件二是放宽时间范围多取几景影像通过中值合成把这些残留区域在时间维上“磨掉”。SCL波段的优势是分类更细可以直接剔除值3云影、值4云、值8云边缘、值9薄卷云、值10高云但有的时候它会把部分高山阴影错误识别为云影。所以在山区做研究时我更倾向于继续用QA60配合qualityMosaic而不是直接用SCL一刀切。6.2 系数转换后NDVI范围异常有人会问我算完NDVI后取值范围是-5000到8000这是怎么回事这个问题九成出在波段系数没统一。如果原始波段是0到10000的整数直接做normalizedDifference结果会是-10000到10000这是正常的但很多后续分析代码默认NDVI范围是-1到1不转换就会出现阈值失效。解决办法就是在计算NDVI前先multiply(0.0001)转成标准反射率或者计算后multiply(0.0001)再clip(-1, 1)。连续性问题的另一种表现是NDVI全图都集中在0.9以上看起来很“假”。这通常说明合成影像里仍有厚云或雪覆盖导致近红外反射率极高。此时要回到掩膜步骤检查mask是否真的生效建议用一次addLayer展示mask波段的二值图。6.3 影像加载为空是什么原因filterDate后s2.size()显示0是新手经常会遇到情况。原因基本有三个时间范围设置太窄、研究区不在数据集覆盖范围内、云量元数据筛选条件过严导致所有影像都被过滤。逐个排查的时候我会先不写filterMetadata只保留filterBounds和filterDate跑通后再逐步加筛选条件。还有一个隐蔽问题是坐标系反了也就是经纬度坐标被写成了经度大于180的数值或者研究区矢量数据本身投影有问题导致filterBounds匹配不到任何影像。在GEE里可以先打印roi.geometry()的坐标值和影像集合的footprint做交集判断就能定位问题。6.4 性能优化减少中间影像、按需加载波段GEE有计算资源限制长时间运行的复杂链式操作可能在浏览器里报“User memory limit exceeded”。性能优化的核心思路是减少中间影像的数量和波段数。比如map去云时可以同时做波段选择只保留后续计算需要的波段比如不需要六十多个波段就在这一步删掉可以省掉很多不必要的内存占用。另一个常见优化点是避免在循环里对ImageCollection做逐个处理能向量化操作的优先用map不要用iterate。iterate效率低还容易出难排查的错误我一般只有在需要累积计算比如逐日累计时才会用它。再一个小技巧如果你只需要统计某个区域的NDVI平均值不需要导出整幅影像可以用reduceRegion直接输出统计值避免生成全分辨率影像占用大量资源。6.5 常见错误速查表症状大概率原因检查/解决方式合成影像有大块白色QA掩膜未正确应用或最薄云系数未覆盖单独可视化mask改用SCL波段NDVI范围异常波段系数未统一先multiply(0.0001)再计算影像集合为空时间/空间/云量筛选过严逐个移除条件排查导出时内存超限影像范围过大或波段过多提前裁剪、select波段、调大maxPixels可视化全灰数据范围设置不对用Map.addLayer中min/max匹配数据范围每年对应时间点结果突变时间窗口内影像数量不足适当放宽时间或补充Landsat数据这张表里的问题都是我处理多个项目时实际遇到的有些坑是在长时间连续分析时才会显现出来。比如“每年对应时间点结果突变”这个问题往往是某一年同期正好赶上连续的阴雨天影像数量不足导致合成影像质量骤降。解决方向就是加入Landsat 8/9做数据补充毕竟两者重访周期加在一起能提供的有效观测天数多不少。6.6 一个提升效率的小习惯建立常用函数库写了几个月的GEE代码后我逐渐意识到一个效率技巧把常用的掩膜函数、指数计算函数、可视化参数统一封装好。每次新项目开始直接加载已有的函数脚本几十行代码就能搞定以前需要写上百行的逻辑。比如把去云函数写成一个公共模块var cloudMask { s2: function(image) { ... }, l8: function(image) { ... }, l9: function(image) { ... } };用的时候直接调用cloudMask.s2(img)代码简洁不少也不容易在复制粘贴中出错。GEE的Scripts面板里可以创建仓库并保存这些公共函数长期积累下来会变成一个很趁手的个人工具箱。再提一个小技巧多用Map.addLayer中以不同palette来检验中间结果。比如mask图层用黑白配色、NDVI用红绿配色这样每次加图层时一眼就能看出哪个环节出了问题而不是等合成结果出来后才发现异常然后再从头排查。写在最后的经验这套“影像加载——去云——波段系数转换”的处理链路表面上只是GEE里的几个API组合但真正支撑它稳定运行的是对每个波段物理含义、每个掩膜位编码、每个系数来源的理解。我在刚开始用GEE的时候也经历过多次“跑通了代码但结果不对”的尴尬后来发现几乎所有问题都出在没有准确理解数据的元信息和波段属性上。如果你也想在自己项目里应用这套流程我从个人经验角度给几个建议第一先在小范围、影像数量少的区域把每一步结果可视化出来确认每一步输出都合理再放大范围批量跑第二不要盲目套用别人代码里的去云阈值和系数先查一查所用数据集的官方说明文档尤其是QA波段的位编码含义第三善用GEE的print函数把中间结果的数据类型、波段名、投影信息都打印出来很多隐藏问题会立刻现形。云端遥感处理给我们的帮助很大但它的核心依然建立在遥感科学的基本原理之上。把基础概念弄扎实了工具切换起来会从容很多。希望这篇文章里记录的思路和坑能帮你少走一些我走过的弯路。
📝

华诺云谱内容团队

资深建站顾问 · 行业研究员

10年+企业数字化服务经验,专注智能建站、SEO优化与品牌营销,持续输出建站技巧、行业洞察与营销干货,已帮助5000+企业实现数字化增长。

你可能需要的服务

订阅华诺云谱资讯周报

每周一封,精选建站技巧、SEO与营销干货,直达邮箱。已有 8,000+ 企业主订阅,助你少走弯路。

↑