基于数值方法的三维激光海面漫反射特性研究

李 恪，姚崇斌，周 勃，尚吉扬

（上海航天电子技术研究所，上海201109）

基于数值方法的三维激光海面漫反射特性研究

李 恪，姚崇斌，周 勃，尚吉扬

（上海航天电子技术研究所，上海201109）

为提高卫星海色遥感技术在实际应用中的效率，对不同海况下不同入射角度的三维激光海面漫反射特性进行了建模仿真计算，并在实验室条件下对激光入射波动水面后产生的散射场能量分布特性进行了研究。根据麦克斯韦尔方程和边界条件分析了散射场在各向分量间的耦合关系，矢量分解后得到了6个标量方程并用矩量法进行离散，给出了离散后的矩阵方程。用稀疏矩阵规则网格法（SMCG）求解矩阵方程，求解时采用共轭梯度法（CGM）计算。在对有限照射区域进行计算时，为减少边界效应产生的误差，用能量强度服从高斯分布的三维锥形波模拟入射激光束。仿真计算了不同粗糙度海面的三维双站散射系数，结果发现：对相同粗糙度海面模型，入射角度越小，粗糙面对激光的漫反射就越大，散射场的能量分布也越均匀。用实验对仿真结果进行验证，结果表明：与基尔霍夫近似（KA）和二维前后向迭代法（FBM）相比，该方法能准确表示波动水面的三维激光漫反射特性，为进一步研究三维激光海面漫反射特性奠定了基础。

粗糙海面；海色；激光；漫反射；三维双站散射系数；稀疏矩阵规则网格法；三维锥形波；快速傅里叶变换

0 引言

海洋占地球表面积的约70%。伴随人类对海洋的开发、利用和保护进程的加快，海洋的资源和环境状况都发生了显著变化。各国对海洋的可持续发展越来越重视，纷纷采用各种高新技术对海洋进行监测和监视。海洋水色卫星主要用于探测海洋水色环境要素，包括叶绿素浓度、悬浮泥沙、海面温度、污染物，以及浅海水深与水下地形等。由于入射激光束经粗糙海面漫反射后，不同偏振状态的激光在各方向上有复杂的耦合关系，因此针对不同海况，获得准确的三维激光海面漫反射特性，对提高卫星海色遥感的准确率有重要意义。因存在耦合关系，为便于计算，以往求解方法多基于一定的近似条件，如KA法［1］。这些近似条件不但对粗糙面的参数有一定的限制，而且会由其固有的物理近似必定带来误差，使计算结果达不到理想效果。FBM虽然是数值计算方法，但其是将三维粗糙界面简化为二维界面后得出的数值模拟，因建立的二维粗糙界面与现实粗糙界面差别较大，故仿真计算结果的准确性较差［2］。为获得准确的激光在波动水面上的漫反射能量分布特性，本文采用SMCG计算激光在波动水面上的双站散射系数（BSC），该算法可减小对计算机内存的需求，而且在每步迭代过程中矩阵与向量的运算采用快速傅里叶变换（FFT）算法，能进一步减少计算时间。

1 激光漫反射原理

假设一束激光入射到波动水面S上，三维激光漫反射的剖面如图1所示。图1中：ki，ks分别为入射激光束和散射激光束的波矢；θi，θs分别为入射角和散射角；ξ（x，y）为波高函数；n为局部法向矢量；k0，k1分别为激光在大气和水中的波数；η0，η1分别为大气和水中的阻抗，且。此处：ε0，μ0和ε1，μ1分别为激光在大气和水中的介电常数与磁导率。

令散射场中的接收点r＝［x y z］，源点r′＝［x′ y′ z′］，则r，r′的距离

式中：Δx，Δy，Δz分别为r，r′在x、y、z向的坐标值之差，Δx＝x－x′，Δy＝y－y′，Δz＝ξ（x，y）－ξ（x′，y′）。则，根据Stratton－Chu公式，以及由麦克斯韦尔方程推导得到散射场的电场和磁场的积分方程与边界条件，对粗糙海面散射场切向和法向进行分析［35］。可得

式中：n′为源点处的局部法向矢量；S′为源点处的局部曲面单元；Hinc（r），Einc（r）分别为入射磁场和入射电场的场量；H，E分别为散射磁场和散射电场的场量；ω为激光在大气中的角频率；i为虚数单位；g0，g1分别为大气和海水中的标量格林函数；为哈密顿算子，且

2 矩阵离散

将光斑在x－y平面上的投影区域Sxy均匀划分成网格Nx×Ny个，该区域的每个采样点位于每个网格的中心处。在x、y方向上，每个网格的边长分别为Δlx＝Lx／Nx，Δly＝Ly／Ny。令m，n分别为该点在网格中的横向和竖向序号，则该区域上的积分顺序如图2所示。由图2可知：采样点的排序是从下到上，再从左至右。

令Fx，Fy分别为n×H在x、y向的分量，Fn为H的法向分量，则

式中：x，y分别为x、y向单位矢量。

同理，令Ix，Iy分别为n×E在x、y向的分量，In为E的法向分量，则

经矢量分解后，由式（2）～（5）可得6个关于标量场的方程。为便于数值计算，本文用矩量法（MOM）离散6个标量方程，离散后的矩阵方程为

式中：p＝1，2，3，…，6，表示6个标量方程的序号；ZAp，…，ZFp分别为第p个方程中不同变量相应的系数矩阵；Up为6个标量方程中的已知量，p＝1， 2，3时，Up＝0，p＝4，5，6时，Up分别为。此处：分别为n× Hinc在x、y向的分量；为Einc的法向分量。

3 稀疏矩阵规则网格法（SMCG）

因源点与接收点的水平距离越近，源点对接收点的贡献就越大，故离接收点近的区域为强区，反之为弱区。SMCG是一种高效率的数值计算方法，它对强区源点的贡献进行详细计算，对弱区源点的贡献则采取近似计算，在保证准确度的同时，极大地节约了计算成本［6］。

令源点与接收点的水平距离Rρ＝为区分强弱区域的临界距离，则Rρ＜rd的区域为强区，Rρ≥rd的区域为弱区。设第（m，n）个点为场点，由图2给出的rd，则图2中的阴影区域为强区，空白区域为弱区。

为便于说明，设存在方程

式中：Z（s），Z（w）分别为强区和弱区矩阵。用稀疏矩阵表示强区矩阵Z（s）；对Z（w）的元素，将其中的格林函数用泰勒级数展开求和的方式近似表示，即

式中：am，bm分别为用泰勒级数对G（R），g（R）进行展开时的第m级系数；下标或上标中的0表示大气中的相应参数，1代表水中的相应参数。

在数值计算中，为获得理想精度，将泰勒级数展开至第6项，M＝5。根据泰勒级数的展开式，可将Z（w）表示为级数求和的形式，即

令式（20）中的第1项为x－y平面的贡献Z（FS），即

则SMCG的迭代过程可简单表述为

定义第n次迭代过程结束后的相对误差

具体求解过程中，本文采用共轭梯度法（CGM）计算［7］。因Z（FS）是一个满矩阵，故CGM在每步的迭代过程中，通过大量计算与向量的乘积。由Z（FS）的表达式可知：这是一个分块托普利兹矩阵，它与向量的乘积满足卷积运算可用FFT实现。令T为Z（FS）的一个分块矩阵，T的第1列为T1，则具体的计算方法为

4 锥形入射波

在对有限照射区域进行数值计算时，为减少边界效应产生的偏差，用具有能量强度服从高斯分布的锥形波模拟入射激光束。设一束激光波数为k0，以方向ki照射在粗糙海面上，则入射电场Einc和入射磁场Hinc可分别表示为

式中：k＝xkx＋yky＋zkz；kx，ky，kz分别为能量在x、y、z向的分量；ETE，ETM分别为光波的水平和垂直偏振状态参数，入射光束为水平偏振时ETE＝1，ETM＝0，入射光束为垂直偏振时ETE＝0，ETM＝1；

其中：kix，kiy分别为ki在x、y向的分量；gx，gy为三维锥形入射波的照射区域控制参数［8］。

入射激光束的总能量

5 数值计算

假设波动水面是一种分布式的随机过程，则采用归一化微分散射截面即散射系数表征粗糙面角度性的散射特征。定义三维激光漫反射在不同偏振状态下的双站散射系数

式中：α，β分别为散射激光和入射激光的偏振状态，α，β＝h表示水平偏振，α，β＝v表示垂直偏振；为不同偏振状态时的散射光强。

通过对双站散射系数进行立体角的积分求和，可得不同偏振状态的入射激光以角度θi入射时产生的散射场的总能量与入射能量的比值，具体可表示为

本文用基于高斯谱的高斯随机过程模拟二维随机粗糙海面［9］。为模拟粗糙海面不同粗糙度，建立了A、B两种海面模型模拟不同粗糙度的波动水面。其中：模型A的均方高度hx＝hy＝0.5cm，相关长度lx＝ly＝6cm；模型B的hx＝hy＝1cm，lx＝ly＝5cm。

为说明三维激光海面漫反射的能量分布特性，以水平偏振为例仿真计算了三维激光海面漫反射的双站散射系数（γ＝γhh）。仿真计算中，分别取θi分别为20°，60°，入射方位角i均为0°。由于激光的能量相当集中，为保证仿真结果的准确性，计算了光斑的全部面积，取参与散射的粗糙面的面积为光斑面积的16倍，则相应的锥形波照射区域的参数和粗糙水面的参数分别为gx1＝gy1＝0.01cm，Lx1＝Ly1＝0.04cm，gx2＝0.02cm，gy2＝0.04cm，Lx2＝0.08cm，Ly2＝0.16cm。令激光波长λ＝0.532μm，为保证激光的电磁特性，在x、y向的采样间隔Δlx＝Δly＝λ／10。

用SMCG进行数值计算，当θi＝20°时，区分强弱区域的临界距离rd1＝40λ；当θi＝60°时，用于区分强弱区域的临界距离rd2＝100λ。为获得散射场能量分布的统计特性，用蒙特卡罗法进行仿真，仿真500次，所得三维激光海面漫反射双站散射系数仿真结果分别如图3、4所示。因Matlab软件中没有直接绘制三维极坐标系中数据图形的命令，直接用θs、散射方位角s作为直角坐标系的变量时，当θs相同时，s＝0°，360°的双站散射系数无法重合。为形象而直观描述散射场三维能量分布的角度特性，将θs，s转换为直角坐标系中的坐标值Px＝sinθs× coss，Py＝sinθssins。

由图3、4可知：当界面的粗糙度相对较小时，散射场的能量分布相对集中，镜面反射方向的峰值能量相对较大，邻近区域的能量相对较小；当界面的粗糙度适度增加时，镜面反射方向的峰值能量减少，邻近区域的能量有所增长，散射场的能量分布更均匀，即散射场能量的覆盖范围变得更大，但在散射场的边缘处，能量衰减速度仍非常快。

随着入射角度的增大，在粗糙面上的光斑面积和粗糙面的截断区域面积都会随之增加，因此就同一粗糙面而言，入射角度不同，粗糙面对激光的漫反射程度亦不同，散射场的能量分布也存在较大差别。比较图3（b）、图4（b）可知：与θi为60°的散射场相比，θi为20°时散射场的能量分布峰值能量相对较弱，总体能量分布更均匀。这说明对相同粗糙度的粗糙面模型，入射角度越小，粗糙面对激光的漫反射程度就会越大，散射场的能量分布就越均匀。

6 实验验证与分析

为验证本文三维激光海面漫反射双站散射系数的数值计算结果，在实验室中进行了实验验证。取θi分别为20°，60°，i为0°，入射高度H＝0.1m；采用功率P＝100mW、光束发散角1mrad、λ＝0.532μm的绿光水平偏振连续激光器。为获得较准确的实验数据，由造波器生成模型A、B两种不同粗糙度的波动水面，用147型激光功率计对散射场中不同方位的采样点的能量进行多次测量，测量结果经处理后与仿真计算结果进行对比。实验中，激光以θi＝60°，H＝0.1m入射在不同粗糙度的波动水面上时，将散射场的能量变化呈现在背景墙上，所得散射场的漫反射能量分布如图5所示。

比较图5（a）、（b）可知：图5（b）的光斑面积明显增大，整个光斑的亮度分布是前向镜面反射方向附近的亮度最高，越靠近边缘，亮度越低，但相对图5（a）而言，光斑中心的亮度明显降低，且在边缘处的能量衰减速度也相对较慢。这说明粗糙度较大的水面使散射场的能量分布更均匀，散射场的覆盖面积明显增大。

由直接观察获得的激光漫反射能量分布特性，与仿真所得激光漫反射能量分布特性是一致的。为获得更客观、量化的激光漫反射能量分布特性，用多次实验测得的数据，对KA，2DFBM，3DSMCG三种方法得到的散射场能量分布特性进行验证、对比和分析。为更清晰地表示散射场能量在前后向上的分布，令s＝0°时的散射角为正，即为前向，s＝180°时的散射角为负，即为后向，则三种不同方法所得的漫反射率以及实验结果分别如图6、7所示。

由图6、7可知：用三维数值计算仿真算得的漫反射率与实验数据更接近。因KA是一种解析算法，只能大致估算散射场的能量分布特性，故是三种方法中最不准确的，甚至在图6（b）中，并未体现出散射场的后向增强现象。二维数值计算的仿真结果大于实验数据。从理论上分析，这主要是因为二维数值计算的仿真结果代表的能量之和是粗糙面上所有源点的贡献之和，即二维散射场的总能量。三维数值计算的仿真结果与真实情况最接近（图6、7中三维漫反射率代表的能量之和只是粗糙面上所有源点产生的贡献的一部分，并非三维散射场的总能量）。虽然三种计算方法的准确度存在差异，但反映的散射场的能量分布特性是相近的。

分析了实验数据与三维双站散射系数存在偏差的原因，认为主要因素是：实验条件不可避免会出现误差；粗糙海面的理论模型与真实水面还有一定差距。

7 结论

本文用数值计算方法对激光在波动水面上的三维漫反射特性进行了研究，并通过实验进行了验证。与其它方法的计算结果进行比较后，认为三维数值计算方法获得的三维BSC是最能准确反映散射场能量分布特性的。研究发现：散射场的能量分布与入射角度和水面粗糙度有密切的关系，区分散射场能量前后向分布大小的入射角度随粗糙度而变，粗糙度越大，入射角的分界点就越大，且较粗糙的水面会使散射场的能量分布更均匀，在散射场的边缘处能量迅速衰减。研究为进一步研究三维激光海面漫反射特性奠定了基础。

［1］ 刘春香，程传福，任晓荣，等.随机表面散射光场的格林函数法与基尔霍夫近似的比较［J］.物理学报，2004，53（2）：427－435.

［2］ 梁华锋，赖建军，彭中良，等.用于气体传感的一维光子晶体窄带红外发射光源［J］.红外与毫米波学报，2009，28（6）：414－417.

［3］ 杨雨川，谭吉春，盛定仪，等.导体板上单极天线对电磁脉冲响应特性的矩量法分析［J］.强激光与粒子束，2008，20（4）：649－652.

［4］ 崔志伟，韩平，李明磊.磁化等离子覆盖导体散射问题的FE／BI方法分析［J］.电波科学学报，2010，（5）：973－978.

［5］ 吴国荣，钟伟芳，吴永东.半空间弹性波强散射双参数非线性反演［J］.华中科技大学学报（自然科学版），2003，31（4）：105－107.

［6］ CHOU H T，JOHNSON J T.A novel acceleration algorithm for the computation of scattering from rough surfaces with the forward－backward method［J］.Radio Sci，1998，33（4）：1277－1287.

［7］ TSANG L，KONG J A.Scattering of electromagnetic waves－numerical simulation［M］.New York：John Wiley，2001.

［8］ 叶红霞，金亚秋.随机粗糙面上介质目标差场散射的快速计算方法［J］.电子学报，2006，34（6）：1074－1079.

［9］ KUANG L，JIN Y.Bistatic scattering from a threedimensional object over a randomly rough surface using the FDTD algorithm［J］.IEEE Transactions on Antennas and Propagation，2007，55（8）：2302－2312.

Study on 3DCharacteristics of Laser Diffusion from Rough Sea Surface Based on Numerical Method

LI Ke，YAO Chong－bin，ZHOU Bo，SHANG Ji－yang
（Shanghai Aerospace Electronic Technology Institute，Shanghai 201109，China）

To improve the application efficiency of ocean color remote sense technology，the characteristics of scattered field of laser from sea surface in different roughness and incident angles were studied.And the energy distribution characteristic of scattered field was analyzed in the lab when laser was incident into sea.The coupling relationships among components of scattering field in different directions were analyzed using Maxwell equation and boundary conditions.Six scalar equations were obtained after vector resolution.The discrete matrix equation was given.The matrix equation was solved by sparse matrix canonical grid method（SMCG）in which conjugate gradient method（CGM）was used in computation.For the calculation of the finite illumination region，the 3Dmodel of incident laser beam was established with 3Dtapered wave to decrease the boundary effect.The 3Dbistatic scattering coefficient（BSC）under different rough sea surface was obtained by simulation.It found that the smaller incident angle，the larger diffusion of the laser and more uniform energy distribution of scattering field for the same rough sea surface.The simulation results were verified by experiment.The result demonstrated that the 3DBSC was enough accurate to express the characteristics of laser diffusion from water surface comparing with Kirchhoff approach and 2Dforward－backward method.The study is useful to the further study about the characteristics of laser diffusion from sea surface.

rough sea surface；ocean color；laser；diffusion；3Dbistatic scattering coefficient；sparse matrix canonical grid method；3Dtapered wave；fast Fourier transform

O451

A

10.19328／j.cnki.1006－1630.2017.01.003

1006－1630（2017）01－0020－07

2016－08－10；

2016－12－20

国家青年科学基金资助（616050110）

李 恪（1980—），男，博士，主要研究方向为光电遥感技术。