基于多时相TanDEM-X极化干涉SAR数据的水稻株高反演

国贤玉，李 坤，邵 芸，Juan M. Lopez-Sanchez，王志勇

1. 山东科技大学测绘科学与工程学院，山东 青岛 266590 2. 中国科学院遥感与数字地球研究所，北京 100101 3. 中国科学院大学资源与环境学院，北京 100049 4. 中科卫星应用德清研究院微波目标特性测量与遥感重点实验室，浙江 德清 313200 5. University of Alicante, Alicante, 99, Spain

引 言

水稻株高是水稻本身以及土壤、水文、气象等因素的综合反映，是水稻长势监测的重要指标。遥感以其宏观性、时效性、周期性等特点，在水稻株高反演中具有不可替代的作用。光学遥感能够直接获取植被水平分布信息，利用光谱反射率、植被指数以及纹理因子等参数，与水稻株高建立定量关系，实现大范围快速反演[1]。但光学遥感对植被垂直分布信息探测难度较大，而且水稻生长季内多云雨天气，又给光学数据获取带来较大的困难，这些原因大大限制了其在水稻株高反演中的应用[2]。合成孔径雷达(synthetic aperture radar，SAR)具有全天时、全天候、穿透性的优势，成为水稻株高反演的重要手段之一。

目前基于SAR的水稻株高反演方法主要有三类：基于辐射传输理论(radiative transfer theory，RTT)的水稻株高反演方法、基于SAR干涉测量(interferometry SAR，InSAR)以及基于极化干涉测量(polarimetric SAR interferometry，PolInSAR)的水稻株高反演方法[3-6]。其中，基于PolInSAR的水稻株高反演方法是通过干涉和极化信息的有效组合实现了观测维度的拓展，将目标竖直方向上不同位置的散射中心分离开来，估计确定具有最高相干性的散射分量，然后利用通道差分法[7]、模型解算法[8-9]，计算水稻株高，该方法既有InSAR对地表植被散射体的空间分布敏感的特性，又具有极化SAR对植被散射体的形状和方向敏感的特性[10]。近年来，PolInSAR已成为植被高度反演研究热点，尤其是在森林高度反演方面得到验证和应用[8，11]。

不同的研究者从不同的角度进行了PolInSAR理论、方法和应用研究。目前基于PolInSAR的植被高度反演方法主要有两类：通道差分法、模型解算法。通道差分法是通过直接选择或计算得到代表植被冠层散射和代表地表散射相位并作差，进而得到水稻株高。该方法计算简单，但选择和计算得到的代表植被冠层散射和代表地表散射分量并不准确，存在较大误差。模型解算法是基于PolInSAR散射模型开展的反演算法，有严密的物理模型的支撑，反演精度较为准确。目前采用的PolInSAR散射模型主要以双层散射模型最为广泛，而双层体散射主要有随机体散射体/地表二层模型RVoG模型[8，12]、方向体积层模型(oriented volume over ground, OVoG)[13]。尽管OVoG模型更贴合真实场景，水稻层内部结构是垂直分布的，消光系数极化状态存在依赖性，但是这种依赖性并不强烈，且OVoG模型应用需要许多假设，涉及参数多，导致对模型的推演和应用较困难。RVoG模型简单，涉及的参数少，反演精度较高，是森林研究中最常见的模型[14]。然而，水稻作为低矮植被的代表，具有其独特的生长环境，下垫面常常为水面，地表散射容易发生镜面散射，传统的RVoG模型并不完全适用[15-16]。

传统的RVoG模型未考虑植被生长环境，反演结果精度不高，对水稻株高反演应用有局限性。因此，针对水稻植株结构与生长环境等特点，利用多时相TanDEM-X数据，对9个时相的TanDEM-X极化干涉SAR数据进行距离向频谱滤波，求取8个复相干系数，进行去相干处理，建立适用于水稻田特性的RVoG模型，构建基于RVoG模型的极化干涉SAR水稻株高反演方法，进而评定反演结果。

1 研究区概况和数据源

研究区位于西班牙塞维利亚西南部，瓜达尔基维尔河的南部，地中海西北沿岸，面积约为30 km2，地理位置位于36.85°—37.38°N，5.91°—6.23°W之间。该区农作物主要为水稻，水稻常年在五月左右种植，十月前后收割，一年一熟。水稻品种大都是普通栽培稻，播种方式为撒播种植，在整个水稻物候期间稻田内有一层水覆盖。图1为Google earth获取QuickBird卫星数据制作的研究区正射影像。

图1 研究区地理位置与范围(右图4个红图框为研究区所选的水稻样区)

由于基于PolInSAR植被株高反演需要较低的时间去相干、大空间基线数据[17]，星载卫星TanDEM-X系统，是双站模式，TanDEM-X卫星发射电磁波，TanDEM-X和TerraSAR-X同时接受回波信号。因此，TanDEM-X系统可以忽略时间去相干和大气的影响, 另外2015年4月到9月，TanDEM-X系统进行一次科学实验，提高到2～3 km的空间基线，是反演农作物的最佳数据。因此，在水稻生育期内获取了9个不同时相的TanDEM-X双极化CoSSC数据，中心波长9.65 GHz，带宽为150 MHz，极化方式为HH和VV极化，该数据为一发双收的双站SAR数据，包括已配准的单视复数据对，详细信息如表1所示。图2给出了TanDEM-X双极化CoSSC数据彩色合成图。

SAR卫星过境的同时，获取了详细的地面数据。9个时相中地面数据的采集大部分为卫星过境当天采集，个别时相样区地面采集数据与卫星过境时间相差1～2 d。地面数据主要来自均匀分布于研究区的4块水稻样区，图1中红色方块为样区。对每一块样区，采集数据包括水稻物候和株高、样区面积(hm2)、测量/播种/收割日期、稻株密度(株·m-2)、稻穗密度(穗·m-2)、产量(kg·hm-2)、灌溉方式、盐碱性、病虫害等信息。其中，利用皮尺测量水稻株高，其精度为cm级。此外，研究区日值气象数据(温度、湿度、风速、风向、降水等)由西班牙农业气候信息系统(SIAR)提供。

表1 TanDEM-XCoSSC信息列表Table 1 List of TanDEM-X CoSSC image information

图2 TanDEM-X双极化CoSSC数据彩色合成图 (获取时间2015年7月7日，R=VV；G=HH-VV；B=HH)(a)：主影像；(b)：辅影像

Fig.2TanDEM-XcolorcompositeimagesofCoSSCdata(time:7July2015,R=VV,G=HH-VV,andB=HH)

(a)：Master image；(b)：Slave image

2 实验部分

采用基于RVoG模型的极化干涉植被高度反演方法，构建复相干系数与水稻株高之间的定量关系。复相干系数与目标场景、传感器特性等多种因素密切相关，可以表示为

(1)

其中，γtemp为时间去相干系数；γgeom为几何去相干系数；γproc为干涉处理过程中造成的去相干系数；γSNR为信噪比去相干系数；γBQ为比特量化去相干系数；为体去相干。在极化干涉SAR水稻株高反演之前，需要剔除水稻场景以外去相干因素对复相干系数的影响。图3给出了极化干涉SAR水稻株高反演基本流程。

图3 极化干涉SAR水稻株高反演基本流程

①首先对9景TanDEM-X CoSSC双极化SAR数据进行距离向频谱滤波补偿[18]，对于TanDEM-X CoSSC数据，由于其空间基线大，传感器空间方位差异对复相干系数的影响较大，γgeom小于0.8，因此必须考虑γgeom的影响。②对距离向频谱滤波后的TanDEM-X数据进行共轭相乘得到干涉图，并利用主辅影像斜距、入射角、中心波长等获取影像时的轨道参数计算平地相位并去除平地效应。③全极化干涉SAR系统中每个散射单元为2×2的散射矩阵，但对于TanDEM-X双极化系统，每个散射单元为2×1的散射矩阵。因此将TanDEM-X数据对应的S矩阵转化为T4矩阵，并利用Boxcar方法对T4矩阵进行滤波处理，窗口大小：21×21，以便减小相干斑点噪声，提高图像的辐射分辨率。④计算8个复相干系数，分别为主辅影像HH，VV，HH-VV，HH+VV复相干系数、两个优化的相干性系数[7](OPT1，OPT2)、PD相位差分法得到的复相干系数(PDHigh，PDLow)。然后进行去相干因素影响的剔除。由于TanDEM-X CoSSC为一发双收的双站模式数据产品，其主辅图像同时获取，不受时间去相干的影响，即γtemp为1。γproc主要是由主辅图像配准引起的，而TanDEM-X CoSSC数据产品具有较高的配准精度，因此γproc可以忽略。对于TanDEM-X and TerraSAR-X数据，考虑到成像区域主要为农田场景，γBQ约为0.965[15]。在干涉测量研究中，信噪比去相干γSNR通常被忽略，因为它只对低散射区造成影响。但对于稻田场景，由于其X波段SAR图像上的后向散射系数通常在-25～-5 dB之间，尤其是幼苗期后向散射系数很低，因此，信噪比去相干γSNR不能忽略。基于TanDEM-X头文件数据计算SNR，采用公式(2)计算信噪比去相干γSNR[4]。

(2)

⑤利用剔除去相干因素影响惑的8个复相干系数，构建一种基于RVoG模型的水稻植株高度反演算法，实现研究区9个时相水稻株高反演。

2.2 基于RVoG模型的水稻植株高度反演算法

RVoG模型包含植被层及其下垫面，其中植被层厚度近似为hv且包含随机取向粒子集合。对于一般植被目标，RVoG主要考虑植被层体散射以及地表面散射。对于水稻而言，由于其生长周期内下垫面多为水面，易发生镜面反射，因此地表面散射可忽略不计，此外水稻与下垫面之间的二次散射较强，需要考虑。因此，对于水稻目标，RVoG模型可以表示为[15]

(3)

其中，kz为垂直波数，φ0为地形相位，hv为植被层厚度也就是我们所要反演的水稻株高，mDB(ω)为二次散射对应的地体幅度比。

γv表示体相干，取决于趴标地物散射体的消光系数和植被层厚度hv

(4)

其中，σ为消光系数，θ为雷达数据入射角。

于是，式(3)和式(4)中包含了4个未知数：hv，mDB，σ，φ0。式(3)可以写成复平面内的线性形式

=eiφ0[γV+L(ω)(γDB-γV)]

(5)

基于上面构建的极化干涉SAR水稻株高反演模型可以看出，共有5个模型参数。这里采用模式搜索法[19]的迭代算法进行水稻株高反演。由于模式搜索法是直接优化方法，不需要求解目标函数的导数。因此，该方法具有方法简单，计算方便快速等优点。总体流程如图4所示。

图4 基于RVoG模型的针对水稻特性的水稻株高反演算法的流程图

①基于剔除各种去相干因素影响的8个不同极化通道的复相干系数，利用最小二乘拟合获得式(5)在复平面内对应的直线。②由于模式搜索法的搜索步长影响搜索速度和全局搜索能力，根据经验知识，确定未知量的搜索步长、搜索范围。给定hv初值、范围、步长，hv∈(0.05, 1.5)，步长Δhv为0.05。然后计算γDB，进而获得地形相位φ0。③设定σ，mDBmax，mDBmin取值范围、初始值、步长。σ∈(0, 10)，步长Δσ为1 dB·m-1；mDBmax∈(-10, 10)，步长ΔmDBmax为0.1 dB；mDBmin∈(-10, 10)，步长ΔmDBmin为0.1 dB。利用设定的hv，σ，mDBmax，mDBmin根据式(5)，计算模拟每一个像元的复相干系数。需要注意的是，消光系数σ在代入模型中要进行转换，其转换系数为0.23。④确定模式搜索法基本参数。确定轴向方向E，加速因子α，缩减率β∈(0, 1)。由于是四个自变量参数，因此，设定轴向方向为E=(0, 0, 0, 1; 0, 0, 0, -1; 0, 0, 1, 0; 0, 0, -1, 0; 0, 1, 0, 0; 0, -1, 0, 0; 1, 0, 0, 0; -1, 0, 0, 0)，加速因子α设为2，缩减率β设为0.5。⑤求出每一个像元8个复相干系数中表真实地表散射的最大复相干系数、最小复相干系数与模拟的复相干系数的距离和的最小值，其具体公式见式(6)。

(6)

基于模式搜索方法求出每一个像元8个复相干系数中表真实地表散射的最大复相干系数、最小复相干系数与模拟的复相干系数的距离的最小值，如果满足2条件(条件1：设定迭代次数为2 000次；条件2：设定函数阈值为0.001)其中之一，终止迭代，提取此时的给定的hv；如果不满足，回到第2步继续运算。⑥遍历每景影像中的每一个像元，从而得到该时相下的水稻株高反演结果。

3 结果与讨论

基于9个时相的TanDEM-X数据，利用第2节所述方法进行研究区水稻株高反演，并利用地块水稻样区的地面数据进行验证与精度评价，如图5所示。

图5 样区1，2，3，4水稻株高反演结果与实测数据对比

The red dotted line represents the measured rice crop height, and the green dots represent the inversion result. What’s more, the error line represents the standard deviation of the inversion result

对样区1来说，在水稻生长前期，第155天(6月4日)、第166天(6月15日)株高实测数据在0.1～0.3 m之间，而反演结果大于0.8 m，明显高于真实值，而且这两个时相反演结果对应的标准差也较大。由此，在水稻生长前期，当株高较小时(<0.3 m)，该方法反演结果较差，株高被过高估计。第177天(6月26日)株高在0.3 m左右，反演结果与实测数据吻合较好，但反演结果标准差较大。在水稻生长中期、成熟期，即第188天(7月7日)、第199天(7月18日)、第210天(7月29日)、第221天(8月9日)、第232天(8月20日)、第243天(8月31日)，株高反演与地面实测值非常接近，绝对差值在0.1 m之内，并且其标准误差也远低于第155天、第166天、第177天的对应值。由此可以看出，当株高大于0.4 m时，该方法反演精度较高，而且稳定性好。

对样区2来说，在第155天、第166天、第188天、第199天，株高实测数据在0.1～0.4 m之间，而这4个时相的反演结果均比其实测结果高0.2 m以上，虽然第177天反演结果较好，但是当株高小于0.4 m时，该方法稳定性较差。第210天、第221天、第232天、第243天株高实测数据在0.4～1.0 m之间，其反演结果与实测株高吻合较好，尤其是第210，221，232和243天，其反演结果与实测值绝对误差在0.1 m以内，而且反演结果标准差也较小。因此，该方法在株高大于0.4 m时，反演结果较好，稳定性较高。

样区3和样区4的反演结果与样区1和样区2类似，当水稻株高小于0.4 m时，反演结果过高估计，当株高高于0.4 m时，反演精度较高(绝对误差<0.15 m)，稳定性也较好。

通过上述4个水稻样区反演结果分析可以看出，当株高低于0.4 m时，反演结果较差，稳定性不高；然而当株高大于0.4 m时，反演结果精度较高而且稳定性较好，与实测值的绝对误差在0.1 m以内，最差不超过0.15 m。图6(a)给出了四个样区9个时相下的水稻株高反演结果。由图6(a)看出，当株高低于0.4 m时，反演结果基本都分布在直线“y=x”上侧，即反演结果被明显高估，绝对误差在0.1～0.8 m之间；当株高大于0.4 m时，反演结果较好，决定系数R2为0.86，均方根误差RMSE为6.79 cm[图6(b)]。

通过反演结果验证与分析发现，基于TanDEM-X数据的水稻株高反演算法，在株高低于0.4 m以下时，效果不好。这主要是由于所使用的TanDEM-X数据，空间基线为2～3 km。根据Cloude[20]提出的株高反演最优范围指数kv=kzhv/2(其中kz为垂直波数，hv为水稻株高)及其最优范围1≤kv≤1.5，结合表1中给出的的kz值，可以得出本文所使用的TanDEM-X数据对应的最佳植被高度反演范围为0.81 m≤hv≤1.21 m。因此，当水稻株高低于0.4 m时，与所使用的TanDEM-X数据对应的最佳植被高度反演范围差距较大，无法较好地反映出水稻植株的较大体散射量，进而导致反演结果较差。

对于株高低于0.4 m的水稻，基于TanDEM-X数据，利用极化干涉的方法进行高度反演需要考虑两方面因素，一是要尽可能的增大基线，提高数据对目标体散射敏感性；二是要保持较高的相干性，不能无限制增大基线。根据频谱偏移超过系统带宽的理论[20]，得出TanDEM-X数据极限空间基线为3.2 km，即当空间基线超过3.2 km时，相干性将为0。根据这一极限条件，TanDEM-X数据，对于高度大于0.33 m的水稻株高反演是有效的。通过kv最优范围计算得出，当TanDEM-X数据基线范围在2.625 km

图6 4个样区9个不同时相水稻株高反演结果与实测水稻株高的对比

4 结 论

利用9个时相的TanDEM-X CoSSC数据，构建了一种基于RVoG模型的极化干涉水稻株高反演算法，实现了多时相的水稻株高反演，结果表明：当株高高于0.4m时，反演精度较高，决定系数R2为0.86，均方根误差RMSE为6.79 cm；然而当株高低于0.4 m时，反演结果被高估而且稳定性差，这主要是因为本研究所使用的TanDEM-X数据空间基线对株高小于0.4 m的水稻敏感性不足。当TanDEM-X数据空间基线范围在2.625 km

致谢：感谢德国航空太空中心(The German Aerospace Center，DLR)(项目编号：NTI-POLI6736)提供本文所用的TanDEM-XCoSSC数据。