高宇潔 陳 健 朱 琳 李雪蓮 李春暉 陳德海

（南京信息工程大學(xué)遙感與測繪工程學(xué)院，江蘇 南京 210044）

0 引言

卷云可以通過散射太陽輻射、吸收地表長波輻射來調(diào)控地-氣系統(tǒng)的溫度。政府間氣候變化專門委員會(huì)（Intergovernmental Panel on Climate Change，IPCC）第四次評估報(bào)告[1]關(guān)注了卷云所具有的氣候效因。傳統(tǒng)方法以實(shí)心勻質(zhì)的球形粒子為前提，利用Mie散射理論研究卷云輻射特性，但實(shí)際卷云多為六角形粒子，因此該理論不適用。近年來，一些學(xué)者開始深入研究非球形粒子的輻射特性。例如王玉文等人利用T-matrix方法建立了對稱非球形冰晶粒子的衰減模型[2]。吳舉秀等人使用離散偶極子近似法（Discrete Dipole Approximation，DDA）研究了取向、形狀和溫度對94 GHz冰晶散射效率的影響[3]。王穎等人研究了非球形粒子的散射矩陣[4]。類成新研究了隨機(jī)取向六角形粒子的散射特征量及敏感程度[5]。探測卷云時(shí)，可見光與近紅外穿透性過弱，微波穿透性過強(qiáng)，而太赫茲波段（30 μm~3 000 μm）[6]和卷云冰晶粒子（50 μm~400 μm）位于同一量級，能夠獲取更全面的回波信息，理論上為探測卷云的最佳波段[7]。該文將利用DDA算法模擬球形粒子在太赫茲頻段下的輻射結(jié)果，并與相同頻率下Mie散射結(jié)果進(jìn)行對比，以評估DDA方法的計(jì)算精度；在該基礎(chǔ)上計(jì)算了340 GHz六角形冰晶粒子的輻射特性，重點(diǎn)分析了不同縱橫比和不同空間取向?qū)椛涮匦缘挠绊憽?/p>

1 DDA方法原理介紹

DDA方法是當(dāng)前被廣泛使用的數(shù)值計(jì)算方法，非常適合用來研究六角形粒子的散射和吸收特性。它可以通過大量離散偶極子模擬計(jì)算不同入射方向、不同入射波長以及不同形狀情況下粒子的散射和吸收系數(shù)[8]。其基本思想如下：散射特性較連續(xù)的物體可近似為可極化點(diǎn)陣。為了從該近似點(diǎn)陣中獲得偶極子，可對該點(diǎn)陣中各個(gè)點(diǎn)對局域電場（入射場及其他點(diǎn)的輻射場都包括在內(nèi)）進(jìn)行研究。在上述方法中，單個(gè)粒子的散射特性是連續(xù)的，因此將這些具有連續(xù)散射特性的單個(gè)粒子轉(zhuǎn)化為多個(gè)離散偶極子陣列，這些偶極子始終具有呈周期性排列及相互作用的特性。在整個(gè)過程中，點(diǎn)陣中每一個(gè)點(diǎn)受局域電場的影響生成偶極矩，最終這些點(diǎn)的輻射總和形成了散射場。

在近似點(diǎn)陣中，已知共有n個(gè)偶極子，假設(shè)第k個(gè)偶極子的坐標(biāo)為rk（k=1，2，…，n），該處電場強(qiáng)度為Ek，偶極矩（即極化強(qiáng)度）為Pk；第j個(gè)偶極子的坐標(biāo)為rj（j=1，2，…，n），該處電場強(qiáng)度為Ej，偶極子極化率為αj，則該點(diǎn)的偶極矩Pj如公式（1）所示。

式中：i為虛數(shù)單位；k為波數(shù)，（λ為入射波長）；rjk為第j個(gè)偶極子與第k個(gè)偶極子間的距離，；I為 3×3 的單位矩陣；為第k個(gè)偶極子到第j個(gè)偶極子的單位向量。

通過公式（3）可解出Pk，由此得出第j個(gè)偶極子電偶極矩為Pj，在該情況下能夠獲得該核殼結(jié)構(gòu)粒子的消光截面Cext、吸收截面Cabs、散射截面Csca以及其他輻射特征量。

該核殼結(jié)構(gòu)粒子的消光截面Cext、消光系數(shù)Qext如公式（4）、公式（5）所示。

式中：E0為原點(diǎn)處的電場強(qiáng)度；E*inc,j為第個(gè)偶極子的入射波電場強(qiáng)度的復(fù)共軛；Qsca為散射系數(shù)。

吸收截面Cabs、吸收系數(shù)Qabs如公式（6）、公式（7）所示。

散射截面Csca、散射系數(shù)Qsca如公式（8）、公式（9）所示。

式中：dΩ為立體角微元；為散射方向單位矢量；aeff為粒子等體積球的等效半徑，；Einc為入射波電場強(qiáng)度；*為該值取復(fù)共軛[9]。

研究結(jié)果表明，離散偶極子近似法在微粒折射率m≈1時(shí)其精度最高，m取值離1越遠(yuǎn)，計(jì)算精度越低。當(dāng)存在較大的折射率|m|或kaeff的情況時(shí)，如果想更精準(zhǔn)地模擬粒子的光散射，就需要更多數(shù)量的偶極子[10]。

2 DDA方法的精度評估

以往DDA的研究多是針對毫米波、可見光的，而針對太赫茲的研究較少，該文首先對大氣透過率較高的340 GHz頻段下DDA的計(jì)算精度進(jìn)行分析。由于Mie散射理論是由麥克斯韋方程組推導(dǎo)出的求解球形粒子散射特性的經(jīng)典算法，因此可以將Mie散射結(jié)果作為標(biāo)準(zhǔn)值，以評估相同條件下DDA方法的計(jì)算結(jié)果精度。鑒于計(jì)算機(jī)性能和DDA算法的限制，該文模擬的球形冰晶粒子等效半徑范圍為50 μm~1 000 μm，步長為10 μm，球形冰晶粒子在340 GHz頻段下的復(fù)折射指數(shù)為1.782+0.0071i（i為虛數(shù)單位），其他參數(shù)選用DDA的默認(rèn)設(shè)置。圖1描述了340 GHz頻段下球形冰晶粒子的輻射特性隨等效半徑變化而變化的規(guī)律。

圖1 球形冰晶粒子輻射結(jié)果

由圖1（a）可知，在340 GHz頻段下由DDA算法和Mie散射理論計(jì)算得到的球形冰晶粒子散射曲線幾乎完全重合。等效半徑小于600 μm時(shí)，散射效率相對誤差穩(wěn)定，在0.5%上下浮動(dòng)；大于該范圍后，隨著等效半徑的增加，相對誤差波動(dòng)上升，最高大約為5%。由圖1（b）可知，在該頻段下2種方法計(jì)算得到的粒子吸收曲線大部分都重合，當(dāng)?shù)刃О霃叫∮?00 μm時(shí)，吸收效率相對誤差較穩(wěn)定，最大不超過2%；大于該范圍后，相對誤差增加，在10%附近上下波動(dòng)且等效半徑為930 μm時(shí)的吸收相對誤差最大。

綜上所述，在340 GHz頻段下隨著等效半徑的增加，散射和吸收效率的相對誤差逐漸增大。對DDA算法來說，這是由于偶極子數(shù)目變少而導(dǎo)致的。無量綱尺度參數(shù)（如公式（10）所示）描述了給定的散射問題，而DDA算法不適用于尺度參數(shù)x非常大的粒子。由粒子大小aeff和偶極子間距d之間的關(guān)系（如公式（11）所示）可知，當(dāng)粒子等效半徑增大時(shí)，偶極子間距也變大，為了滿足約束條件，偶極子數(shù)目將會(huì)減少，計(jì)算結(jié)果的精度也隨之下降。因此計(jì)算高頻率下的大尺度粒子時(shí)得到的結(jié)果誤差較大。

比較圖1（a）和圖1（b）還可以發(fā)現(xiàn)，冰晶粒子的吸收效率遠(yuǎn)小于散射效率，因?yàn)榱Ｗ颖旧砦蛰^小，此時(shí)偶極子數(shù)量減少所帶來的相對誤差就比較明顯。綜合圖1可以得到以下結(jié)論：當(dāng)有效半徑小于600 μm時(shí)，無論是散射和吸收，其相對誤差均小于2%，因此使用DDA方法仍然會(huì)具有較高精度。當(dāng)有效半徑大于600 μm時(shí)，誤差隨半徑增加而增加，但由于吸收效率的絕對值較小，因此仍然可以使用DDA算法，但其可靠性會(huì)降低。

3 六角形冰晶粒子散射和吸收特性研究

與球形粒子相比，六角形粒子的縱橫比、空間取向不同，對應(yīng)電磁波的散射和吸收特性也不同，因此有必要研究不同縱橫比、不同空間取向以及不同粒子大小情況下冰晶粒子的散射和吸收特性，為準(zhǔn)確地進(jìn)行卷云參數(shù)遙感反演提供依據(jù)。該文通過DDA方法建立六角形冰晶粒子的物理模型，計(jì)算在340 GHz頻率下不同縱橫比、不同空間取向下不同有效半徑的粒子散射和吸收效率，進(jìn)而分析縱橫比和空間取向的影響。其中縱橫比分別取4、2、1、1/6和1/20；空間取向的極化角φ、β都取0°，θ角取值范圍為0°~90°且步長為1°。

3.1 球形粒子與六角形粒子的比較

假設(shè)卷云中存在的冰晶粒子均為球形，此處選取極化角θ為0°、縱橫比為1的球形冰晶粒子，即塊柱狀六角形冰晶粒子為代表。如圖2（a）所示，描述了2種冰晶粒子的散射特性隨有效半徑增大而波動(dòng)變化的特點(diǎn)，當(dāng)有效半徑較大時(shí)，塊柱狀六角形冰晶粒子的散射特性比球形冰晶粒子的波動(dòng)幅度更大。如圖2（b）所示，描述了2種冰晶粒子的吸收特性隨有效半徑增大而波動(dòng)增加的變化特點(diǎn)，但球形冰晶粒子與六角形冰晶粒子的吸收特性在粒子大小不同的情況下，其數(shù)值并不一致。綜上所述，球形與六角形的冰晶粒子在340 GHz頻段下會(huì)產(chǎn)生不同的散射和吸收特性，進(jìn)一步證明太赫茲為探測卷云的最佳波段，如果將卷云粒子假設(shè)為球形粒子會(huì)帶來較大的誤差，因此該假設(shè)不成立，無法使用球形冰晶粒子代替六角形冰晶粒子。

圖2 球形與六角形冰晶粒子輻射結(jié)果

3.2 不同縱橫比的影響

六角形粒子的形狀可以由高L和邊長a來描述，縱橫比定義為L/2a，當(dāng)縱橫比分別為4、2、1、1/6和1/20時(shí)，六角形粒子可以分別稱為長柱、厚柱、塊柱、厚板和薄板。該文研究了在340 GHz頻段下不同縱橫比粒子的散射和吸收特性，如圖3所示。

圖3 六角形冰晶粒子輻射結(jié)果

由圖3可見，在340 GHz頻段下，縱橫比為4的六角形冰晶粒子散射曲線波動(dòng)很大，存在3個(gè)峰值，分別位于360 μm、560 μm以及830 μm處；縱橫比為2、1的六角形冰晶粒子散射特性曲線波動(dòng)相對明顯，存在2個(gè)峰值，且隨著有效半徑的增加，其峰處對應(yīng)數(shù)值變?。欢v橫比為1/6、1/20的粒子散射特性總體趨勢隨有效半徑的增加而增加，在有效半徑大于400 μm的部分散射數(shù)值大于另外3個(gè)散射數(shù)值。相比于散射特性，不同縱橫比下的六角形冰晶粒子的吸收特性曲線變化規(guī)律更一致，在所有縱橫比情況下都呈現(xiàn)波動(dòng)上升的趨勢，縱橫比越大，曲線的峰值越多，且相同半徑下不同縱橫比的粒子吸收數(shù)值相差不大。

綜上可得，在340 GHz頻段下相同空間取向的冰晶粒子縱橫比越大，輻射特性曲線越陡峭、峰值越多且散射曲線變化越復(fù)雜，而吸收曲線變化規(guī)律相同且都為上升趨勢。吳舉秀等人[3]對六角形冰晶的94 GHz毫米波后向散射特性進(jìn)行模擬、計(jì)算和分析，在冰晶后向散射效率對形狀的敏感性方面得出了相似的結(jié)論。當(dāng)粒子的縱橫比遠(yuǎn)離1時(shí)，它的形狀偏離球形，相同尺度下粒子的最大截面積也越大，因此厚板和薄板粒子的后向散射也就越大。

3.3 不同空間取向的影響

當(dāng)粒子的形狀為非球形時(shí)，其散射和吸收特性參數(shù)與光束入射角有關(guān)，相同粒子在固定的坐標(biāo)系中會(huì)因空間取向不同而引起不同的光學(xué)特性，因此在研究非球形粒子的散射和吸收特性時(shí)應(yīng)該考慮其空間取向?？臻g取向指非球形粒子在空間所擺放的位置，建立如圖4所示的坐標(biāo)系：系統(tǒng)坐標(biāo)系和固定于粒子的當(dāng)?shù)刈鴺?biāo)系a1-a2-a3[11]，取向主要由極化角θ、φ和β來決定，其中θ為a1與軸之間的夾角，φ為a1與平面之間的夾角，β為a2與Xa1平面之間的夾角。該文選取340 GHz頻段下、縱橫比為2的厚柱和縱橫比為1/6的厚板狀六角形冰晶粒子在空間中從橫放到豎放的變化過程，模擬正常情況下不同入射角時(shí)粒子空間取向?qū)ζ漭椛涮匦缘挠绊懀碚撋蠘O化角θ與入射角的大小相同，因此θ從0°到90°以1°為步長各取1個(gè)空間取向進(jìn)行分析，其他極化角為0°。

圖4 非球形粒子在系統(tǒng)坐標(biāo)中的取向

當(dāng)縱橫比為1/6時(shí)，在340 GHz頻段下不同空間取向的厚板狀六角形冰晶粒子的散射和吸收特性如圖5所示（以空間取向?yàn)?°、15°、30°、45°、60°、75°以及90°的輻射特性為代表）。由圖5（a）可知，當(dāng)粒子等效半徑范圍為50 μm~1 000 μm時(shí)，由于空間取向的變化，隨著粒子等效半徑的增大，散射特性曲線的變化趨勢各不相同，整體上曲線先波動(dòng)上升，增加到一個(gè)峰值后再波動(dòng)下降。并且隨著空間取向的增加，散射特性曲線變化越來越劇烈、峰值處對應(yīng)的等效半徑越來越小。在橫縱比相同的情況下，當(dāng)粒子尺度較小時(shí)，空間取向值大的粒子的散射效率大于空間取向小的粒子；隨著粒子尺度的增加，空間取向值大的粒子的散射效率逐漸小于空間取向值小的粒子。

由圖5（b）可知，吸收特性曲線整體上呈現(xiàn)出上升的趨勢；隨著空間取向的增加，曲線越來越陡峭，且在粒子尺度較大的部分，波動(dòng)越來越劇烈。在粒子尺度小于500 μm 的部分，空間取向越大的粒子的吸收效率就越大；而后半部分的粒子吸收效率隨著空間取向的增大先增加后減少。

圖5 六角板的散射和吸收特性

當(dāng)縱橫比為2時(shí)，在340 GHz頻段下不同空間取向的厚柱狀六角形冰晶粒子的散射和吸收特性如圖6所示。由圖6（a）可知，隨著空間取向的增大，散射特性曲線由先上升到峰值后下降又繼續(xù)上升的多峰值波動(dòng)趨勢，變?yōu)檎w波動(dòng)先上升后下降的趨勢，即曲線的整體趨勢從兩端凸起、中間下凹變成中間凸起。并且隨著空間取向的增大，散射峰值先減少后增加，曲線變化越來越平緩。

由圖6（b）可知，吸收特性曲線波動(dòng)上升，隨著空間取向的增大，曲線波動(dòng)由較劇烈變?yōu)槠骄彙Ｇ以谙嗤臻g取向的情況下，隨著粒子尺度的增加，曲線波動(dòng)越來越明顯。

圖6 六角柱的散射和吸收特性

綜上所述，在縱橫比相同的情況下，散射和吸收特性隨粒子等效半徑變化而變化的趨勢與空間取向有關(guān)。當(dāng)粒子尺度較小時(shí)，不同空間取向?qū)α切伪ЯＷ虞椛涮匦缘挠绊戄^小，整體上隨等效半徑的增加而增加；當(dāng)粒子尺度較大時(shí)，不同取向引起的散射截面差異被放大。這是因大尺度粒子具有較大的散射截面，微小的空間取向差異會(huì)使其散射截面產(chǎn)生較大的變化而導(dǎo)致的。對比柱狀和板狀的六角形冰晶粒子對應(yīng)的散射和吸收特性曲線與空間取向的關(guān)系可以發(fā)現(xiàn)，六角板冰晶粒子受空間取向的影響較大，其原因可能為柱狀冰晶比板狀冰晶更接近球形，其光學(xué)特性受空間取向影響的敏感性與粒子形狀有關(guān)。粒子縱橫比偏離1越多，其形狀越偏離球形，不同取向下散射截面的差異也就越大，光學(xué)特性受空間取向影響越敏感[3]；粒子縱橫比越接近1，其形狀越接近球形，輻射特性就越近似于各向同性，光學(xué)特性受空間取向的影響較小。

4 結(jié)論

通過對比驗(yàn)算得出在DDA算法可用尺度內(nèi)[12]，粒子輻射的計(jì)算具有較高的精度。因此該文在340 GHz頻段使用該方法計(jì)算等效半徑為50 μm~1 000 μm的六角形冰晶粒子輻射特性，著重分析了不同縱橫比、不同空間取向?qū)碓票ЯＷ虞椛涮匦缘挠绊?。結(jié)果表明：1） 在340 GHz頻段下，通過比較Mie散射理論和DDA計(jì)算得到的球形粒子輻射結(jié)果可以發(fā)現(xiàn)，2種方法得到的輻射曲線幾乎一致，但隨著等效半徑的增加，相對誤差整體增大，且散射誤差最大不超過5%，吸收誤差最大約為10%。因此在340 GHz頻段使用DDA方法計(jì)算常見冰晶粒子（尺度在1 000 μm以內(nèi)）的輻射特性具有較高的精度。2） 空間取向一定時(shí)，不同縱橫比的六角形冰晶粒子輻射特性不同?？v橫比越大，粒子的散射特性曲線更陡峭、波動(dòng)更明顯，而吸收特性曲線整體呈上升的趨勢。粒子縱橫比偏離1越多，相同尺度下粒子的最大截面積就越大，厚板和薄板冰晶粒子的后向散射也就越大。3） 縱橫比一定時(shí)，空間取向差別較大的六角形冰晶粒子輻射特性不同。隨著等效半徑的增大，不同取向引起的散射截面差異被放大。粒子縱橫比偏離1越多，微小的空間取向差異會(huì)使其散射截面產(chǎn)生較大的變化，光學(xué)特性受空間取向影響越敏感；反之，輻射特性就越近似于各向同性，光學(xué)特性受空間取向的影響較小。