《電子技術(shù)應(yīng)用》
您所在的位置:首頁 > 嵌入式技術(shù) > 設(shè)計應(yīng)用 > 基于改進(jìn)OMP的超寬帶穿墻雷達(dá)稀疏成像方法
基于改進(jìn)OMP的超寬帶穿墻雷達(dá)稀疏成像方法
2015年電子技術(shù)應(yīng)用第11期
晉良念1,,2,,錢玉彬2,申文婷2,,劉 琦2,,張 燕2
(1.廣西無線寬帶通信與信號處理重點實驗室,廣西 桂林541004,; 2.桂林電子科技大學(xué) 信息與通信學(xué)院,,廣西 桂林541004)
摘要: 在超寬帶穿墻雷達(dá)成像應(yīng)用中,壓縮感知理論應(yīng)用可以獲得高分辨率稀疏成像,。但這必須建立在目標(biāo)恰好落在預(yù)設(shè)網(wǎng)格點的前提條件下,,一旦目標(biāo)偏離預(yù)設(shè)網(wǎng)格點,目標(biāo)像會發(fā)生偏移,,甚至還會產(chǎn)生虛假像,。本文提出一種基于梯度優(yōu)化的貝葉斯假設(shè)檢驗正交匹配追蹤(GBTOMP)稀疏成像方法。該方法以傳統(tǒng)正交匹配追蹤(OMP)為基礎(chǔ),,從預(yù)設(shè)網(wǎng)格點空間位置出發(fā),,以梯度優(yōu)化的最速上升方法搜索到目標(biāo)真實空間位置,并由此修正模型中的感知矩陣,再利用修正后的感知矩陣恢復(fù)目標(biāo)散射系數(shù),??紤]到正交匹配追蹤會帶來冗余下標(biāo),利用貝葉斯假設(shè)檢驗設(shè)置一個合適的門限去冗余以保證目標(biāo)真實像的準(zhǔn)確恢復(fù),。仿真和實驗結(jié)果表明,,該方法校正了目標(biāo)偏離預(yù)設(shè)網(wǎng)格所帶來的模型誤差,稀疏成像效果明顯改善,。
中圖分類號: TN957.51
文獻(xiàn)標(biāo)識碼: A
DOI:10.16157/j.issn.0258-7998.2015.11.038

中文引用格式: 晉良念,,錢玉彬,申文婷,,等. 基于改進(jìn)OMP的超寬帶穿墻雷達(dá)稀疏成像方法[J].電子技術(shù)應(yīng)用,,2015,41(11):135-139.
英文引用格式: Jin Liangnian,,Qian Yubin,Shen Wenting,,et al. Sparse imaging for ultra-wideband through-the-wall radar based on modified OMP algorithm[J].Application of Electronic Technique,,2015,41(11):135-139.
Sparse imaging for ultra-wideband through-the-wall radar based on modified OMP algorithm
Jin Liangnian1,,2,,Qian Yubin2,Shen Wenting2,,Liu Qi2,,Zhang Yan2
1.Guangxi Key Laboratory of Wideband Communication & Signal Processing,Guilin 541004,,China,; 2.School of Information and Communication Engineering,Guilin University of Electronic Technology,,Guilin 541004,,China
Abstract: Compressive sensing theory is widely applied in ultra-wideband through-the-wall radar imaging, and its performance is very good. However, the target is just falling on the preset grid point. Once the target is off the grid point, the target image will be away from the true position and the ghost will be generated. This paper proposes a Bayesian hypothesis testing orthogonal matching pursuit sparse imaging method based on gradient optimization(GBTOMP). The method is based on the conventional Orthogonal Matching Pursuit(OMP) algorithm. It searches the target true position from preset-grid-point coordinate through the steepest ascent with gradient optimization and corrects the sensing matrix in the model. Then the target scattering coefficient is reconstructed with corrected sensing matrix. Moreover, this paper uses the Bayesian hypothesis testing to set up a suitable threshold to reduce the redundancy, and guarantee the accurate recovery of the target image. The simulation and experiment results show that the proposed method can correct the model errors caused by off-grid target, effectively improve the imaging quality and is better than conventional sparse reconstruction method.
Key words : through-the-wall radar imaging;off-the-grid,;orthogonal matching pursuit,;Bayesian hypothesis testing

 

0 引言

  穿墻雷達(dá)成像(Through the Wall Radar Imageing,TWRI)使用電磁波穿透墻體和內(nèi)部封閉的建筑物進(jìn)行監(jiān)測,,以確定建筑結(jié)構(gòu)布局,,區(qū)分建筑物內(nèi)部的活動情況,檢測,、識別和跟蹤運動目標(biāo),。TWRI技術(shù)在搜救、墻后目標(biāo)檢測和城區(qū)環(huán)境的監(jiān)視等方面?zhèn)涫荜P(guān)注[1-4],。TWRI系統(tǒng)通常利用超寬帶信號和大孔徑天線陣列提供高的距離向和方位向分辨率,,能夠在較短的采樣時間內(nèi)利用較小的存儲空間獲得高質(zhì)量圖像,,與之對應(yīng)的后向投影一類算法需要大的空間和時間(或頻率)采樣數(shù)據(jù)才能獲得高質(zhì)量成像。

  近年來,,壓縮感知理論的興起,,國內(nèi)外學(xué)者將該理論很好地應(yīng)用到TWRI中,提出了一些高分辨稀疏成像方法[5-6],。它們首先將成像區(qū)域劃分為有限個網(wǎng)格,,當(dāng)目標(biāo)恰好落在網(wǎng)格點上,稀疏成像效果很好,。一旦目標(biāo)偏離了網(wǎng)格點,,基矩陣不匹配,從而引起成像模糊,,甚至散焦,。為此,文獻(xiàn)[5]通過建立off-gird稀疏表示模型,,提出了一種基于L21范數(shù)懲罰項的FISTA算法,,該算法雖能較準(zhǔn)確地恢復(fù)目標(biāo)散射系數(shù),但是算法的性能受人工參數(shù)設(shè)置的影響,,在實際應(yīng)用中很難正確地選擇,。文獻(xiàn)[6]針對穿墻雷達(dá)偏離網(wǎng)格目標(biāo)成像中提出了一種貝葉斯推理方法,克服了FISTA算法性能受人工參數(shù)設(shè)置的影響,。但是,,該方法在成像區(qū)域劃分網(wǎng)格數(shù)過大時處理速度會很慢,所以需要尋求新的方法提高處理速度,。

  本文提出一種改進(jìn)的正交匹配追蹤方法(GBTOMP算法)實現(xiàn)偏離網(wǎng)格目標(biāo)超寬帶穿墻雷達(dá)稀疏成像,。該方法以正交匹配追蹤(Orthogonal Matching Pursuit,OMP)算法為基礎(chǔ),,先通過最速上升梯度優(yōu)化方法搜索到目標(biāo)的真實位置,,逐漸減小網(wǎng)格偏移量,修正模型中的感知矩陣,,再由修正的感知矩陣通過最小二乘方法計算目標(biāo)散射系數(shù),,但在此計算過程中,會引入冗余下標(biāo),。本文利用貝葉斯假設(shè)檢驗去冗余,,這樣既可以減小重構(gòu)過程的運算量,還可以增強重構(gòu)的精度和抗噪性能,。仿真和實驗結(jié)果驗證了本文方法的有效性和可行性,。

1 稀疏表示模型


001.jpg

  一般的穿墻傳播模型可以等效為如圖1所示的兩層傳播模型[7]。假設(shè)均勻介質(zhì)墻體的厚度為d,相對介電常數(shù)為?著,。將成像區(qū)域在方位向和距離向上離散成Nx Ny個網(wǎng)格點,,對應(yīng)的像素值用?滓k,l,,k=1,,…,Nx,,l=1,,…,Ny表示,,當(dāng)網(wǎng)格點有目標(biāo)時k,,l不為0,反之k,,l等于0,。為了方便表示,將k,,l按照順序重新排列成Nx Ny維的目標(biāo)散射系數(shù)列向量,,令所有元素的網(wǎng)格位置用矩陣表示,其中的元素用p=[xp,,yp]T表示,對應(yīng)的像素值用p表示,,則第n個天線回波信號的離散數(shù)據(jù)模型為:

  1.png

  式中,,yn=[yn(t1),…,,yn(tM)]T,,An為M×Nx Ny的矩陣,第m行的元素為[s(tn,,m-n,,1),…,,s(tn,,m)],這里n,,p為第n個天線到第p個網(wǎng)格點之間的傳播時延,。

  根據(jù)圖1所示的傳播模型,n,,p可以近似表示為:

  2.png

  式中,,xn、yn分別表示天線的坐標(biāo)位置,n,,p為傳播的入射角,,當(dāng)成像區(qū)域距離天線較遠(yuǎn)時可以近似為:

  3.png

  將N個天線到成像區(qū)域中所有像素點的傳播時延按照式(2)和式(3)來計算,并把所有天線的接收數(shù)據(jù)進(jìn)行堆疊,,由此所構(gòu)造的離散數(shù)據(jù)模型為:

  4.png

  然而,,實際情況中目標(biāo)很可能偏離預(yù)設(shè)網(wǎng)格,如圖1所示,,這種偏移會影響A′的精確計算,。定義成像空間所有網(wǎng)格的偏移量矩陣為,A′關(guān)于的函數(shù)為A′,。另外,,對于穿墻成像場景來說滿足稀疏特性,因此可以利用較少的測量值恢復(fù),。令表示從y中隨機選擇選擇Q1個天線單元和每個天線單元隨機選擇Q2個數(shù)據(jù)得到的Q1Q2個數(shù)據(jù)矢量,,并定義這種隨機選擇形式構(gòu)成的測量矩陣用表示[8],由此得到大小為Q1Q2×1的測量數(shù)據(jù)(或稀疏成像)模型(為了書寫方便將A為:

  5.png

  根據(jù)式(5)的稀疏成像模型可以看出,,感知矩陣A有關(guān),,因此求解式(4)反演需要事先準(zhǔn)確估計?專參數(shù)。當(dāng)且僅當(dāng)A滿足約束等距性質(zhì)時,,可以通過求解如下優(yōu)化問題將高維信號從低維的測量值中準(zhǔn)確重構(gòu)出來,。

  6.jpg

  由于最小0范數(shù)問題是一個非確定性多項式難問題,目前的重構(gòu)算法都是基于0范數(shù)的變換或逼近處理,,其中貪婪算法以其運算速度快,,實現(xiàn)簡單的優(yōu)勢被用于穿墻雷達(dá)成像中。

2 GBTOMP的稀疏成像方法

  2.1 GBTOMP算法

  本文提出的GBTOMP算法是在正交匹配追蹤算法[9]的迭代過程中采用梯度優(yōu)化最速上升算法校正網(wǎng)格偏移量,,然后采用貝葉斯假設(shè)檢驗[10]剔除冗余下標(biāo),。該算法整體上可以分為兩大步:第一,根據(jù)殘差與感知矩陣的內(nèi)積最大值找到最匹配原子位置(預(yù)設(shè)網(wǎng)格位置),。由于目標(biāo)可能會偏離預(yù)設(shè)的網(wǎng)格位置,,因此在該位置需要進(jìn)行梯度優(yōu)化的最速上升,以搜索到目標(biāo)的真實位置,,從而以目標(biāo)真實位置建立新的感知矩陣,,并由該矩陣?yán)米钚《朔ㄇ蟮蒙⑸湎禂?shù)向量的近似值,同時通過交替迭代更新下標(biāo)集求解來得到目標(biāo)散射系數(shù),。第二,,OMP算法在有噪聲和沒有噪聲的情況下都會存在冗余下標(biāo)問題。為解決此問題,,GBTOMP算法最后采用貝葉斯假設(shè)檢驗理論給出判決準(zhǔn)則對第一步輸出的預(yù)選下標(biāo)集中原子進(jìn)行篩選,,得到等于或者大于信號真實下標(biāo)集估計后,,再用最小二乘算法進(jìn)行重構(gòu),以改善算法的性能,,增強其抗噪能力,。GBTOMP算法流程如下。

  初始化:構(gòu)造測量向量,,網(wǎng)格偏移矩陣,,感知矩陣A。

  (1)采用OMP算法找到匹配原子對應(yīng)的空間網(wǎng)格位置,,利用最速上升梯度優(yōu)化方法對偏離網(wǎng)格進(jìn)行估計,,再利用最小二乘法對目標(biāo)散射系數(shù)進(jìn)行粗估計,輸出預(yù)選下標(biāo)集P和粗重構(gòu)目標(biāo)散射系數(shù)X,。

  (2)從預(yù)選下標(biāo)集P出發(fā),,通過貝葉斯假設(shè)檢驗準(zhǔn)則設(shè)置一個適當(dāng)門限Thj對P中的元素進(jìn)行篩選,將大于門限的元素都保留在最終的下標(biāo)集F中,。

  輸出:根據(jù)最終下標(biāo)集F,,利用最小二乘法得到目標(biāo)散射系數(shù)最終的估計值。

  2.2 網(wǎng)格偏移量的估計

  在GBTOMP算法的迭代過程中需要估計網(wǎng)格偏移量,,其代價函數(shù)的構(gòu)造是關(guān)鍵,。OMP算法是以測量數(shù)據(jù)和感知矩陣A構(gòu)造代價函數(shù)的。為了能夠使得估計準(zhǔn)確,,并使處理過程更加穩(wěn)定,,本文使用OMP算法相關(guān)函數(shù)的平方開方項來表示[11]:

  7.png

  首先初始化迭代更新次數(shù)k=1,殘差向量r,。計算感知矩陣A所有列和殘差向量r的內(nèi)積值,,找到內(nèi)積值最大時所對應(yīng)的列下標(biāo)pk,接著進(jìn)入搜索目標(biāo)真實位置過程,,滿足搜索條件|Ji-Ji-1|<a結(jié)束,這里的i代表搜索過程中的迭代次數(shù),,a代表很小的一個正數(shù),。在此過程中,從?仔中取出第pk列賦給?仔0,,以該點坐標(biāo)(xp,,yp)計算出與隨機選擇Q1個天線組成新集合的第q個索引號天線位置之間的雙程傳播時延為:

  89.png

  為求解代價函數(shù)J關(guān)于所選擇坐標(biāo)(xp,yp)的梯度,,需要事先求解?子q,,p關(guān)于(xp,yp)的偏導(dǎo),,分兩種情況進(jìn)行:

  (1)當(dāng)xq=xp時,,即q,,p=0,q,,p簡化為:

  1011.png

  (2)當(dāng)xq≠xp時,,即q,p≠0,,q,,p關(guān)于xp和yp的偏導(dǎo)數(shù)為:

  1213.jpg

  假設(shè)輻射脈沖波形為高斯脈沖s(t),則第q個天線和所選目標(biāo)坐標(biāo)(xp,,yp)對應(yīng)的感知矩陣Aq,、B是大小為Q2×1的矩陣,它們的第q2行分別為:

  1415.jpg

  式中是當(dāng)前(xp,,yp)對應(yīng)的感知矩陣,,僅僅是A的一個子集,這樣處理是為了提高運行速度,,減小存儲空間,。B1x和B1y表示對矩陣分別關(guān)于xp和yp的偏導(dǎo)數(shù),它們的具體表達(dá)式與s(t)有關(guān),。以梯度搜索方法更新所選空間的坐標(biāo)為:

  16.png

  式中是搜索步長,,本文中以網(wǎng)格尺寸為參考的很小的正數(shù)。由此得到更新后的重新計算代價函數(shù)Ji,,搜索迭代次數(shù)i=i+1,,一旦代價函數(shù)滿足條件|Ji-Ji-1|<a就終止搜索過程。終止后,,將所選的網(wǎng)格空間位置?仔1所對應(yīng)的感知矩陣保存在A2矩陣中,,同時將所選的感知矩陣列下標(biāo)pk添加到預(yù)選下標(biāo)集P中,由A2通過最小二乘法計算出第k次粗估計的散射系數(shù)xpre,,并更新殘差向量,。

  17.png

  繼續(xù)搜索下標(biāo)并校正網(wǎng)格偏移量,如此循環(huán)迭代,,k=k+1,,直到殘差向量r滿足條件,粗估計結(jié)束,。將k次循環(huán)求得的散射系數(shù)xpre對應(yīng)其下標(biāo)集P賦給Nx Ny維列向量X,。

  X{P}=xpre(18)

  2.3 貝葉斯假設(shè)檢驗去冗余

  為剔除OMP算法輸出的預(yù)選集P中的冗余下標(biāo),根據(jù)貝葉斯假設(shè)檢驗?zāi)P蚚10]構(gòu)造似然函數(shù),,其似然比檢驗公式為:

  19.jpg


  根據(jù)判決準(zhǔn)則可知,,當(dāng)(zj-mj)2>Thj時,下標(biāo)j對應(yīng)的散射系數(shù)?滓j的值為非零,,更新此時的j為最終下標(biāo),,反之,,下標(biāo)j對應(yīng)的散射系數(shù)?滓j為零,則從預(yù)選下標(biāo)集中剔除此時的j值,。對預(yù)選下標(biāo)集中元素進(jìn)行篩選,,最終剩下的真實下標(biāo)集為F,接下來利用最小二乘法進(jìn)行目標(biāo)散射系數(shù)重構(gòu)完成整個GBTOMP算法,。

3 仿真與實驗結(jié)果分析

  3.1 仿真結(jié)果分析

  利用GPRMAX產(chǎn)生仿真數(shù)據(jù)[12],,發(fā)射信號采用高斯脈沖,脈沖中心頻率為1 GHz,。墻體厚度為0.20 m,,相對介電常數(shù)為6,電導(dǎo)率為0.03 S/m,。

  場景1:天線陣列由16個收發(fā)共置單元組成,,隨機選擇其中的10個天線單元,陣元間距為0.1 m,。成像區(qū)域是[-1 m,,1 m]×[0.3 m,1.3 m],,距離向和方位向的網(wǎng)格尺寸均為0.05 m,。假設(shè)成像區(qū)域中有三個點目標(biāo),其坐標(biāo)分別為(-0.414 m,,0.514 m),,(0.025 m,0.725 m),,(0.325 m,,0.965 m)。圖2為OMP方法和本文GBTOMP方法成像結(jié)果,。從圖2(a)可以看出,,OMP方法成像出現(xiàn)嚴(yán)重散焦,同時帶來虛假像,,從圖2(b)可以看到,,估計的目標(biāo)位置和預(yù)設(shè)目標(biāo)位置完全吻合,同時目標(biāo)像聚焦程度高,,很明顯本文方法性能優(yōu)于傳統(tǒng)OMP方法,。

002.jpg

  場景2:天線陣列由45個收發(fā)共置單元組成,,隨機選擇其中的10個天線單元,,陣元間距為0.04 m,與墻體之間的距離為0.2 m,。圖3給出了2維GPRMAX模型,,圖中左邊的長方形目標(biāo)中心位于(0.7 m,,3.25 m),長為0.4 m,,寬為0.3 m,,右邊的圓形目標(biāo)體的中心位置位于(1.25 m,3.05 m),,半徑為0.2 m,。GPRMAX網(wǎng)格單元尺寸為0.005 m,采樣時窗為30 ns,。假設(shè)回波信號由有目標(biāo)時減去無目標(biāo)時的信號進(jìn)行模擬,。

  圖4給出了兩個擴展目標(biāo)的成像結(jié)果,成像區(qū)域的方位向和距離向的網(wǎng)格尺寸設(shè)為0.03 m,。圖中的矩形和圓形虛線框代表兩目標(biāo)的真實位置,。從圖4可以看出,傳統(tǒng)OMP方法的成像效果較差,,雖然能夠大致識別目標(biāo)位置,,但是引入了很多的虛假目標(biāo)像,而本文方法GBTOMP成像效果明顯改善,,目標(biāo)像清晰,。

  3.2 實驗結(jié)果分析


005.jpg

  使用美國GSSI公司的探地雷達(dá)SIR-20搭建穿墻實驗場景,如圖5所示,。實驗墻體為實心磚墻,,墻體厚度0.20 m,相對介電常數(shù)6.4,。選用1 GHz的喇叭天線,,架高1.2 m,貼著墻體,。水平移動喇叭天線共掃描21個測點,,數(shù)據(jù)處理時隨機選擇10個測點,測點間距為0.05 m,,合成孔徑長度為1 m,,要求在每個天線孔徑測試點處測量2次,包括有人和無人的場景,。人體目標(biāo)高度為1.80 m,,體型寬度約0.40 m,站在墻后1 m處,。SIR-20系統(tǒng)數(shù)據(jù)采集的參數(shù)設(shè)置:每道采樣1 024點,,每秒60道,時間窗15 ns,。

006.jpg

  首先將SIR-20系統(tǒng)采集的數(shù)據(jù)波形做平均處理,,然后經(jīng)過去噪和雜波相消,、自動增益控制等信號處理得到較好的目標(biāo)回波,所有成像都是以天線陣列中心為坐標(biāo)原點,。圖6給出實驗數(shù)據(jù)成像結(jié)果比較,,圖中的矩形虛線框代表目標(biāo)的真實位置。從圖中可以看出,,傳統(tǒng)OMP方法目標(biāo)成像散焦現(xiàn)象較嚴(yán)重,,目標(biāo)區(qū)域不能確定,且周圍還有較多的虛假像,。而本文GBTOMP方法成像則不同,,雜波相對較少,且目標(biāo)位置清晰可見,,由此看來通過校正目標(biāo)偏離網(wǎng)格能夠很好解決傳統(tǒng)稀疏成像所出現(xiàn)的散焦現(xiàn)象,。

4 結(jié)論

  本文提出的GBTOMP稀疏成像方法,交替迭代得到目標(biāo)真實空間位置和目標(biāo)散射系數(shù),,目標(biāo)位置的校正保證了基矩陣的匹配,,實現(xiàn)目標(biāo)的高質(zhì)量稀疏成像。仿真和實測數(shù)據(jù)處理結(jié)果表明:與OMP成像相比,,本文GBTOMP方法成像效果在目標(biāo)位置準(zhǔn)確度還是目標(biāo)聚焦程度有明顯的提高,。但是,由于以點稀疏模型對擴展目標(biāo)進(jìn)行成像,,存在圖像不夠平滑,,也不連續(xù)等問題,下一步工作將對擴展目標(biāo)的稀疏成像方法展開研究,。

參考文獻(xiàn)

  [1] TIVIVE F H C,,BOUZERDOUM A,AMIN M G.A subspaceprojection approach for wall clutter mitigation in through--wall radar imaging[J].IEEE Transactions on Geoscienceand Remote Sensing,,2015,,53(4):2108-2122.

  [2] LI G,BURKHOLDER R J.Hybrid matching pursuit for distributed through-wall radar imaging[J].IEEE Transactionson Antennas and Propagation,,2015,,63(4):1701-1711.

  [3] 蔣留兵,韋洪浪,,許騰飛,,等.EEMD生命探測雷達(dá)人體數(shù)量識別技術(shù)[J].電子技術(shù)應(yīng)用,2014,,40(5):122-125.

  [4] LIU J,,CUI G,JIA Y,et al.Sidewall detection using mul-tipath in through-wall radar moving target tracking[J].IEEEGeoscience and Remote Sensing Letters,,2015,12(6): 1372-1376.

  [5] XIA S,,LIU F.Off-grid compressive sensing through-the-wall radar imaging[C].SPIE Defense+Security.International Society for Optics and Photonics,,2014:90771F-90771F-8.

  [6] 晉良念,錢玉彬,,劉慶華,,等.超寬帶穿墻雷達(dá)偏離網(wǎng)格目標(biāo)稀疏成像方法[J].儀器儀表學(xué)報,2015,,36(4).

  [7] JIN T,,CHEN B,ZHOU Z.Imaging-domain estimation of wall parameters for autofocusing of through-the-wall SAR imagery[J].IEEE Transactions on Geoscience and Remote Sensing,,2013,,51(3):1836-1843.

  [8] LAGUNAS E,AMIN M G,,AHMAD F,,et al.Joint wall mitigation and compressive sensing for indoor image reconstruction[J].IEEE Transactions on Geoscience and Remote Sensing,2013,,51(2):891-906.

  [9] 李少東,,裴文炯,楊軍,,等.貝葉斯模型下的OMP重構(gòu)算法及應(yīng)用[J].系統(tǒng)工程與電子技術(shù),,2015,37(2).

  [10] 甘偉,,許錄平,,蘇哲,等.基于貝葉斯假設(shè)檢驗的壓縮感知重構(gòu)[J].電子與信息學(xué)報,,2011,,33(11):2640-2646.

  [11] GURBUZ A C,TEKE O,,ARIKAN O.Sparse ground-penetrating radar imaging method for off-the-grid target problem[J].Journal of Electronic Imaging,2013,,22(2):021007-021007.

  [12] XIE J L,XU J L.Ground penetrating radar-based experimental simulation and signal interpretation on roadway roof separation detection[J].Arabian Journal of Geosciences,,2014:1-8.


此內(nèi)容為AET網(wǎng)站原創(chuàng),,未經(jīng)授權(quán)禁止轉(zhuǎn)載。