目前国内外利用GNSS-R技术对海上溢油的研究尚处于起步阶段，大多研究使用星载方式，并且进展多停留于仿真阶段，基于海岸平台的研究则相对较少。西班牙****Valencia等利用溢油后的大面积黏稠油膜对海面均方斜率(Mean Square Slopes，MSS)的影响，通过仿真生成包含溢油区海面的时延-多普勒图像(Delay-Doppler Maps，DDM)，将DDM进行反卷积处理同时配合空间滤波的方法消除定位模糊，反演出海面散射系数的分布图进而凭此区分溢油区域，并将结果与欧洲航天局(ESA)的Envisat ASAR得到的数据进行比对，得到其分辨率(约2 km)与SAR基本一致，平均误差低于10%，但反卷积处理与振铃效应造成图像的某些地方误差增加至35%[11-12]。Li和Huang则在此基础上进行了更加深入地研究，将原有仿真中卫星高度角为90°的特殊场景扩展为一般场景，并引入了空间积分法(Spatial Integration Approach，SIA)[13]来替代雅可比矩阵，同时应用双天线方式[14]，降低了在DDM去模糊过程中由约束最小二乘(Constrained Least Square，CLS)滤波造成的误差[15-16]。吴军等则在假设反射表面光滑的前提下，提出了在岸基条件下的双天线溢油探测模型[17]。
由于在海岸平台上进行GNSS-R的接收，信号的空间衰减较星载和机载平台更小，设备安装也更易操作。本文针对石油与纯净海水介电性能的固有差异，提出了岸基条件下三天线GNSS-R观测模式的海上溢油探测方法。本文首先论述了GNSS直射与反射信号的数学表达方法，以及信号相关值和反射率的计算方法，分析并结合电磁波多极化现象建立了利用GNSS-R左右旋分量探测水域相对介电常数的反演模型与信号处理流程。接着介绍了基于此溢油探测模型的实验，并对实验数据的处理结果进行了分析与讨论。最后对所提出的探测方法、实验结果及不足之处进行了总结。
1 岸基GNSS-R海上溢油探测原理 纯净海水及几种常见石油产品的相对介电常数范围如表 1所示，虽然不同油品的介电特性不尽相同，但都与海水存在较大差异，因而可以通过反演目标海域的相对介电常数来判断海面的溢油状况。
表 1 海水和不同油品相对介电常数 Table 1 Relative dielectric constant of sea and different kinds of oil
| 物质名称 | 海水 | 原油 | 重油 | 汽油 | 煤油 | 柴油 |
| 相对介电常数 | 68~84 | 2~2.4 | 2.6~3.0 | 1.9 | 2.8 | 2.1 |
表选项
1.1 GNSS信号数学描述 GNSS信号可视为准单色、相位调制球面波信号[18]，并且对于某一固定接收装置而言大多数反射信号经由围绕着镜面点(Specular Point，SP)的第一菲涅尔区反射而来。如图 1所示，θ为卫星高度角，RHCP表示入射右旋极化波，LHCP表示反射信号中左旋分量，在接收点处天线接收到的直射信号ud可以表示为
 |
| 图 1 GNSS-R几何关系 Fig. 1 Geometrical relationship of GNSS-R |
| 图选项 |
 | (1) |
式中：Ad为接收到直射信号的幅度；RT, R为卫星与接收点之间的距离；D(·)为导航电文数据码信号；C(·)为二进制扩频序列；c为光速；i2=-1；t为时间；k为卫星和接收点之间的载波数；fc为信号的载波频率；fD1为由发射机和接收机相对运动引起的多普勒频移。
同样，天线接收到反射信号ur可以表示为
 | (2) |
式中：R为镜面反射点处的菲涅尔反射系数；RT, SP和RSP, R分别为卫星到镜面反射点、镜面反射点到接收点之间的距离；fD2为由发射机和接收机相对运动以及由反射元相对运动引起的多普勒频移。
1.2 直射与反射信号的相关值与反射率的计算 由于GNSS信号属扩频体制，在处理时需要进行解调解扩才可实现信号同步。在某个特定的多普勒频移fD0下，任意时刻t0的本地PRN复制码C(·)与接收天线输出的信号ud, r(·)的一维相关函数定义为
 | (3) |
式中：Ti为相干积分时间；τ为接收信号码片与本地码之间的时间延迟；
 | (4) |
式中：Γpq为p极化入射波对应q极化反射波时的反射率；Δτ为由反射信号相对于直射信号的路径差所产生的时延(岸基平台所接收的反射信号主要来自于镜面反射点附近反射面的反射，反射信号相对于直射信号的时延近似为0，二者的多普勒频移近似相等)；Yr, q(t0+Δτ, 0)为接收为q极化时反射信号的相关峰值；Yd, p(t0, 0)为接收为p极化时直射信号的相关峰值；〈·〉为时间平均；Ar为反射信号幅值；Pr和Pd分别为反射信号功率和直射信号功率。
1.3 多极化现象 通常情况下，GNSS信号入射至地表时反射的电磁波能量可以由散射矩阵确定。而散射矩阵可以通过菲涅尔反射系数表示为如下形式[19]：
 | (5) |
 | (6) |
 | (7) |
式中：菲涅尔反射系数
 |
| 图 2 石油和海水的圆极化反射系数与卫星高度角的关系 Fig. 2 Relationship between reflection coefficient of circular polarization (oil and sea) and satellite elevation |
| 图选项 |
由此可见GNSS信号经地表反射后会产生去极化现象，部分RHCP极性反转；并且LHCP随着卫星高度角的增加而增大，右旋分量则随卫星高度角的增加而减小。
2 岸基GNSS-R海上溢油探测模型 进行海上溢油探测，无疑需要考虑到粗糙度对电磁波的影响。对于海面而言，可以认为下层媒介是近似均匀的，仅仅在海面上发生散射。而海浪是一种十分复杂的现象，在空间和时间上都具有不确定性。海水从海面风场的压力和摩擦作用中吸收能量同时产生波动传递能量，当能量的增加和消耗达到平衡——海浪达到成熟状态时，根据随机过程理论海浪的统计规律已经基本不变，可以将此时的海浪视为一种具有各态历经性的零均值平稳随机过程[20]，因而可以用平稳随机过程理论来分析海浪。自Longuet-Higgins于1952年在波面位移为正态过程和窄谱的假定下推导出海浪波高分布函数以来[21]，众多****在此假定下又提出了多种推导方法，并得到了一致的结论即波高服从瑞利(Rayleigh)分布
 | (8) |
式中：f(ξ)为无因次波高ξ的分布密度函数[22]。
根据Beckmann和Spizzichino[19]于1963年得到的电磁波经由平稳随机过程所产生的粗糙面表面散射后，散射系数与粗糙面高度Z的特征函数χ(vz)的关系为
 | (9) |
式中：
假设E1为入射场，将式(9)两边同乘以
 | (10) |
即
 | (11) |
结合本文1.3节叙述的GNSS-R信号中同时存在左旋和右旋分量的多极化现象，由于无论何种信号分量在同一时刻都经过了相同的反射面，可认为二者都在一定程度上包含了反射面粗糙度特性的信息。这意味着若通过左旋圆极化和右旋圆极化天线分别接收不同极化方式的反射波，分别与直射信号进行功率比得到左旋反射率Γrl和右旋反射率Γrr，再结合式(11)将二者相比，得
 | (12) |
这样便消去了χ(vz)所带来的影响。结合式(5)~式(7)，可将不同极化方式下反射率的比值定义为相对介电常数和卫星高度角的函数f(εr, θ)，即
 | (13) |
因此，在卫星高度角已知的情况下，可以通过测量GNSS-R反射率的极化比来反演目标海域的相对介电常数，进而依据石油和海水相对介电常数的不同探测海上石油的泄漏。
综上所述，岸基条件下通过GNSS-R三天线观测模式探测海上溢油的方法总结如图 3所示。
 |
| 图 3 GNSS直射、反射信号接收处理与海上溢油反演过程 Fig. 3 Process of GNSS direct and reflected signals for retrieval of marine oil spill |
| 图选项 |
通过在海岸平台上的接收端架设用于接收GNSS信号的右旋圆极化天线，以及接收GNSS-R信号的右旋和左旋圆极化天线，将信号分为3个独立通道进行处理。直射信号处理通道完成直射信号的捕获和跟踪，精确获取直射信号的码延时和载波多普勒，使本地信号与直射信号精确同步，并进行导航电文解算，得到卫星高度角θ和一维直射信号相关值。2个反射信号处理通道分别完成反射信号左、右旋分量的开环跟踪，得到一维反射信号相关值。通过相关峰值与反射率之间的关系，得到左旋反射率Γrl和右旋反射率Γrr，并结合天线增益方向图对其进行修正；再根据左、右旋反射率之比与反射面相对介电常数的映射关系，获得目标海域的相对介电常数值，凭此进行溢油情况的判断。
3 实验及分析 3.1 实验 为了探索本文提出的基于GNSS-R技术的海上溢油探测方法是否可行，课题组在2015年9月7日于山东省东营市一处污油处理池(北纬37°52′，东经119°2′)，以及2016年11月28日于东营市刁口镇的一处潮间带(北纬38°8′，东经118°42′)进行了实验。污油池实验中油池长约20 m，宽约5 m，油罐车向池内持续倾倒污油，同时从另一出口流出注入下一级处理池；潮间带实验则选择在涨潮与落潮期间进行GNSS直射、反射信号的采集。
2次实验均采用三天线观测模式，分别于油池短边中间点处和潮间带处架设用于接收GNSS直射信号的RHCP天线、接收反射信号的LHCP和RHCP天线，天线架设高度分别为1.5 m和3.2 m，主瓣方位角分别为180°和167°，其中直射天线斜向上45°放置，反射天线下倾45°放置，保证彼此间无相互遮挡，使用馈线将天线与四通道采集卡相连，再将采集卡连接至笔记本电脑。在设备安装完毕并进行测试后，进行持续的数据观测采集(见图 4)。
 |
| 图 4 实验场地 Fig. 4 Experiment scenario |
| 图选项 |
3.2 结果与分析
3.2.1 污油池实验 结合本文所述的海上溢油探测模型，对采集到的导航卫星直射、反射信号数据进行了相关的处理；同时为保证GNSS-R信号镜面反射点轨迹落于油池内部，结合天线架高、朝向以及采集时段内每颗卫星的高度角、方位角等信息，综合分析了卫星数据的质量，经过进一步筛选后，选取北斗1、3、6、7号星(12:20~12:59时段)，以及北斗1、3、8号星(17:32~18:17时段)的数据。反演得到相对介电常数随观测时间的变化情况如图 5所示。
 |
| 图 5 污油池实验相对介电常数的反演结果 Fig. 5 Inversion results of relative dielectric constant of heavy-oil basin experiment |
| 图选项 |
图 5(a)和(b)中，横坐标的每一个时间点都分别对应了采集时段内的每一分钟。可以看出，在水面存在浮油的情况下利用反射率极化比方式得到的相对介电常数值其范围为2~10。
反演结果的均值及标准差如表 2所示。
表 2 污油池实验反演结果的均值与标准差 Table 2 Mean and standard deviation of inversion results in heavy-oil basin experiment
| 卫星 | 均值 | 标准差 |
| GEO卫星 | 2.34 | 1.28 |
| IGSO卫星 | 2.13 | 1.70 |
| GEO+IGSO卫星 | 2.14 | 1.44 |
表选项
从表 2可以看出，极化比方式下，利用实验采集时段内可用北斗地球同步轨道(GEO)卫星导航信号所反演出的相对介电常数均值为2.34，标准差为1.28；北斗倾斜地球同步轨道(IGSO)卫星信号反演出的相对介电常数均值为2.13，标准差为1.70；所有卫星反演结果的均值为2.14，标准差为1.44。从本次实验结果来看，反演值与重油相对介电常数的理论值(由于污油池中油质主要是原油在提取了汽油、柴油后的剩余重质油，因而本实验采用分子量大、黏度高的重油的相对介电常数作为理论值进行比对)2.6~3.0较为一致，远小于海水的相对介电常数。并且，由于GEO卫星相对地球表面静止的特性，反演结果覆盖了整个采集区间，利用GEO卫星信号反演的相对介电常数标准差比利用IGSO卫星信号反演的标准差小，这表明GEO卫星更有利于实现对固定海域溢油状况的长期持续监测。相比于文献[17]中以假定反射表面光滑为建模前提并直接利用信号左旋反射率来反演相对介电常数，得到的均值为3.6、标准差为2.13的结果而言，本文提出的通过反射率极化比的方法进行处理所得到的结果在一定程度上更加接近石油相对介电常数的真实范围，同时在数值上更加稳定。
此外，从图 5(a)中可以看出反演结果在0~4 min的时段内随着时间的增加整体呈上升趋势，在5~30 min内逐渐下降并趋于平缓；PRN1的数值从7开始微微上升后下降至2，PRN3和PRN6则从开始的3分别增加至5和10左右之后开始减小，在反演结果降低至2后基本保持不变，而PRN7的反演值则在3~4之间略微起伏。从图 5(b)可知，在采集开始时段PRN1、PRN3和PRN8号星的反演结果在1~2之间，直至22 min后(即15:54时左右)呈平缓上升的趋势，并在分别上升至4、8和5后开始略微下降。
分析曲线整体随时间上升及下降的原因：①如图 5(a)所示的采集时段，在初始时刻，污油池内油层由薄变厚，探测表面反射能量逐渐减弱，反射率极化比也随之减小，导致反演结果呈现出略微下降的趋势；图 5(b)所示同理，随着后来油膜厚度逐渐增加，反演结果呈现出随时间增加而逐渐上升的现象。②污油池中的乳状液并不是纯粹的油(由于浮油下层为含有杂质的水，使得油层内也含也有一定量的水)，根据文献[23]提出的油水混合物的介电常数公式:
 | (14) |
式中：εo和εw分别为油和水的介电常数。
可知，油水混合物的介电常数εm会根据混合物中水的体积分数w的不同而有所波动，即含水量越高，混合物介电常数越高，反之则越低。
3.2.2 潮间带实验 图 6为潮间带实验的反演结果, 从图中反演的相对介电常数来看，在潮汐的作用下，伴随着涨潮这一过程，海滩的相对介电常数由采集之初的13左右逐渐增长至55，在维持一段时间(伴随有一定程度的起伏)之后，随着海水的退去，被水淹没的沙滩又逐渐重新露出水面，此时段的反演结果又逐渐减小，回归至26以下。需要说明的是，在卫星高度角一定时(PRN1为GEO卫星)，由于天线架设较高，使得信号镜面反射点离天线架较远，在涨潮时段内天线所接收到的总能量中包含了一小部分高潮线以上的沙滩所反射的功率，因此相对应的反演数值介于沙滩和海水之间，与纯海水的相对介电常数理论值存在一定差异。
 |
| 图 6 潮间带实验相对介电常数的反演结果(PRN1) Fig. 6 Inversion results of relative dielectric constant of intertidal zone (PRN1) |
| 图选项 |
综上，利用本文提出的海上溢油探测方法，污油池实验的反演结果与石油相对介电常数的理论值基本相一致，潮间带实验的反演结果与海水相对介电常数的理论值较为接近；并且在水面存在油膜时反演得到的相对介电常数值远小于潮间带在潮水处于平潮期时段所反演出的相对介电常数值，表明本文所提出的探测方法具有一定的可行性和有效性。
4 结论 本文在分析GNSS-R信号特点与多极化现象以及考虑海面粗糙度因素的情况下，建立了GNSS-R左、右旋信号分量的反射率与相对介电常数的映射关系，提出了一种利用GNSS-R技术在海岸平台下对小规模海域进行溢油探测的方法；并结合本文所建立的模型，通过在污油处理池和潮间带进行实地实验，利用接收设备进行信号的采集与处理，反演得到的相对介电常数分别与石油和海水相对介电常数的理论值基本吻合，同时二者存在明显区别，表明在岸基条件下该方法在溢油监测领域具有一定的可行性。不过，由于受到实验环境及条件的限制，污油池内重油表面的粗糙度并不完全吻合海洋表面的粗糙度，同时在潮间带实验中风力与海浪的等级也较低，今后将在高风速下进行海上溢油探测实验，深入研究以进一步优化反演模型与算法。
可以预见，随着研究的不断深入以及对海上溢油探测需求的逐渐增大，具有高时空分辨率、全时段覆盖、设备功耗低、抗干扰等诸多优点的新兴无源遥感技术GNSS-R，在探测海上石油泄漏这一领域的应用前景十分广阔。同时，伴随卫星导航系统的飞速发展与日益完善，尤其是GEO卫星的加入，使得对固定海域溢油情况进行长期稳定地监测成为可能。
参考文献
| [1] | HALL C D, CORDEY R A. Multistatic scatterometry[C]//Preceeding of the IEEE Intemational Geoscience and Remote Sensing Symposium. Piscataway, NJ: IEEE Press, 1988: 561-562.http://ieeexplore.ieee.org/xpls/abs_all.jsp?arnumber=570200 |
| [2] | ALONSO-ARROYO A, CAMPS A, PARK H, et al. Retrieval of significant wave height and mean sea surface level using the GNSS-R interference pattern technique:Results from a three-month field campaign[J].IEEE Transactions on Geoscience & Remote Sensing, 2015, 53(6): 3198–3209. |
| [3] | RODRIGUEZ-ALVAREZ N, AKOS D M, ZAVOROTNY V U, et al. Airborne GNSS-R wind retrievals using delay-Doppler maps[J].IEEE Transactions on Geoscience & Remote Sensing, 2013, 51(1): 626–641. |
| [4] | CARDELLACH E, FABRA F, NOGUéS-CORREIG O, et al. GNSS-R ground-based and airborne campaigns for ocean, land, ice, and snow techniques:Application to the GOLD-RTR data sets[J].Radio Science, 2011, 46(6): 3604–3629. |
| [5] | FABRA F, CARDELLACH E, RIUS A, et al. Phase altimetry with dual polarization GNSS-R over sea ice[J].IEEE Transactions on Geoscience & Remote Sensing, 2012, 50(99): 1–10. |
| [6] | SABIA R, CAPARRINI M, RUFFINI G. Potential synergetic use of GNSS-R signals to improve the sea-state correction in the sea surface salinity estimation:Application to the SMOS mission[J].IEEE Transactions on Geoscience & Remote Sensing, 2007, 45(7): 2088–2097. |
| [7] | JACOBSON M D. Dielectric-covered ground reflectors in GPS multipath reception-Theory and measurement[J].IEEE Geoscience & Remote Sensing Letters, 2008, 5(3): 396–399. |
| [8] | RODRIGUEZ-ALVAREZ N, BOSCH-LLUIS X, CAMPS A, et al. Soil moisture retrieval using GNSS-R techniques:Experimental results over a bare soil field[J].IEEE Transactions on Geoscience & Remote Sensing, 2009, 47(11): 3616–3624. |
| [9] | JI W, XIU C, LI W, et al. Ocean surface target detection and positioning using the spaceborne GNSS-R delay-Doppler maps[C]//IEEE International Symposium on Geoscience and Remote Sensing. Piscataway, NJ: IEEE Press, 2014: 3806-3809.http://ieeexplore.ieee.org/xpl/articleDetails.jsp?reload=true&arnumber=6947313 |
| [10] | 夏威. 基于遥感技术的海上溢油应急响应[J].中国水运, 2015(9): 20–23. XIA W. Response of marine oil spill emergency based on remote sensing technique[J].China Water Transport, 2015(9): 20–23.(in Chinese) |
| [11] | VALENCIA E, CAMPS A, MARCHAN-HERNANDEZ J F, et al. Ocean surface's scattering coefficient retrieval by delay-Doppler map inversion[J].IEEE Geoscience & Remote Sensing Letters, 2011, 8(4): 750–754. |
| [12] | VALENCIA E, CAMPS A, RODRIGUEZ-ALVAREZ N, et al. Using GNSS-R imaging of the ocean surface for oil slick detection[J].IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2013, 6(1): 217–223. |
| [13] | LI C, HUANG W. Sea surface oil slick detection from GNSS-R delay-Doppler maps using the spatial integration approach[C]//Radar Conference. Piscaway, NJ: IEEE Press, 2013: 1-6.http://ieeexplore.ieee.org/xpls/abs_all.jsp?arnumber=6585990 |
| [14] | PARK H, VALENCIA E, RODRIGUEZ-ALVAREZ N, et al. New approach to sea surface wind retrieval from GNSS-R measurements[C]//Geoscience and Remote Sensing Symposium. Piscataway, NJ: IEEE Press, 2011: 1469-1472.http://ieeexplore.ieee.org/xpls/icp.jsp?arnumber=6049344 |
| [15] | LI C, HUANG W. Simulating GNSS-R delay-Doppler map of oil slicked sea surfaces under general scenarios[J].Progress in Electromagnetics Research B, 2013, 48(48): 61–76. |
| [16] | LI C, HUANG W, GLEASON S. Dual antenna space-based GNSS-R ocean surface mapping:Oil slick and tropical cyclone sensing[J].IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2015, 8(1): 425–435. |
| [17] | 吴军, 张波, 洪学宝, 等. 基于北斗卫星反射信号的海面溢油探测方法及试验[J].北京航空航天大学学报, 2017, 43(1): 193–199. WU J, ZHANG B, HONG X B, et al. Sea surface oil spill detection method using BeiDou satellite reflected signal and its experiment[J].Journal of Beijing University of Aeronautics and Astronautics, 2017, 43(1): 193–199.(in Chinese) |
| [18] | 杨东凯, 张其善. GNSS反射信号处理基础与实践[M].北京: 电子工业出版社, 2012: 76. YANG D K, ZHANG Q S. GNSS reflected signal processing:Fundamentals and applications[M].Beijing: Publishing House of Electronics Industry, 2012: 76.(in Chinese) |
| [19] | BECKMANN P, SPIZZICHINO A. The scattering of electromagnetic waves from rough surfaces[M].New York: Pergamon Press, 1963: 90-200. |
| [20] | 文圣常. 海浪原理[M].济南: 山东人民出版社, 1962: 148, 167-172. WEN S C. The principle of the waves[M].Jinan: Shandong People Press, 1962: 148, 167-172.(in Chinese) |
| [21] | LONGUET-HIGGINS M S. On the statistical distribution of the heights of sea waves[J].Journal of Marine Research, 1952, 11: 245–266. |
| [22] | 孙孚. 三维海浪要素的统计分布[J].中国科学:数学, 1988, 31(5): 501–508. SUN F. Statistical distribution of three dimension ocean wave[J].Science in China:Mathematics, 1988, 31(5): 501–508.(in Chinese) |
| [23] | 罗海清. 基于MSP430的原油含水率测量[D]. 大连: 大连理工大学, 2003: 6-7. LUO H Q. Measurement of water content of oil based on MSP430[D]. Dalian: Dalian University of Technology, 2003: 6-7(in Chinese).http://cdmd.cnki.com.cn/article/cdmd-10141-2003061578.htm |
