岳榮剛 王世濤 王虎妹 尹可 金挺 藏潔 晏廷飛 李歡 唐紹凡
基于慣性參考基準(zhǔn)的像移測量與圖像復(fù)原
岳榮剛1王世濤1王虎妹1尹可1金挺1藏潔2晏廷飛3李歡4唐紹凡4
(1 遙感衛(wèi)星總體部,北京 100094)(2 北京空間飛行器總體設(shè)計部,北京 100094)(3 北京衛(wèi)星環(huán)境工程研究所,北京 100094)(4 北京空間機(jī)電研究所,北京 100094)
衛(wèi)星上存在高、中、低頻率的微振動源,引起光學(xué)載荷光軸抖動,造成遙感圖像成像品質(zhì)下降。文章提出了一種基于慣性參考基準(zhǔn)的像移測量與圖像復(fù)原方法,介紹了該方法的原理和涉及的關(guān)鍵部件。為驗證該方法的可行性和性能,研制了一套地面試驗驗證系統(tǒng),包括驗證相機(jī)(含慣性參考基準(zhǔn)單元、哈特曼傳感器)、積分球、靶標(biāo)、平行光管、隔振平臺、激振器、六自由度支撐臺,以及相應(yīng)的圖像復(fù)原軟件等。測量了試驗環(huán)境的背景噪聲,開展了像移測量精度試驗,在此基礎(chǔ)上完成了多種工況下的圖像恢復(fù)試驗。結(jié)果表明,基于慣性參考基準(zhǔn)的像移測量誤差均方根值不超過0.12像元,基于測得的像移數(shù)據(jù)可顯著提升像質(zhì),圖像的調(diào)制傳遞函數(shù)一般可提升至原圖的1.38~1.52倍。該方法獲得的像移數(shù)據(jù)還可以反饋給穩(wěn)像鏡,實時補(bǔ)償星上微振動,有助于直接獲取高品質(zhì)圖像。
慣性參考 微振動 像移 圖像復(fù)原 遙感衛(wèi)星
隨著空間應(yīng)用能力的提高,對光學(xué)遙感衛(wèi)星、激光通信衛(wèi)星等的性能要求也越來越高,涌現(xiàn)出一系列高分辨率光學(xué)成像衛(wèi)星,如:GeoEye-1[1]、WorldView-4[2]、“高分四號”衛(wèi)星[3],以及未來的“巡天號”光學(xué)艙等,鏡頭口徑越來越大,角分辨率越來越高。
然而,衛(wèi)星上普遍存在多種頻率的擾振源,如:太陽翼驅(qū)動機(jī)構(gòu)、天線驅(qū)動機(jī)構(gòu)、控制力矩陀螺、推進(jìn)器、制冷機(jī)等等[4],加上衛(wèi)星的姿態(tài)不穩(wěn)定因素,綜合起來會引起星上顫振,其覆蓋極低頻到上千赫茲高頻的寬頻帶[5-8]。這些因素會引起光學(xué)載荷的光軸抖動,導(dǎo)致成像相對焦平面產(chǎn)生移動或者轉(zhuǎn)動,造成遙感圖像成像品質(zhì)下降[9-10]。
美國宇航局戈達(dá)德航天中心的研究表明,顫振的能量主要集中在中低頻區(qū)[11]。為測量顫振對光軸的影響,研究人員提出了多種像移測量方法。應(yīng)用比較多的是基于光學(xué)聯(lián)合變換相關(guān)器(Joint Transform Correlator,JTC)進(jìn)行空間相機(jī)像移測量的方法[12-17]。該方法需要的探測設(shè)備包括高速CCD和光學(xué)相關(guān)器,高速CCD測得的相鄰兩幀圖像經(jīng)JTC處理后得到相應(yīng)的像移量,所得像移量被反饋到穩(wěn)像鏡,確保了成像與CCD間的相對穩(wěn)定,可得到優(yōu)良的成像。這種方法將測量位置置于焦面處,像移測量精度高、速度快。
文獻(xiàn)[18-22]則采用了基于遙感圖像的像移探測方法,在相機(jī)積分過程中,用高速相機(jī)記錄下多幀連續(xù)圖像,通過配準(zhǔn)算法得到相鄰幀間的位移量,從而實現(xiàn)像移測量,基于像移構(gòu)建圖像復(fù)原所需的點擴(kuò)散函數(shù)(Point Spread Function,PSF),通過后期處理的方式進(jìn)行圖像復(fù)原。該方法能夠測量成像位置像移的真實狀態(tài),測量精度較高,但是盲點頻率處的像移無法測量,且存在滯后性。
無論是基于JTC的像移測量方法,還是基于遙感圖像的像移探測方法,均要求被拍攝對象的圖像層次豐富,且對信噪比要求較高[17]。一旦對遇到遙感圖像較平滑或信噪比低的情況,精度難以保證。針對上述問題,本文提出了一種基于慣性參考基準(zhǔn)的高精度像移測量方法,該方法不依賴于地面景物特點,也不依賴于光照條件,以自身提供的參考基準(zhǔn)保證像移測量系統(tǒng)可全天時、全天候有效工作。
基于慣性參考基準(zhǔn)的像移測量方法是一套帶慣性傳感器的參考激光系統(tǒng)工作。將激光與成像光共光路引入相機(jī)系統(tǒng),直接測量像移信息,用于驅(qū)動穩(wěn)像鏡實時進(jìn)行補(bǔ)償,或者用像移信息后期恢復(fù)退化的圖像,減小顫振對成像品質(zhì)的影響。測量方法原理如圖1所示。慣性參考基準(zhǔn)單元發(fā)出一束參考激光束并通過角錐棱鏡引入光學(xué)系統(tǒng),參考光束和成像光通過相同的光路(相機(jī)主鏡、相機(jī)次鏡、穩(wěn)像鏡等)后到達(dá)位置探測器(位置探測器位于相機(jī)焦面的邊緣視場,并與相機(jī)焦面剛性固定在一起);位置探測器探測到參考激光束的抖動信息;該信息攜帶了參考光束相對于焦面的抖動,再融合慣性參考基準(zhǔn)單元測量到的相機(jī)姿態(tài)變化信息;該融合信息通過閉環(huán)反饋控制穩(wěn)像鏡,穩(wěn)像鏡通過快速運動補(bǔ)償微振動對成像的影響,從而得到優(yōu)良的遙感圖像?;蛘呋谙褚菩畔?,用軟件恢復(fù)退化的圖像,也可以減小星上顫振對成像的影響。
圖1 基于慣性參考基準(zhǔn)的像移測量方法原理
慣性參考基準(zhǔn)單元主要包括慣性傳感器和激光源(如圖2所示)。兩軸慣性傳感器的敏感軸互相垂直,用于高精度測量(積分時間內(nèi))兩個垂直方向上衛(wèi)星的姿態(tài)變化和微振動;激光器用于發(fā)出一束參考光束。由于激光器與慣性傳感器剛性連接,其光束在慣性空間的絕對指向可精確測量。慣性傳感器長時間工作后的累積誤差可用星敏感器在軌標(biāo)校。
圖2 慣性參考基準(zhǔn)單元模型
像移測量試驗系統(tǒng)原理如圖3所示。靶標(biāo)置于平行光管的焦面上以模擬無窮遠(yuǎn)的目標(biāo),靶標(biāo)被積分球照亮后,其像經(jīng)過相機(jī)光學(xué)系統(tǒng)在CCD上成像。平行光管和靶標(biāo)安裝在隔振平臺上,隔離地面微振動的影響,保證目標(biāo)光的穩(wěn)定性。相機(jī)系統(tǒng)(含慣性參考基準(zhǔn)單元)安裝到六自由度支撐臺上,模擬相機(jī)在軌自由狀態(tài),激振器產(chǎn)生模擬的微振動,施加到相機(jī)上。
試驗中所用的相機(jī)主鏡光學(xué)口徑330mm,焦距3 300mm,CCD分辨率1 024像元×1 024像元,單像元尺寸12μm。相機(jī)上集成了慣性參考基準(zhǔn)單元、角錐棱鏡、位置探測器(Position Sensitive Detector, PSD)等,如圖4所示。為了驗證像移測量精度,另外集成了1臺哈特曼傳感器,該傳感器通過對子孔徑圖像進(jìn)行高速成像,計算每幀圖像的偏移量,實現(xiàn)像移的直接測量。經(jīng)驗證,該哈特曼傳感器對試驗相機(jī)的像移測量精度優(yōu)于0.1像元(均方根值)。
圖3 基于慣性參考基準(zhǔn)的像移測量試驗系統(tǒng)原理
圖4 驗證相機(jī)三維模型
像移測量試驗依賴于精密的測量系統(tǒng),為了減小測量誤差,需保證平行光管和靶標(biāo)系統(tǒng)的穩(wěn)定,因此將其置于隔振臺上用于隔離地面的微振動。試驗開始前用高精度角位移傳感器測量了隔振臺工作面繞方向和方向的角位移噪聲(坐標(biāo)定義如圖4所示),測量環(huán)境條件與試驗時相同:打開積分球風(fēng)扇,停止周邊其余振動試驗,關(guān)閉試驗區(qū)的通風(fēng)裝置,遮光簾處于關(guān)閉狀態(tài)等。
在以上條件下測得兩個方向的微振動曲線如圖5所示,相應(yīng)的噪聲(均方根值)如表1所示??梢姼粽衽_在1~200Hz范圍內(nèi)的背景噪聲引起的繞、軸的角位移最大值為0.037″,小于1/20像元。而在3~200Hz頻段上,隔振臺背景噪聲引起的繞、軸的角位移最大值為0.014″,小于1/50像元??蓾M足本試驗精度需求。
圖5 隔振臺背景噪聲微振動曲線
表1 隔振臺背景噪聲引起的角位移均方根值
Tab.1 The Root mean square value of angular displacement caused by the background noise of the vibration isolator
為了驗證基于慣性參考基準(zhǔn)的像移測量精度,將其測量數(shù)據(jù)與哈特曼傳感器所測數(shù)據(jù)進(jìn)行了比較,該哈特曼傳感器對相機(jī)的像移測量精度優(yōu)于0.1像元,作為精度驗證的參考。試驗過程如下:
1)用激振器對相機(jī)施加沿軸方向的激勵;
2)分別用激光陀螺、PSD和哈特曼傳感器記錄相應(yīng)的測量數(shù)據(jù);
3)將激光陀螺和PSD的測量數(shù)據(jù)進(jìn)行融合處理,得到基于慣性參考基準(zhǔn)的像移測量曲線;
4)對哈特曼傳感器所測數(shù)據(jù)進(jìn)行處理,得到另一條像移曲線。
比較兩條曲線(如圖6、圖7所示),可得到以下結(jié)論:
圖6 20Hz定頻正弦激勵條件下像移測量曲線對比
圖7 110Hz定頻正弦激勵條件下像移測量曲線對比
1)在圖6、圖7所示的各種工況下,兩條曲線像移測量誤差均方根值不超過0.12像元??梢娀趹T性參考基準(zhǔn)的像移測量方法也能實現(xiàn)較高的測量精度。
2)激光陀螺和PSD的采樣率均為10kHz,遠(yuǎn)高于哈特曼傳感器的550Hz,所以其測量曲線包含了像移的高頻信息,能夠更精細(xì)地描述像移過程。
為進(jìn)一步驗證本方法的像移測量性能和應(yīng)用效果,開展了圖像恢復(fù)試驗。用測得的像移信息對退化的圖像進(jìn)行恢復(fù),得到復(fù)原后的圖像,并分別計算退化圖像和復(fù)原圖像的MTF值,定量評估圖像恢復(fù)效果。
基于慣性參考基準(zhǔn)的圖像恢復(fù)試驗方法和流程如圖8所示,說明如下:
1)用相機(jī)拍攝一張退化的圖像,同時用慣性參考基準(zhǔn)系統(tǒng)測得光學(xué)系統(tǒng)的像移信息;
2)基于像移信息對像移模型進(jìn)行重構(gòu),得到積分時間內(nèi)像移綜合退化模型;
3)基于綜合退化模型,用自研的圖像恢復(fù)軟件對退化的圖像進(jìn)行恢復(fù),得到清晰圖像;
4)分別計算退化圖像和復(fù)原圖像在奈奎斯特頻率處的MTF值,進(jìn)行定量判斷。
圖8 基于慣性參考基準(zhǔn)的圖像恢復(fù)試驗方法和流程
圖像復(fù)原方面,非盲復(fù)原算法取得了長足的發(fā)展。幾種具有代表性的圖像非盲復(fù)原算法如維納濾波算法[23]、總變分(Total Variation,TV)正則化圖像復(fù)原算法[24]等,都在各自適用的領(lǐng)域獲得了理想的圖像復(fù)原效果。另外,文獻(xiàn)[25]提出了一種分段局部正則化RL(Richardson-Lucy)方法,構(gòu)建了一種新型的正則項,使得復(fù)原圖像中的噪聲和邊緣振鈴效應(yīng)得到有效控制。
本文研制的圖像恢復(fù)軟件提供了三種算法,分別是TV(Total Variation,總變分)方法、分段局部正則化RL方法和一種盲恢復(fù)方法。在本文的圖像恢復(fù)過程中,比較了三種圖像恢復(fù)算法的效果,結(jié)論如下:TV方法對噪聲和振鈴的抑制效果明顯,但圖像細(xì)節(jié)損失較大;盲恢復(fù)方法在一定程度上能抑制振鈴,但會產(chǎn)生較大噪聲,且會引入偽信息;分段局部正則化RL方法較真實地保留了圖像細(xì)節(jié),并可較好地抑制振鈴,產(chǎn)生的噪聲也不明顯,是圖像恢復(fù)效果最好的一種算法。因此,本文中的圖像恢復(fù)均采用了分段局部正則化RL方法。
針對建筑靶標(biāo)進(jìn)行了圖像恢復(fù)試驗,試驗原圖如圖9所示。試驗結(jié)果如圖10和圖11所示,分別給出了不同振動激勵工況下,退化圖像與復(fù)原圖像的對比圖:圖10試驗工況條件為20Hz定頻激勵,抖動量3.5像元,積分時間80ms;圖11試驗工況條件為110Hz定頻激勵,抖動量2像元,積分時間80ms。從視覺上定性判斷,恢復(fù)后的圖像無論是對比度還是銳度都得到了明顯提升,無明顯振鈴,噪聲也在可控范圍內(nèi)。兩種工況下,圖像的MTF值分別可提升至原圖的1.52和1.38倍??梢姡?dāng)像元的抖動量在3.5個像元以內(nèi),基于慣性參考基準(zhǔn)的像移測量方法可明顯提升圖像品質(zhì)。
圖9 試驗用靶標(biāo)原始圖像
圖10 MTF從0.044提升至0.067
圖11 MTF從0.037提升至0.051
本文針對提出的基于慣性參考基準(zhǔn)的像移測量方法,研制了一套相應(yīng)的試驗驗證系統(tǒng)?;谠撓到y(tǒng)開展了多種工況的像移測量試驗,并利用獲得的像移數(shù)據(jù)對退化的圖像進(jìn)行了復(fù)原,結(jié)果表明:基于慣性參考基準(zhǔn)的像移測量方法實現(xiàn)了較高的測量精度,與哈特曼傳感器的測量數(shù)據(jù)相比,相對誤差不超過0.12像元(RMS);當(dāng)像元的抖動量在3.5個像元以內(nèi),基于慣性參考基準(zhǔn)的像移測量方法可明顯提升成像品質(zhì),圖像的MTF可提升至原圖的1.38~1.52倍;獲得的像移數(shù)據(jù)還可以反饋給穩(wěn)像鏡,實時補(bǔ)償星上微振動對成像品質(zhì)的影響,有效提升光學(xué)遙感衛(wèi)星圖像品質(zhì)。該方法具有測量精度高、速度快、不受地面光照條件影響、可全天時工作的優(yōu)點,可用于其余有高精度光軸抖動測量需求的場合。
[1] Wikipedia.GeoEye-1[EB/OL]. [2020-10-17]. https://en.wikipedia.org/wiki/GeoEye-1.
[2] Wikipedia.WorldView-4[EB/OL]. [2020-10-17]. https://en.wikipedia.org/wiki/ WorldView-4.
[3] 王殿中, 何紅艷. “高分四號”衛(wèi)星觀測能力與應(yīng)用前景分析[J]. 航天返回與遙感, 2017, 38(1): 98-106. WANG Dianzhong, HE Hongyan. Observation Capability and Application Prospect of GF-4 Satellite[J]. Spacecraft Recovery & Remote Sensing, 2017, 38(1): 98-106. (in Chinese)
[4] 李賢, 胡渝. 星間光通信中振動抑制的研究[J]. 宇航學(xué)報, 2002, 23(3): 77-80. LI Xian, HU Yu. Vibration Suppression in Optical Inter-satellite Communications[J]. Journal of Astronautics, 2002, 23(3): 77-80. (in Chinese)
[5] PARK G, LEE D O, HAN J H. Development of Multi-degree-of-freedom Micro Vibration Emulator for Efficient Jitter Test of Spacecraft[J]. Journal of Intelligent Material Systems and Structures, 2014, 25(9): 1069-1081.
[6] TOYOSHIMA M, TAKAYAMA Y, KUNIMORI H, et al. In-orbit Measurements of Spacecraft Micro Vibrations for Satellite Laser Communication Links[J]. Optical Engineering, 2010, 49(8): 083604.
[7] HE S, XU Z, WANG X, et al. Design and Testing of a Parallel Manipulator for Space Micro-vibration Simulation[C]∥Conference Towards Autonomous Robotic Systems, July 19-21, 2017, Guildford, UK. Springer, 2017: 86-100.
[8] PRASHANT A R, MADHESWARAN M, KARTIK V, et al. System Development for Micro Vibration Measurements on Spacecrafts[C]∥International Conference on Control, Instrumentation, Communication and Computational Technologies (ICCICCT), December 16-17, 2016, Kumaracoil, India. New York: IEEE, 2016: 98-103.
[9] 葛任偉, 吳清文, 王運, 等. 基于聯(lián)合變換相關(guān)器的像移測量方法研究[J]. 計算機(jī)仿真, 2010, 27(3): 215-219. GE Renwei, WU Qingwen, WANG Yun, et al. The Study of Image Motion Measurement Based on Joint Transform Correlator[J]. Computer Simulation, 2010, 27(3): 215-219. (in Chinese)
[10] 劉海秋, 徐抒巖, 王棟. 基于空間相機(jī)時間延遲積分傳感器拼接區(qū)圖像的像移測量[J]. 光學(xué)學(xué)報, 2014, 34(2): 108-114. LIU Haiqiu, XU Shuyan, WANG Dong. Space Camera Image Motion Measurement Based on Images from Time Delayed Integration Sensors Overlapped Area[J]. Acta Optica Sinica, 2014, 34(2): 108-114. (in Chinese)
[11] HAYDEN W L, MCCULLOUGH T, RETH A. Wideband Precision Two-axis Beam Steer Tracking Servo Design and Test Results[J]. Proceedings of SPIE, 1993, 1866: 271-279.
[12] JANSCHER K, TCHERNYKH V. Optical Correlator for Image Motion Compensation in the Focal Plane of a Satellite Camera[J]. Space Technology, 2001, 21(4): 127-132.
[13] VALERIJ T, SERGEI D, KLAUS J, et al. Smartsan Hardware Test Results for Smart Optoelectronic Image Correction for Push Broom Cameras[J]. Proceeding of SPIE, 2002, 4814: 264-272.
[14] TCHERNYKH V, BECK M, JANSCHEK K. Optical Correlator Based Optical Flow Processor for Real Time Visual Navigation[M]//Goro Obinatn, Ashish Dutta. Vision Systems: Applications. I-Tech Education and Publishing, 2007.
[15] YI H, ZHAO H, LI Y, et al. Improved Digital Processing Method Used for Image Motion Measurement Based on Hybrid Opto-digtial Joint Transform Correlator[J]. Chinese Optics Letters, 2010, 8(10): 989-992.
[16] NAYAR S K, BEN-EZRA M. Motion-based Motion Deblurring[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2004, 26(6): 689-698.
[17] 樊超, 李英才, 傅洪亮, 等. 光學(xué)相關(guān)法測量空間相機(jī)像移的性能研究[J]. 光學(xué)學(xué)報, 2011, 3(7): 1-5. FAN Chao, LI Yingcai, FU Hongliang, et al. Research on Measurement Method of Image Motion of Space Camera Based on Optical Correlator[J]. Acta Optica Sinica, 2011, 3(7): 1-5. (in Chinese)
[18] MCEWEN A, BANKS M, BAUGH N, et al. The High Resolution Imaging Science Experiment (Hirise) During MRO’S Primary Science Phase (PSP)[J]. Icarus, 2010, 205(1): 2-37.
[19] BELY P Y, LUPIE O L, HERSHEY J L. Line-of-sight Jitter of the Hubble Space Telescope[J]. Proceedings of SPIE, 1993, 1945: 55-61.
[20] MATTSON S, BOYD A, KIRK R L. HiJACK: Correcting Spacecraft Jitter in HiRISE Images of Mars[J]. Health Management Technology, 2009, 33(5): A162.
[21] WANG M, FAN C, PAN J, et al. Image Jitter Detection and Compensation Using a High-frequency Angular Displacement Method for Yaogan-26remote Sensing Satellite[J]. ISPRS Journal of Photogrammetry & Remote Sensing, 2017, 130: 32-43.
[22] TONG X H, XU Y S, YE Z, et al. Attitude Oscillation Detection of the ZY-3 Satellite by Using Multispectral Parallax Images[J]. IEEE Transactions on Geoscience and Remote Sensing, 2015, 53(6): 3522-3534.
[23] GONZALEZ R C, WOODS R E. Digital Image Processing[M]. 2nd ed. Publishing House of Electronics Industry, 2002.
[24] RODRIGUEZ P, WOHLBERG B. Efficient Minimization Method for a Generalized Total Variation Functional[J]. IEEE Transactions on Image Processing: A Publication of the IEEE Signal Processing Society, 2009, 18: 322-332.
[25] 董文德. 基于光纖陀螺顫振探測的圖像復(fù)原技術(shù)研究[D]. 杭州: 浙江大學(xué), 2013. DONG Wende. Research on Image Restoration Based on Vibration Detection Using Fiber Optic Gyroscope[D]. Hangzhou: Zhejiang University, 2013. (in Chinese)
Image Motion Measurement and Image Restoration Based on the Inertial Reference
YUE Ronggang1WANG Shitao1WANG Humei1YIN Ke1JIN Ting1ZANG Jie2YAN Tingfei3LI Huan4TANG Shaofan4
(1 Institute of Remote Sensing Satellite, Beijing 100094, China)(2 Beijing Institute of Spacecraft System Engineering, Beijing 100094, China)(3 Beijing Institute of Spacecraft Environment Engineering, Beijing 100094, China)(4 Beijing Institute of Space Mechanics & Electricity, Beijing 100094, China)
There are lots of high, medium and low frequency micro vibration sources on the satellite, which cause optical axis jitter of the optical load, resulting in the degradation of remote sensing image quality. Aimed to that question, a method of image motion detection and image restoration was proposed based on an inertial reference, and the principle and the key components of the method was introduced. To verify the performance of the method, a test system was developed, including a space camera, an inertial reference unit, a Hartmann sensor, an integrating sphere, a simulated image target, a parallel light pope, a vibration isolation platform with the vibration generator, a six degrees of freedom platform, and the software of image restoration. The background noise of the test environment was measured, and the image motion measurement accuracy test was carried out. Based on that work, the verification tests of image restoration under various working conditions were completed. The test results show that the error of image motion detection is less than 0.12 pixels (Root Mean Square, RMS). Using the image motion data to improve the image quality, MTF (Modulation Transfer Function) of the restored image was improved to 1.38-1.52 times of the original image MTF. The image motion data can be used as a feedback to the fast steering mirror to compensate the satellite jitter in real time, which helps to get high quality images directly.
inertial reference; micro vibration; image motion; image restoration; remote sensing satellite
V447+.1
A
1009-8518(2021)01-0125-10
10.3969/j.issn.1009-8518.2021.01.015
岳榮剛,男,1980年生,2010年獲北京航空航天大學(xué)機(jī)械電子工程專業(yè)博士學(xué)位,高級工程師。主要研究方向為遙感衛(wèi)星、移動機(jī)器人。E-mail:beijing2008-v@126.com。
2020-09-18
國家重點研發(fā)計劃(2016YFB0500502)
岳榮剛, 王世濤, 王虎妹, 等. 基于慣性參考基準(zhǔn)的像移測量與圖像復(fù)原[J]. 航天返回與遙感, 2021, 42(1): 125-134.
YUE Ronggang, WANG Shitao, WANG Humei, et al. Image Motion Measurement and Image Restoration Based on the Inertial Reference[J]. Spacecraft Recovery & Remote Sensing, 2021, 42(1): 125-134. (in Chinese)
(編輯:龐冰)