丁輝,牛生杰,張澤嬌

(南京信息工程大學 中國氣象局氣溶膠與云降水重點開放實驗室,江蘇 南京 210044)

?

利用微脈沖激光雷達探測大氣混合層高度和氣溶膠垂直結構的初步研究

丁輝,牛生杰,張澤嬌

(南京信息工程大學 中國氣象局氣溶膠與云降水重點開放實驗室,江蘇 南京 210044)

采用2008年5月14日—12月28日安徽壽縣地區(qū)(116°46′55.02″E,32°33′30.18″N)MPL(micro pulse lidar,微脈沖激光雷達)資料,分析了當?shù)鼗旌蠈痈叨?mixed layer height,MLH)的日變化特征。利用梯度法、拐點法及歸一化梯度法反演了當?shù)鼗旌蠈痈叨?并將3種方法的反演結果與SONDE探空確定的MLH做相關性分析,結果表明:由梯度法反演得到的MLH與探空資料確定的MLH高度相關(相關系數(shù)為0.976),是3種方法中最適合確定壽縣地區(qū)MLH的反演算法。利用梯度法反演了晴天無云天氣條件下MLH的日變化,研究了當?shù)貧庀髼l件對MLH的影響,結果發(fā)現(xiàn):MLH增長率與前1 h地表溫度增長率有很好的相關性;地表風速的變化也會導致MLH的急劇變化;MLH的發(fā)展滯后于氣溶膠光學厚度0.5～3 h。

微脈沖激光雷達;混合層高度;歸一化相對后向散射;反演方法

0 引言

大氣邊界層指最靠近下墊表面的對流層底層,厚度從數(shù)百米到一兩千米,是大氣與下墊面直接發(fā)生相互作用的層次,是地球大氣之間物質和能量交換的橋梁,和人類關系最為密切。大氣邊界層的科學進展是氣象發(fā)展的重要標志之一(熊超超等,2010)。由于大氣邊界層結構,如混合層高度、卷夾層厚度等,與氣象預報、氣候預測、環(huán)境保護等密切相關,大氣邊界層的研究一直是大氣科學研究的熱點問題之一(Emeis et al.,2004)。

目前,通過實測資料確定混合層高度主要是利用探空資料獲得的位溫廓線(陸春松等,2011),這種方法切實可行且有一定的可靠性,但其觀測手段難以進行連續(xù)觀測。大氣激光雷達精細的空間分辨率、時間上可連續(xù)測量以及較高的探測精度,可實現(xiàn)對大氣邊界層結構的實時監(jiān)測。MPL(micro pulse lidar,微脈沖激光雷達)是一種彈性后向散射激光雷達,該系統(tǒng)結構緊湊,發(fā)射能量對人眼安全,其可靠性已在國內外很多大型試驗中得到了驗證,并部署在世界各處進行大氣氣溶膠和云層的長期自動監(jiān)測。MPL因其測量準確度高、時空分辨率高及探測范圍廣,是當前一種重要的主動遙感工具,被廣泛應用于大氣氣溶膠以及云層的連續(xù)、自動觀測。MPL區(qū)別于傳統(tǒng)激光雷達的收發(fā)共享光學路徑,也會引發(fā)其他的問題,諸如探測器駐留脈沖訂正和近端填充函數(shù)確定等,這些問題在數(shù)據(jù)處理過程中必須仔細考慮。

Chen et al.(2001)在日本駐波市進行了一年的微脈沖激光雷達連續(xù)觀測,發(fā)現(xiàn)混合層高度(mixed layer height,MLH)的變化趨勢與氣溶膠光學厚度(aerosol optical thickness,AOT)一致,只是時間上略有滯后。對MLH日變化及年變化的分析發(fā)現(xiàn),城市混合層的變化與氣象條件的季節(jié)變化有很好的相關性。Welton et al.(2002)在印度洋上空進行了氣溶膠垂直廓線及光學特性的觀測研究,發(fā)現(xiàn)在此區(qū)域海洋性邊界層的高度低于1 km,并提出了一個污染性海洋大氣氣溶膠垂直消光的新模式。郭本軍等(2008)利用微脈沖激光雷達觀測了大連地區(qū)2005年11月6日的一次沙塵天氣過程,發(fā)現(xiàn)PM10的時間變化與由MPL反演的近地面消光系數(shù)相當一致,它們之間的轉換系數(shù)在1.94到6.50 mg·m-3·km-1之間,說明MPL可以實現(xiàn)近地面氣溶膠消光系數(shù)的有效觀測。Voss et al.(2001)在1999年1月14日—2月8日進行的Aerosols99移動觀測計劃中,通過對微脈沖激光雷達和太陽光度計同時觀測得到的氣溶膠光學厚度的比較,發(fā)現(xiàn)兩者的觀測結果具有較好的一致性。賀千山和毛節(jié)泰(2005)利用北京地區(qū)的微脈沖激光雷達遙感資料及325 m鐵塔的觀測資料,基于導數(shù)法提出了一種自動判別邊界層高度的歸一化梯度法,并探討了利用MPL的遙感結果提高對流混合層發(fā)展模式預報精度的可行性。Campbell et al.(2008)利用改進的微脈沖激光雷達數(shù)據(jù),提出了一個基于閾值法計算云及氣溶膠層高度的新算法。Huang et al.(2008)利用2007年4—5月共55 d的微脈沖激光雷達數(shù)據(jù),結合太陽光度計資料,分析了黃土高原上空氣溶膠的垂直分布特征。

1 MPL基本原理及資料來源

MPL的探測原理可以用米散射激光雷達方程說明,其形式如下(Campbell et al.,2002):

nr(r)=[(Oc(r)CEβ(r)T2/r2)+nb(r)+ nap(r)]/D[n(r)]。

(1)

BNR(r)=Cβ(r)T2(r)。

(2)

可見NRB信號只受大氣總后向散射系數(shù)或消光系數(shù)的影響,表明該信號的時空分布能夠在一定程度上反映氣溶膠濃度的變化趨勢,它隨高度的變化表征了實際大氣中氣溶膠的垂直分布特征,是計算邊界層高度和相關大氣光學參數(shù)的基本量,可用于計算云滴散射截面、云的光學厚度、行星邊界層高度、氣溶膠消光系數(shù)廓線及光學厚度等(Campbell et al.,2002)。而線性去偏振比δL可作為區(qū)分氣溶膠粒子類型及形狀的一個潛在工具(Flynn et al.,2007)。研究發(fā)現(xiàn):液態(tài)球形粒子對線性去偏振光的后向散射信號是完全去極化的δL=0;純分子大氣δL典型值為1.4%;大氣氣溶膠δL量級小于10%;非球形顆粒物的δL值在40%～70%之間;簡單冰晶δL值為50%,當有液體顆粒物存在時這個值會顯著減小(Mishra et al.,2010)。結合NRB信號與線性去偏振比可以有效地分析邊界層結構的日變化趨勢。

隨著雷達探測靈敏度的提高,使得人們有可能在晴空條件下獲得更多的大氣結構細節(jié)(魏鳴等,2007)。本文根據(jù)安徽壽縣地區(qū)2008年5月14日—12月28日晴天無云條件下的MPL遙感探測

資料,結合SONDE探空資料、常規(guī)氣象資料以及Micaps資料分析當?shù)鼗旌蠈拥淖兓卣鳌?/p>

圖1 由NRB信號反演得到的2008年11月29日17時(a)和12月9日13時(b)1 h平均氣溶膠消光系數(shù)的垂直廓線(單位：km-1)

2 氣溶膠消光系數(shù)的反演

2.1 Fernald反演方法

根據(jù)Fernald方法(Fernald et al.,1972),如果已知某一高度zc處氣溶膠粒子消光系數(shù),則zc高度以下各高度上氣溶膠粒子的消光系數(shù)的后向積分公式為:

(3)

其中:X(i)為NRB信號強度;l為MPL的距離分辨率(0.015 km);Sa=50(單位:sr);Sm=σm(z)/βm(z)=8π/3(單位:sr)(袁松等,2005)。σm(z)利用SONDE探空資料根據(jù)大氣分子瑞利散射相關理論計算。參考高度zc選取近乎不含氣溶膠粒子的清潔大氣層所在的高度(邱金桓等,2003),在此高度上X(i)/βm(z)值應最小,一般在對流層頂附近,受MPL信噪比的限制,本文zc的選取如下:白天(08—18時)zc在3～6 km內尋找,夜晚zc在5～10 km內尋找。參考高度處的氣溶膠后向散射系數(shù)通過氣溶膠散射比:R(zc)=βa(zc)/βm(zc)=1.02來確定。

2.2 大氣分子光學參數(shù)的確定

按照瑞利散射理論,考慮到非各向同性,在同性球形小粒子的散射系數(shù)公式的基礎上增加一個訂正項(賀千山,2006)后,得到根據(jù)SONDE探空資料反演大氣分子消光系數(shù)的公式為:

T0=288 K,

P0=1 013.25 hPa,λ=532 nm。

(4)

其中:r是粒子半徑;λ是激光波長。圖1給出了由NRB信號反演得到的2008年11月29日17時和12月9日13時的1 h平均氣溶膠消光系數(shù)垂直廓線。相應時刻的AOT分別為0.39和0.57。11月29日17時(圖1a)對流層下部有一充分發(fā)展的較厚的氣溶膠層即混合層,厚度約為2 km。大氣氣溶膠主要集中混合層中,混合層以上大氣比較清潔,氣溶膠消光系數(shù)垂直廓線呈指數(shù)遞減,其變化趨勢接近于瑞利散射。12月9日(圖1b)氣溶膠消光系數(shù)隨高度的增加降低較緩慢,污染物的垂直延伸范圍較大。氣溶膠消光系數(shù)垂直分布是一個明顯的多層結構,1.3 km、2.8～4.0 km處存在消光系數(shù)較大的抬升氣溶膠層,圖2可以看出12月9日15時之前在1.5～4.0 km持續(xù)存在著線性去偏振比較大的抬升氣溶膠層,此氣溶膠層有可能是外地移入本地或者前一天殘留的污染物。

圖2 2008年12月9日歸一化相對后向散射信號(a)及線性去偏振比(b)的時間—高度剖面

3 混合層高度的反演

3.1 混合層高度反演方法

混合層與自由大氣層相比更潮濕、氣溶膠含量更大,對激光束的散射量更大,因此,激光雷達很容易探測到混合層與自由大氣之間的邊界(Melfi et al.,1985)。Flamant et al.(1997)將激光雷達距離訂正信號對高度的一階導數(shù)最小值所在的高度定義為混合層頂高,即為梯度法。Menut et al.(1999)提出的拐點法將距離訂正信號對高度的二階導數(shù)最小值所在高度定義為混合層頂高。賀千山和毛節(jié)泰(2005)提出了歸一化梯度法,去除了信號強度隨高度變化對反演結果的影響。

本文首先對NRB信號做6 min平均,此平均時間與邊界層中大氣湍流的特征時間尺度一致,分別利用上述三種方法反演MLH,并將反演結果與探空資料進行比較。SONDE探空得到的溫度,位溫、相對濕度及風速等特征量在混合層頂均會表現(xiàn)出較大的梯度,可驗證MLH反演方法的正確性(李麗平等,2013)。圖3a、b、c是上述3種方法得到的11月12日13時25分的MLH,分別為0.877、0.832及0.892 km,與圖3d、e、f、g中強逆溫出現(xiàn)的高度(0.867 km)比較一致,拐點法的結果稍低于梯度法,與Dupont et al.(1994)的研究結果一致。溫度及位溫垂直廓線在0.867 km處開始出現(xiàn)強逆溫,此高度以上為頂蓋逆溫層,該層之下水汽幾乎均勻混合,同時刻的NRB垂直廓線也有相似的變化趨勢,表明大氣氣溶膠在混合層中的分布類似于水汽分布較均勻。其中,位溫廓線是利用探空資料中溫度與氣壓垂直廓線根據(jù)位溫定義的計算值而非實測值。

圖3 梯度法、拐點法及歸一化梯度法反演的混合層高度與探空廓線的對比 a.梯度法由NRB信號(實線)得到的一階導數(shù)(虛線);b拐點法由NRB信號(實線)得到的二階導數(shù)(虛線);c.歸一化梯度法由NRB距離訂正信號(實線)得到的一階導數(shù)(虛線);d.溫度(K)的垂直廓線;e.位溫(K)的垂直廓線;f.相對濕度(%)的垂直廓線;g.風速(m/s)的垂直廓線

選擇10—12月共31 d的數(shù)據(jù),分別對由上述3種方法反演得到的MLH與SONDE探空確定的MLH做相關性分析(圖4),結果發(fā)現(xiàn):梯度法反演得到的MLH與探空資料確定的MLH相關性最好,r=0.976;其次為歸一化梯度法,r=0.835;拐點法的結果最差。其中梯度法與歸一化梯度法確定的MLH之間的相關系數(shù)r=0.851,說明梯度法和歸一化梯度法均可以很好地確定壽縣地區(qū)的MLH。本文選擇相關性最好的梯度法反演MLH的日變化,由于MPL的低信噪比,有時這種反演方法會產(chǎn)生錯的MLH值。另外,其他氣象原因,如:云層、霾等也會偶爾阻礙MLH的反演。為了去除上述錯誤值,檢查整個MLH數(shù)據(jù)列去除其中被認為是錯誤的數(shù)據(jù)點,最后對整個MLH數(shù)據(jù)列做5點滑動平均得到MLH的日變化。

圖4 由梯度法(a)、拐點法(b)及歸一化梯度法(c)反演得到的MLH與SONDE探空廓線確定的MLH的相關性分析

3.2 MLH日變化

圖5給出了NRB信號及由3種方法反演得到的MLH的日變化特征:日出前,邊界層比較穩(wěn)定,混合層較淺;08時之后逆溫層被打破混合層開始緩慢加厚,10時左右混合層開始迅速發(fā)展且殘余層被完全破壞,邊界層結構從早晨的分層結構轉變?yōu)槌浞只旌系幕旌蠈?12時左右MLH達到最大值并一直維持至20時左右。從圖5中可以看出:3種反演方法得到的MLH在白天比較一致,當有殘余層存在時,梯度法追蹤的是沉降層的邊界,而歸一化梯度法追蹤的是殘余層頂?shù)母叨?拐點法的結果介于梯度法與歸一化梯度法之間。圖6a是12月6日MLH與地表溫度的日變化,白天(10—18時)MLH的變化趨勢與地表溫度的變化趨勢比較一致,只是MLH的增長略滯后于地表溫度1～2 h。利用11月11—13日的數(shù)據(jù)分析了MLH增長率與前1 h地表溫度之間的相關性(圖6b),兩者的線性相關系數(shù)r=0.787,與Stull(1985)的研究結果一致,說明邊界層中熱力強迫是混合層發(fā)展的主導動力。

圖5 2008年11月12日NRB時間—高度剖面以及由梯度法、歸一化梯度法及拐點法反演得到的MLH的日變化(單位：km)

圖6 MLH(km)和地表溫度(K)的日變化(a)以及MLH增長率與前1 h地表溫度增長率的相關性分析(b)

圖7 2008年12月6日地表風速(m/s)及風向的日變化

在分析的過程中發(fā)現(xiàn)地表風速是控制MLH發(fā)展的另一個重要因素,地表風速、風向的變化也會導致MLH的劇烈變化。從圖7可以看出,08時30分之后風向由偏南風轉為北風,且風速急劇增加,導致相應時刻MLH的急劇增長。11—16時風速保持在3 m/s左右,MLH也持續(xù)在1.2 km高度以上,比較穩(wěn)定,此后風速與MLH表現(xiàn)出了相似的波動性。

圖8給出了對流層中氣溶膠光學厚度AOT的日變化,MLH與AOT的變化趨勢基本一致,只是MLH的發(fā)展在時間上略有滯后。通過對長期資料的分析研究發(fā)現(xiàn),MLH略滯后于AOT的發(fā)展0.5～3 h,MLH這種滯后性的可能原因是清晨近地面霾的擴散速度比邊界層中熱力對流快。

3.3 MLH的逐月變化

根據(jù)Micaps云量數(shù)據(jù)選擇MPL觀測期間晴天無云條件下的MPL數(shù)據(jù)(其中,7月沒有符合要求的觀測數(shù)據(jù)),利用梯度法反演每日的MLH并計算月平均值,分析MLH的逐月變化特征。圖9給出了2008年5—12月壽縣地區(qū)MLH及地表溫度的逐月變化,發(fā)現(xiàn)MLH開始急劇增加的時間由從5月的06時左右逐漸推遲到了12月的09時左右。最大混合層高度也逐月推遲,10月月平均MLH最大值為1.23 km,出現(xiàn)于13時;11月月平均MLH最大值為1.24 km,出現(xiàn)于14時;12月月平均MLH的最大值為0.834 km,出現(xiàn)于14時。地表溫度月平均值的增長提前于MLH月平均值1～2 h,與日平均值得變化趨勢一致。12月混合層發(fā)展最弱,平均MLH低于0.834 km。5—9月由于數(shù)據(jù)樣本資料較少,不能夠反映真實的月變化趨勢。

圖8 2008年12月6日(a)和9日(b)MLH(km)和AOT的日變化

圖9 2008年5月(a)、6月(b)、8月(c)、9月(d)、10月(e)、11月(f)、12月(g)月平均的MLH(km)和地表溫度(K)的日變化(7月無符合要求的觀測數(shù)據(jù))

4 結論

1)線性去偏振比作為氣溶膠粒子類型及形狀的指示劑,可結合NRB時間高度分布來分析大氣氣溶膠及邊界層的結構特征。

2)用梯度法、拐點法及歸一化梯度法反演了MLH,選擇晴天無云天氣條件下的MLH與探空資料確定的MLH做相關性分析,結果表明梯度法是3種方法中最適合壽縣地區(qū)的MLH反演算法。

3)通過對長期資料的分析發(fā)現(xiàn),MLH與AOT的日變化趨勢基本一致,只是MLH的發(fā)展滯后于AOT約0.5～3 h;地表風速、風向的變化也會導致MLH的劇烈變化。

4)結合地表溫度、風速、能見度等常規(guī)氣象資料分析MLH的日變化,結果發(fā)現(xiàn)相同氣象條件下MLH的日變化趨勢有可能會表現(xiàn)出較大的差異,即MLH的日變化與引發(fā)它的氣象條件并不一一對應。

由于有的月數(shù)據(jù)樣本較少,MLH月平均值的結果缺少代表性,因此,MLH逐月變化特征還有待進一步研究。

郭本軍,賀千山,李成才,等.2008.大連一次沙塵天氣的遙感研究[J].大氣科學學報,31(5):624-632.

賀千山.2006.微脈沖激光雷達研究對流層氣溶膠光學屬性和分布特征[D].北京:北京大學.

賀千山,毛節(jié)泰.2005.北京城市大氣混合層與氣溶膠垂直分布觀測研究[J].氣象學報,63(3):374-384.

李麗平,成麗萍,靳莉莉,等.2013.北太平洋混合層深度異常及其與中國夏季降水的關系[J].大氣科學學報,36(1):20-28.

陸春松,牛生杰,岳平,等.2011.南京冬季霧多發(fā)期邊界層結構觀測分析[J].大氣科學學報,34(1):58-65.

邱金桓,鄭斯平,黃其榮,等.2003.北京地區(qū)對流層中上部云和氣溶膠的激光雷達探測[J].大氣科學,27(1):1-7.

魏鳴,秦學,王嘯華,等.2007.南京地區(qū)大氣邊界層晴空回波研究[J].大氣科學學報,30(6):736-744.

熊超超,謝麗萍,吳維.2010.大氣邊界層探測技術初步分析[J].氣象水文海洋儀器,3:77-80.

袁松,辛雨,周軍.2005.合肥市郊低層大氣的激光雷達研究[J].大氣科學,29(3):387-395.

Campbell J R,Hlavka D L,Welton E J,et al.2002.Full-time,eye-safe cloud and aerosol lidar observation at atmospheric radiation measurement program sites:Instrument and data processing[J].J Atmos Oceanic Technol,19:431-442.

Campbell J R,Sassen K,Welton E J,et al.2008.Elevated cloud and aerosol layer retrievals from micropulse lidar signal profiles[J].American Meteorological Society,25:685-700.

Chen W,Kuze H,Uchiyama A,et al.2001.One-year observation of urban mixed layer characteristics at Tsukuba Japan using a micro pulse lidar[J].Atmos Environ,35:4273-4280.

Dupont E,Pelon J,Flamant C,et al.1994.Study of the moist convective boundary layer structure by backscatter lidar[J].Bound-layer Meteor,69:1-25.

Emeis S,Munkel C,Vogt S,et al.2004.Atmospheric boundary-layer structure from simultaneous sodar,rass,and ceilometer measurements[J].Atmos Environ,38(2):273-286.

Fernald F G,Herman B H,Reagan J A,et al.1972.Determination of aerosol height distribution by lidar[J].J Appl Meteor,11:482-489.

Flamant C,Pelon J,Flamant P H,et al.1997.Lidar determination of the entranment zone thickness at the top of the unstable marine atmospheric boundary layer[J].Bound-layer Meteor,83:247-284.

Flynn C J,Mendoza A,Zheng Y H,et al.2007.Novel polarization-sensitive micropulse lidar measurement technique[J].Optics Express,15(6):2785-2790.

Huang J P,Huang Z W,Bi J R,et al.2008.Micro-pulse lidar measurements of aerosol vertical structure over the Loess Plateau [J].Atmospheric and Oceanic Science Letters,1(1):8-11.

Melfi S H,Sphinhirne J D,Chou S H,et al.1985.Lidar observations of the vertically organized convection in the planetary boundary layer over the ocean[J].Climate Appl Meteor,24:806-821.

Menut L,Flamant C,Pelon J,et al.1999.Urban boundary layer height determination form lidar measurements over the Paris area[J].Applied Optics,38(6):945-954.

Mishra M K,Rajeev K,Bijoy V,et al.2010.Micropulse lidar observations of mineral dust layer in the lower troposphere over the southwest coast of Peninsular India during the Asian summer monsoon season[J].Journal of Atmospheric and Solar-terrestrial Physics,72(17):1251-1259.

Stull R B.1985.An introduction to boundary layer meteorology[M].Netherlands:Kluwer Academic Publisher.

Voss K J,Welton E J,Quinn P K,et al.2001.Lidar measurements during aerosols99 [J].J Geophys Res,106(D18):20821-20831.

Welton E J,Voss K J,Quinn P K,et al.2002.Measurements of aerosol vertical profiles and optical properties during INDOEX 1999 using micropulse lidars[J].J Geophys Res,107(D19):1029-1049.

(責任編輯：張福穎)

Measurements of mixed layer height and aerosol vertical profiles using a micro pulse lidar

DING Hui,NIU Sheng-jie,ZHANG Ze-jiao

(Key Laboratory for Aerosol-Cloud-Precipitation of China Meteorological Administration,NUIST,Nanjing 210044,China)

Based on the MPL(micro pulse lidar) data over Shouxian area(32°33′30.18″N,116°46′55.02″E),Anhui Province from 14 May to 28 December 2008,this paper analyzes the diurnal variation of mixed layer height(MLH).Comparisons are made between MLHs retrieved by MPL retrieval methods(the gradient method,the inflection point method and the normalized gradient method) and MLH determined by the sounding data measured by SONDE.Results show that MLH retrieved from the lidar data using the gradient method are highly correlated with MLH measured by SONDE(the correlation coefficient is 0.976),and the gradient method is the most suitable retrieval algorithm in Shouxian area.The diurnal cycle of MLH on sunny days are retrieved by the gradient method,and the effect of meteorological conditions on MLH are also studied.Results show that the growth rate of MLH agrees well with the hourly increase of surface temperature one hour before.Surface wind speed and direction can also lead to rapid change of MLH.The development of MLH exhibits temporal behavior similar to that of the aerosol optical thickness,but with 0.5—3 h lag time.

micro pulse lidar;mixed layer height;normalized relative backscatter;retrieval method

2012-01-12；改回日期：2014-11-06

江蘇省青藍工程“云霧降水與氣溶膠研究”創(chuàng)新團隊項目資助;江蘇省高校優(yōu)勢學科建設工程資助(PAPD);公益性行業(yè)(氣象)科研專項(GYHY201006036)

牛生杰,教授,博士生導師,研究方向為云霧降水物理學,niusj@nuist.ecu.cn.

10.13878/j.cnki.dqkxxb.20120112001.

1674-7097(2015)01-0085-08

P407.5

A

10.13878/j.cnki.dqkxxb.20120112001

丁輝,牛生杰,張澤嬌.2015.利用微脈沖激光雷達探測大氣混合層高度和氣溶膠垂直結構的初步研究[J].大氣科學學報,38(1):85-92.

Ding Hui,Niu Sheng-jie,Zhang Ze-jiao.2015.Measurements of mixed layer height and aerosol vertical profiles using a micro pulse lidar[J].Trans Atmos Sci,38(1):85-92.(in Chinese)