简介面向遥感数据处理人员与相关专业学生这份操作流程文档系统梳理了使用ENVI软件完成Sentinel-2数据大气校正的完整方案。内容在ENVI 5.3与5.4版本中实测可用先讲大气校正原理与辐射率转换公式再说明L1C辐射率定标工具的两种安装方法包括应用商店安装与手动拷贝文件同时解释十米四波段、二十米九波段、六十米十一波段三种输出方案的波段构成和参数含义并补充FLAASH大气校正的关键设置与常见问题提示可帮助读者快速完成从原始影像到大气校正结果的流程减少参数选择的盲目性。资源包仅有一个docx文档大小约515KB图文形式呈现。目前已有九千七百四十二人浏览学习对需要独立完成哨兵二号影像预处理的在校学生和科研人员是一份实用的参考材料。1. 为什么 Sentinel-2 的大气校正一定要在 ENVI 里重走一遍做 Sentinel-2 数据处理的人应该都有过这种经历下载下来的 L1C 数据在 QGIS 里看挺好植被绿油油的水体也清亮但一拿去算 NDVI、做分类或者跟历史地表反射率产品对比数值就全对不上。原因很简单L1C 是表观反射率大气散射和吸收把真实的地表信号搅浑了尤其蓝光和红光波段气溶胶影响能占到 30% 以上。ENVI 的 FLAASH 模块是大家用得最多、资料最全、结果也最容易验证的一套大气校正方案网上能搜到的流程帖不少但很多只写到点一下 FLAASH 就行真正跑过的人会告诉你中间藏着大量细节坑——波段顺序、单位换算、中心波长、气溶胶模型每一样不对出来的图都是一张废图。这篇流程整理把我亲测能跑通的完整步骤和踩过的坑都写清楚适合正在用 ENVI 处理哨兵 2 数据、需要拿到可靠地表反射率的从业者照着复现。2. FLAASH 为什么是首选先搞懂辐射传输模型再动手2.1 FLAASH 和 Sen2Cor 的差异什么时候必须用 ENVI处理 Sentinel-2 大气校正大多数人第一个想到的是 Sen2Cor——欧空局官方提供的 L2A 处理器命令行一跑就能输出地表反射率。那为什么还要用 ENVI 的 FLAASH 再走一遍我的理由很实际Sen2Cor 是个黑匣子它输出的 L2A 产品质量不错但你很难干预中间过程比如气溶胶厚度反演失败时它直接给你一个默认值你不知道它什么时候反演失败、为什么失败。FLAASH 用的是 MODTRAN 辐射传输模型参数全部暴露在面板上每个设置都能控制出了异常结果你知道往哪个方向找原因。还有一类场景是必须用 ENVI 的当你需要把 Sentinel-2 跟其他传感器数据如 GF-5 高光谱、无人机多光谱做交叉定标或联合反演时各数据的处理链路必须一致不能一个用 Sen2Cor、一个用 FLAASH。ENVI 里能统一用 FLAASH 走完保证大气校正的算法内核和参数体系一致数据之间的可比性才有保障。如果你只是想要一张能用的 L2A 图不关心过程和参数Sen2Cor 够用但如果你要处理批量的哨兵 2 数据并且要做质量控制和参数调优FLAASH 更可控。2.2 从 L1C 到地表反射率要过哪几道换算这里先把大气校正的物理链路说清楚因为我发现很多人卡在中间某一环其实是没搞明白每一步在干什么。Sentinel-2 L1C 产品里存的原始值是 DN 值数字量化值它跟辐射亮度之间的关系是线性的辐射亮度 L DN × 缩放因子。这个缩放因子通常印在 .jp2 文件的元数据里哨兵 2 一般是 0.0001意味着原始 DN 值 10000 对应辐射亮度 1 W/m²·sr·μm。但如果直接用 ENVI 打开 L1C 的 .jp2 文件ENVI 会把这个缩放因子读进去界面上显示的数值已经是辐射亮度了——这是第一个容易混乱的地方因为很多教程说先做辐射定标再大气校正那是针对 Landsat 或其他卫星数据的老流程Sentinel-2 的 L1C 其实已经是定标过的辐亮度数据。大气校正要做的是把辐亮度数据里的大气贡献剥离掉得到一个地表反射率。FLAASH 做的事可以简化为一个公式地表反射率 (观测辐亮度 - 大气路径辐射) / (大气透过率 × 太阳辐射照度)。这里面每一项参数都由 MODTRAN 根据你设定的大气模型、气溶胶模型、观测几何来计算。所以 FLAASH 面板里的每一项都不是摆设——传感器高度、飞行时间、中心点经纬度、大气模型、气溶胶模型它们共同决定了一次完整的辐射传输模拟。理解了这一层你就知道为什么参数填错一点结果图看起来都是花的。3. 用 ENVI 读入并预处理 Sentinel-2 L1C辐射定标与格式转换3.1 读入 jp2 的正确姿势别直接 File Open最常见的翻车操作就是直接File Open把 L1C 的 .jp2 文件拖进 ENVI然后发现波段全挤在一起、没有坐标系信息、显示成灰蒙蒙一片。这不算错但会给后面 FLAASH 挖坑。正确做法是用 ENVI 自带的 Sentinel-2 专用读取器。在较新版本的 ENVI 里File Open弹出的文件选择窗口左下角有个 Open 按钮旁边的下拉箭头选择 Open Sentinel-2如果你用的是旧版 ENVI 5.x需要先安装对应的 Sentinel-2 读入插件或者直接用File Open As Optical Sensors Sentinel-2。选择 L1C 数据包里的MTD_MSIL1C.xml文件ENVI 会自动把 13 个波段按顺序读进来并带上对应的中心波长和坐标系。这里有个关键点ENVI 读入 Sentinel-2 时会把 B1海岸气溶胶波段和 B10卷云波段也读进来但 FLAASH 大气校正里它们基本不参与。FLAASH 要求输入数据是辐亮度单位且波长覆盖范围要在 0.4~2.5μm 之间所以后面做波段选择时我们一般保留 B2-B8A、B11、B12 这 9 个常用波段把 B1、B9、B10 去掉。B9水汽波段理论上可用于水汽反演但 FLAASH 对波段位置有严格匹配要求实操中经常因为它导致报错不如直接去掉省心。3.2 辐射定标与 Scale Factor把 DN 值变成辐亮度如果你读入 .jp2 后看到的数值范围在 0~10000 左右说明 ENVI 没有自动应用缩放因子此时显示的是原始 DN 值必须手动做辐射定标。如果看到的数值是 0~1 左右的小数说明缩放因子已经被读进去了可以直接进 FLAASH。判断标准很简单在 ENVI 的 Layer Manager 里右键点击数据层选 View Metadata找到 Data Ignore Value 和 Scaling Factor 两项。Sentinel-2 L1C 的 Scaling Factor 记录在MTD_TL.xml或MTD_MSIL1C.xml里的Radiometric_Scale_Value通常为 0.0001。ENVI 的专用读取器在多数版本里会自动应用这个值但如果你是通过.jp2单文件读入的这个缩放因子可能不生效。保险起见我一般直接在 FLAASH 之前先做一次统一的辐射定标操作用 ENVI 的Radiometric Calibration工具Toolbox 里搜 Radiometric Calibration选择读入的哨兵 2 数据把 Scale Factor 设置为 0.0001或者从元数据里查到的值输出为浮点型Float得到辐射亮度数据。这一步做完了Layer Manager 里的数据单位应该显示为 µW/(cm²·sr·nm) 或者 W/(m²·sr·μm)具体单位取决于 ENVI 设置不统一也没关系FLAASH 会自动做单位换算。# 如果你偏好用 ENVI Modeler 或 IDL 批处理辐射定标这一步可以写成这样 ; 假设你已经在 ENVI 中打开了读取后的 Sentinel-2 数据获取其 raster raster E.OpenRaster(D:/S2_L1C/MTD_MSIL1C.xml) ; 创建辐射定标任务 task ENVITask(RadiometricCalibration) task.INPUT_RASTER raster task.CALIBRATION_TYPE Top-of-Atmosphere Reflectance task.OUTPUT_RASTER_URI D:/S2_L1C/calibrated.dat task.Execute ; 注意这里的 CALIBRATION_TYPE 如果选 TOA Reflectance得到的是反射率而非辐亮度 ; FLAASH 要求输入辐亮度因此应使用 ENVITask(ConvertToRadiance) 或直接让 FLAASH 处理。实际经验是FLAASH 面板里有一个 Use Raster Resolution 和 Set Scale Factors 的选项在 FLAASH 弹出的 Input Radiance Image 选择窗口里会自动弹出 Radiance Scale Factors 对话框默认是 1.0。如果你之前已经做过辐射定标并生成辐亮度数据这里保持 1.0 不变如果你是直接从 L1C 读入、数值范围 0~10000这里要填 0.001 还是 0.0001看你的数据单位。最简单的办法是先做辐射定标再进 FLAASH避免在 FLAASH 面板里反复试错。3.3 波段顺序、中心波长和 FWHMFLAASH 能不能跑通的命门FLAASH 对输入数据的波段顺序没有硬性要求但要求你告诉它每个波段的中心波长和半高宽FWHM而且波长单位必须是 nm微米也行但要跟面板里选的一致。这一步是整条流程里最容易出鬼的地方因为 Sentinel-2 有多个分辨率级别10m/20m/60m 波段的像素尺寸不同直接混合输入 FLAASH 会出现网格不对齐的问题。我的做法是在读入 .jp2 之后用 ENVI 的Layer Stacking工具把所有需要的波段叠成一个多波段文件分辨率统一设为 20mB5、B6、B7、B8A、B11、B12 本来就是 20mB2-B4、B8 是 10m重采样到 20m 信息损失在接受范围内。这样做的好处是 FLAASH 不会因为像素尺寸不一致而报错处理速度也快不少。注意Layer Stacking里一定要勾选 Output Metadata 或保持默认让中心波长信息跟着走。中心波长和 FWHM 在 FLAASH 面板里怎么填取决于你读入的数据有没有把 metadata 里的波长信息带进来。推荐的做法是在 FLAASH 的 Wavelength 设置里选择 Read array from file指向你从元数据导出的哨兵 2 波段列表文件。这个文件格式是 ENVI 标准波长文件两列中心波长nm和满宽半高值nm9 个波段就是 9 行。或者直接在 FLAASH 打开的波长设置表格里手动敲进去表格长这样波段序列中心波长nmFWHMnmB2蓝49065B3绿56035B4红66530B5红边 170515B6红边 274015B7红边 378320B8近红外842115B8A窄近红外86520B11短波红外 1161090B12短波红外 22190180这里有个很多人忽略的坑FLAASH 提供的通用传感器类型列表里有一个 Sentinel-2 的预设但你选了这个预设之后它自动填充的波段列表可能跟你手头数据的波段顺序不一致——尤其是你把 10m 和 20m 波段混合叠加后波段顺序往往跟预设对不上。每次改波段顺序FLAASH 面板里的波长表不会自动同步更新很容易出现错位。我习惯每次新建 FLAASH 任务时都手动核对一遍波长表确认中心波长跟波段一一对应这件事值得花两分钟检查后面省两小时的返工时间。4. FLAASH 大气校正参数逐个设置照着填就能出图4.1 必填参数中心点经纬度、传感器类型、飞行时间FLAASH 面板打开后Toolbox Radiometric Correction FLAASH Atmospheric Correction第一步是选择输入文件就是上一步生成的辐亮度多波段文件。选完之后弹出 Radiance Scale Factors按前面说的保持 1.0。接着进入 FLAASH 主面板有几个参数必须填对。Scene Center Location场景中心经纬度可以从 Sentinel-2 元数据里读打开MTD_MSIL1C.xml搜Geoposition里面有LAT和LON注意单位是十进制度数WGS84 坐标系。直接填进去就行。这个参数决定了 MODTRAN 计算大气透过率用的太阳-目标几何填错了整幅图的校正结果都偏。Sensor Type传感器类型选择 Sentinel-2 的预设如果版本里有没有的话选Generic Multispectral并把 Sensor Altitude 填成 786飞行时间填数据获取日期。这里要注意FLAASH 的 Sensor Altitude 单位是千米哨兵 2 的轨道高度是 786km填 786 就行不是 786000。Ground Elevation地面平均高程对地形起伏大的区域很重要。FLAASH 默认是 0.5km如果你处理的是青藏高原或者云贵高原的数据把它改成实际的区域平均海拔单位 km比如青藏高原填 4.0 左右。这个参数直接影响大气路径辐射的计算平原地区填 0~0.5 都无所谓高海拔地区不填对蓝光波段会校正过头。4.2 气溶胶模型与大气模型怎么选按季节和纬度查表FLAASH 面板里的 Atmospheric Model 和 Aerosol Model 是两个下拉框很多人随便选其实它们决定了 MODTRAN 使用的标准大气廓线和气溶胶类型。选错的结果是植被光谱曲线形状看着对但绝对值偏差 15%~30%。大气模型按纬度和季节查表这是 FLAASH 官方推荐的做法我做了个简化表格纬度范围夏季5-10月冬季11-4月0-20°热带TropicalTropical20-45°亚热带Mid-Latitude SummerMid-Latitude Winter45-70°中高纬Mid-Latitude SummerSub-Arctic Winter70° 以上Sub-Arctic SummerSub-Arctic Winter对于中国大部分区域华北、华东、华南夏季选 Mid-Latitude Summer冬季选 Mid-Latitude Winter东北北部和新疆北部冬季可以选 Sub-Arctic Winter。这个按月份边界是可微调的我的经验是 4 月和 10 月这种换季月份如果当天温度明显偏高或偏低可以手动微调但影响不大按表格选就行。气溶胶模型分四种Rural乡村、Urban城市、Maritime海洋、Tropospheric对流层。城市区域选 Urban郊区或植被覆盖高选 Rural沿海城市有时候选 Maritime 反而更准因为海盐气溶胶占主导。这个选择影响气溶胶光学厚度的反演如果你不确定选 Rural 是最保守的结果不会太离谱。初始能见度Initial Visibility这个参数也值得讲一下。FLAASH 用初始能见度作为气溶胶反演的起点默认是 40km。如果你处理的场景当天有霾或者轻微污染把它降到 15~20km如果雨后空气特别干净可以填 50km 以上。这个参数选多少不是关键关键是你得给 FLAASH 一个合理的起点它自己会通过暗像元算法迭代修正。填 40 一般都能收敛。4.3 水汽反演与 865nm 波段设置了不等于能用水汽反演Water Retrieval是 FLAASH 里一个高级选项它利用 820nm 附近的水汽吸收带特征来反演大气水汽含量。问题是 Sentinel-2 的波段设置里没有标准的 820nm 水汽吸收通道只有一个 865nm 的窄带B8A它在水汽吸收窗口的边缘信噪比不够FLAASH 用它做水汽反演容易失败或者反演出离谱的负值。我在实际操作中直接关掉水汽反演Water Retrieval 选 No然后让 FLAASH 使用默认的 1.5g/cm² 水汽含量这个值在 MODTRAN 的 Mid-Latitude 模型中会结合大气模型自动调整。对照过用探空数据做的辐射传输模拟结果地表反射率在可见光波段的偏差小于 5%完全在可接受范围内。如果你非要开水汽反演必须保证输入波段里有 865nm 且传感器类型预设正确然后仔细检查反演输出的水汽分布图——如果出现斑块状异常说明反演挂了赶紧切回 No。K-T 下垫面Kaufman-Tanre的 Aerosol Retrieval 也建议选上它使用暗像元算法反演气溶胶光学厚度。但请注意如果你的场景里水体太多比如处理湖区或沿海地区暗像元选择会出错这时候把 K-T 下垫面关掉改用固定能见度。还有一种情况是冬季北方植被枯黄没有足够多的暗像元也会失败同样的处理关 K-T给初始能见度 25km 左右出来的结果可接受。所有参数填完后Output Reflectance 勾上分块处理Tiles建议设置成 1 个或按内存大小自动分块。点击 Apply 开始跑一个标准的 10980×10980 像素 L1C 场景10m 全分辨率在 32GB 内存、i7 处理器上要跑 20~40 分钟属于正常速度。5. Sentinel-2 大气校正避坑实录现象、原因、解决5.1 波谱曲线全是锯齿中心波长没有按哨兵 2 实际波段填现象大气校正跑完了用 Z Profile 看一个植被像元的光谱曲线形状大体对但每个波段之间出现不自然的跳变曲线像锯齿状尤其红边和近红外过渡处。原因FLAASH 面板里的波长列表和实际输入文件的波段顺序错位。最常见的情况是把 B5、B6、B7红边波段的顺序填反了或者波段文件里实际排布和波长表差了 1~2 个位置。解决回到 3.3 节的波长表对照输入文件波段顺序在 FLAASH 面板里逐行核对而不是只看一眼就关掉。如果你用 ENVI 的 Layer Stacking 叠过波段Metadata 里的波长信息可能没同步进来建议直接用Edit Edit Metadata打开数据的波长列表手动粘贴 9 个波长的值。这个检查花 3 分钟能省掉一上午的重跑。5.2 大气校正后植被近红外反射率反而变低现象校正前 L1C 的近红外反射率看起来有 0.45FLAASH 校正完变成 0.30比植被正常值还低而且 NDVI 也跟着掉。原因大概率是大气模型选错了常见的是冬天选了 Summer 模型。还有一种情况是气溶胶模型选了 Maritime但场景是内陆农田海盐气溶胶的光学特性跟实际的大陆型气溶胶差异大导致 MODTRAN 计算的大气透过率偏大校正时把地表信号压过头了。解决按 4.2 的表格重新选大气模型气溶胶模型改 Rural 或 Urban重跑一次。这里提醒一句FLAASH 校正出来的反射率应该比 L1C 的表观反射率在近红外波段偏高一点因为去掉了大气吸收如果结果反而偏低 15% 以上说明参数有问题不用怀疑数据本身。5.3 边缘波段出现负值或黑白点现象校正完的 B2蓝光和 B12SWIR 2波段的图像上大片像元的数值是负的或等于 0图像显示成黑色和白色雪花点。原因有两个常见来源。第一是元数据里的坏像元比如传感器边缘的探测器异常没被掩膜掉DN 值本身是异常的大气校正后变成了负反射率。第二是 FLAASH 里输入数据的 Scale Factor 填错辐亮度被放大或缩小了 10 倍超过了 MODTRAN 的计算范围。解决先检查原始 L1C 数据在对应位置是不是有 0 值或异常值如果有用 ENVI 的Mask工具先做掩膜或者干脆只做有效区域的子集。检查 FLAASH 面板里 Data Ignore Value 是否读入正确哨兵 2 的无值像元一般是 0这个要勾选 Apply data ignore value。如果原始数据没问题而校正后依然有负值大概率是 Scale Factor 填大了 10 倍——改小后重跑。5.4 K 分布渲染报错内存与分块设置现象FLAASH 跑了几分钟后弹出一个错误框提示 Invalid Tilename 或者 Tile size is too large 之类的信息进程直接中断。原因FLAASH 的 K 分布方法对内存占用比较敏感当你把 Tiles 数设成 1全图一块处理而数据又特别大时内存不够就会崩。另外如果输入文件路径里有中文字符或空格某些版本会出现 Tile 文件命名失败的问题。解决把 Tiles 改成 4 或 6让 FLAASH 分块处理。输入文件路径改成纯英文且不带空格如D:/S2_data/corrected.dat。如果还是崩把输入文件裁剪成一个测试子区500×500 像元左右先跑通验证参数没问题再全图重跑。5.5 结果图出现横向条纹现象大气校正完的反射率影像在显示时出现横向的明暗条纹跟在原始影像上叠加了一层百叶窗一样。这在单张无人机影像去条纹的案例里也常见原理类似。原因哨兵 2 的 L1C 产品在探测器拼接处存在系统性辐射不均匀性大气校正本身不会新产生条纹但不均匀性在辐亮度转反射率的过程中会被放大。另一个原因是 FLAASH 的分块处理在不同块之间使用了略有差异的初始参数导致块边界出现带状差异。解决如果条纹存在于块边界调整 Tiles 设置或勾选 Spectral Polishing 和 Spatial Smoothing这两个选项在 FLAASH 面板上可勾选默认关闭它们会对光谱维和空间维做轻微平滑能够有效压制块状差异。如果条纹存在于固定位置属于传感器自身的条纹噪声需要在 FLAASH 之后用 ENVI 的Striping Removal工具做后期修复Tools Image Processing Striping Removal参数一般用默认的 Moving Window 就能处理干净。6. 结果验证与批量处理肉眼过关不算过关6.1 用典型地物光谱剖面验证大气校正跑完后先别急着出图。一个能不用依靠外部实测数据、快速判断校正结果是否有效的办法是看典型地物的光谱形态植被的红边效应B5-B7 之间反射率陡升、水体在近红外波段的低反射B8 应该接近 0、裸土在 SWIR 波段的平坦光谱特征。在 ENVI 里打开校正结果用Spectral Profile工具分别点几个明显是植被、水体、裸土的区域对照这些特征做判断。还有一个量化指标是 NDVI 的范围。校正后的植被 NDVI 一般在 0.6~0.9 之间水体在 -0.3~0.1 之间裸土在 0.1~0.3 之间。如果整幅 NDVI 的方差很小比如全部集中在 0.1 左右说明校正效果不好回看参数设置是不是出了问题。另外有条件可以和 ESA 的官方 L2A 产品做对比——从 Copernicus Open Access Hub 下载同一个时间、同一景的 L2A 数据在 ENVI 里重采样到相同分辨率逐像元计算差值。正常情况下 FLAASH 和 Sen2Cor 的结果差异应该在 ±0.05 反射率以内如果差异超过 0.1肯定是某个参数填错了。6.2 批量大气校正的脚本骨架处理一个季度的哨兵 2 数据一景一景地打开 ENVI 点 FLAASH 面板操作量很大。用 ENVI Modeler 把流程串起来是更省力的做法。这里给一个批处理脚本的框架基于 ENVI 的ENVITask接口原理上就是把前面手动操作的步骤自动化# 这是 ENVI 的 Python API 批处理框架在 ENVI 的 Python 环境里运行 import os from envy import envipy from ENVI import ENVITask # 设置输入和输出路径 input_folder rD:\S2_L1C\scenes output_folder rD:\S2_L1C\corrected # 遍历所有 L1C 数据包 for scene_folder in os.listdir(input_folder): xml_file os.path.join(input_folder, scene_folder, MTD_MSIL1C.xml) # Step 1: 读入数据 raster envi.open_raster(xml_file) # Step 2: 辐射定标如果数据还是 DN 值 cal_task ENVITask(RadiometricCalibration) cal_task.input_raster raster cal_task.scale_factor 0.0001 cal_task.output_raster_uri os.path.join(output_folder, scene_folder _rad.dat) cal_task.execute() # Step 3: 大气校正FLAASH flaash_task ENVITask(FLAASH) flaash_task.input_raster cal_task.output_raster flaash_task.scene_center_location [30.5, 114.3] # 从 MTD 文件中读取 flaash_task.sensor_type Sentinel-2 flaash_task.sensor_altitude 786.0 flaash_task.atmospheric_model Mid-Latitude Summer flaash_task.aerosol_model Rural flaash_task.initial_visibility 40.0 flaash_task.water_retrieval False flaash_task.output_raster_uri os.path.join(output_folder, scene_folder _sr.dat) flaash_task.execute()脚本里的参数需要根据实际数据逐景修改中心点经纬度、飞行时间、大气模型是必须动态更新的这三项错了整批数据全部作废。我的习惯是把场景信息先用一个 CSV 文件管理起来每一景一行记录场景名、中心经纬度、获取日期、对应的大气模型脚本逐行读取再跑这样既不会漏参出问题也好定位是哪一景的参数写错了。批量跑完之后和 L2A 官方产品做一次抽检对比是值得投入的时间——抽一景比较植被/水体典型像元的反射率差值。如果差异在可接受范围就可以放心地批量使用这套处理结果作为后续分类、反演或变化检测的数据源。这也是我踩过不少坑之后固化下来的习惯新拿到一批数据不要急着全量处理先挑一景跑通、验证确认参数合理后再放量跑。前面省下的二十分钟后面会以几倍的时间还回来。这个流程我跑了两年最大的体会是大气校正别迷信任何官方默认参数每个场景都要自己过一遍参数表再好的工具也需要人确认它算出来的东西是不是符合物理规律。希望这篇流程整理能帮你少走一些我走过的弯路。本文还有配套的精品资源点击获取