楊 鑫,魏 兵
(1.西安電子科技大學物理與光電工程學院,陜西西安 710071; 2.西安電子科技大學信息感知協(xié)同創(chuàng)新中心,陜西西安 710071)
一種射線跟蹤建模的新方法
楊 鑫1,2,魏 兵1,2
(1.西安電子科技大學物理與光電工程學院,陜西西安 710071; 2.西安電子科技大學信息感知協(xié)同創(chuàng)新中心,陜西西安 710071)
在射線跟蹤的研究中,針對計算區(qū)域進行網(wǎng)格劃分問題,以單位球內(nèi)接二十面體為基礎(chǔ),提出了一種基于降維度計算的空間曲面三角形網(wǎng)格劃分方案.通過對球面的三角形劃分,介紹了該方法的具體實現(xiàn)過程,給出了兩種網(wǎng)格的編號方案;結(jié)合電磁輻射的方向圖對射線進行篩選;給出了曲面三角形劃分以及射線選擇的算例.數(shù)值算例證明了該算法的有效性.
射線跟蹤;降維方法;網(wǎng)格劃分;射線選擇
在射線跟蹤法(Ray Tracing,RT)的計算中,首先要解決的問題就是計算域的網(wǎng)格劃分問題,即將輻射區(qū)域劃分為一定數(shù)量的射線管束.然后,再對射線進行跟蹤計算便可得到特定的目標區(qū)域的電磁響應(yīng)參量.而射線的跟蹤計算通??煞譃檎蚋櫽嬎愫头聪蚋櫽嬎銉深怺1-6],其中,正向跟蹤計算的實現(xiàn)方法是:先對劃分所得的所有射線做一遍跟蹤計算,然后將有貢獻的射線的計算結(jié)果累加至目標區(qū)域.這種方法的優(yōu)點是,方便簡潔且易于程序的實現(xiàn),但同時也需要花費大量的時間計算對目標區(qū)域沒有貢獻的射線,因此效率較低;而反向射線跟蹤計算的方法為:首先篩選出對目標區(qū)域有貢獻的射線,然后再依次對這些篩選出來的射線進行跟蹤計算,便可得到目標區(qū)域的電磁響應(yīng)參量.這種方法克服了正向跟蹤計算效率較低的缺點,但是具有明顯的數(shù)學以及程序?qū)崿F(xiàn)的復雜性,特別是在多徑問題的跟蹤計算中.
射線跟蹤法常將經(jīng)過輻射球面的射線管束劃分為三角形、四邊形以及六邊形網(wǎng)格形式,但是三角形網(wǎng)格是更為基礎(chǔ)和常用的形式,在其基礎(chǔ)上就可方便地實現(xiàn)其他網(wǎng)格形式的劃分.因此,該問題的實質(zhì)是曲面的三角形網(wǎng)格劃分的問題.而在曲面的三角形網(wǎng)格劃分方面,前人已做了許多的工作[7-14].文獻[7]在重心一致Delaunay三角化方法基礎(chǔ)上得到了平面的三角網(wǎng)格劃分.文獻[8]在最優(yōu)Delaunay三角化方法的基礎(chǔ)上將曲面劃分為局部網(wǎng)格和全局網(wǎng)格,并通過減小插值誤差的方法實現(xiàn)了網(wǎng)格的平滑優(yōu)化.文獻[9]采用epsilon-net方法和泊松表面重構(gòu)法實現(xiàn)了對巖石表面的網(wǎng)格劃分及近似.文獻[10]采用調(diào)和映射的方法對三角網(wǎng)格曲面進行參數(shù)化,提出一種三角網(wǎng)格曲面等參數(shù)線刀軌生成算法.文獻[11]在將曲面劃分為可展曲面和不可展曲面的基礎(chǔ)上,利用Watson算法實現(xiàn)了曲面的三角網(wǎng)格劃分.上述的網(wǎng)格劃分方法普適性好,但其數(shù)學過程和程序?qū)崿F(xiàn)一般比較復雜,在具體的工程應(yīng)用中實現(xiàn)起來難度較大.因此,針對上述問題以及為滿足射線跟蹤法的需要,筆者以球內(nèi)接二十面體為基礎(chǔ),提出了一種基于降維度計算的空間曲面的三角形網(wǎng)格劃分方案.該方案的思路是:將三維空間曲面的網(wǎng)格劃分轉(zhuǎn)化為三維空間平面的網(wǎng)格劃分,然后將三維空間平面轉(zhuǎn)化為二維平面的網(wǎng)格劃分,接著把二維空間平面的劃分轉(zhuǎn)化為一維線段的劃分,這樣就可得到空間曲面的三角形網(wǎng)格劃分.而對于球面特定部分需要劃分的場景,可對二十面體進行選擇性的劃分,以此便可在一定程度上達到提高計算效率的目的.
射線跟蹤法中通常需要建立全向模型,但不需要進行全向計算,而只需要計算特定方向上的空間范圍.比如在計算天線輻射場時,就可根據(jù)天線的輻射方向圖對空間進行劃分,在有輻射場方向上對空間進行計算,而對沒有輻射場或是輻射場很弱的空間范圍則無需計算.鑒于此,筆者設(shè)計了一種基于正向跟蹤的網(wǎng)格選擇方案以及網(wǎng)格編號方法.該方案首先對球面劃分所得的三角網(wǎng)格進行編號,然后根據(jù)輻射方向圖對特定方向上的射線進行選擇跟蹤,從而提高了跟蹤計算的效率.
球面三角形網(wǎng)格劃分在其內(nèi)接正二十面體的基礎(chǔ)上實現(xiàn)[6,15].以單位球面網(wǎng)格劃分為例,首先,在待劃分球面內(nèi),內(nèi)接正二十面體(如圖1所示);然后,可對二十面體的各個頂點進行編號,以便于后續(xù)網(wǎng)格計算.在球坐標系下,可方便地求出二十面體的各頂點的坐標值[15]:A(1,0,0)、B(1,π,0)、C(1,a,0)、D(1,a,a)、E(1,a,2a)、F(1,a,3a)、G(1,a,4a)、H(1,π-a,2.5a)、I(1,π-a,3.5a)、J(1,π-a,4.5a)、K(1,π-a,0.5a)、L(1,π-a,1.5a),其中,a=1.107 149.這樣,對球面的網(wǎng)格劃分就可轉(zhuǎn)化為對二十面體的網(wǎng)格劃分,二十面體的網(wǎng)格劃分結(jié)果只需要徑向投影至球面就可得到球面的網(wǎng)格劃分.而這轉(zhuǎn)化過程只需用直角坐標與球坐標系之間的互換關(guān)系就可以實現(xiàn).
圖1 單位球內(nèi)接正二十面體
單位球面內(nèi)接二十面體的三角面編號為:△ACG為T1,△ADC為T2,△AED為T3,△AFE為T4,△AGF為T5,△GCJ為T6,△CKJ為T7,△CDK為T8,△DLK為T9,△DEL為T10,△EH L為T11,△EFH為T12,△FIH為T13,△FGI為T14,△GJI為T15,△BJK為T16,△BKL為T17,△BHL為T18,△BIH為T19,△BJI為T20.
圖2 母三角形3階網(wǎng)格劃分及編號結(jié)果
以上述20個三角形中的任意一個三角形為母三角形,采用對其邊取中點的方式可得到高階劃分的子三角形網(wǎng)格(如圖2所示).例如,對母三角形的每條邊取中點,就可得到母三角形的1階網(wǎng)格劃分,從而形成4個子三角形網(wǎng)格;然后,再對所得的4個子三角形以取每邊中點的方式可得到2階網(wǎng)格劃分,形成16個子三角形;以此類推,就可得到母三角形的更高階網(wǎng)格劃分.這樣,原二十面體的每個母三角形在進行n階網(wǎng)格劃分后,其所含的子三角形的行數(shù)為2n,每行含2m-1個子三角形(m為該子三角形所在行的行編號),而每個母三角形所含子三角形總數(shù)為22n.于是每個子三角形可取(i,j,k)編號形式,其中,i為母三角形編號(i≤20),j為子三角形在母三角形中的行編號(j≤2n),k為該子三角形在其所在行內(nèi)的位置編號(k≤2m-1).對于每個格點的編號,可采用如下方法:以順時針方向?qū)γ總€子三角形的3個點排序,每個點可以編號為(i,j,k,l),其中,l為取值1、2或3,分別代表三角網(wǎng)格的3個組成點.另外,還可采取直接用點編號方式代替對三角網(wǎng)格的編號.對于進行n階網(wǎng)格劃分后的每個母三角形中的點可編號為(i,j,k),其中,i為該點所在母三角形的編號(i≤20),j為該編號點所在的行編號(j≤2n+1),k為其所在行內(nèi)的位置編號(k≤m,m為該點所在的行編號即行數(shù)).
在文中的計算中,采用上述兩種編號結(jié)合的方法,即分別給三角網(wǎng)格點和網(wǎng)格面編號.首先,由網(wǎng)格劃分方法計算出每個網(wǎng)格點的坐標值;然后,再由網(wǎng)格點計算相應(yīng)的三角網(wǎng)格的面積;最后,得到平面及曲面的所有網(wǎng)格劃分點和網(wǎng)格面積.
對空間曲面進行網(wǎng)格劃分,可先將曲面劃分為若干個對應(yīng)的平面,然后可對這些平面進行網(wǎng)格劃分,最后將劃分好的平面投影到曲面之上,即可得到曲面的網(wǎng)格劃分.由此可見,對曲面網(wǎng)格劃分的關(guān)鍵是如何得到平面的網(wǎng)格劃分.下面給出文中算法的具體實現(xiàn)方案.
圖3 網(wǎng)格劃分的降維示意圖
以空間平面三角形網(wǎng)格劃分為例,如圖3所示,在直角坐標系Oxyz中,對三角形ABC進行n階三角網(wǎng)格劃分,則可先將其投影到任意兩個坐標面(圖3中的投影面選為x Oz和x Oy面)形成二維投影三角形(△ABCxz和△ABCxy),然后對兩投影三角形進行n階網(wǎng)格劃分,最后將劃分所得的網(wǎng)格投影坐標點合成,即可得到空間平面三角形的n階網(wǎng)格劃分點.而對于二維的投影三角形,其n階網(wǎng)格點劃分則可轉(zhuǎn)化為一維線段的n階劃分:首先,對三角形的3邊進行劃分,計算其n階網(wǎng)格劃分點;然后,由3邊計算所得的劃分點再計算該三角形內(nèi)部所需的劃分點.其中,三角形的內(nèi)部所需的劃分點可以分為3個類型,即分別與3邊平行的線段上的劃分點;最后,由這3類線段和三角形3邊的n階網(wǎng)格劃分點,合成即可得到坐標投影平面內(nèi)的二維三角形的n階網(wǎng)格劃分點.圖4給出了三角形頂點分別位于點A(-50,0,50)、點B(0,70,0)以及點C(50,0,50)處的3階網(wǎng)格劃分點的分布情況.
圖4 三角形的3階網(wǎng)格劃分點分布
這樣,利用上述的降維方法,就得到空間三角形的網(wǎng)格劃分.然后將該方法依次作用于球面內(nèi)接二十面體的各個表面之上,就可得到二十面體的n階網(wǎng)格劃分點.
在進行高階劃分后的二十面體的網(wǎng)格點,可利用坐標變換公式
將計算所得的劃分點投影到球面上,即可得球面的內(nèi)接高階多面體及相應(yīng)的球面三角網(wǎng)格劃分.圖5給出了4階網(wǎng)格劃分下的球面網(wǎng)格點分布圖.
對于其他類型的曲面,同樣可先通過將其劃分為一定數(shù)量的空間平面母三角形,然后再利用文中所使用的方法對這些平面母三角形網(wǎng)格劃分,最后再投影回原曲面上,即可得到該曲面的任意階網(wǎng)格劃分.圖6給出了曲面z=x2+x y3的5階網(wǎng)格點分布情況.
圖5 球面的4階網(wǎng)格劃分點分布
圖6 曲面的5階網(wǎng)格劃分點分布
在射線跟蹤法中,為提高跟蹤計算的效率,常將跟蹤射線分為有效射線和無效射線.嚴格上來講,有效射線即是對接收點有能量貢獻的射線;反之,沒有貢獻的射線即稱為無效射線.而廣義上來說,即便射線對接收點有能量貢獻,但若該貢獻小于一定的閾值,則依然可將其視為無效射線.但總的說來,可用設(shè)定閾值的方法來劃分有效射線和無效射線.把射線進行如此劃分之后,就可對射線進行選擇性的跟蹤,以此就可達到提高計算效率的目的,即對有效射線進行跟蹤計算而放棄對無效射線的跟蹤,這樣,就需要解決如何選擇跟蹤射線或者網(wǎng)格點的問題了.
在無線電波的射線跟蹤問題中,電波一般可視為天線所發(fā)出的,因此,電波的分布可用方向圖來描述.而方向圖一般為(θ,φ)的函數(shù),這時就可通過控制θ和φ變量,實現(xiàn)格點的選擇.圖7是矩形波導天線的輻射方向圖,在上述6階網(wǎng)格劃分的情形下會生成10 242根射線.利用文中算法可自由地對上述射線根據(jù)強度選擇,從而提高射線跟蹤的計算效率.不同閾值下有效射線數(shù)量、原有射線跟蹤所需時間與選擇后射線跟蹤所需時間的比值(加速比)如表1所示.由表1可見,隨著閾值的降低,需要跟蹤的射線條數(shù)有所增加,然而相比于全向跟蹤算法選擇后的計算時間明顯降低.
圖7 矩形方孔輻射方向圖
表1 控制閾值與加速比觀察表
在某些電磁輻射問題的射線選擇中,可能還存在方向函數(shù)難以計算的問題,比如非常復雜的方向函數(shù)或者難以尋找甚至是不存在的方向函數(shù).對于這些問題,由于球面是在內(nèi)接二十面體的基礎(chǔ)上實現(xiàn)網(wǎng)格劃分的,所以,對于上述問題中的射線選擇就可以轉(zhuǎn)化為對二十面體的面的選擇,即計算有輻射母三角面而忽略沒有輻射的母三角面的計算,以此也可達到提高計算效率的目的.
最后,對于已知場點的情況下尋找有效輻射面元的問題,則可根據(jù)該場點坐標值(R,θ,φ)找到電磁波的運動軌跡,然后以其軌跡與球面求交點的方式即可找到該點的有效輻射面元.而對于直射問題,面元的尋找就顯得相對容易,只需找出與場點有相同θ和φ分量的輻射面元即可.
為解決電波射線的空間劃分以及提高計算效率等問題,文中給出了一種根據(jù)射線跟蹤的空間曲面的三角形網(wǎng)格劃分方法,在該方法的基礎(chǔ)上也可容易地實現(xiàn)曲面的四邊形以及六邊形網(wǎng)格劃分;同時,在已知電波的輻射特性的條件下,比如輻射方向圖或者電波的大致覆蓋范圍等,文中還設(shè)計了一種根據(jù)正向射線跟蹤法建立的射線選擇方案,該方案有效剔除了無需計算的冗余射線,從而提高了射線跟蹤的計算效率.最后,整個方案的數(shù)值結(jié)果表明了文中算法的正確有效性,并且由該方法的實現(xiàn)原理和結(jié)果亦可表明,其具有較好的普適性.總的來說,該方法具有簡單高效、易于工程中程序?qū)崿F(xiàn)的特點.
[1]張忠波.基于射線跟蹤技術(shù)的室內(nèi)電波傳播預測研究[D].西安:西安電子科技大學,2012.
[2]de ADANA F S,BLANCO O G,DIEGO I G,et al.Propagation Model Based on Ray Tracing for the Design of Personal Communication Systems in Indoor Environments[J].IEEE Transactions on Vehicular Technology,2000,49(6): 2105-2112.
[3]劉忠玉,郭立新,種稚萌,等.城市微蜂窩環(huán)境下一種改進的射線跟蹤預測模型[J].西安電子科技大學學報,2014,41 (2):137-143. LIU Zhongyu,GUO Lixin,ZHONG Zhimeng,et al.Improved Ray Tracing Prediction Model in Urban Microcellular Environments[J].Journal of Xidian University,2014,41(2):137-143.
[4]李祥震,韓香娥.均勻橢球粒子的彩虹角分析[J].西安電子科技大學學報,2010,37(4):731-736. LI Xiangzhen,HAN Xiange.Rainbow Angle Analysis of a Homogeneous Spheroid[J].Journal of Xidian University,2010,37(4):731-736.
[5]楊紹華,張福順,焦永昌.天線罩電磁特性的仿真分析[J].西安電子科技大學學報,2004,31(6):873-876. YANG Shaohua,ZHANG Fushun,JIAO Yongchang.An Analysis of the Parameters for Antennas with the Dome[J]. Journal of Xidian University,2004,31(6):873-876.
[6]袁正午.移動通信系統(tǒng)終端射線跟蹤定位理論與方法[M].北京:電子工業(yè)出版社,2007.
[7]WANGA B,KHOOC B C,XIE Z Q,et al.Fast Centroidal Voronoi Delaunay Triangulation for Unstructured Mesh Generation[J].Journal of Computational and Applied Mathematics,2014,280:158-173.
[8]CHEN L,HOLST M.Efficient Mesh Optimization Schemes Based on Optimal Delaunay Triangulations[J].Computer Methods in Applied Mechanics and Engineering,2011,200:967-984.
[9]LAI P,SAMSON C,BOSE P.Surface Roughness of Rock Faces Through the Curvature of Traingulated Meshes[J]. Computers and Geosciences,2014,70:229-237.
[10]陳曉兵,廖文和,戴寧.三角網(wǎng)格曲面等參數(shù)線刀軌生成算法[J].中國機械工程,2013,24(8):1047-1051. CHEN Xiaobing,LIAO Wenhe,DAI Ning.Algorithm for ISO-parametric Tool Path Generation for Triangular Mesh Surface Machining[J].Chinese Mechanical Engineering,2013,24(8):1047-1051.
[11]陳永府,張華,陳興.任意曲面的三角形網(wǎng)格劃分[J].計算機輔助設(shè)計與圖形學學報,1997,9(5):396-401. CHEN Yongfu,ZHANG Hua,CHEN Xin.Triangle Meshing for Surface[J].Journal of Computer Aided Design and Computer Graphics,1997,9(5):396-401.
[12]胡曉娟,盧兆林,葛德彪,等.基于三角面元的涂層目標FDTD共形網(wǎng)格生成技術(shù)[J].系統(tǒng)工程與電子技術(shù),2010,32 (9):1884-1888. HU Xiaojuan,LU Zhaolin,GE Debiao,et al.Conformal FDTD Mesh-generating Scheme for Coated Targets Based on Triangle-patch[J].Systems Engineering and Electronics,2010,32(9):1884-1888.
[13]ZHONG Z C,LIANG S,MAO J,et al.Anisotropic Surface Meshing with Conformal Embedding[J].Graphical Models,2014,76(5):468-483.
[14]CHEN X,CHEN L,SHI M D.A Highly Solid Model Boundary Preserving Method for Large-scale Parallel 3D Delaunay Meshing on Parallel Computers[J].Computer-Aided Design,2015,58:73-83.
[15]王毅.極低頻電磁波時域方法的理論研究及其應(yīng)用[D].南京:南京航空航天大學,2012.
(編輯:齊淑娟)
New modeling method for ray tracing
YANG Xin1,2,WEI Bing1,2
(1.School of Physics and Optoelectronic Engineering,Xidian Univ.,Xi’an 710071,China;2.Collaborative Innovation Center of Information Sensing and Understanding,Xidian Univ.,Xi’an 710071,China)
It is often necessary to divide the spatial domain into grids in the ray tracing method.Based on the unit sphere inscribed icosahedron in this paper,a method of triangulation of the space curved surfaces by dimension reduction calculation is proposed.The process of this method is presented in detail by the triangulation of a unit sphere and two grid numbering schemes are given simultaneously.Then,combined with the direction of the electromagnetic radiation figure,the rays are selected.Finally,the examples of surface triangulation and tube selection are presented,which show the validation and effectiveness of this algorithm.
ray tracing;dimension reduction;surface meshing;ray selection
TN011
A
1001-2400(2016)03-0073-05
10.3969/j.issn.1001-2400.2016.03.013
2015-02-20
時間:2015-07-27
國家自然科學基金資助項目(61231003,61571348,61401344)
楊 鑫(1988-),男,西安電子科技大學博士研究生,E-mail:e?yangx@126.com.
魏 兵(1970-),男,教授,E-mail:bwei@xidian.edu.cn.
http://www.cnki.net/kcms/detail/61.1076.TN.20150727.1952.013.html