叶 秋 李 钢 张仁斌

空中复杂结构目标杂散辐射快速计算算法

叶 秋 李 钢 张仁斌

(合肥工业大学计算机与信息学院 安徽 合肥 230009)

针对在自然光环境下天地背景对空中目标的辐射作用，引入光学系统中的杂散辐射信息的计算，结合有限元算法采用并改进传统光线追迹计算模型求解目标散射信息。并且利用Kd-tree算法对面元结构目标建立空间分割描述结构实现加速计算，同时借助Modtran大气传输模型模拟计算自然环境下目标表面的入射光信息，从而利用BRDF模型计算获得复杂目标在天地背景中自然光辐射的散射信息。仿真计算结果表明，该方法相对于有限元求解算法具有相同的计算准确度，并且提升了6到10倍传统杂散辐射计算方法的计算速度。对目标杂散辐射的计算结果表明该方法对目标光学探测、识别和导弹制导等方面的研究具有很高的应用参考价值。

空中目标复杂结构 背景辐射 杂散辐射 Kd-tree算法

0 引 言

空中目标的光学特性分析作为空中光学遥感、目标探测和识别的重要理论和数据基础[1,2]，为现有空中目标检测方面以及目标探测与识别系统的波段选择优化以及雷达等探测器的探测能力分析提供依据，同时也对导弹成像制导、目标的探测、跟踪与识别等空间突防技术具有十分重要的工程应用价值[3-5]。空中目标散射信息的辐射源包括天空大气和地面反射以及太阳的辐射，天地背景的辐射作用对典型目标的复杂外形如挂载武器、进风口凹腔、旋转翼等处易产生一种非直接散射即杂散辐射，杂散辐射对于光学系统而言，会直接影响成像面的对比度，降低像面的信噪对比度，并且可能在成像时产生光斑，导致成像质量下降，影响导弹成像制导、空间监控系统的探测识别。因此研究空中目标的杂散辐射信息对研究空中目标光学特性具有重要的应用参考价值。

国内外对空中目标光学散射特性进行了大量研究。其中国内研究主要针对空中目标进行建模求解[6-8]，提出天地背景辐射近似计算，以及应用计算几何面元网格化思想[9]，提出并实现基于面元网格化的空间目标光学特性计算方法。总结国内研究应用的主要计算模型有理论计算模型[11,12]、目标有限元计算模型[6,7]，但这两种算法只考虑单次入射而未对目标体面间的多次散射予以计算，所以得到的结果不适合精确散射信息分析。相对国内而言，国外较早对目标散射特性进行广泛研究，提出了多种计算思想，如微扰法计算模型，基尔霍夫切平面法，蒙特卡罗法等[10]，以上思想在具体应用时主要分为[13,14]：蒙特卡洛追迹方法、能流密度追迹法和传统光线追迹法。在上述方法中充分考虑了光线在目标体面间多次散射作用，并且对目标的杂散辐射进行了准确计算，但计算模型构建以及计算过程相对复杂，尤其前两种方法针对性较强，一般主要用于目标典型特性与边界观测条件下的局部计算，并不适用于一般散射信息计算情况。

综合上述杂辐散射计算方法的优缺点，传统光线追迹法因为计算模型相对简单而且相对其它计算杂散辐射的方法更具实时性，因此符合实时应用需求。然而为了保证计算准确度，由大量光线入射造成的大计算量成为其应用的瓶颈，本文针对该算法的计算瓶颈，提出采取构建Kd-tree结构对空中复杂目标进行空间划分，从而对复杂结构目标散射信息的计算进行加速。

1 复杂目标对背景辐射的散射的计算

1.1 空中典型目标的杂散辐射源分析

在背景辐射下典型目标产生的散射信息中，不仅包括直接散射信息也包括间接散射信息。其中间接散射在可见光通道称为杂散光，若在红外通道则称作杂散辐射。如图1所示，光学探测中除了成像光线外，同样存在作用于探测器外表面上的非成像光线辐射，以及通过非正常散射光线轨迹抵达探测器的成像光线辐射。一般对空中目标的辐射作用的辐射源可分为内部作用和外部作用，所以空中典型目标的杂散辐射源按照来源可以分为三种[15]：第一种为目标外部的背景辐射源，如太阳光照射、地面的散射光和漫射光以及大气漫射光等，在辐射能量到达作用目标后，如图2所示，经过典型目标内部构件的多次反射、折射或衍射最终到达探测器，此类辐射源称为外部杂散辐射或外杂光，这也是本文研究的空中背景环境对典型目标辐射的主要辐射源；第二种是典型目标内部辐射源，像目标与空气摩擦产生的热源、目标发动机运行时的热能以及温度较高的大功率部件等产生的红外波段辐射，在经过目标表面的反射、折射以及衍射后抵达探测器，这部分光能称之为内部杂散辐射或内杂光；第三种是成像光线经非正常传播而抵达探测器的辐射能量，我们称之为成像杂散光。

图1 背景辐射产生杂散辐射示意

图2 面元间散射作用

1.2 构建复杂结构散射信息计算模型

光线追迹法是解决光线入射后产生的杂散信息的计算的常用方法，在计算目标光学散射信息时，往往为保持结果准确度而采用大量光线密集入射，由此造成的几何级倍数的计算成本的增加是制约其广泛应用的瓶颈。针对光线追迹优化加速的改进关键点有两点：一是减少非兴趣光线入射的数目，二是加快光线与面元求交的速度。对于第一点可首先采用适应网格划分算法迭代缩小兴趣光线范围，可以有效地减小非兴趣入射光数量，另外针对最后出射面元光线可利用有限元方法采用一般光线替代，只针对二次及以上反射作用的光线作单独追迹处理；第二点关键是使用加速数据结构来减少求交次数，一旦建立好加速数据结构，光线追迹会更加有效率。

在空中目标天地背景辐射源确定的情况下，假定探测接收方向上的追迹光线为平面光波，根据Multi-Grid算法[15,16]对空中目标模型在接收方向上的投影进行网格划分。如图3所示，将入射平面波均匀划分成射线管网格，并且保证射线管之间互相平行。其中每根射线管包括5个光线矢量，并假定射线管中心的射线包含全部光线能量,此能量将用来判断是否终止继续追迹。此外为保证计算结果的准确，在划分网格时要确定所有入射光线的射线管截面的面积总和要大于所测目标模型的最大包围盒在探测接收方向上的投影面积。

图3 典型目标射线Multi-Grid分割前后示意

针对采用Multi-Grid算法精简首次入射光线数量还是会带来大数量的光线处理的状况，本文采用最后出射光线一般化处理，只针对二次反射及以上的光线采用单独光线追迹，从而大大减少计算量。如图4所示。

图4 a:有限元方法处理示意；b:本文算法处理示意

图4中最后出射光线由于是平面光光波，有限元方法采用面元消隐法即将光线L1替换所有同等角度方向的出射光线，并由面元S1面积与光线密度乘积计算最终总的出射信息，本文利用同样方法计算最后出射光线，与此同时针对二次即以上反射的出射光线如图4中光线L2利用自定义光线密度参数采样追迹，然后利用Kd-tree算法进行后续相交求解。本文算法通过以上操作保证了与有限元法具有相同的计算精度，而且由于杂辐散射的引入，使得到的散射信息更符合实际光线作用结果。

Kd-tree算法是由Harvan于2001年提出的[17]，它有效解决了一些传统追迹算法中对目标划分时间长、划分的空网格多、查找网格节点效率低等问题。如图5所示，Kd-tree是一种自上而下递归的对空间进行剖分的数据结构，同时是一个非均匀的网格剖分算法。根据目标空间的面元分布利用与坐标轴垂直的平面进行剖分，并且当每步递归剖分构建的时候，将一个节点处理成内部节点或叶节点。如果该节点包含的面元数比自定义的阈值小或者节点深度大于预设的深度阈值时则停止划分该节点，另外当子结点中面片数较多，但是本身已经达到规定分割大小的时候，也要停止进行划分[18]。在内部节点中，最优划分位置是基于光线追踪成本预估值，该成本包括内部节点遍历成本和叶节点光线与面元的求交成本。成本预估采用著名的Surface Area Heuristic(SAH)算法[19]。其中内部节点N的遍历成本公式CN为：

CN=Ct+[nlCiP(Nl|N)+nrCiP(Nr|N)]

(1)

其中nl、nr是内部节点N左右子树Nl、Nr所含面片的数目，Ct为遍历花费，Ci为单次光线相交计算花费。面片分割示意见图6所示。

图5 左：一个二维Kd-tree.内部节点其最近分割面标记;右：图形表示左边的Kd-tree

图6 坐标轴上划分候补，每个面元拥有两个待分割位置，如a1、a2是面元A的待分割点

1.3 复杂结构目标散射信息的计算流程

根据文献[6]在使用MODTRAN4软件进行背景的辐射亮度计算时，以天顶角135°的地面背景辐射近似代替下半球空间的地面或海面背景平均辐射；以天顶角45°的天空背景辐射近似代替上半球空间的天空背景平均辐射；天空太阳辐射模型因为相对简单则单独计算后叠加。

本文采用五参数半经验统计模型来进行参数优化建模[20,21]，该统计模型是针对Torrance-Sparrow模型进行修正，是将粗糙表面起伏分布和斜率分布拓展到非高斯分布的一般情况。目标的背景辐射如图7所示，目标几何模型在地球坐标系中建立,即设地球坐标系和目标坐标系重合。图中面元本地法向量z′(n)向上，Ki,sky和Ki,earth分别为地表和天空辐射，Ks为某方向上散射光波，计算公式如下：

图7 面元本地入射散射角

(2)

式(2)中λ1和λ2分别代表波谱起始值和终止值，Li在计算时分别带入天空和地表辐射，对方位角积分求解，最后对波谱积分求解即可。

利用该模型计算最终获得典型目标散射信息的角信息图。另外计算过程中参照的太阳方位角以正北为准，角方向以顺时针方向为正；观测角同样以正北为准，逆时针为正角方向。目标散射信息计算的主要计算流程如图8所示。

图8 目标散射信息计算流程图

1.4 计算实例与分析

下面以表面为铝蜂窝材料的某带挂载战机目标为例，计算目标对太阳辐射和天地背景辐射的散射，对比使用本文算法和有限元算法的计算结果。计算波段为：3～5 μm，8～12 μm。设飞机目标位于东经117°，北纬31°，高度2 km，时间为2010年7月10日上午10时；大气模式为：1976美国标准大气；地表温度为:290 K。BRDF计算模型采用五参数经验模型。计算设备为：CPU为Intel(R)Dual-Core，主频2.20 GHz，程序计算了围绕目标360°，每1度设1个数据点，共360个接收方位的角分布曲线，接收方位与目标姿态如图9所示，图中箭头曲线为接收方位变化路径，圆点为起始位置。计算结果如图10-图15所示。

图9 某战机在地球坐标系下不同观察角度

从图10中得出，在3～5μm波段战机姿态1时,太阳辐射相对天地背景辐射对目标辐射能量较大,目标战机对天地背景辐射的散射亮度在机背近330°上出现峰值。由此可见,太阳辐射作用较天地背景辐射更为明显即贡献值更大，而且散射亮度出现的峰值方位与太阳入射方向有关，由图10和图11也可以看出，在目标的机背被太阳辐射的可见部分的散射亮度较大，并且由图10-图12散射结果可以得出本文算法与有限元算法结果差异不大；在8～12 μm波段，由于太阳辐射与天地辐射的此消彼长而且由于天地背景辐射都基本上均匀分布、强度相当,并且飞机自身具有对称性，目标的散射亮度基本呈对称分布。由图13-图15散射信息结果可以看出，本文算法结果在目标的凹腔，机翼挂载等处观察到散射亮度与有限元算法计算结果出现明显差异，这是由于有限元算法因为采取面元消隐，只考虑单次辐射作用因而得到的结果较平滑，而本文算法因为考虑到天地背景辐射下，周围光线在入射以上部位后进行多次作用产生的杂散辐射的缘故，在上述部位所得结果出现峰值。相对前者，本文算法结果更符合实际情况从而更准确。

图10 3～5 μm波段战机在姿态1两种算法散射亮度结果对比图

图11 3～5 μm波段战机在姿态2两种算法散射亮度结果对比图

图12 3～5 μm波段战机在姿态3两种算法散射亮度结果对比图

图13 8～12 μm波段战机在姿态1两种算法散射亮度结果对比图

图14 8～12 μm波段战机在姿态2两种算法散射亮度结果对比图

图15 8～12 μm波段战机在姿态3两种算法散射亮度结果对比图

综合以上实验信息可以得出在研究目标对天地背景辐射的散射特性时,在3～5 μm波段,对目标辐射作用贡献中，太阳辐射贡献突出而天地背景辐射贡献可忽略不计。而且目标对太阳辐射散射的峰值角度与太阳入射角度相当，并且由于典型目标背部较平滑，几乎没有杂散辐射的产生；而在8～12 μm波段,此时天地背景辐射增强而太阳辐射贡献相对减弱，并且由于天地背景辐射的均匀分布以及空中目标本身外形具有对称性,从而导致目标的散射亮度呈现对称分布，同时在凹腔、机翼挂载等部位，由于辐射光线在目标内部体面间的多次作用产生的杂散辐射信息较其它部位明显增强，而且散射亮度的峰值和方位与目标外形结构密切相关。由此得出：目标的几何外形、观测目标的方位等因素相互作用,共同决定着空中目标的散射特性。

表1是针对单个面元，两个面元，8463个面元分别采用本文改进追迹算法和Octree光学追迹在采用OPENMP优化下的计算的时间比较，计算花费的时间均为多次测量平均值。从表中可以看出，在模型面片较少时采用本文方法计算速度约为采用Octree方法计算速度的2倍左右，随着面元数增加，两者计算时间都随之增加，而采用本文方法的光线追迹算法较采用Octree算法计算速度提高到3倍左右，体现了改进算法在计算大型模型时的优势。

表1 计算SH-2直升机散射信息花费时间

2 结 语

本文通过对传统光学追迹方法进行改进和加速，实现了对空中目标在天地背景辐射下散射信息的加速计算，并且与主流的有限元算法计算结果对比，本文算法提高了计算准确度，为获取和分析空中目标在红外和可见波段的散射特性提供了有效的手段，与此同时在目标和环境红外、可见光成像制导、空中监控以及实时场景仿真等方面具有应用价值。但是本文模型的精确度依旧存在着一些误差，主要是没有考虑内部辐射作用和对光线密度选取造成的，后续工作主要集中在计算目标内部辐射作用和对入射光线密度与结果精确度的关系进行建模选取最佳数值，综合实时性和精确性进而改善空中目标散射信息计算模型。

[1] 李道勇,王增玉,张艳群.空间目标可见光散射特性与可视条件研究[J].光学与光电技术,2005,3(5):183-187.

[2] 张伟,汪洪源,王治乐,等.空间目标可见光散射特性建模方法研究[J].光子学报,2008,37(12):2462-2467.

[3] 高国旺,刘上乾,秦翰林,等.复杂背景下的红外目标自动跟踪算法[J].光电工程,2010,37(6):78-83.

[4] 刘文哲,张婉怡,董会,等.光学相关红外目标识别算法研究[J].仪器仪表学报,2011,32(4):850-855.

[5] 许晓航,肖刚,云霄,等.复杂背景及遮挡条件下的运动目标跟踪[J].光电工程,2013,40(1):23-30.

[6] 李良超,牛武斌,吴振森.空中复杂目标对背景红外辐射的散射的并行计算[J].系统工程与电子技术,2011,33(12):2573-2576.

[7] Wu Zhensen,Liu Anan.Scattering of solar and atmospheric background radiation from a targets[J].Int JIR and Mill Waves,2002,3(6):903-917.

[8] 鲍文卓,丛明煜,张伟,等.基于面元网格化的空间目标光学特性计算方法[J].哈尔滨工业大学学报,2010,42(5):710-715.

[9] Moldosanov Kamil A,Kashirin Victor A.Black coating forstray light and thermal control application[C]//SPIE,2001,4458:87-94.

[10] Rask J D.Modeling of Diffuse Photometric Signatures of Sattellites for Space Objects Identification[D].Ohio State,USA:Air Force Institute of Technology,1983:20-119.

[11] Zissis G J.The Infrared & Electro-optical Systems Handbook vol. 1 Sources of Radiation[M].Bellingham,Wash:SPIE Optical Engineeing Press,1993:151-157.

[12] Zhao N,Xue Y,Wang J.Analysis of stray radiation from infrared optical system with monte-carlo method[J].Chinese Journal of Optics and Applied Optics,2010,3(6):665-670.

[13] 黄强.空间光学系统的杂散光分析[J].红外,2006,27(1):26-33.

[14] 原育凯.光学系统杂散光的消除方法[J].大气与环境光学学报,2007,2(1):6-10.

[15] Suk S,Seo T I,Park H S,et al.Multiresolution grid algorithm in the SBR and its application to the RCS calculation[J].Microwave and Optical Technology Letters,2001,29(6):394-397.

[16] Bang J K,Kim B C,Suk S H,et al.Time consumption reduction of ray tracing for RCS prediction using efficient grid division and space division algorithms[J].Journal of Electromagnetic Waves and Applications,2007,21(6):829-840.

[17] Roccia J P,Paulin M,Coustet C.Hybrid CPU/GPU KD-Tree Construction for Versatile Ray Tracing[C]//Eurographics (Short Papers),2012:13-16.

[18] Jiang Huiyan,Zhao Yudong,Li Ning.The Study of 3D Reconstruction Method Based on Dynamic Threshold Method and Improved Ray Casting Algorithm[C]//Washington: IEEE Computer Society. Proceedings of the 2008 International Conference on Internet Information Hiding and Multimeda Signal Processing,2008:402-405.

[19] MacDonald J D,Booth K S.Heuristics for ray tracing using space subdivision[J].The Visual Computer,1990,6(3):153-166.

[20] 张涵璐,吴振森,张昌民,等.BRDF的遗传算法和遗传模拟退火算法建模及比较[J].系统工程与电子技术,2010,32(7):183-185.

[21] 吴振森,谢品华.粗糙表面激光散射统计建模的遗传算法[J].光学学报,2002,22(8):897-901.

FAST ALGORITHM OF CALCULATING STRAY RADIATION OF AERIAL TARGETS WITH COMPLEX STRUCTURE

Ye Qiu Li Gang Zhang Renbin

(SchoolofComputerandInformation,HefeiUniversityofTechnology,Hefei230009,Anhui,China)

Aiming at the radiation effect of sky and earth background on aerial targets in natural light environment, we introduced the stray radiation information computation in optical system, and calculated the target scattering information in combination with the finite element algorithm as well as using and improving traditional ray tracing calculation model. Moreover, we adopted the Kd-tree algorithm to establish the spatial segmentation description structure for bin structural targets to achieve the acceleration of calculation, at the same time, with the help of Modtran atmospheric transmission model we simulated the computation of the incident light information of target surface in natural environment, so that used the computation of BRDF model to obtain the scattering information of natural light radiation of complex targets in sky and earth background. Simulation computation results showed that the method in the paper had the same computational accuracy relative to finite element solution algorithm, and raised the calculation speed by six to ten times of that of traditional stray radiation calculation method. Result of the target stray radiation calculation showed that this method had very high applied reference value for researches in regard to target optical detection, recognition and missile guidance.

Aerial target complex structure Background radiation Stray radiation Kd-tree method

2014-10-20。国家自然科学基金项目(61271121)；中央高校基本科研业务费专项资金项目(J2014HGBZ0129)。叶秋，硕士，主研领域：系统仿真，目标探测。李钢，博士。张仁斌，博士。

TP391.9

A

10.3969/j.issn.1000-386x.2016.04.014