王振武， 卜異亞， 馬鍵強

(中國礦業(yè)大學(xué)(北京) 機電與信息工程學(xué)院，北京 100083)

?

地質(zhì)雷達數(shù)據(jù)中克里金插值采樣數(shù)據(jù)選擇算法

王振武， 卜異亞， 馬鍵強

(中國礦業(yè)大學(xué)(北京) 機電與信息工程學(xué)院，北京 100083)

采樣數(shù)據(jù)選擇的合理性是克里金算法插值結(jié)果是否準(zhǔn)確的前提。針對如何從地質(zhì)雷達數(shù)據(jù)中選擇恰當(dāng)?shù)目死锝鸩逯挡蓸訑?shù)據(jù)的問題，本文提出了兩種克里金插值采樣數(shù)據(jù)選擇算法，即動態(tài)三維球體覆蓋(D3DGC)選點算法和雙次反距離選點(DRD)算法。兩種算法在基于粒子群優(yōu)化的普通克里金算法下驗證了采樣數(shù)據(jù)選擇的合理性，實驗結(jié)果分析表明，D3DGC算法在最佳選點個數(shù)下誤差率比傳統(tǒng)選點方法降低了6.7%，在其他選點個數(shù)下誤差率(<10%)也明顯低于傳統(tǒng)選點方法，DRD算法在參數(shù)e和c為4以及最佳選點個數(shù)的條件下誤差率也達到了6.2%，在其他情況下誤差率也低于傳統(tǒng)選點方法。

克里金算法；插值；地質(zhì)雷達；采樣數(shù)據(jù)選擇算法；動態(tài)三維球體覆蓋選點算法；雙次反距離選點算法；粒子群優(yōu)化；誤差率

地質(zhì)雷達(ground penetrating radar，GPR)技術(shù)是一種基于電磁波反射原理，用于淺層地質(zhì)構(gòu)造探測和工程治理檢測的地球物理方法[1]。為了深入勘查或檢測探測區(qū)域的地質(zhì)構(gòu)造或地下結(jié)構(gòu)，需要對測線數(shù)據(jù)做插值處理來獲取探測區(qū)域的整體數(shù)據(jù)。在常見的空間插值方法中以克里金(Kriging)插值方法為代表的地質(zhì)統(tǒng)計方法提供了揭示空間變量結(jié)構(gòu)相關(guān)性與隨機性的理論，是地質(zhì)領(lǐng)域應(yīng)用最廣泛的插值方法之一[2]。在對克里金算法的研究中，有些學(xué)者根據(jù)不同的應(yīng)用情況，將不同類型的克里金算法組合使用。牛文杰等將同位置協(xié)同克里金(Collocated CoKriging)和貝葉斯克里金方法相結(jié)合，推導(dǎo)出一種新的貝葉斯同位置協(xié)同克里金估值方法，綜合了上述兩種方法的優(yōu)點[3]。李曉軍等針對不連續(xù)地層，提出了一種采用指示克里金和普通克里金相結(jié)合的方法，即先用指示克里金法估計地層分布范圍，再用普通克里金方法估計分布范圍內(nèi)的地層厚度，更多的學(xué)者將克里金方法和其他方法結(jié)合使用[4]。王輝贊等針對變異函數(shù)難以刻畫實際數(shù)據(jù)分布的不足，采用最小二乘支持向量機從實際資料場擬合重構(gòu)變異函數(shù)，表現(xiàn)出較好的針對性和客觀性[5]。賈雨等將粒子群算法和克里金算法相結(jié)合，在粒子群優(yōu)化過程中，通過高斯變異、樣本點權(quán)重系數(shù)設(shè)定、搜索范圍約束等方式提高了插值精度[6]。顏華等針對三維溫度場重建問題，首先采用最小二乘法獲得一個由少量像素描述的溫度場，然后運用克里金內(nèi)插和外推運算獲得整個被測溫度場的細致描述，這種方法可實現(xiàn)復(fù)雜的三維溫度場的高精度重建[7]。徐馳等將高光譜反演與協(xié)同克里金方法相結(jié)合，將土壤表層含水率反演結(jié)果作為協(xié)同克里金插值的協(xié)變量來輔助進行協(xié)同克里金插值，提高了預(yù)測效率[8]。ZHANG等先采用克里金算法進行大型球磨機的功能適配，然后利用遺傳算法優(yōu)化了齒輪傳送的可靠性問題[9]。FANG等將序列克里金算法和人工蜂群算法結(jié)合來優(yōu)化點焊接頭以改善其疲勞壽命，該方法用克里金模型來近似表達點焊接頭設(shè)計模型的目標(biāo)函數(shù)，然后采用人工蜂群算法進行迭代優(yōu)化[10]。

在上述的研究中，對克里金插值方法的研究多集中在對不同類型克里金算法的組合[3-4]以及與其他算法理論的結(jié)合[5-11]中，鮮有學(xué)者研究克里金插值方法中采樣數(shù)據(jù)的選擇問題。一般情況下，按照空間相關(guān)性原理，地質(zhì)數(shù)據(jù)應(yīng)該選取距離插值點較近的數(shù)據(jù)點作為采樣數(shù)據(jù)。文獻[12]針對某儲層數(shù)據(jù)，研究了空間相關(guān)分析的主要影響因素與克里金估計結(jié)構(gòu)的關(guān)系，指出影響插值結(jié)構(gòu)的因素和參數(shù)包括采樣點的位置分布特征、變差函數(shù)模型的塊金效應(yīng)、變差函數(shù)模型的變程、變差函數(shù)模型的類型、空間結(jié)構(gòu)的各向異性以及空間搜索范圍等。文獻[11]進一步指出，距離估計點較近的采樣數(shù)據(jù)對估計結(jié)構(gòu)的影響較大，且大多數(shù)情況下應(yīng)該選用球狀模型或指數(shù)模型。但文獻[11]也只是給出了一般性的選點原則，對具體數(shù)據(jù)而言，應(yīng)該如何選擇采樣數(shù)據(jù)并沒有進行深入分析。文獻[2]針對傳統(tǒng)的克里金方法采用歐式距離測度描述空間結(jié)構(gòu)相關(guān)性的缺點，基于最小方差準(zhǔn)則，提出了綜合曼氏、歐式和切氏三種距離測度的多測度加權(quán)克里金方法，并通過實驗表明該方法提高了插值精度和可靠性，但是該方法也僅僅是對距離測度公式做了擴展，沒有深入討論采樣數(shù)據(jù)如何選擇的問題，而選擇合適的采樣數(shù)據(jù)是提高克里金插值準(zhǔn)確性的前提。論文針對這一問題提出了動態(tài)三維球體覆蓋(dynamic three-dimensional globe cover，D3DGC)選點算法和雙次反距離(double reverse distance，DRD)選點算法，并基于粒子群優(yōu)化的普通克里金算法驗證了兩種算法的有效性。

1 克里金插值采樣數(shù)據(jù)選擇算法描述

1.1 基本概念

本文在已獲取地質(zhì)雷達測線數(shù)據(jù)的基礎(chǔ)上研究克里金插值的采樣數(shù)據(jù)選擇問題。為了對所提出的算法進行描述，首先給出相關(guān)概念。

定義1 地質(zhì)數(shù)據(jù)點δ=〈ζ，γ〉：δ是一個二元組，其中ζ=〈x，y，z〉表示三維坐標(biāo)值，γ表示地質(zhì)點的屬性值，如物質(zhì)對電磁波的反射值等。

(1)

式中：δij表示測線數(shù)據(jù)上的第i行j列的地質(zhì)數(shù)據(jù)點，且1≤i≤m，1≤j≤n。為表述方便，設(shè)一共p條測線數(shù)據(jù)，第l條測線表示為l，將測線上第h行的地質(zhì)數(shù)據(jù)點集合定義為δh，h∈[1，m]，l∈[1，p]。

定義3 測線數(shù)據(jù)集Γ：Γ表示對同一探測區(qū)域的同一次探測活動中獲取的測線數(shù)據(jù)總和，即

Γ={1，2，…，lp}

(2)

定義4 測線數(shù)據(jù)凸包Δ：Δ表示當(dāng)前測線數(shù)據(jù)集所描述的最大三維探測區(qū)域。

Δ={δ|δ.ζ.x∈[Xmin，Xmax]，δ.ζ.y∈[Ymin，Ymax]，
δ.ζ.z∈[Zmin，Zmax]}

Xmin=min{k.δij.ζ.x}，Xmax=max{k.δij.ζ.x}
Ymin=min{k.δij.ζ.y}，Ymax=max{k.δij.ζ.y}
Zmin=min{k.δij.ζ.z}，Zmax=max{k.δij.ζ.z}

(3)

式中：k∈[1，p]，i∈[1，mk]，j∈[1，nk]，Xmin和Xmax分別表示x坐標(biāo)的最大和最小值，min{ }和max{ }表示取集合的最小值和最大值操作。

定義5 待插值數(shù)據(jù)點集ψ：ψ=Δ-Γ，克里金插值的目的就是從Γ中選擇合適的采樣點對ψ進行插值。

定義6 凸包水平分割層η：η表示Δ按照測線平行方向的層次劃分，即

(4)

Ψ={φ1，φ2，…,φu}

(5)

式中：φr表示第r層的待插值數(shù)據(jù)點的集合， 1≤r≤u且存在關(guān)系測線l中的δi滿足δi?l.δi，l∈[1，p]，i∈[1，m]。

定義7 待插值點的采樣數(shù)據(jù)點集φr[ε]→π[μ]：假設(shè)選擇μ個采樣點進行克里金插值，則φr中第ε個插值點的μ個采樣點的集合表示為φr[ε]→π[μ]。

定義8 距離計算公式τ：論文采用歐幾里得距離作為距離測度公式。設(shè)ψ中的某待插值點為θ(a，b，c)，測線上某點l中δ的ζ值為(x，y，z)，則距離τ定義為

(6)

1.2 采樣數(shù)據(jù)選擇算法原理

地質(zhì)數(shù)據(jù)具有空間相關(guān)性和隨機性的特點，針對不同數(shù)據(jù)的特點，本文提出了D3DGC選點算法和DRD選點算法。

1.2.1 D3DGC 算法

根據(jù)地質(zhì)數(shù)據(jù)的相關(guān)性原理，一般而言選取距離插值點較近的數(shù)據(jù)點作為采樣點，針對地質(zhì)雷達數(shù)據(jù)特點，論文提出了D3DGC選點算法，算法原理如圖1所示。

圖1 D3DGC算法原理圖Fig.1 Theory graph of D3DGC

定義9 球與直線交點方程：設(shè)待插值點為θ(a，b，c)， 球的半徑為R，l上的δh與球的交點φ(x，y，z)計算公式為

(7)

設(shè)測線垂方向采樣點步長為d，每個待插值點選取μ個采樣點進行插值?；谌缟隙x，D3DGC算法的流程如圖2所示。

圖2 D3DGC算法流程圖Fig.2 The flow graph of D3DGC

對于每一個待插值點θ(a，b，c)，球的半徑以d為單位遞增，計算l的δh與球體交截后包含在球體內(nèi)的地質(zhì)數(shù)據(jù)點個數(shù)kl，這里直到條件滿足時為止。通過逐層計算ηi上的φi，直到所有待插值點的采樣點選擇完畢為止。

如果θ(a，b，c)的插值采樣點只從ηi中的δi選擇，其中δi屬于l，1≤≤p則D3DGC算法就簡化為二維模式，如圖3所示，稱之為動態(tài)圓覆蓋(dynamic circle cover，DCC)選點算法。

準(zhǔn)時制（JIT)，可概括為“任何時間、任何地點、任何事情”，都可以及時解決和處理，以達到零時間、零距離的高效作業(yè)模式。審批平臺可以引入此理論在審批業(yè)務(wù)流程的設(shè)計和模型的架構(gòu)當(dāng)中，以實現(xiàn)審批的高效、便捷。具體JIT模式的優(yōu)勢見圖1。

1.2.2 DRD算法

圖3 DCC算法原理圖Fig.3 Theory graph of DCC algorithm

定義10 點到直線的距離?。設(shè)φi中當(dāng)前的待插值點為θ(a，b，c)，l上的δh可由式(7)的方程2描述，若將方程2簡記為則θ(a，b，c)到l上的δh的距離可定義為

(8)

設(shè)待插值點選取μ個樣本點進行插值，.δh的點數(shù)為n(為插值準(zhǔn)確性和效率考慮，通常μ?n)，按照式(8)計算距離后，分配到第l條測線上的采樣點數(shù)由式(9)計算：

(9)

(10)

圖4 DRD算法原理圖Fig.4 Theory graph of DRD

圖5 DRD算法流程圖Fig.5 The flow graph of DRD

2 算法實驗分析

由于主觀設(shè)置的變差函數(shù)模型參數(shù)對克里金算法的插值結(jié)果影響較大，為了客觀準(zhǔn)確地驗證兩種采樣點選擇算法，論文采用了粒子群優(yōu)化(particle swarm optimization，PSO)算法對普通克里金算法的變差函數(shù)模型參數(shù)進行優(yōu)化。

PSO算法將優(yōu)化問題的候選解假設(shè)為一個無體積和質(zhì)量的粒子，在d維空間中粒子通過自身和其他粒子的經(jīng)驗來調(diào)整自己的位置。設(shè)有n個粒子組成的粒子群，其中第i個粒子在d維空間中的速度為vi=[vi1vi2…vid]T，位置為xi=[xi1xi2…xid]T，粒子i的當(dāng)前個體最優(yōu)位置為pid，整個粒子群的當(dāng)前全局最優(yōu)位置為gd，粒子i通過式(11)的迭代來更新自己的速度和位置。

(11)

定義11 適應(yīng)度函數(shù)F:本文給定的PSO算法的適應(yīng)度函數(shù)表示為

(12)

式中：F(δo) 為第o個粒子(地質(zhì)數(shù)據(jù)點)的適應(yīng)度函數(shù)值，V*(δo)表示粒子的克里金估計量，V(δo)表示例子的實際值。論文采用普通Kriging算法進行實驗結(jié)果的驗證，變差函數(shù)采用球狀模型，拱高f、塊金值g、變程s，參數(shù)采用PSO算法進行優(yōu)化，粒子的速度和位置用向量〈f,g,s〉表示，粒子群的最大迭代次數(shù)α=1 000，適應(yīng)度函數(shù)閾值β=10-6。PSO優(yōu)化的克里金(PSO-Kriging)算法的流程如圖6所示。

圖6 PSO-Kriging算法流程圖Fig.6 The flow graph of PSO-Kriging

本文基于地質(zhì)雷達探測數(shù)據(jù)來驗證所提出的兩種采樣點選擇算法的有效性。地質(zhì)雷達儀為中國礦業(yè)大學(xué)(北京)自主研發(fā)的MTGR-4F車載地質(zhì)雷達，雷達參數(shù)為：探測深度5m，時間窗口100ns，采樣點數(shù)512，天線主頻200MHz。實驗數(shù)據(jù)探測路段為北京市朝陽區(qū)東三環(huán)中路，探測目的有3個：地底的分層結(jié)構(gòu)、有無地下管線以及有無空洞。實驗數(shù)據(jù)取值范圍為[-23 679，15 421]；測線條數(shù)為3，采樣道數(shù)分別為18 632、18 088和18 838，實驗所用數(shù)據(jù)的采集工作照如圖7所示。

圖7 實驗數(shù)據(jù)的獲取Fig.7 The acquirement of experiment data

實驗硬件和軟件環(huán)境信息如下：CPU為Pentium(R)Dual-CoreT4300 2.10GHz，內(nèi)存為RAM2GB，操作系統(tǒng)為Windows7 32位，實驗平臺為VS.NET2012，開發(fā)語言為C++。

傳統(tǒng)的采樣點選擇方式是對同一ηi，1≤i≤u上所有l(wèi)的δi數(shù)據(jù)按照式(6)進行排序，選擇距離最近的某些地質(zhì)數(shù)據(jù)點作為克里金插值的采樣點。本文將所提出的兩種算法和這種傳統(tǒng)的選點方法進行了比較，為準(zhǔn)確起見，誤差率統(tǒng)計的為N0=2 000個點的平均誤差率，誤差率ρ由式(13)計算：

(13)

實驗探測數(shù)據(jù)的最優(yōu)采樣數(shù)據(jù)點數(shù)與采樣道數(shù)、采樣點數(shù)等具體的地質(zhì)雷達儀參數(shù)有關(guān)?？紤]到插值算法的運算效率和插值結(jié)果精度等因素，本實驗選取的數(shù)據(jù)點個數(shù)μ∈[160，300]，依次按10個點的數(shù)量遞增。如圖8所示，在不同的選點個數(shù)下，D3DGC算法的誤差率明顯低于傳統(tǒng)選點算法。當(dāng)選點個數(shù)為210時，兩種算法的誤差率都為最小，D3DGC算法的誤差率僅為1.3%，而傳統(tǒng)選點算法的誤差率卻為8%，而且當(dāng)選點個數(shù)在180～270時，D3DGC算法的誤差率都小于10%，而傳統(tǒng)選點算法只有在200～210誤差才小于10%，且都大于8%，在其他范圍內(nèi)的最大誤差達到了31.1%。與傳統(tǒng)選點算法相比，D3DGC算法效果好的原因在于，該算法在三維空間中考慮相關(guān)性，按照球狀半徑遞增的原理逐步尋找采樣數(shù)據(jù)點，比傳統(tǒng)選點算法更為準(zhǔn)確合理。

圖8 D3DGC算法實驗結(jié)果Fig.8 The experiment results of D3DGC

表1 DRD算法實驗結(jié)果

3 結(jié)論

1)D3DGC算法是在三維空間中按照球狀半徑遞增原理選取的采樣點，因此其克里金插值誤差率明顯地比傳統(tǒng)選點方法好；2)針對數(shù)據(jù)量大且相鄰數(shù)據(jù)值較為接近的地質(zhì)雷達數(shù)據(jù)，DRD算法綜合考慮了地質(zhì)數(shù)據(jù)的空間相關(guān)性和隨機性，在參數(shù)e和c取值為4的條件下也取得了不錯的效果。

D3DGC和DRD算法均已在道路塌陷的檢測問題中得到了應(yīng)用，使得克里金插值的精度得到了提高。對于DRD算法，針對不同的地質(zhì)雷達數(shù)據(jù)，如何實現(xiàn)參數(shù)e和c的自適應(yīng)取值是將來需要進一步研究的問題。

[1]李晉平，邵丕彥，谷牧.地質(zhì)雷達在鐵路隧道工程質(zhì)量檢測中的應(yīng)用[J].中國鐵道科學(xué)， 2006， 27(2)： 56-59.

LI Jinping， SHAO Piyan， GU Mu. Application and analysis of GPR in railway tunnel engineering quality inspection[J]. China railway science， 2006，27(2)： 56-59.

[2]劉志平，何秀鳳，張淑輝.多測度加權(quán)克里金法在高邊坡變形穩(wěn)定性分析中的應(yīng)用[J].水利學(xué)報， 2009， 40(6)： 709-715.

LIU Zhiping， HE Xiufeng， ZHANG Shuhui. Multi-distance measures weighted Kriging method for deformation stability analysis of steep slopes[J]. Shuili xuebao， 2009，40(6)： 709-715.

[3]牛文杰，朱大培，陳其明.貝葉斯同位置協(xié)同克里金(Bayesian Collocated CoKriging)估值法的研究[J]計算機輔助設(shè)計與圖形學(xué)學(xué)報， 2002， 14(4)： 343-347.

NIU Wenjie， ZHU Dapei， CHEN Qiming. Bayesian collocated CoKriging[J]. Journal of computer-aided design & computer graphics， 2002， 14(4)： 343-347.

[4]李曉軍，張振遠.基于指示和普通克里金的不連續(xù)地層厚度估計方法[J].巖土力學(xué)， 2014， 35(10)： 2881-2887.

LI Xiaojun， ZHANG Zhenyuan. A combined indicator-ordinary Kriging method for estimating thicknesses of discontinuous geological strata[J]. Rock and soil mechanics， 2014， 35(10)： 2881-2887.

[5]王輝贊，張韌，劉巍，等.支持向量機優(yōu)化的克里金插值算法及其海洋資料對比試驗[J].大氣科學(xué)學(xué)報， 2011， 34(5)： 567-573.

WANG Huizan， ZHANG Ren， LIU Wei， et al. Kriging interpolation method optimized by support vector machine and its application in oceanic data[J]. Transactions of atmospheric sciences， 2011， 34(5)： 567-573.

[6]賈雨，鄧世武，姚興苗，等.基于約束粒子群優(yōu)化的克里金插值算法[J].成都理工大學(xué)學(xué)報：自然科學(xué)版， 2015， 42(1)： 104-109.

JIA Yu， DENG Shiwu， YAO Xingmiao， et al. Kriging interpolation algorithm based on constraint particle swarm optimization[J]. Journal of Chengdu University of Technology：science & technology edition， 2015， 42(1)： 104-109.

[7]顏華，李欣，王善輝.基于最小二乘法和克里金插值的三維溫度場重建[J].沈陽工業(yè)大學(xué)學(xué)報， 2014， 36(3)： 303-307.

YAN Hua， LI Xin， WANG Shanhui. Reconstruction of three-dimensional temperature field based on least-square method and kriging interpolation[J]. Journal of Shenyang University of Technology， 2014， 36(3)： 303-307.

[8]徐馳，曾文治，黃介生，等.基于高光譜與系統(tǒng)克里金的土壤耕作層含水率反演[J].農(nóng)業(yè)工程學(xué)報， 2014， 30(13)： 94-103.

XU Chi， ZENG Wenzhi， HUANG Jiesheng， et al. Modelling soil water content in plow layer using co-kriging interpolation based on hyperspectral data [J]. Transactions of the Chinese society of agricultural engineering， 2014，30(13)： 94-103.

[9]ZHANG Guanyu， WANG Guoqiang， LI Xuefei， et al. Global optimization of reliability design for large ball mill gear transmission based on the Kriging model and genetic algorithm[J]. Mechanism and machine theory， 2013， 69:321-336.

[10]FANG Jianguang， GAO Yunkai， SUN Guangyong， et al. Optimization of spot-welded joints combined artificial bee Colony algorithm with sequential sriging optimization[J]. Advances in mechanical engineering， 2014: 1-10.

[11]DENG Mingguo， LI Wenchang， LI Bo， et al. Application of log kriging on estimated reserves of the 10-9 ore body of lutangba in the gejiu tin deposits[J]. Journal of China university of mining & technology， 2007， 17(2)： 286-289.

[12]劉永社，印興耀，賀維勝.空間相關(guān)分析因素對儲層建模中克里金估計結(jié)果的影響[J].石油大學(xué)學(xué)報：自然科學(xué)版， 2004， 28(2)： 24-29.

LIU Yongshe， YIN Xingyao， HE Weisheng. Effects of the factors for spatial correlation analysis on Kriging estimation in reservoir modeling[J]. Journal of the University of Petroleum： science & technology edition， 2004， 28(2)： 24-29.

本文引用格式：

王振武， 卜異亞， 馬鍵強. 地質(zhì)雷達數(shù)據(jù)中克里金插值采樣數(shù)據(jù)選擇算法[J]. 哈爾濱工程大學(xué)學(xué)報， 2017， 38(5)： 784-790.

WANG Zhenwu，BU Yiya，MA Jianqiang. Sampling data selection algorithm for Kriging interpolation based on ground penetrating radar data[J]. Journal of Harbin Engineering University， 2017， 38(5)： 784-790.

Sampling data selection algorithm for Kriging interpolation based on ground penetrating radar data

WANG Zhenwu，BU Yiya，MA Jianqiang

(School of Mechanical Electronic & Information Engineering， China University of Mining & Technology (Beijing)， Beijing 100083， China)

The rationality of the data sampling selection method is a prerequisite for judging whether the Kriging interpolation results are accurate. Considering the problem of how to properly choose sampling data from the probing data of ground penetrating radar， in this paper， we propose two sampling selection algorithms for Kriging interpolation: the dynamic three-dimensional globe cover (D3DGC) algorithm and the double reverse distance (DRD) algorithm. We verified the rationality of these two methods with respect to the ordinary Kriging algorithm based on particle swarm optimization. The experimental results show the D3DGC algorithm to reduce the error rate to 6.7% compared with the traditional method for an optimal number of sampling points. With respect to other sampling numbers， D3DGC also demonstrated obvious advantages (error rate<10%). The DRD algorithm error rate was 6.2% for the optimal number of sampling points when the parameter e and c values were 4， respectively. For other sampling numbers， the DRD algorithm also obtained better results than the traditional method.

Kriging algorithm; interpolation; ground penetrating radar; sampling selection algorithm; dynamic three-dimensional globe cover; double reverse distance; particle swarm optimization(PSO); error

2016-01-07.

日期：2017-04-26.

國家高技術(shù)研究發(fā)展計劃重大專項(2012AA12A308)；核設(shè)施退役及放射性廢物治理科研項目(FZ1402-08)；國家自然科學(xué)基金項目(61302157)；北京市高等學(xué)校青年英才計劃項目(YETP0939).

王振武(1978-)， 男， 副教授，博士.

王振武，E-mail: wangzhenwu@126.com.

10.11990/jheu.201601022

TP39

A

1006-7043(2017)05-0784-07

網(wǎng)絡(luò)出版地址：http://www.cnki.net/kcms/detail/23.1390.u.20170426.1152.046.html