Sentinel-1 SAR数据入门到DInSAR实战：SNAP处理全流程与避坑指南

1. 哨兵1号SAR卫星到底能做什么第一次接触Sentinel-1数据的人十有八九是被它那种“黑白颗粒感”的影像劝退的。没有光学卫星那种色彩鲜艳、一眼就能看懂的视觉效果满屏都是灰白噪点看起来像老式电视机没信号的画面。但恰恰是这种“不好看”的数据在洪水监测、地表形变分析、海冰识别、农作物估产这些领域里干着光学卫星根本干不了的活。哨兵1号Sentinel-1是欧洲空间局哥白尼计划里的一颗C波段合成孔径雷达卫星分1A和1B两颗1B在2021年12月因电源故障退役目前主要靠1A和后续的1C接力。它最大的特点就是主动微波成像——自己发射雷达波再接收地面反射回来的信号。这意味着它不依赖太阳光照白天黑夜都能拍微波能穿透云层阴雨天照样出图。对于常年多云多雨的南方地区或者需要应急响应的洪涝灾害场景这个能力就是刚需。这颗卫星能提供的数据产品大致分几类原始RAW数据、SLC单视复数数据、GRD地距探测数据还有OCN海洋产品。SLC保留了相位信息是做InSAR形变监测的核心数据GRD已经做了多视和投影处理适合做强度影像分析和分类OCN则专门服务于海洋风场、波浪谱这些应用。不同产品对应不同处理链路选错了数据类型后面所有工作都白搭。这篇文章面向的是想上手Sentinel-1数据、但被各种参数和处理流程绕晕的人。不管你是遥感专业的学生还是做地灾监测的工程师或者只是对SAR影像好奇的技术爱好者我都会从数据获取、参数理解、处理流程到常见坑点把这条链路讲清楚。尤其是InSAR处理这块我会以SNAP软件为操作平台把DInSAR的实战流程拆到每一步都能照着做的程度。2. 搞懂Sentinel-1数据之前必须弄明白的几个核心概念2.1 波段、极化与成像模式选数据的三个关键维度Sentinel-1工作在C波段波长大约5.6厘米。这个波段的选择是有讲究的——波长太长比如L波段穿透力强但空间分辨率难做高波长太短比如X波段分辨率高但对植被和大气太敏感。C波段算是一个平衡点既能保证一定的穿透能力对植被冠层有一定穿透但不至于完全穿透又能做到米级到十米级的分辨率。极化方式是另一个必须搞清楚的参数。Sentinel-1支持单极化HH或VV和双极化HHHV或VVVH。简单说HH是水平发射水平接收VV是垂直发射垂直接收HV和VH是交叉极化。不同地物对不同极化组合的响应完全不同VV对地表粗糙度和含水量敏感适合做土壤湿度和洪水监测HH对垂直结构敏感适合做海冰和城市区域分析交叉极化HV或VH则对植被体积散射特别敏感做森林生物量估算时经常用。成像模式直接决定了你拿到的数据覆盖范围和分辨率。Sentinel-1有四种主要模式成像模式幅宽分辨率极化方式典型用途SM条带80 km5×5 mHHHV/VVVH小区域精细观测IW干涉宽幅250 km5×20 mHHHV/VVVH陆地默认模式InSAR主力EW超宽幅400 km20×40 mHHHV/VVVH海洋和极地大范围监测WV波浪20×20 km5×5 mVV海浪谱和风场IW模式是绝大多数陆地应用的首选也是InSAR处理最常用的数据源。它用TOPS渐进扫描合成孔径雷达技术通过天线在方位向扫描来获取宽幅影像代价是每个子条带的方位向带宽变窄对配准精度要求极高。这也是为什么做InSAR时IW数据的配准必须用专门的Enhanced Spectral Diversity方法普通配准根本不够用。2.2 SLC、GRD、OCN三种产品到底该选哪个拿到Sentinel-1数据下载页面时新手最容易懵的就是产品类型选择。我见过太多人下载了GRD数据想做InSAR结果发现根本没有相位信息白忙一场。SLCSingle Look Complex是单视复数产品每个像素包含实部和虚部保留了完整的相位和幅度信息。这是做InSAR、PS-InSAR、SBAS这些形变监测应用的唯一选择。SLC数据没有做多视处理所以看起来会有很多斑点噪声空间分辨率是最高的但数据量也最大。一个IW模式的SLC数据大概4-8 GB处理起来对内存和算力都有要求。GRDGround Range Detected是地距探测产品已经做了多视处理和地距投影只保留幅度信息没有相位。数据量小很多大概1-2 GB适合做强度影像分析、变化检测、分类这些不需要相位的应用。GRD又分高分辨率和中分辨率两种高分辨率保留更多细节中分辨率数据量更小。OCN是海洋产品专门为海洋应用设计包含风场、波浪谱、洋流等信息。如果你不做海洋相关的研究基本用不到。选数据的原则很简单要做形变监测和干涉处理必须下SLC只做强度分析和目视解译GRD就够了海洋应用才考虑OCN。2.3 轨道参数与重访周期时间基线怎么控制Sentinel-1的重访周期是12天单颗卫星两颗卫星同时运行时可以缩短到6天。但实际做InSAR时时间基线两景影像的获取时间差直接决定了干涉相干性的好坏。时间基线越长地表变化越大相干性越差干涉条纹越模糊。一般来说做DInSAR形变监测时间基线最好控制在12-24天以内最长不要超过48天。超过这个范围植被生长、土壤湿度变化、甚至大气条件差异都会导致相干性严重下降。做PS-InSAR这种时序分析时因为要选大量影像做长时间序列时间基线可以放宽到几年但需要足够多的影像数量来保证可靠性。轨道参数方面Sentinel-1提供精密轨道星历POD和恢复轨道星历RES。POD精度更高厘米级但发布时间有延迟大概20天左右RES精度稍低分米级但可以快速获取。做高精度形变监测时一定要等POD轨道出来再处理否则轨道误差会直接混入形变信号里。3. Sentinel-1数据获取与预处理实操3.1 从哪下载数据几种常用渠道对比下载Sentinel-1数据有几个主要渠道各有优劣。最官方的是欧空局的Copernicus Open Access Hub数据最全但下载速度有时候不太稳定尤其是国内访问时。另一个常用的是ASFAlaska Satellite Facility的数据门户下载速度相对好一些界面也更友好。还有PEPS平台主要提供欧洲区域的数据。我个人的习惯是先用ASF的Vertex界面搜索和预览确定好数据ID后再去下载。ASF支持按区域、时间、模式、极化方式等多条件筛选还能直接在地图上画框选区域比欧空局那个老界面好用不少。下载之前有几个参数一定要确认清楚产品类型SLC还是GRD、成像模式IW还是EW、极化方式VV还是VH、时间范围、以及轨道方向升轨还是降轨。升轨和降轨的几何不同对同一地区的形变敏感方向也不同做形变分析时经常需要两个方向的数据配合使用。3.2 用SNAP做基础预处理轨道校正、辐射定标、斑点滤波SNAPSentinel Application Platform是欧空局官方推出的遥感处理软件免费开源对Sentinel系列数据支持最好。虽然界面看起来有点老旧但功能确实强大尤其是InSAR处理这块基本是标配工具。基础预处理的流程一般是Apply Orbit File应用精密轨道→ Radiometric Calibration辐射定标→ Speckle Filtering斑点滤波→ Terrain Correction地形校正。如果是做InSAR中间还要加配准和干涉步骤这个后面单独讲。Apply Orbit File这一步很多人会跳过觉得没必要。但实测下来用精密轨道和用预测轨道干涉图的相位质量差别很明显。尤其是做长时间序列分析时轨道误差会累积最后形变结果里混入的轨道残差可能比真实形变还大。辐射定标是把DN值转成后向散射系数σ0或γ0这一步对强度分析是必须的对InSAR处理其实可以跳过因为干涉用的是相位信息。但如果你后面要做相干性计算或者幅度图显示定标还是做一下比较好。斑点滤波是SAR处理里绕不开的话题。SAR影像天生就有斑点噪声这是相干成像系统的固有特性。滤波方法有很多Lee滤波、Frost滤波、Gamma MAP滤波、Refined Lee滤波等等。做InSAR时一般用Refined Lee滤波它在保持边缘和抑制噪声之间平衡得比较好。但要注意滤波要在干涉图生成之后做不能提前对SLC数据滤波否则会破坏相位信息。地形校正用Range-Doppler Terrain Correction需要DEM数据。SNAP自带SRTM DEM也可以自己提供更高精度的DEM。地形校正后影像就投影到地理坐标系下了可以和光学影像、GIS数据叠加分析。3.3 数据导入与格式转换的注意事项SNAP支持直接打开Sentinel-1的.zip压缩包不需要提前解压。但如果你要用其他软件处理比如GAMMA、SARscape可能需要先转成ENVI格式或者GeoTIFF。SNAP的Export功能支持多种格式输出但要注意导出SLC数据时一定要选复数格式Complex否则相位信息就丢了。还有一个坑是数据命名。Sentinel-1的产品命名规则很长包含了卫星平台、模式、产品类型、时间、轨道号、极化方式等信息。比如S1A_IW_SLC__1SDV_20230101T120000_20230101T120030_046000_058000_ABCD这里面S1A是卫星IW是模式SLC是产品类型1SDV表示双极化VVVH后面是时间戳和轨道号。搞清楚命名规则批量处理时不容易搞混。4. DInSAR处理全流程从SLC到形变图4.1 配准InSAR处理里最容易翻车的一步配准是InSAR处理里最基础也最关键的一步。配准精度不够干涉相位里就会混入地形相位残差最后形变图上一片乱七八糟的条纹根本没法解释。Sentinel-1的IW模式用TOPS技术每个子条带在方位向是扫描成像的不同子条带之间的配准精度要求极高——方位向配准误差要控制在千分之一像素以内。普通配准方法根本达不到这个精度必须用SNAP里的Enhanced Spectral DiversityESD方法。具体操作是先用Apply Orbit File把两景影像的轨道信息更新然后用Back Geocoding做粗配准接着用Enhanced Spectral Diversity做精配准。ESD的原理是利用重叠区域的干涉相位来估计方位向偏移量通过迭代优化把配准精度提到千分之一像素级别。这一步做完可以看一下配准报告里的方位向和距离向偏移量正常情况下应该在0.001像素量级。配准没做好后面全白搞。我见过太多人干涉图上一堆平行条纹以为是形变信号其实是配准误差导致的地形相位残差。4.2 干涉图生成与去平地效应配准完成后下一步是生成干涉图。SNAP里的Interferogram Formation工具会同时完成干涉图生成和去平地效应。干涉图就是两景影像的相位差平地效应是因为地球曲率导致的相位变化去平地效应就是把这块相位去掉让干涉条纹只反映地形和形变信息。干涉图生成后会同时输出相干性图。相干性是衡量干涉质量的重要指标值域在0到1之间。相干性高于0.3的区域干涉相位才比较可靠低于0.3的区域基本就是噪声做解缠时也会被掩膜掉。植被覆盖区、水体、农田这些地方相干性通常很低这是正常的。去平地效应之后干涉图里还会包含地形相位。如果只是做差分干涉DInSAR需要用外部DEM模拟地形相位然后从干涉图里减掉。SNAP里的Topographic Phase Removal工具就是干这个的需要提供DEM数据。SRTM 1弧秒DEM是常用的选择覆盖全球精度也够用。4.3 滤波与解缠让相位变得可用干涉图里的噪声主要来自两个方面一是斑点噪声二是失相干噪声。滤波的目的就是压制这些噪声让干涉条纹更清晰。SNAP里常用的滤波方法是Goldstein滤波和Boxcar滤波。Goldstein滤波对条纹的增强效果更好但参数需要调Boxcar滤波简单粗暴效果一般但稳定。滤波之后就是相位解缠。解缠是把缠绕在[-π, π]之间的相位恢复成连续相位的过程。这是InSAR处理里算法最复杂的一步也是误差最容易传播的一步。SNAP里提供了几种解缠方法Snaphu、SNAPHU外部调用、以及自己的解缠算法。Snaphu是最常用的效果好但速度慢SNAP自带的解缠速度快但精度稍差。解缠前一定要做相干性掩膜把相干性低于阈值的区域去掉。阈值一般设0.2-0.3具体看数据质量。掩膜没做好解缠误差会像瘟疫一样蔓延到整个干涉图。4.4 相位转形变与地理编码解缠后的相位是雷达视线向LOS的相位差要转成形变量需要做相位转形变。公式是形变量 (相位 × 波长) / (4π)。C波段波长5.6厘米所以一个2π的相位变化对应大约2.8厘米的LOS向形变。这里要注意DInSAR得到的是LOS向形变不是水平或垂直形变。LOS向是卫星视线方向跟垂直方向有一个夹角入射角所以LOS向形变是垂直形变和水平形变在视线方向的投影。如果要得到垂直形变需要假设水平形变为零或者用升降轨数据联合解算。最后一步是地理编码把雷达坐标系下的形变图投影到地理坐标系。SNAP里的Range-Doppler Terrain Correction可以同时完成地理编码和地形校正。地理编码后形变图就可以和GIS数据叠加做进一步分析了。5. 常见问题与排查技巧实录5.1 干涉图全是噪声怎么办这是最常见的问题原因可能有很多。先看相干性图如果整体相干性都很低那可能是时间基线太长或者地表变化太大换一景时间更近的影像试试。如果相干性局部低那可能是那个区域本身是植被或水体这是正常的。如果相干性看起来还行但干涉图还是乱检查配准精度。配准报告里的偏移量如果超过0.001像素重新做ESD配准。还有可能是轨道误差试试用POD精密轨道重新处理。5.2 解缠结果出现大面积跳变解缠跳变通常是因为相干性掩膜没做好或者解缠参数设置不当。先检查相干性掩膜阈值太低会导致噪声区域参与解缠太高会丢失有效信号。Snaphu解缠时可以调整--mtf参数最小相干性阈值一般设0.2-0.3。还有一种可能是相位不连续比如形变梯度太大一个条纹周期内相位变化超过π解缠算法就懵了。这种情况需要先做相位解缠的预处理比如用最小费用流算法或者增加解缠的迭代次数。5.3 形变结果里出现规律性条纹规律性条纹大概率是轨道误差或者DEM误差。轨道误差导致的条纹通常是长波长的覆盖整个影像DEM误差导致的条纹跟地形相关在山区特别明显。区分方法如果条纹跟地形等高线平行那就是DEM误差如果条纹是全局性的那就是轨道误差。轨道误差可以用轨道精化Orbit Refinement来校正SNAP里有对应的工具。DEM误差可以通过迭代去趋势来减弱或者用更高精度的DEM。5.4 处理速度太慢怎么优化SNAP处理大区域SLC数据时确实很吃资源。几个优化方向一是用子区域处理只处理感兴趣的区域不要整景处理二是调整SNAP的JVM内存参数在snap.conf里把-Xmx调大建议至少8GB三是用命令行模式gpt批量处理比GUI快不少四是考虑用云平台或者高性能计算资源。常见问题可能原因排查方法解决思路干涉图全噪声时间基线长/失相干查看相干性图换短时间基线数据解缠跳变掩膜不当/相位不连续检查掩膜阈值调整阈值或换解缠算法规律条纹轨道误差/DEM误差对比地形等高线轨道精化/换DEM处理慢内存不足/数据量大查看内存占用子区域处理/调JVM参数6. 光学与SAR协同云去除和变化检测的进阶玩法6.1 光学影像云去除SAR来补位光学卫星最大的软肋就是云。南方地区一年到头能用的光学影像没几景这时候SAR的优势就体现出来了。但SAR影像的目视解译性不如光学影像所以一个常见的做法是用SAR数据来填补光学影像被云遮挡的区域。具体思路是先用光学影像做分类或变化检测识别出被云遮挡的区域然后用同一时间的SAR影像通过强度分析或极化分解提取被遮挡区域的地物信息最后把SAR提取的信息融合到光学影像的分类结果里。这种方法在土地利用监测和灾害评估里很实用。6.2 变化检测SAR强度影像的用武之地SAR强度影像做变化检测最常用的方法是比值法和差值法。比值法是把两景影像的后向散射系数相除取对数后得到变化图差值法是直接相减。比值法对辐射定标误差更鲁棒所以更常用。变化检测的关键是阈值确定。阈值太高会漏检太低会误检。常用的阈值确定方法有双峰法假设变化和未变化区域在直方图上形成两个峰、Otsu法自动阈值分割、以及基于统计的K-S检验。实际应用中往往需要结合目视解译来调整阈值。6.3 时序InSAR从DInSAR到PS-InSAR和SBASDInSAR只能得到两景影像之间的形变而且受大气延迟和失相干影响很大。要做长时间序列的形变监测就需要PS-InSAR或SBAS这些时序InSAR方法。PS-InSAR永久散射体InSAR的思路是在多景影像里找出那些相位稳定的点比如建筑物、岩石这些硬目标然后对这些点做时序分析分离出形变信号、大气延迟、轨道误差等成分。优点是精度高毫米级缺点是依赖人工地物在植被覆盖区效果不好。SBAS小基线集InSAR的思路是把多景影像按时间基线和空间基线分组组内做干涉然后用最小二乘或奇异值分解方法解算形变时间序列。优点是能处理自然地表缺点是精度比PS-InSAR稍低。SNAP里做时序InSAR需要配合外部工具比如StaMPS或者MintPy。SNAP负责预处理和干涉图生成StaMPS或MintPy负责时序解算。这个流程比较复杂涉及大量参数调优建议先从DInSAR练手熟悉了再上时序方法。7. 几个容易踩坑的细节和我的实操心得7.1 数据下载前的检查清单下载数据前我习惯列一个检查清单卫星平台1A还是1C、成像模式IW还是EW、产品类型SLC还是GRD、极化方式VV还是VH、时间范围、轨道方向升轨还是降轨、以及数据质量标记。尤其是数据质量标记有些数据会有标注“degraded”或者“partial”这种数据尽量别用处理到一半发现有问题更浪费时间。7.2 SNAP内存配置的坑SNAP默认的JVM内存是2GB处理IW模式的SLC数据根本不够用跑着跑着就卡死或者报OutOfMemoryError。解决办法是修改snap.conf文件把default_options里的-Xmx参数调大。我一般设成8GB到16GB具体看机器配置。但也不要设太大超过物理内存反而会频繁触发swap更慢。7.3 解缠参数的经验值Snaphu解缠有几个关键参数--mtf最小相干性阈值、--init-method初始化解缠方法、--cost-mode代价模式。根据我的经验--mtf设0.25左右比较平衡--init-method用MCF最小费用流效果最好--cost-mode用DEFO形变模式适合形变监测场景。当然这些参数不是绝对的具体数据具体调。7.4 形变结果验证的重要性DInSAR得到的形变结果一定要验证。验证方法有几种一是跟GNSS观测数据对比如果有同区域的GNSS站点直接对比形变时间序列二是跟已有研究成果对比看看量级和趋势是否一致三是做内部一致性检查比如用升降轨数据交叉验证。没有验证的形变结果只能算半成品。7.5 批量处理的脚本化思路如果要做大量数据的处理GUI操作效率太低了。SNAP提供了命令行工具gpt可以把处理流程写成XML格式的图谱文件然后用gpt批量执行。我一般会把整个DInSAR流程写成一个图谱文件然后用shell脚本循环调用一晚上能跑几十景数据。这个思路在做时序InSAR时尤其重要因为动辄要处理上百景影像。最后分享一个小技巧SNAP的图谱文件可以用GUI操作一遍后自动生成在Tools菜单里选“Create Graph”就能把当前操作流程导出成XML。然后在这个基础上改参数、加循环比从头写快得多。这个内容后续还可以这样扩展一是结合Cesium做卫星轨道和波束可视化的三维展示把SAR成像几何直观地呈现出来二是用Python调用SNAP的snappy接口把处理流程嵌入到自动化工作流里三是探索Sentinel-1C的新数据特性看看跟1A相比有哪些改进和变化。