基于I/Q支路相干積分觀測濾波的GPS接收機信號跟蹤方法

沈 鋒 李偉東李 強

(哈爾濱工程大學(xué)自動化學(xué)院 哈爾濱 150001)

基于I/Q支路相干積分觀測濾波的GPS接收機信號跟蹤方法

沈 鋒 李偉東*李 強

(哈爾濱工程大學(xué)自動化學(xué)院 哈爾濱 150001)

針對傳統(tǒng)GPS接收機在弱信號環(huán)境下跟蹤誤差大，收斂速度慢的缺點，該文提出一種基于I/Q支路相干積分觀測濾波的GPS信號跟蹤方法。將接收機I/Q支路相干積分輸出為觀測量，應(yīng)用無跡卡爾曼濾波(Unscented Kalman Filter, UKF)算法構(gòu)建卡爾曼濾波器，得到基帶數(shù)字信號處理濾波模型，閉合接收機跟蹤環(huán)路。該方法能夠有效減小傳統(tǒng)GPS跟蹤環(huán)路中信號參數(shù)的估計誤差，提高接收機抗干擾能力和弱信號環(huán)境下環(huán)路跟蹤性能。仿真對比結(jié)果表明，不同載噪比環(huán)境下相比傳統(tǒng)GPS信號跟蹤的方法，基于I/Q支路相干積分觀測的信號跟蹤算法能夠提高跟蹤精度，加快跟蹤收斂速度。

全球定位系統(tǒng)；弱信號；跟蹤環(huán)路；相干積分

1 引言

全球定位系統(tǒng)(GPS)定位誤差不隨時間積累，全天候提供高精度的定位和導(dǎo)航服務(wù)已經(jīng)實現(xiàn)全球覆蓋，在軍民領(lǐng)域均具有廣闊的用途和發(fā)展前景[1,2]。近年來，全球用戶對于GPS應(yīng)用提出了越來越高的要求，應(yīng)用范圍也越來越廣泛。

傳統(tǒng)GPS接收機在高載噪比(C/N0)強信號環(huán)境下跟蹤性能良好，但是在弱信號、多路徑以及強電磁干擾環(huán)境下，接收機跟蹤性能受到嚴(yán)重影響，對弱信號處理一直是接收機技術(shù)研究的難點[3,4]。針對復(fù)雜環(huán)境下的GPS信號跟蹤問題，文獻(xiàn)[5]介紹了用載波跟蹤環(huán)和碼跟蹤環(huán)結(jié)構(gòu)分別對接收GPS信號中的碼相位和多普勒頻移進行同步跟蹤，確保本地復(fù)現(xiàn)信號與輸入信號的相關(guān)積分值達(dá)到最大；文獻(xiàn)[6]提出一種位同步方法(K-P法)，能夠在低載噪比條件下檢測電文數(shù)據(jù)位邊緣，但前提條件是沒有載波相位和頻率誤差；文獻(xiàn)[7]利用擴展卡爾曼濾波(Extended Kalman Filter, EKF)消除信號跟蹤時的多路徑效應(yīng)，其設(shè)計模型較為簡單，但不能有效地跟蹤微弱信號；文獻(xiàn)[8]提出一種將EKF濾波算法和貝葉斯估計結(jié)合的方法，實現(xiàn)對弱信號載波和偽碼相位的跟蹤；文獻(xiàn)[9,10]給出了基于自適應(yīng)卡爾曼濾波的信號跟蹤方法，通過歷史信息估計噪聲協(xié)方差矩陣，跟蹤變化中的系統(tǒng)模型，由于動態(tài)協(xié)方差矩陣的非正定性，在微弱信號的環(huán)境下容易失鎖；文獻(xiàn)[11]在分析EKF濾波算法的基礎(chǔ)上，將無跡KF(UKF)濾波算法引入到接收機跟蹤環(huán)路中，完成對載波相位和碼相位的估算，但未在弱信號跟蹤問題上做深入研究；文獻(xiàn)[12]給出了UKF濾波算法在接收機碼環(huán)中的應(yīng)用，并通過采用仿真噪聲干擾下的微弱GPS偽碼中頻信號對該算法進行了性能分析，但文中選取的狀態(tài)量維數(shù)低且觀測量信息少，未能更好地反應(yīng)跟蹤環(huán)路信息。

為解決弱信號環(huán)境下的信號跟蹤，本文在分析傳統(tǒng)GPS接收機載波和碼跟蹤環(huán)路工作的基礎(chǔ)上，基于I/Q支路相干積分觀測提出了利用卡爾曼濾波器實現(xiàn)GPS信號跟蹤方法，該方法利用接收機I/Q支路相干積分輸出為觀測量，克服傳統(tǒng)鑒別器輸出時帶來的誤差，對載波環(huán)和碼環(huán)的特征量進行聯(lián)合估計并應(yīng)用UKF濾波構(gòu)建濾波器模型，通過對非線性函數(shù)概率密度函數(shù)的估計來獲取狀態(tài)估計，得到基帶信號處理濾波模型，能夠有效地解決微弱信號下的信號跟蹤問題。文中基于UKF濾波算法建立Kalman濾波器模型，采用有限數(shù)量的Sigma點集χ對載波頻率誤差、載波相位誤差和碼相位誤差等參數(shù)進行估計，并對不同載噪比環(huán)境下的跟蹤性能進行了仿真對比分析，驗證了所提方法的正確性與有效性。

2 GPS信號跟蹤模型

GPS衛(wèi)星首先利用偽碼(C/A碼)對數(shù)據(jù)碼進行擴頻調(diào)制，再將偽碼與數(shù)據(jù)碼的組合碼通過雙相移位鍵控機制(BPSK)對載波(L1)進行調(diào)制。在不考慮噪聲的情況下，由接收機天線接收到的第i顆衛(wèi)星中頻信號可寫成

信號跟蹤階段，信號通道從捕獲階段獲得的對當(dāng)前衛(wèi)星信號載波頻率和碼相位的粗略估計值出發(fā)，通過跟蹤環(huán)路逐步精細(xì)對這兩個信號參量的估計[1]。根據(jù)鑒別器的不同，載波跟蹤環(huán)路通常有相位鎖定環(huán)路(Phase Lock Loop, PLL)和頻率鎖定環(huán)路(Frequency Lock Loop FLL)兩種形式，碼跟蹤環(huán)路采用延遲鎖定環(huán)路(Delay Lock Loop, DLL)。傳統(tǒng)的跟蹤環(huán)路如圖1所示。

圖1 傳統(tǒng)GPS跟蹤環(huán)路總體結(jié)構(gòu)圖

在歷元k時刻，載波環(huán)鑒別器使用即時碼(P)相關(guān)器的輸出值IP和QP，得到載波頻率差異δf和相位差異δφ, FLL和PLL分別采用叉積鑒頻法和二象限反正切函數(shù)鑒相器[13]，載波頻率誤差和相位誤差分別為

式(2)中，叉積Pcross與點積Pdot分別等于

延遲鎖定環(huán)路(DLL)通過超前碼(E)、滯后碼(L)相關(guān)器的輸出來計算碼相位誤差，采用非相干超前減滯后功率法，碼相位誤差為

式(4)中，IE,QE,IL和QL為超前支路與滯后支路相干積分值。

環(huán)路濾波器通常是一個低通濾波器，其目的在于降低環(huán)路中的噪聲，載波和碼鑒別器的輸出經(jīng)過環(huán)路濾波后，將濾波結(jié)果反饋給數(shù)控振蕩器(Numerical Controlling Oscillator, NCO)進而調(diào)整本地載波頻率、載波相位和偽碼相位，載波環(huán)路盡力使得本地復(fù)現(xiàn)載波信號與衛(wèi)星載波信號的相位或頻率保持一致[2]。碼NCO復(fù)制出3份不同相位的E、P和L復(fù)制C/A碼，與接收到的C/A碼做相關(guān)運算，推算出C/A碼自相關(guān)函數(shù)的主峰位置，從而檢測出碼相位差異。

在弱信號環(huán)境下，傳統(tǒng)GPS接收機接收到的信號容易受到內(nèi)部熱噪聲和動態(tài)應(yīng)力誤差的影響，使GPS接收機跟蹤環(huán)路中信號參量估計誤差加大，且鑒別器在鑒別相位時存在著非線性因素，因此跟蹤環(huán)路容易產(chǎn)生失鎖現(xiàn)象。由于傳統(tǒng)跟蹤方法處理微弱信號方面的先天性不足，本文提出利用接收機內(nèi)部跟蹤環(huán)路I/Q支路相干積分值作為觀測量，應(yīng)用UKF濾波方法構(gòu)建濾波器模型代替?zhèn)鹘y(tǒng)算法中載波跟蹤和碼跟蹤的鑒別器。

3 基于I/Q支路相干積分觀測信號跟蹤算法

3.1 基于相干積分觀測濾波的信號跟蹤模型

圖2為本文所提出的基于I/Q支路相干積分觀測濾波的GPS信號跟蹤方法，其核心思想是采用Kalman濾波器代替?zhèn)鹘y(tǒng)GPS跟蹤環(huán)路中的鑒別器。Kalman濾波是一種利用最小方差準(zhǔn)則的最優(yōu)估計算法，可以對系統(tǒng)狀態(tài)參量進行精確估計，且可以更為有效地對系統(tǒng)進行控制[14,15]。為了能夠準(zhǔn)確地跟蹤輸入信號的頻率和相位，選擇濾波器狀態(tài)量為

圖2 Kalman濾波器跟蹤環(huán)路結(jié)構(gòu)圖

式(5)中，A為歸一化信號幅值，δφ為載波相位差，δw 為載波頻率差，δa為載波頻率變化率，δτ為碼相位差。卡爾曼濾波器的狀態(tài)方程如式(6)所示。

式(6)中，λL為L1載波波長，λCA為C/A碼波長， Qk為過程噪聲，記為：零均值白噪聲序列。

將I/Q支路相干積分輸出值作為濾波器的觀測量，建立與狀態(tài)量相關(guān)的非線性方程為

式(7)中，δ為本地C/A碼超前滯后間隔，R(εi)為C/A碼自相關(guān)函數(shù)。Vk為量測噪聲，噪聲方差陣為

式(8)中，σN是相關(guān)信號I, Q處理后的噪聲強度[11]。跟蹤過程中，觀測量與狀態(tài)量存在非線性關(guān)系，因此采用UKF濾波算法進行環(huán)路濾波處理，得到估計出的載波相位差δφ、載波頻率差δw和碼相位差δτ。

3.2 基于UKF的卡爾曼濾波器

UKF濾波算法是采用UT變換和傳統(tǒng)卡爾曼濾波架構(gòu)的有效結(jié)合，通過對非線性函數(shù)概率密度函數(shù)的估計來獲取狀態(tài)估計，UT變換是UKF算法的核心和基礎(chǔ)。在確保采樣均值和協(xié)方差為x和px的前提下，選擇2n+1個Sigma點集，將非線性變換應(yīng)用于采樣的每個Sigma點，得到非線性轉(zhuǎn)換后的點集和py是變換后Sigma點集的統(tǒng)計量。

卡爾曼濾波器離散化后的非線性系統(tǒng)模型為

狀態(tài)變量XK初始分布均值為，方差矩陣為P。UT變換產(chǎn)生Sigma點χ和對稱采樣的加權(quán)序列w，表示為

時間更新

測量更新

UKF濾波算法是對非線性函數(shù)的概率密度分布進行近似，采用一系列確定的Sigma點集來逼近狀態(tài)的后驗概率密度，避免了EKF濾波算法求導(dǎo)計算雅克比矩陣帶來的大量計算，克服了EKF濾波算法線性化忽略高階項而帶來的跟蹤精度問題。UKF濾波算法的計算量大小取決于UT變換選擇Sigma點集平方根分解運算，當(dāng)狀態(tài)變量維數(shù)較高時，2n+1個Sigma點增多，UKF濾波算法計算量也會相應(yīng)地增加，優(yōu)化平方根分解方式可以減小計算量，同時隨著計算機技術(shù)的不斷發(fā)展，計算速度越來越快，濾波算法復(fù)雜度帶來的滯后問題將會較好地解決。

4 仿真與分析

為了驗證基于I/Q支路相干積分值為觀測信息的GPS信號跟蹤性能與正確性，首先利用傳統(tǒng)二階跟蹤環(huán)路對GPS中頻模擬信號進行跟蹤仿真，通過鑒別器輸出載波相位差和碼相位差。然后利用本文提出的Kalman濾波器模型代替?zhèn)鹘y(tǒng)的跟蹤環(huán)路，利用I/Q支路相干積分值為觀測量建立系統(tǒng)模型，并輸出跟蹤信號的誤差信息。仿真過程中，數(shù)字中頻信號頻率和采樣頻率分別設(shè)置為1.5 MHz和5 MHz，碼速率1.023 Mcps，載波和碼環(huán)路帶寬分別設(shè)為8 Hz 和4 Hz，預(yù)檢測積分時間設(shè)為1 ms。

圖3為載噪比42 dBHz下的載波相位和碼相位跟蹤誤差對比圖，由圖可以看出傳統(tǒng)二階跟蹤環(huán)路在強信號環(huán)境下能夠較好地跟蹤GPS衛(wèi)星信號，載波跟蹤環(huán)更新360次左右達(dá)到穩(wěn)態(tài)，碼環(huán)更新580次左右達(dá)到穩(wěn)態(tài)，穩(wěn)態(tài)時載波和碼環(huán)跟蹤誤差分別控制在10°和0.1個碼片內(nèi)。基于I/Q支路相干積分值為觀測量的Kalman濾波模型下的信號跟蹤性能更為優(yōu)越，載波和碼跟蹤初始階段均有一個較小的超調(diào)量，在環(huán)路更新240次左右就能達(dá)到穩(wěn)態(tài)值，快速性得到充分體現(xiàn)，且跟蹤精度更加準(zhǔn)確，能夠?qū)崟r、快速地跟蹤接收到的衛(wèi)星信號。

圖4是在32 dBHz下所做的載波相位和碼相位跟蹤誤差對比圖，由圖可以看出，相對42 dBHz傳統(tǒng)二階跟蹤環(huán)路信號跟蹤性能有所下降，穩(wěn)態(tài)時載波相位和碼相位跟蹤誤差分別擴大到16°和0.2個碼片范圍內(nèi)，基于I/Q支路相干積分值為觀測量的Kalman濾波模型超調(diào)量有所增加，但是其值變化不大。碼跟蹤環(huán)路更新310次左右達(dá)到穩(wěn)態(tài)值，較傳統(tǒng)二階碼跟蹤環(huán)路630次更快，快速性仍然能夠得到體現(xiàn)，且跟蹤精度較傳統(tǒng)二階跟蹤環(huán)路更高，誤差波動范圍更小。

表1和表2給出了傳統(tǒng)和Kalman濾波模式下跟蹤數(shù)據(jù)對比結(jié)果，由表可以看出，載噪比42 dBHz 和32 dBHz環(huán)境下，卡爾曼濾波跟蹤模式載波和碼相位誤差的均方差均比傳統(tǒng)跟蹤模式小，且波動范圍較后者要小。在載噪比32 dBHz環(huán)境下，傳統(tǒng)跟蹤模式載波相位均方差9.90°, Kalman濾波模式均方差2.44°，后者要比前者要小很多。傳統(tǒng)跟蹤模式載波相位誤差最大值達(dá)到近30°，而Kalman濾波模式最大值僅為6.42°，遠(yuǎn)小于前者。由此也可以得出，Kalman濾波模式的跟蹤性能更優(yōu)。

圖5為傳統(tǒng)二階環(huán)路和Kalman濾波器模式在低載噪比24 dBHz下的I/Q支路輸出信息對比圖，由圖5可以看出傳統(tǒng)二階環(huán)路在低載噪比環(huán)境下，輸出的信息是一些發(fā)散的數(shù)據(jù)點，起始階段跟蹤過程就不理想，隨著環(huán)路鑒別誤差逐漸增加，最終導(dǎo)致環(huán)路跟蹤失鎖。而Kalman濾波器模式I支路信息正確解調(diào)，Q支路輸出的只是一些噪聲。Kalman濾波器跟蹤模式下系統(tǒng)利用接收機最原始的I/Q支路相干積分值作為觀測量，有效減小了使用鑒別器觀測時所帶來的誤差量，在低載噪比環(huán)境下仍能正確解調(diào)出跟蹤的衛(wèi)星信號導(dǎo)航數(shù)據(jù)。

圖4 傳統(tǒng)和Kalman濾波模式載噪 比32 dBHz下跟蹤誤差對比

圖5 低載噪比下二者I/Q 支路信息對比圖

表1 傳統(tǒng)和Kalman濾波模式42 dBHz下跟蹤數(shù)據(jù)對比

表2 傳統(tǒng)和Kalman濾波模式32 dBHz下跟蹤數(shù)據(jù)對比

為了進一步分析濾波方法的時效性，對傳統(tǒng)二階環(huán)路和Kalman濾波器模式的耗時情況進行仿真與分析。

圖6為傳統(tǒng)二階環(huán)路和Kalman濾波器模式在環(huán)路更新過程中的時效性對比圖，由圖可以看出，跟蹤環(huán)路更新500次，兩者每次環(huán)路更新耗時數(shù)量級都在毫秒級，時間消耗上差別很小，在濾波器中引入UKF濾波算法，不但在跟蹤性能上得到改善，而且算法復(fù)雜度引起的時間滯后問題也隨著計算機技術(shù)的發(fā)展而得到改善。

5 結(jié)束語

圖6 傳統(tǒng)和Kalman濾波模式下時效性對比

在對傳統(tǒng)跟蹤環(huán)路通過鑒別器輸出誤差信息分析的基礎(chǔ)上，提出基于I/Q支路相干積分觀測的Kalman濾波模型及應(yīng)用UKF濾波算法進行仿真。通過Kalman濾波處理與傳統(tǒng)跟蹤環(huán)路在不同載噪比環(huán)境下的性能仿真對比，可以得出：應(yīng)用UKF濾波算法處理非線性信號快速收斂、對動態(tài)參量能夠精確估計，使其具有更高的跟蹤精度，更快的收斂速度，在低載噪比的環(huán)境下能夠較好地跟蹤衛(wèi)星信號，正確地解調(diào)出衛(wèi)星信號中的數(shù)據(jù)碼信息，比傳統(tǒng)跟蹤環(huán)路在弱信號環(huán)境下具有更為優(yōu)越的跟蹤性能。

[1] 田甜, 安建平, 王愛華. 高動態(tài)環(huán)境下無數(shù)據(jù)輔助的擴展Kalman濾波載波跟蹤環(huán)[J]. 電子與信息學(xué)報, 2013, 35(1): 63-67. Tian Tian, An Jian-ping, and Wang Ai-hua. Non-data-aided extended Kalman filter based carrier tracking loop in high dynamic environment[J]. Journal of Electronics & Information Technology, 2013, 35(1): 63-67.

[2] 謝鋼. GPS原理和接收機設(shè)計[M]. 北京: 電子工業(yè)出版社,2009: 1-2.

[3] Lashley M, Bevly D M, and Hung J Y. Performance analysis of vector tracking algorithms for weak GPS signals in high dynamics[J]. IEEE Journal of Selected Topics in Signal Processing, 2009, 3(4): 661-673.

[4] 趙琳, 丁繼成, 孫明, 等. 極弱信號環(huán)境下GPS位同步和載波跟蹤技術(shù)[J]. 航空學(xué)報, 2010, 31(6): 1205-1212. Zhao Lin, Ding Ji-cheng, Sun Ming, et al.. Bit synchronization and carrier tracking for very weak GPS signals[J]. Acta Aeronautica et Astronautica Sinica, 2010, 31(6): 1205-1212.

[5] Wu Ling-juan and Cui Ying-ying. Signal acquisition and tracking for software GPS receivers[C]. Proceedings of the Asia Pacific Conference on Postgraduate Research in Microelectronics & Electronics, Shanghai, China, 2009: 384-387.

[6] Kokkonen M and Pietila S. A new bit synchronization method for a GPS receiver[C]. Proceedings of the IEEE Position Location and Navigation Symposium, Palm Springs, CA, 2002: 85-90.

[7] Jee G, Kim H S, and Lee Y J. A GPS C/A code tracking loop based on extended Kalman filter with multipath mitigation [C]. Proceedings of the Institute of Navigation GPS, Portland, USA, 2002: 446-451.

[8] Psiaki M L and Jung H. Extended Kalman filter methods for tracking weak GPS signals[C]. Proceedings of the Institute of Navigation GPS 2002, Portland, 2002: 2539-2553.

[9] E Shi. An improved real-time adaptive Kalman filter for low-cost integrated GPS/INS navigation[C]. Proceedings of the 2012 International Conference on Measurement, Information and Control (MIC), Harbin, China, 2012: 1093-1098.

[10] Janos G and Andrei C. A simplified adaptive Kalman filter algorithm for carrier recovery of M-QAM signals[C]. Proceedings of the 17th International Conference on Methods and Models in Automation and Robotics (MMAR), Szczecin, Poland, 2012: 303-307.

[11] 張敏虎, 任章, 華春紅. UKF在深組合GPS/INS導(dǎo)航系統(tǒng)中的應(yīng)用[J]. 中國慣性技術(shù)學(xué)報, 2009, 17(6): 697-700. Zhang Min-hu, Ren Zhang, and Hua Chun-hong. Application of UKF in deeply coupled GPS/INS navigation system[J]. Journal of Chinese Inertial Technology, 2009, 17(6): 697-700.

[12] 曾慶喜, 唐琳琳, 王慶, 等. 基于UKF的GPS L1 C/A碼微弱信號跟蹤算法[J]. 艦船電子工程, 2014, 35(1): 58-62. Zeng Qing-xi, Tang Lin-lin, Wang Qing, et al.. Weak signal tracking algorithm base on UKF for L1 C/A code[J]. Ship Electronic Engineering, 2014, 35(1): 58-62.

[13] Li Hui and Yang Jing-shu. Analysis and simulation of vector tracking algorithms for weak GPS signals[C]. Proceedings of the International Asia Conference on Informatics in Control, Automation and Robotics, Hefei, China, 2010: 215-218.

[14] Luo Yong, Ravindra Babu, Wu Wen-qi, et al.. Double-filter model with modified Kalman filter for baseband signal pre-processing with application to ultra-tight GPS/INS integration[J]. GPS Solutions, 2012, 16(4): 463-476.

[15] 謝非, 劉建業(yè), 李榮冰, 等. 基于環(huán)路相關(guān)積分觀測的SINS/ GPS 深組合導(dǎo)航算法[J]. 中國慣性技術(shù)學(xué)報, 2013, 21(4): 472-477. Xie Fei, Liu Jian-ye, Li Rong-bing, et al.. Deeply SINS/GPS integrated navigation algorithm based on tracking loop correlation measurements[J]. Journal of Chinese Inertial Technology, 2013, 21(4): 472-477.

[16] 王璐, 李光春, 喬相偉, 等. 基于極大似然準(zhǔn)則和最大期望算法的自適應(yīng)UKF算法[J]. 自動化學(xué)報, 2012, 38(7): 1200-1210. Wang Lu, Li Guang-chun, Qiao Xiang-wei, et al.. An adaptive UKF algorithm based on maximum likelihood principle and expectation maximization algorithm[J]. Acta Automatica Sinica, 2012, 38(7): 1200-1210.

沈 鋒： 男，1981年生，副教授，研究方向為信號檢測、衛(wèi)星導(dǎo)航技術(shù)研究.

李偉東： 男，1989年生，碩士生，研究方向為衛(wèi)星導(dǎo)航技術(shù)研究、SINS/GPS深組合導(dǎo)航.

李 強： 男，1987年生，博士，研究方向為衛(wèi)星導(dǎo)航、GPS/INS組合導(dǎo)航.

GPS Receiver Signal Tracking Method Based on I/Q Branch Coherent Integration Measurements Filter

Shen Feng Li Wei-dong Li Qiang
(College of Automation, Harbin Engineering University, Harbin 150001, China)

Due to the disadvantages of a traditional GPS receiver in the environment of weak signal are that slow convergence rate and serious tracking error, a GPS signal tracking algorithm is proposed, in which the coherent integration of I/Q branch is considered as measurement of Unscented Kalman Filter (UKF). The filter model of the baseband signal processing is constructed and tracking loop is closed by UKF. Then the signal parameters of GPS tracking loop can be accurately estimated. Moreover, the anti-jamming capability and tracking ability of receiver are improved in the presence of weak signals. Simulation results demonstrate that the tracking ability and convergence rate of proposed algorithm outperform traditional method in the case of different Carrier to Noise ratio (C/N).

GPS; Weak signal; Tracking loop; Coherent integration

TN966

A

1009-5896(2015)01-0037-06

10.11999/JEIT140314

2014-03-10收到，2014-06-23改回

國家自然科學(xué)基金(61102107, 61374208)和中央高校基本科研業(yè)務(wù)費專項基金(HEUCFX41310)資助課題

*通信作者：李偉東 lwd52531258@163.com