InSAR相位拆解与星载/地基形变监测工程实践

1. 相位里到底藏了什么：把InSAR的观测量拆开看

项目评审会上被追问"你这个毫米级形变是怎么算出来的"，很多人答不上来——不是因为不会点软件按钮，而是从没把干涉相位这一项一项拆开过。合成孔径雷达干涉测量（InSAR）看起来是"两景影像一相减"这么朴素的操作，实际上那条彩色条纹图里叠着地形、形变、大气、轨道误差、噪声至少五六种成分，你把哪一项当成了信号、把哪一项当成了噪声，直接决定了最后出来的数字能不能用。DInSAR形变监测、DEM重建、时序分析这三条技术路线，本质上都是在回答同一个问题：怎么把不属于目标的那几项从相位里剔干净。

先把最核心的一条关系摆出来。雷达测的是斜距，两景影像在同一个像元上的相位差，对应的是两次观测到目标距离之差：

φ = -4π/λ · ΔR

λ 是波长，ΔR 是斜距差。这个负号容易被忽略，它意味着目标朝雷达方向靠近（距离变短）时相位是增加的。把 ΔR 拆成视线向形变 d（规定远离雷达为正），就得到从业者最熟悉的那条换算：

d = -λ/(4π) · φ_def

拿常见波段代进去感受一下量级：C波段波长约5.6 cm，每弧度相位对应视线向约4.4 mm；L波段波长约23.6 cm，每弧度对应约18.8 mm。反过来，一个完整的2π条纹对应 λ/2 的视线向变化，C波段就是2.8 cm，L波段是11.8 cm。这个数字非常关键——它同时是InSAR的灵敏度来源，也是它的模糊度来源。相位是被"缠绕"在 [-π, π] 里的，真实相位可能是 φ、φ+2π、φ+4π，你解缠时少算或多算一个2π，形变量就整整齐齐地偏了2.8 cm，而且这个错误在图上看起来一点都不突兀，只是相邻区块之间多了一条不明显的界线。

完整写出来，一景干涉图的相位大致是这么几块的和：

这张表是我认为初学者最该背下来的东西。后面不管讲星载还是地基、讲DEM重建还是DInSAR，都是在跟这六行打交道。而相干性——两个像元复数向量的归一化互相关——则是判断这张干涉图"有没有资格进入下一步"的第一道闸门。相干性低于0.3的区域，相位基本是随机数，你硬解缠出来的东西就是好看但没意义的色块。

还有一个容易被新手忽略的点：InSAR测的是视线向（LOS）的一维投影。卫星侧视观测，入射角通常20°到45°，一个纯水平的东西向位移，投影到视线向上的分量可能只有三到五成；而南北向位移因为卫星飞行方向接近极轨，几乎投影不上去。所以拿到一张"沉降速率图"就宣布某地每年下沉多少毫米，是很危险的动作——你必须说清楚这是LOS向速率，若要分解成垂直和水平分量，得靠升降轨联合或者外部约束。

2. 星载与地基：两套系统在工程上根本不是一回事

标题里把"星载/地基系统"并列，但做过这两类项目的人都知道，它们从硬件到数据处理思路几乎没有可复用的部分，唯一的共同点是都用了合成孔径和干涉相位。把它们混为一谈，最常见的后果是拿星载的精度预期去承诺地基的活，或者反过来用星载的处理惯性去处理地基数据，结果卡在配准和相位闭合上出不来。

星载系统的特点是重复轨道、被动等待、覆盖广。同一颗卫星按重访周期回到近似同一轨道位置，垂直基线由轨道控制精度决定，一般在几十米到几百米之间波动。它的优势是区域覆盖和长期历史数据积累，一个县、一个流域甚至一个省范围的地面沉降普查，星载是唯一现实的选择。劣势也很明确：重访周期从几天到几十天不等，短时突发形变（比如滑坡加速阶段）可能直接被漏掉；侧视几何决定它对南北向位移不敏感；C波段在植被区穿透能力有限，夏秋季节相干性掉得很厉害。

地基系统（常说的GB-SAR一类）则是固定平台、主动连续观测。设备架在边坡对面或者桥墩侧面，几秒到几分钟就能完成一次扫描，测的是距离向和方位向的二维图像，采样密度极高，形变精度在条件好的时候可以做到亚毫米。它特别适合三类场景：一是矿山边坡、水库岸坡这类需要预警的高风险点；二是桥梁、大坝、塔架这类结构物的局部微变形；三是作为星载结果的地面真值校验。代价是覆盖范围小（通常几百米到几公里）、受通视条件限制、设备成本和维护成本高，雨雪天气对信号也有影响。

选型时我一般按这个顺序问三个问题：要的是面还是点、要的是历史还是现在、要的是普查还是预警。要面、要历史、要普查，星载；要点、要现在、要预警，地基。预算够的项目，两者搭配是最稳的组合——地基设备负责关键点的高频监测和预警阈值触发，星载时序负责把点上的结论外推到整个区域，互相校验。

至于机载系统，它在DEM重建这个方向上有独特位置：飞行高度低、基线可以设计得很短、分辨率高，做小区域高精度地形比星载更容易控制基线误差。但机载的轨迹控制精度和姿态稳定性对干涉的影响比星载更苛刻，运动补偿做不好，干涉图直接就废了。

提醒一句：地基雷达的"毫米级"和星载时序的"毫米级"不是同一个统计口径。前者通常指单次观测的精度，后者往往是长时间序列拟合速率的精度，是除以了时间跨度的。拿两个数字直接对比，很容易在汇报时说错话。

3. 干涉处理链的六个关口：从原始回波到干涉图

数据处理这条链，我不喜欢把它讲成"点几下就出图"的流程，因为真正决定成败的往往是被跳过的那几步判断。按我的习惯，把它拆成六个关口，每个关口都有一个"过不去就得回头"的硬指标。

3.1 第一关：干涉对的几何设计

干涉对选错了，后面所有工序都是在给错误数据做美颜。核心控制量是三个基线：垂直基线、时间基线、多普勒质心差。

垂直基线决定地形相位的灵敏度，也决定了相干性衰减。经验上，垂直基线超过临界基线的三分之一到二分之一，相干性就开始明显下降；超过临界基线，理论上完全没有干涉条纹。临界基线近似为：

B_c ≈ λ · R · tanθ / (2 · ρ_r)

其中 R 是斜距、θ 是入射角、ρ_r 是斜距分辨率。实际项目里我不会去精确算这个值，而是把它当作一个"相对尺度"——同一颗卫星、同一种成像模式，历史数据里相干性还不错的那些对的基线范围，就是你的安全区。

时间基线要看地表类型。城区、裸岩、沙漠，几年时间基线还能保持不错的相干性；农田、草地、林地，一个生长季过去就基本失相干。做形变监测时，我通常先按地表覆盖把研究区分层，再给不同层设不同的时间基线上限，而不是一刀切。

多普勒质心差这一项最容易被忽略。两颗卫星的视角如果差得太多，方位向频谱不重叠，干涉相位会被严重污染。处理时要做频谱平移（common band filtering），把重叠部分提出来再干涉，否则相干性会莫名其妙地低。

3.2 第二关：配准精度

配准没到亚像素，干涉条纹就是一片模糊的斜纹。要求是方位向和距离向都要优于十分之一像元，做时序分析时甚至要更严。

我的做法是两级：先用轨道参数做粗配准，把偏移量压到几个像元以内；再用相干性最大化或者条纹频率估计做精配准，在距离向和方位向上分别拟合偏移多项式。这里有个实操细节——精配准的多项式阶数不要贪高，二次到三次足够，阶数过高会在数据边缘产生虚假的偏移量，把好好的边缘区域拉花。另外要检查偏移量残差的分布，如果某个角落的残差明显偏离整体，多半是那一片地形起伏大或者已经失相干，得把它从拟合样本里剔掉。

3.3 第三关：干涉图生成与平地相位扣除

这一步的产物是缠绕干涉图和相干图。平地相位有两种扣法：用轨道几何直接算，或者用多项式拟合条纹频率。前者物理意义清晰，后者对轨道误差更宽容。我的习惯是先用几何法，然后看残余条纹——如果还剩下明显的线性或二次趋势，说明轨道有残差，再用多项式补一刀。

3.4 第四关：滤波的窗口大小是个取舍

滤波的目的是压制噪声、让条纹连续，但窗口一大，细节就没了，形变梯度的边界会被抹平。常用的自适应滤波（按条纹频率调整方向）比固定窗口好，但方向估计本身也依赖噪声水平。

窗口怎么定？我一般按"目标条纹间距的三到五分之一"来估。如果预估研究区最大形变梯度会产生50个像元一个条纹，窗口取10到17比较合适。窗口给到31以上，基本只能看到大尺度趋势，小范围形变坑会被平滑掉。滤波之后一定要回头再看一次相干图，滤波能让噪声区看起来"有条纹了"，那是假象，不能作为解缠可信的依据。

3.5 第五关：解缠

解缠是整个链条里最容易埋雷的一步，因为它的错误表现形式非常"温和"——只是某些区域整体平移了整数个2π。

主流方法分两类：一类是路径积分（枝切法、质量引导），靠质量图决定积分路径，把残差点连起来绕开；另一类是最小二乘加网络流（最小费用流），把解缠转成全局优化问题。前者快、适合质量较好的数据，后者更稳、对噪声区更鲁棒，但计算量大。

我的经验是：解缠之前先划掩膜。相干性低于阈值的区域、水体、阴影区，直接掩掉不参与解缠，宁可在图上留空洞，也不要用错误的积分路径穿过去。掩膜阈值别设得太教条，城区0.3、植被区0.2，按数据实际情况调。解缠之后必须做一致性检查：把解缠相位重新缠绕回去，和原始缠绕相位逐像元比对，残差超过一定比例的区域要标记出来重新处理。

3.6 第六关：相位到物理量的转换

到这一步才是"出结果"。转换本身只是一个线性缩放，但前提工作全在前面：轨道要精确、入射角要逐像元算（不能用一个平均值糊弄整个幅宽）、坐标系要统一。斜距-地距转换、地理编码，都要用到外部高程数据和精密轨道，这一步的迭代次数和插值方法会影响最终配准精度，做DEM生产时尤其明显。

4. DEM重建：把缠绕相位变成高程的几个必答题

用干涉相位反演地形，是InSAR最早也最成熟的用途。原理一句话说完：地形相位与垂直基线、入射角、斜距成比例，比例系数知道了，相位就能换成高程。但工程上要答的题远不止这一道。

4.1 高程模糊度决定了你能做多精细的地形

地形相位对高程的敏感度，用高程模糊度描述：

h_amb = λ · R · sinθ / (2 · B_⊥)

含义是"相位变化一个2π对应多少米的高程变化"。垂直基线越大，模糊度越小，对地形越敏感；但基线越大，相干性越低、解缠越难、大气和轨道误差的放大倍数也越大。这是一个典型的此消彼长。做DEM时我倾向选择中等偏小的基线，宁可牺牲一点灵敏度，换解缠的稳定性和误差的可控性——因为解缠错一个条纹，在高模糊度下就是几米甚至十几米的高程跳变，而且这种跳变在晕渲图上不一定看得出来。

4.2 基线误差是长波误差的元凶

DEM产品最常见的质量问题是"整体倾斜"或者"碗状起伏"，这类长波误差十有八九来自基线估计不准。基线在几百米的量级上，几个厘米的误差，乘以几千倍的放大系数，就能在几十公里范围内造成数米的系统性偏差。

控制办法有三条：一是用精密轨道（精密定轨产品，厘米级）而不是快速轨道；二是用地面控制点或者已有高程数据做基线精化，把残余的线性趋势拟合掉；三是干涉对尽量选同轨道、同成像模式、时间接近的，减少轨道状态差异。

4.3 叠掩、阴影和透视收缩：几何决定的盲区

侧视成像在陡峭地形上必然出现几何畸变。迎坡面被压缩（透视收缩），陡到一定程度发生叠掩，两个不同位置的地面点落到同一个像元里，相位直接失效；背坡面则可能出现阴影，完全没有回波。这些区域不是"处理得不好"，而是"物理上测不到"，只能在成果里明确标注为无效区。

处理策略是多视角联合。升轨和降轨的几何互补，一个方向的叠掩区在另一个方向可能是良好观测区。做山区DEM时，如果只有单侧视角数据，山脊背坡的漏洞是补不上的，这一点在交付前必须跟使用方讲清楚。

4.4 植被穿透带来的系统性高程偏差

这是一个特别容易被忽视的坑。雷达信号会穿透植被冠层，在树冠和地面之间形成体散射，干涉相位反映的是某种"等效散射中心"的高度，通常落在冠层内部某个位置。结果就是：林区的高程系统性偏高，偏高量取决于波长、树高、密度。C波段穿透浅，偏差相对小；L波段穿透深，偏差可能到几米。

如果拿这样的DEM和光学摄影测量或者激光点云对比，会在林区出现整片"正偏差"，很多人第一反应是算法有问题，其实是物理机制差异。做精度验证时必须按地表覆盖分类统计，林区、裸地、城区分开看，否则那个总体RMSE数字毫无意义。

4.5 地理编码与多视：分辨率和噪声的平衡

多视处理是拿空间分辨率换信噪比。多视数越大，相位噪声越低，但地形细节越糊。做1:50000比例尺的地形产品，多视比大概在距离向2到4、方位向4到8这个区间，具体要看原始分辨率和最终成图比例尺。多视之后像元不再是正方形，地理编码时要重采样回规则格网，这一步的插值方法（最近邻、双线性、双三次）对陡峭地形的边缘影响明显，我一般用双线性，兼顾平滑和细节。

地理编码需要迭代：用初始高程把斜距坐标投到地图坐标，插值得到该位置的高程，再重新计算投影，反复两三次直到收敛。轨道精度差的时候，这个迭代可能不收敛，表现是成果图整体偏移几十米，这时候要回去检查轨道文件版本和坐标系定义，而不是继续调插值参数。

5. DInSAR差分：把地形相位剔干净之后剩下的东西

DInSAR是差分干涉，逻辑上是"再减一次"：在干涉相位里把地形相位减掉，剩下的就是形变。听起来简单，但减得干不干净，决定了这个结果是毫米级还是"仅供参考"。

5.1 三种差分路线的取舍

实际操作中二轨法用得最多，因为外部高程数据易得，流程可控。但它有个前提：外部DEM的精度要优于你要测的形变量级对应的等效高程误差。用30米级、高程精度几米的公开DEM去测厘米级形变，残余地形相位会随着垂直基线放大，形成与地形高度相关的条纹，看起来像形变，实际是DEM误差。

识别这类假形变有个简单办法：把残余相位和高程做回归，如果相关系数很高、空间分布与地形高度吻合，基本就是DEM误差而不是形变。垂直基线越大，这个效应越明显，所以差分干涉对也要控制垂直基线，不能太大。

5.2 差分相位里还剩什么

减掉地形之后，差分相位大致是：

φ_diff = φ_def + φ_atm + φ_orbit + φ_dem_error + φ_noise

• 形变相位：目标信号，通常空间上局部化，与地质构造、地下水抽取、地下开采、边坡活动对应。
• 大气相位：分湍流和分层两部分。湍流是高频随机，分层与地形高度强相关，在山区可能达到几厘米，是单对DInSAR最大的误差源。
• 轨道误差：表现为全场的线性或二次趋势，用远离形变区的稳定区域拟合一个多项式扣掉即可。
• DEM误差：与垂直基线成正比，与地形相关。
• 噪声：失相干区域的随机相位。

反演形变时，相位转视线向位移仍然是那个换算，但要注意解缠基准。差分相位解缠后整体带一个未知常数，必须选一个假定不动的参考点把常数定下来。参考点选在形变区里，整个图就整体偏移，这是最常见的低级错误。

5.3 形变梯度和一个条纹的天花板

InSAR有个硬限制：相邻像元的相位差不能超过π，否则就无法正确解缠。换算成形变梯度，就是一个像元距离内形变量不能超过 λ/4。C波段约1.4 cm，如果像元地面尺寸是20米，那形变梯度上限大约是0.7毫米/米。超过这个梯度，比如地裂缝两侧、矿山塌陷边缘、滑坡后缘，条纹会密集到解缠失败。

应对办法有两个：一是用高分辨率数据，像元变小，同样梯度下相位差减小；二是用L波段，波长长，容忍的梯度大一倍多。所以矿区、滑坡这类大梯度场景，L波段数据几乎是首选。另外，用外部先验形变模型做相位模拟和解缠引导（相位解缠的模型辅助），也是常用手段。

5.4 升降轨联合分解

单轨道只能得到视线向投影，要分离垂直和水平分量，需要升降轨联合。设升轨视线向单位矢量为 (e, n, u)，降轨为 (e', n', u')，两个观测方程解三个未知量，还需要额外约束（比如假定南北向位移为零，或者引入外部GNSS观测）。

实际操作中，我通常只做垂直向和东西向的分解，南北向留空或者用外部数据约束——因为卫星几何对南北向本来就近乎不敏感，强行解出来的结果误差极大。分解前要确保两个轨道的参考时间一致、参考点一致，否则解出来的分量会带系统性偏差。

6. 时序InSAR：从一堆干涉图里抠出毫米级趋势

单对DInSAR受大气和失相干限制，精度往往停在厘米级。要把精度推到毫米级，就必须用时序方法——用多景影像、多组干涉对，把大气这种"每次观测都不一样的随机项"在时间维上平均掉，同时用永久散射体或者分布式散射体把稳定点找出来。

6.1 PS和SBAS是两种思路，不是两个软件

PS思路：选一景主影像，所有影像都和它干涉，靠幅度离差指数挑出那些长期稳定的点（建筑物、裸岩、人工设施），要求这些点在所有干涉对里都保持高相干。它的优势是精度高、能捕捉单点非线性形变；劣势是选点依赖人工地物，自然地表覆盖区点太稀，而且对主影像的选择敏感，影像数量少于一定数量（经验上25景以上）时解算不稳定。

SBAS思路：不设唯一主影像，而是组合很多短时空基线的干涉对，形成干涉对网络，用奇异值分解解决秩亏问题，得到每个时间段的形变增量。它能在自然地表上获得较密的观测，适合大范围缓变形变；但对单个干涉对的解缠误差更敏感，需要做相位闭合检查。

我的实际选择是看场景：城区沉降普查用PS类方法，成果点位清晰、便于和建筑对照；山区、农田、大范围缓变形变用SBAS类方法，覆盖更完整；复杂场景（城区+城郊过渡带）用两者结合，先分布式选点再联合解算。

6.2 处理链上的几个关键决策

干涉对组合：SBAS类方法里，时间基线阈值和空间基线阈值是核心参数。阈值给太松，失相干的干涉对混进来，污染网络；给太紧，网络连通性不够，分成很多小连通分量，解算时只能各自定基准，拼起来会出现台阶。我的做法是先画基线图，看数据分布，把阈值设在能保证全网连通的最小值附近，再人工剔掉明显低相干的边。

相位闭合：三个干涉对构成一个三角形，三个相位差之和理论上应该是零（或者2π整数倍）。实际不为零的部分就是解缠误差和噪声的度量。做时序时我会统计所有闭合环的残差，残差大的环对应的干涉对要重点检查，甚至直接剔除。这一步能捞回不少隐蔽的解缠错误。

大气校正：时序方法本身能通过时间平均削弱湍流大气，但分层大气和地形相关，时间平均不掉。常用做法是相位-高程回归（假设同一时刻大气相位是高程的函数），或者引入外部气象再分析数据做改正。在山区，不做大气校正的时序结果经常出现"形变跟着等高线走"的假象。

参考点与基准：时序结果整体带一个常数，参考点一般选在研究区外、远离形变、相干性长期良好的稳定区域。选点前先看一遍区域地质资料和已有监测数据，别把参考点选在正在沉降的城区边上。基准问题也一样——不同处理方法给出的形变速率可能差几毫米每年，交叉验证时必须换算到同一参考基准。

6.3 速率精度能到多少

一个粗略的精度估算关系：

σ_v ≈ λ · σ_φ / (4π · T · √N)

λ 是波长，σ_φ 是单幅干涉图的相位标准差（弧度），T 是时间跨度，N 是有效干涉图数量。从这个式子能看出两件事：时间跨度越长、干涉图越多，速率精度越高；波长越短，同样相位噪声下的位移精度越高（但相干性也更容易掉）。这也是为什么C波段做城区沉降能到毫米级，而植被区往往只能给出厘米级的原因——不是算法不行，是相干性决定了有效观测数量。

7. 我踩过的坑与精度验证的土办法

技术流程讲完，真正拉开水平差距的是排查经验。下面几条都是实际项目中反复遇到的，写出来供对照。

相位跳变整块偏移。成果图上某个区块和周边差了一个明显的台阶，边界是一条不规则曲线——十有八九是解缠时多算或少算了一个2π。判断方法：把该区块的形变量除以 λ/2，如果接近整数，基本可以确认。修复不能靠"加回去"手动修补，必须回到解缠那一步，检查掩膜和质量图，找到错误积分路径的入口。

轨道文件版本不一致。升轨和降轨分别处理时用了不同版本的精密轨道，结果联合分解出来的垂直分量出现系统性倾斜。这个问题极其隐蔽，因为单看每个轨道的成果都正常。我的固定动作是：所有数据统一用同一批次、同一版本的轨道产品，处理前记录版本号，交付时写进元数据。

参考点选在形变区。整个研究区速率整体偏移，但相对形变形态是对的，所以内部检查发现不了。避免办法是选参考点时用外部数据交叉看一眼，或者干脆选多个候选参考点，比较不同选择下速率图的整体平移量，平移量大的说明参考点本身在动。

季节相干性被忽略。夏天做的干涉对植被区全糊，秋天做的同样基线却很好。如果不分季节混在一起做时序，网络连通性会时好时坏，解算权重失衡。按季节分层组合干涉对，是比较土但很有效的办法。

把LOS速率当垂直速率汇报。这是最危险的一条，因为它涉及结论的正确性。一个纯水平位移占主导的滑坡，LOS速率可能只有实际位移的三到五成，按垂直向解释会严重低估风险。汇报时必须写明"视线向速率"，需要垂直分量时说明分解方法和假设条件。

至于验证，除了常规的野外点位比对，我常用三个"土办法"：

1. 看速率图有没有条带。时序成果出现沿方位向或距离向的条带，通常是解缠块拼接或者轨道残差没扣干净，这类图案自然界不会长成那样。
2. 看时间序列曲线的形状。真实的沉降曲线通常是渐变的、有季节波动的，如果某个点的曲线是阶梯状跳变，多半是解缠在某个时刻出错。
3. 多方法交叉。同一区域用不同波段、不同轨道、不同处理方法各做一遍，取交集看结论是否一致。一致的部分可以放心用，不一致的部分必须逐个查原因，不能取平均糊弄过去。

最后分享一个参数习惯：我在每个项目开始前会先做一个"最小可行试验"——挑两三个干涉对，把整条链快速跑一遍，看相干性分布、条纹密度、解缠残差，用这个结果反推参数范围，再铺开做全量。这比一上来就批量处理几十景、跑到一半发现基线选错要省时间得多。干这行久了会发现，InSAR的功夫不在软件里，而在你决定"这一对能不能用""这个条纹是真是假"的那几秒钟判断上。