0 引言

合成孔径雷达（synthetic aperture radar，SAR）成像具备全天时、全天候、分辨率高、穿透性强等优势，因此被广泛应用于冰情监测和海洋遥感等领域^[1]。电磁波极化对海冰的介电特性和几何结构非常敏感。极化SAR（polarimetric SAR，PolSAR）通过在两个正交极化通道上分时发射和同时接收全极化电磁波，能够有效获取海冰的完整极化散射信息，这使得对海冰的精细监测成为可能。目前国际上已有多套可用于海冰遥感的星载PolSAR系统，例如加拿大的RADARSAT-2、日本的ALOS-PALSAR-1/2、中国的高分三号（GF-3）和陆地探测一号（LT-1）等。这些系统提供了大量的海冰PolSAR图像数据，全极化海冰散射信息解译也成为了当前的研究热点之一。

PolSAR获取的海冰数据往往是一个3×3的相干矩阵，需要采用极化目标分解等方法从相干矩阵中有效解译出可用于识别分类的海冰几何信息及散射信息。基于模型的极化目标分解方法是将复杂目标的后向散射分解为若干标准散射分量的加性组合，该方法物理意义明确，已成功应用于海冰监测中^[2-13]。SHARMA等^[2]提出将冰川冰的后向散射分解为风致雪表面散射分量和体散射分量的非相干叠加，构造了一个含有13个未知数的散射平衡方程组，并通过一系列近似和简化求解散射功率。但该分解中总存在一个可被最小化到L2范数、但又无法消除的残余分量。ZHANG等^[7]将海冰后向散射分解为表面散射、二面散射和体散射等分量之和，但该分解中仍不可避免地存在残余协方差矩阵分量。PARRELLA等^[8]将冰盖后向散射简化为X-Bragg表面散射和体散射的非相干叠加，但该分解中也存在类似的残余分量。因此，基于模型的海冰分解方法依然存在极化信息利用不充分的问题。此外，为简化起见，当前针对海冰的极化分解仅采用单一的体散射模型，即默认海冰粒子只存在一种定向分布。例如文献^[2]在分解中仅使用了随机定向体散射模型，而文献^[7]只考虑了粒子呈现垂直定向的情况。事实上，垂直、水平以及随机定向分布海冰粒子都有可能出现在冰层中，这要求分解方法能根据实际冰层中海冰粒子的定向方式，自适应地匹配恰当的体散射模型。

为了实现对极化信息的完整利用，本文提出基于模型的星载PolSAR海冰自适应极化分解法，该算法包括体散射分量自适应分解以及表面和二面散射分量自适应分解两个方面。在海冰体散射自适应分解方面，构建覆盖海冰粒子水平、垂直和随机定向分布的体散射模型，并基于最大功率拟合的自适应模型选择策略来实现体散射模型的自适应匹配。在表面和二面散射分量自适应分解方面，结合标准Huynen二分法，根据实际冰面的散射机制自适应选取分解方法，并保证分解结果不退会化。

1 海冰分解的散射模型

1.1 表面散射和二面散射模型

海冰目标的散射过程通常可以利用表面散射、二面散射和体散射等不同的散射机制来表征^[7，11]，海冰典型散射机制示意图如图1所示。表面散射描述了电磁波在微粗糙表面的散射过程；二面散射描述了电磁波在两个不同表面的散射过程；体散射则描述了电磁波多次散射过程的叠加。

体现不同散射的海冰冰况图如图2所示。

浮冰表面比较粗糙，表现为表面散射；冰脊受到洋流的挤压呈现二面角结构，表现为二面散射；海冰的水平切面图表明海冰中夹杂着较多椭球状的盐水包裹体，电磁波在海冰中存在多次散射，表现为体散射。

与文献^[7]类似，本文也将海冰目标散射相干矩阵T表示为表面散射分量、二面散射分量以及体散射分量相干矩阵的非相干叠加，即

式中：f_S，f_D，f_V为表面散射、二面散射和体散射分量在总后向散射中所占的权重；T_S，T_D，T_V为表面散射、二面散射和体散射分量的相干矩阵。下面对分解过程中使用的散射模型进行介绍。

假设空气-海冰界面为微粗糙表面，由于雷达回波受表面粗糙度影响较大，如果海冰的表面比较粗糙（例如破碎冰、光滑浮冰以及开阔水域等），就会发生表面散射。在对三分量模型的Freeman-Durden分解（Freeman-Durden decompo-sition，FDD）中使用一阶Bragg表面散射模型来描述微粗糙面的表面散射，表面散射的相干矩阵T_S的表达式^[14]为

其中

式中：β为表面散射复参数；|·|为复数取模运算符；*为复共轭运算符；R_h，R_v为水平极化和垂直极化的Bragg反射系数。在实际场景中，移动的浮冰在风和洋流的作用下会发生碰撞和汇聚，从而导致冰表面发生变形，形成冰脊，进而呈现二面角结构。当极化电磁波与由两个相互垂直的散射面构成的二面角散射体相互作用时，所产生的散射过程为二面散射。在FDD中，二面散射的相干矩阵T_D的表达式^[14]为

其中

式中：α为二面散射复参数；γ_h，γ_v为水平和垂直极化通道的相位；R_Th，R_Tv为平整海冰表面的水平极化和垂直极化菲涅尔反射系数；R_Gh，R_Gv为变形海冰表面的水平极化和垂直极化菲涅尔反射系数。

1.2 体散射模型

如图1所示，海冰的体散射分量主要源于冰层中夹杂的椭球状盐水包裹体。对海冰的体散射建模即是对一定体积内大量相互独立的椭球状盐水包裹体的非相干散射过程的建模。海冰的次表层结构模型为形状和大小相同的均匀粒子云，粒子云由大量的椭球状盐水包裹体构成。文献^[15]给出了三维均匀定向椭球状粒子云的广义体散射解析模型，鉴于极化体散射响应会受到粒子云定向分布的影响，本节基于海冰的结构特性对该解析模型进行简化，得到最终的体散射模型。

与文献^[2，7，16]类似，本节根据几何特征将椭球状盐水包裹体简化成独立偶极子，根据方位向特征将海冰次表层结构建模成满足反射对称性的偶极子云，即海冰粒子的方位角φ在[-π，π]范围内均匀分布。

海冰冰层中盐水包裹体在俯仰向上不一定总是以俯仰角ψ=π/2垂直排列，而是存在垂直、水平以及随机等不同定向分布。盐水包裹体的不同定向分布如图3所示，其中ψ₀为中心俯仰角，Δψ为俯仰角的分布宽度。图中水平定向分布的盐水包裹体的中心俯仰角ψ₀为零，俯仰角分布宽度为π/4；垂直定向分布的盐水包裹体的中心俯仰角ψ₀为π/2，俯仰角分布宽度为π/4；随机定向分布的盐水包裹体的中心俯仰角ψ₀为零，俯仰角分布宽度为π/2。

垂直、水平以及随机定向分布粒子都有可能出现在冰层中。不失一般性，本节选定水平定向分布粒子云的中心俯仰角ψ_h0=0、垂直定向分布粒子云的中心俯仰角ψ_v0=π/2进行分析，并赋予相同的俯仰角分布宽度，即Δψ_h=Δψ_v=π/4。

根据上述条件，对体散射解析模型进行简化，构建可完全涵盖盐水包裹体水平、垂直和随机等不同定向分布的体散射模型的相干矩阵${\mathit{T}}_{\mathrm{V}}^{\mathrm{h}}$，${\mathit{T}}_{\mathrm{v}}^{\mathrm{v}}$，${\mathit{T}}_{\mathrm{V}}^{\mathrm{r}}$，其表达式分别为

其中

式中：θ为入射角。

2 自适应极化分解

2.1 体散射分量的自适应提取

散射建模是为了更好地描述构成海冰散射的潜在分量，而用于目标识别和分类的关键参数（如散射功率等）则需通过极化分解获得。基于模型的海冰分解方法普遍存在极化信息利用不充分的问题^[2，7，11]，除此之外，FDD^[14]、四分量分解^[17]等基于模型的极化分解方法还同时存在负散射功率问题。为了避免负散射功率的出现，文献^[18]提出了一般意义上的非负特征值分解方法，即完全模型分解法（complete model-based decompo-sition，CMD）。CMD首先采用广义特征分解来求解体散射功率，准确地将相干矩阵的秩降为1；然后通过模型适配或特征分解对剩余相干矩阵进行二分量分解。CMD不仅解决了基于模型的分解方法中的负散射功率问题，而且还充分利用了极化信息。因此，针对海冰分解方法中存在的极化信息利用不充分的问题，借鉴CMD中采用的广义特征分解来求解体散射功率，推导出使相干矩阵精确降秩且不违反物理可实现条件的体散射功率的解析表达式。根据T_S，T_D，T_V的表达式，可得

式中：rank（·）为矩阵取秩函数。可知，矩阵T_S，T_D的秩小于矩阵T_V，T的秩。为了充分利用T中包含的极化信息，f_VT_V应能精确地使T的秩减少1，即

式中：T_N为剩余相干矩阵。

根据广义特征分解原理^[18]，可直接得到

式中：λ（·）为矩阵特征值的取值函数；-1为矩阵求逆运算符。从降秩的角度来看，1/f_V可以取矩阵T^-1T_V三个正特征值中的任意一个。然而，权重f_S和f_D的物理非负性要求f_V应该尽可能小，从而使得秩为2的剩余相干矩阵T_N满足半正定性，此时可唯一确定

式中：λ_max（·），λ_min（·）为矩阵最大特征值和最小特征值的取值函数。则体散射功率分量

式中：trace（·）为矩阵的迹的取值函数。由式（6）~式（10）可分别推导出与三个体散射模型对应的水平、垂直和随机等不同定向分布的体散射功率分量$P_{\mathrm{V}}^{\mathrm{h}}，P_{\mathrm{V}}^{\mathrm{v}}，P_{\mathrm{V}}^{\mathrm{r}}$。为了从$P_{\mathrm{V}}^{\mathrm{h}}，P_{\mathrm{V}}^{\mathrm{v}}，P_{\mathrm{V}}^{\mathrm{r}}$中选择一个分量作为最终的PV，本节基于最大功率拟合准则进行体散射模型的自适应选取，即在提取体散射分量时选择$P_{\mathrm{V}}^{\mathrm{h}}，P_{\mathrm{V}}^{\mathrm{v}}，P_{\mathrm{V}}^{\mathrm{r}}$中取值最大的作为最终的体散射功率，有

式中：max（·），min（·）为最大值和最小值的取值函数；λ_l（·）为矩阵第l个特征值的取值函数。

在${\mathit{T}}_{\mathrm{V}}^{\mathrm{h}}，{\mathit{T}}_{\mathrm{V}}^{\mathrm{v}}，{\mathit{T}}_{\mathrm{V}}^{\mathrm{r}}$中选取能获得最大功率的模型作为T_V的最终模型，即

2.2 表面和二面散射分量的自适应提取

CUI等^[18]提出了模型适配法（model fitting，MF），在完成体散射分量提取后，采用MF从剩余相干矩阵T_N中提取了表面散射分量和二面散射分量。针对模型适配法对二面散射分量功率和表面散射分量功率的求解过程复杂且物理意义不明确的问题，AN等^[19]以理想情况下的表面散射分量和二面散射分量的交叉极化项等于零为前提，提出了改进的CMD（improved CMD，ICMD），实现对表面散射和二面散射分量提取过程的修正。CMD和ICMD在确保散射功率非负的同时解决了信息利用不完全的问题，但经分析发现，CMD和ICMD仍无法实现完整的三分量分解，在很多情况下分解会退化成二分量分解。为了解决这个问题，本文提出了改进的模型适配法。

在CMD和ICMD中，去取向后的剩余相干矩阵T_N可表示为

其中

式中：T_SorD，T″分别为表面散射和二面散射分量的剩余相干矩阵。可以发现，式（18）的分解实际上采用的是标准Huynen二分法（canonical Huynen dichotomy，CHD）。LI等^[20]在2015年提出了统一Huynen二分法（unified Huynen dichotomy，UHD），该方法是对Huynen分解（Huynen decomposition，HD）的拓展。HUYNEN^[21]提出了HD，采用唯象学理论来提取雷达目标的物理结构特性，并指出由于现实世界具有对称性和规则性，因此HD的结果总会偏好目标的表面散射分量，无法用于不对称和不规则散射体的分析。文献^[20]对HD进行了拓展，提出了CHD。CHD保留了与目标非对称性和非规则性相关的参数，提供了三种适用于相干矩阵的二分法，在分解过程中可自适应地选择偏好于表面散射分量、二面散射分量或体散射分量。文献^[19]采用了CHD中偏好体散射的分解方法对剩余相干矩阵T_N进行二分量分解，实现对表面散射和二面散射分量的提取。在提取表面散射和二面散射分量之前，已经求得了体散射分量，因此采用模型适配法得到的两个单目标散射体并不总是对应表面散射分量和二面散射分量，即分解结果中存在二面散射分量或表面散射分量为零的情况。这导致基于模型的三分量分解退化成二分量分解，使得对地物后向散射机制的分析不够全面。为了解决这一问题，本文提出了一种改进的模型适配法，利用CHD的另外两种二分法，即偏好目标表面散射分量的分解方法和偏好目标二面散射分量的分解方法，来分解剩余相干矩阵T_N。

本文所提的改进的模型适配法首先利用CHD中偏好目标表面散射分量的方法对T_N进行二分量分解，此时的偏好矢量q₁=[1，0，0]^T[20]，其中T为矩阵转置运算符。与q₁对应的归一化目标矢量k₁和目标的相干矩阵T₁的表达式为

式中：H为矩阵共轭转置运算符；T₁为偏好表面散射分量的表面散射相干矩阵。进一步可给出对T_N进行偏好表面散射分量的二分量分解后的结果，即

其中

式中：T′₁为偏好表面散射分量的二面散射相干矩阵。然后利用CHD中偏好目标二面散射分量的分解方法对T_N进行二分量分解，此时偏好矢量q₂=[0，1，0]^T[20]，对应的T_N二分量分解结果为

其中

式中：T₂为偏好二面散射分量的二面散射相干矩阵；T′₂为偏好二面散射分量的表面散射相干矩阵。

这里通过比较T_N11和T_N22的大小，根据实际海冰的散射机制，自适应地选取式（22）或式（24）作为分解方法，来提取表面散射和二面散射分量。由于T_N满足半正定性，其二阶主子式应当是非负的，因此可得到不等式组

当T_N11≥T_N22时，结合式（26），有

进一步推导，可得

如果T_N11≥T_N22，根据文献^[19]的判定准则，T_N对应的散射机制中表面散射占主导，所以此时选取偏好表面散射分量的二分法，即式（22）。根据式（23），求解表面散射功率P_S和二面散射功率P_D，其表达式为

如果T_N11＜T_N22，则T_N对应的散射机制中二面散射占主导，此时选取偏好二面散射分量的二分法，即式（24）。根据式（25）求解二面散射功率P_D和表面散射功率P_S，其表达式为

2.3 算法小结

综上所述，可给出基于模型的海冰自适应极化分解算法的流程，如图4所示。

本文提出的基于模型的海冰自适应极化分解法涉及体散射分量自适应分解以及表面和二面散射分量自适应分解两方面。

在海冰体散射自适应分解方面，本文推导了使相干矩阵精确降秩且不违反物理可实现条件的体散射功率的解析表达式；针对冰层中呈现不同定向的盐水包裹体，给出了水平、垂直和随机三种体散射模型；提出了基于最大功率拟合的自适应模型选择策略，即从水平、垂直和随机体散射分量中选取功率最大的作为最终的体散射分量，该分量所对应的散射模型为最终的体散射模型，从而根据冰层粒子云的定向方式来自适应匹配体散射模型。

在表面和二面散射分量自适应分解方面，本文先分别利用CHD中偏好表面散射分量和二面散射分量的方法对剩余相干矩阵进行降秩分解；再根据剩余相干矩阵中的元素T_N11，T_N22的相对大小关系，判断冰面的实际散射机制，进而自适应地选取对应的分解方法，得到表面散射和二面散射分量功率；最终在实现极化散射信息完整利用的同时，保证分解结果不会退化。

3 实验验证

3.1 普利兹湾海冰数据集

2016年发射的GF-3卫星是中国首颗多极化C波段高分辨率SAR卫星。本文海冰分解实验所使用的PolSAR数据由GF-3卫星于2019年4月29日在南极洲普利兹湾东北部以全极化模式获取，PolSAR的分辨率为8m。此数据的中心地理坐标为东经73°6′，南纬68°1′，所在区域毗邻的中国南极中山站的地理坐标为东经76°22′14.28″，南纬69°22′24.76″。普利兹湾位于印度洋，区域海冰受海面温度、风和上升气流的影响，具有冻结期（3月至9月）和融化期（10月至次年2月）交替变化的特点。由于地处高纬度海域，普利兹湾的海冰变化在全球气候变化研究方面发挥着重要作用。

GF-3卫星获取的全极化数据为Level-1A单视复数产品数据，因此，本文利用欧洲航天局（ESA）的PolSARpro软件进行多视和Refined Lee滤波来抑制斑点噪声。图5为研究区域的Pauli分解RGB合成图像。图6为美国国家冰雪中心于2019年5月2日获取的南极局部海冰冰况图，图中不同颜色对应不同类型的海冰区域，蓝色为无冰区，紫色为初冰区，黄色为一年冰区，绿色为一年薄冰区，棕色为多年冰区，灰色为冰脊区。红星标记处为研究区域所在位置。可以看出，研究区域的海冰以初冰为主。根据世界气象组织（world meteorological organization，WMO）按照海冰形成和发展阶段对海冰类型的划定，初冰是尼罗冰和一年冰之间的过渡，厚度从10cm到30cm不等，初冰冰层中包含许多椭球状的盐水包裹体^[16]。某些季节的初冰可能被一层薄雪覆盖，由于雪层盐度低，通常被视为透明介质。因此对于初冰而言，由粗糙表面和雪冰界面所引发的表面散射占主导地位，次要的体散射则主要由冰层内的盐水包裹体贡献。

3.2 基于模型的三分量分解结果

基于本文所提的海冰自适应极化分解法得到的三分量分解结果如图7所示。

设M×N像素大小的图像中像素点（i，j）的散射功率为P_qij，由分解结果可生成散射分量功率矩阵P_q=[P_qij]，其中P_qij为像素点（i，j）的表面散射（S）、二面散射（D）和体散射（V）分量功率，q∈{S，D，V}，i=1，2，···，M，j=1，2，···，N，与数据集大小匹配的M=1 820，N=1 824。

由图7可以直观地看出，研究区域海冰的主导散射机制为表面散射，其次是体散射，二面散射分量很少。为了定量分析研究区域海冰的散射机制，需计算分解结果中各像素点散射分量的功率与其总后向散射功率之比。像素点（i，j）的总后向散射功率

则像素点（i，j）各散射分量的功率与其总后向散射功率的比定义为P_qij/P_Aij×100%，q∈{V，S，D}。计算结果与图7所示的分解结果一致：表面散射分量功率占比最高，占总功率的75.73%；体散射分量为次要散射机制，其功率占总功率的21.42%；二面散射分量功率只占总功率的2.85%。结合图7和计算结果可知研究区域海冰的主散射机制是表面散射。这是因为雷达后向散射易受海冰表面粗糙度的影响，浮冰和破碎冰的表面粗糙度较高。分解结果中的次要散射机制是由冰层中的盐水包裹体引发的体散射，结合图6所示冰况，可知初冰冰层中含有大量的盐水包裹体，这与文献^[16]中的结论一致。由于初冰表面通常相对平坦，没有大规模的二面角结构（例如多年冰表面的雪丘以及冰脊等），所以分解结果中二面散射分量占比最少。

3.3 体散射分量提取结果

本节详细分析基于最大功率拟合选取策略所得到的体散射分量分解结果。设M×N像素大小的图像中像素点（i，j）对应的体散射功率为P_Vij，进而可生成体散射功率矩阵P_V=[P_Vij]，其中i=1，2，···，M，j=1，2，···，N。从功率的角度来看，在${\mathit{T}}_{\mathrm{V}}^{\mathrm{h}}，{\mathit{T}}_{\mathrm{V}}^{\mathrm{v}}，{\mathit{T}}_{\mathrm{V}}^{\mathrm{r}}$中选取T_V意味着对P_V进行三分量分解，即

其中

式中：${\mathit{P}}_{\mathrm{V}}^{\mathrm{h}}，{\mathit{P}}_{\mathrm{V}}^{\mathrm{V}}，{\mathit{P}}_{\mathrm{V}}^{\mathrm{r}}$为水平、垂直和随机定向分布的体散射功率分量矩阵；P^k_Vij为像素点（i，j）对应的水平、垂直和随机定向分布的体散射功率分量。如果不同定向分布的体散射相干矩阵${\mathit{T}}_{\mathrm{V}}^k$被确定为像素点（i，j）所匹配的体散射相干矩阵T_V，则

体散射功率分量提取结果如图8所示。与相干矩阵匹配的像素点数占比统计示意如图9所示，图中红色、绿色、蓝色分别对应与${\mathit{T}}_{\mathrm{V}}^{\mathrm{h}}，{\mathit{T}}_{\mathrm{V}}^{\mathrm{v}}，{\mathit{T}}_{\mathrm{V}}^{\mathrm{r}}$匹配的像素点。

由图9可以看出，P_V的大部分被分解为垂直分布分量${\mathit{P}}_{\mathrm{V}}^{\mathrm{v}}$，这在图8中也有反映，垂直分布的相干矩阵${\mathit{T}}_{\mathrm{V}}^{\mathrm{v}}$与约70%的像素点匹配，即研究区域冰层中大部分盐水包裹体在俯仰向上偏好于垂直定向，这一结果与对海冰结构的现有分析结论是一致的^[16]。海冰是柱状多晶体，其主轴（即c轴）几乎平行于水平面，偏差只有几度。这样的几何结构使冰层中的冰板垂直于水平面，因此，夹在冰板之间的椭球状盐水包裹体往往优先呈现垂直定向。柱状多晶体的c轴相对水平面排列的一致程度主要取决于深度，海冰上层多晶体的c轴方向实际上是随机的，但随着深度的增加，多晶体的c轴逐渐呈现水平排列。多晶体c轴的定向从随机到水平的过渡过程被称为几何选择。除受深度的影响外，多晶体c轴的方向还可能受到底层洋流的影响^[22]，底层洋流会使原本垂直排布的盐水包裹体变成随机排布甚至水平排布。

由上述分析可知，研究区域内海冰冰层内盐水包裹体的取向以垂直定向为主导，也存在次级随机定向和不可忽略的水平定向。因此，本文建议采用三种不同定向的体散射模型来描述海冰的体散射机制，而不是采用单一的垂直定向的体散射模型，以实现对盐水包裹体水平、垂直和随机定向分布分析的完全覆盖。

3.4 表面和二面散射分量的提取算法对比

为了验证本文提出的表面散射分量和二面散射分量功率自适应提取算法三分量分解的完整性，基于实验数据，统计了本文算法和ICMD提取的表面散射功率或二面散射功率为零的像素点的个数。统计结果表明，本文算法的分解结果中不存在像素点表面或二面散射功率为零的情况。而在基于最大功率拟合的体散射模型选择策略下，ICMD的分解结果中有3 109 893个像素点的分解退化为二分量分解，退化的像素点数占总像素点数的92.76%，大量像素点的二面散射分量功率为零。可知，ICMD虽然充分利用了极化自由度，但无法精准分离出研究区域海冰目标的二面散射分量信息，而本文算法在极化信息完全利用的前提下实现了完整的三分量分解，这对充分分析海冰后向散射机制的散射分量有重要意义。

4 结论

本文首先推导了使相干矩阵精确降秩且不违反物理可实现条件的体散射功率的解析表达式，提出了利用标准Huynen分解根据冰面实际散射机制自适应提取表面和二面散射分量的算法，在实现极化散射信息完整利用的同时保证分解结果不退化。然后，基于最大功率拟合的自适应模型选择策略，解决了现有海冰模型分解法无法实现对冰层粒子定向方式自适应匹配的问题。最后，基于GF-3南极洲普利兹湾海冰数据开展了体散射模型和分解算法有效性的实验验证。实验结果首次从极化分解角度揭示了冰层中盐水包裹体的取向以垂直定向为主导、同时存在次级随机定向以及不可忽视的水平定向的定向特征；验证了海冰多晶体的c轴随深度增加逐渐对准、呈现水平排列的几何选择机制。

参考文献