曹輝 宋有建 于佳禾 師浩森 胡明列 王清月
(天津大學精密儀器與光電子工程學院,光電信息技術教育部重點實驗室,天津 300072)
長度是國際單位制中的七個基本單位之一,同時也是科學研究與工業生產中最常用的物理量之一.近年來,衛星編隊飛行、大型工件形貌測量等領域對大量程、高精度的絕對距離測量提出了很高的要求.
傳統的光學方法測量絕對距離主要分為干涉法[1]和飛行時間法[2].干涉法使用單色光測距,其分辨率受限于單色光的波長穩定性,分辨率可以達到納米量級.但是光干涉法的絕對距離測量范圍僅為測量光的半波長,對于更長的待測距離只能通過增量式測量獲取距離值.飛行時間法通過測定激光脈沖的發射時刻和其回波返回時刻之間的時間差推算距離信息.該方法具有測量范圍大、發射功率高以及抗干擾能力強等優點.但是,該方法的時間分辨率受限于光電探測器的響應帶寬,僅能達到皮秒量級,即毫米級的距離分辨率.
基于飛秒激光的飛行時間測距,有效地克服了干涉法以及傳統飛行時間法的上述缺點,有效地提升了絕對距離測量精度,在航空航天、遙感和表面形貌測量等領域有著很好的應用前景.2000年,Minoshima和Matsumoto[3]首次提出基于飛秒激光的絕對距離測量,實現了240 m長的高精度飛行時間測量.2009年,Coddington等[4]使用兩臺相干的光纖激光頻率梳(雙光梳)在約1.5 m的測試距離下實現了5 nm的測試精度,并同時實現了對運動物體的追蹤鎖定.2010年,Lee等[5]在飛秒激光測距系統中使用平衡光學互相關(BOC)技術,通過等效采樣原理成功地避免了采用光電探測器對光脈沖直接進行光電接收引入的時間分辨率損失,將1 s內采集到的測距信號進行平均,在0.7 km的距離下達到了7 nm的絕對測量精度.2012年,秦鵬等[6]基于平衡互相關技術,采用摻Yb鎖模光纖激光器作為測距光源,在52 m的待測距離上,獲得了12 nm的絕對距離測量精度.2013年,邢書劍等[7]采用基于改進的邁克耳孫干涉原理的任意絕對距離測量系統,使用飛秒光學頻率梳在0.6 m時獲得了0.5μm的測距精度.2016年,張曉聲等[8]使用飛秒激光模間拍頻法實現了優于3μm的絕對測量精度.
在已有的飛秒激光測距方法中,雙光梳方法由于測量精度高、更新速率快,并可以同時跟蹤多個目標,引起了廣泛的關注.尤其是通過使用兩臺自由運轉鎖模激光器(簡化的雙光梳),實驗裝置可以極大地簡化,同時能夠保證微米級的測量精度,對工業測量具有極大的吸引力.但是,雙梳測距極容易受到激光器中量子噪聲的影響,從而增加測量結果的不確定度,使測量精度下降[9].盡管可以使用滑動平均等方法減輕上述問題,可是必須以降低系統的更新速率為代價.另外,滑動平均會在距離信號本身帶有高頻成分或者噪聲源更加復雜的情況下失效.
奇異譜分析(SSA)是一種從含噪信號中提取量化信息的方法.它基于Takens的嵌入理論,在1986年由Broomhead和King共同建立[10].SSA方法是一種無參數頻譜估計技術,它用于分析無先驗知識的動態系統采集到的時間序列,并從中提取若干有實際意義的成分[11].同時,SSA是一種無模型方法,可以應用于線性和非線性系統中.SSA廣泛應用于地理學[12]、經濟學[13]和社會學[14]等數據較為復雜的學科中.此外有報道SSA方法用于季節性GPS信號的分析[15],顯示了SSA對于提取GPS中的調制振蕩信號是一種完備的方法.在高光譜遙感中,SSA也表現出有效的信號提取功能[16?19].
本文應用SSA方法分析由高更新速率雙光梳測距系統獲得的飛行時間絕對測距數據,證明了在無先驗隨機噪聲模型的前提下,SSA方法可以用來提取含噪一維階梯結構信號.進一步在實驗中利用200 Hz更新速率的簡化的雙光梳測距系統,對移動的合作靶的位置進行了測定,通過SSA方法對獲得的信號進行分析,經由后處理所提取出的信號標準差達到0.15μm.

圖1 自由運轉雙梳激光測距系統,插圖為經“時間展寬”的脈沖回波信號Fig.1. Dual-comb free-running laser ranging system working at 1550 nm. The inset illustrates the equivalently-sampled Target/Reference-re fl ected pulses in the “stretched” time frame.HWP(half wave plate);LO(local oscillator);LUT(laser under test);PBS(polarization beam splitter);PD(photo detector);PPKTP(periodically poled KTiOPO4);QWP(quarter wave plate).
雙光梳測距系統的測距原理是基于光學等效采樣方法將超快的周期性脈沖信號在時間上“展寬”,從而可以應用現有的電子學手段提取并解算脈沖飛行時間,進而計算得出飛秒激光源和目標之間的絕對距離.系統結構如圖1所示,主要由飛秒激光源、光學互相關系統和數據采集處理系統組成.實驗中采用兩臺相互獨立的自由運轉被動鎖模飛秒激光器,分別作為本地振蕩器(local oscillator,LO)和信號激光源(laser under test,LUT).兩臺激光器均為全保偏光纖激光器,使用碳納米管(carbon nanotube,CNT)作為可飽和吸收材料,工作波長約為1550 nm,輸出平均功率分別為10 mW和13 mW,重復頻率約為74 MHz,且兩臺激光器之間有約為2 kHz的重復頻率差.信號光經過準直器和半波片(half wave plate,HWP)后,經由偏振分光棱鏡(polarization beam splitter,PBS)分束,分得的兩束光分別經過參考鏡和目標鏡反射,經反射的參考脈沖序列和目標脈沖序列分別與LO脈沖序列合束,經過基于PPKTP(periodically poled KTiOPO4)晶體的光學互相關系統采樣,并使用光電探測器(photo detector,PD)接收采樣信號.由于本地振蕩器和信號激光源具有一定的重復頻率差?fr,該值遠小于激光器本身的重復頻率fr,因此產生等效采樣效果,即飛行時間被展寬,其展寬倍數為N=fr/?fr.在實際應用中,并不需要測定放大率N,只需要測得到達探測器獲得的互相關信號脈沖峰值的時間即可由下式計算絕對距離D:

其中c表示真空光速,ng表示空氣群折射率,tref1和ttar分別表示一個周期內的參考脈沖和目標脈沖到達探測器的時刻,tref2表示下一個周期的參考脈沖到達探測器的時刻.激光實驗裝置的詳細設定參照文獻[9].
在沒有運動模型和隨機過程模型的前提下,為了實現信號與噪聲的分離,本文用SSA方法對飛秒激光距離信號進行處理.SSA方法將觀測到的時間序列信號分解為一系列獨立分量的總和,并賦予每個分量相應的物理意義(例如,慢變信號、周期振蕩信號和隨機噪聲).經典SSA方法分為以下四個步驟:數據內嵌、奇異值分解(singular value decomposition,SVD)、對角平均和分組重建.
首先,通過等間隔采樣獲得包含N個采樣點的初始數據集YN,

其中yi表示第i個采樣點.單個采樣點包含的數據個數稱為該觀測數據的維度(例如,對于同時采集目標位置與速度的系統,每個采樣點包含位置和速度兩個數據,即觀測數據維度為2).數據采集時需要保證采樣點數N遠遠大于觀測數據的維度.選定時間嵌入窗口,其長度為整數L,滿足2≤L≤N/2,對初始數據集YN進行延遲坐標映射,即選取L個連續采樣點,以獲得時滯向量xi,

因此可以由初始數據集YN獲得N?L+1個時滯向量.這些時滯向量又可以組成軌跡矩陣X,

在第二步奇異值分解中,對軌跡矩陣X進行奇異值分解,即得到

其中矩陣U和V均為列歸一化矩陣,矩陣Σ為對角矩陣.U的列向量組{ui}和V的列向量組{vi}分別為矩陣XXT和XTX的正交歸一化本征矢.列向量組{ui}和{vi}又被分別稱為矩陣X的經驗正交函數(EOF)和主成分(PC).每一個經驗正交函數和主成分對{ui,vi}都對應一個矩陣X的奇異值σi,其大小可衡量相應經驗正交函數和主成分在原始信號中的占比,該值大小等于矩陣XXT中與相應經驗正交函數ui所對應的本征值的平方根.
考慮到經過SVD分解后所得的矩陣數量d=rankX,軌跡矩陣X還可以表示為

其中Xi表示通過SVD分解獲得的第i個矩陣.對于每個矩陣Xi,本征矢和本征值對(ui,σi,vi)按照奇異值降序進行排列,即σ1≥σ2≥···≥σd≥0.
第三步,對角平均.由于上一步獲得了分解后的若干軌跡矩陣分量,這里需要將其還原為原始信號的數據結構.因此還原后的信號為對應矩陣的副對角線元素的平均值.設還原后的信號重建成分(RC)Y(i)的各元素和相應矩陣中元素的對應關系如下:

最后,將得到的還原后的RC進行分組,進而可以獲得重建信號和噪聲的分離.通常按照趨勢、振蕩和噪聲進行分類.實際應用中根據經驗,一般前幾個RC歸為趨勢類,即重建信號本身.歸類后的RC需要按照類別進行累加,以獲得相應分類的重建信號.由于本文所應用場景不包含振蕩信息,因此信號本身僅僅包含趨勢和噪聲成分.
首先使用仿真數據驗證SSA對任意有色噪聲的提取能力.仿真了雙光梳測距系統對一維物體的掃描以得到其表面形貌的應用場景.待測一維物體為1μm高的臺階,如圖2青綠線所示.仿真測量以200 Hz更新速率持續2 s,獲得了400個采樣點,測量結果如圖2黑線所示.這里測量信號采用仿真紫色噪聲進行劣化,噪聲標準差為0.3μm.

圖2 原始一維仿真測距數據(青綠線);被0.3μm標準差的紫色噪聲劣化后的仿真數據(黑線).Fig.2.The original 1-D signal with a 1μm step(cyan curve)and the same signal embedded in computergenerated violet noise with a standard deviation 0.3μm(black curve).

圖3 SSA方法應用于信號提取及噪聲分離 (a)一維含噪聲形貌信號及重建信號;(b)提取的紫色噪聲時域波形;(c)噪聲功率譜密度Fig.3.Signal extraction simulation by using SSA method:(a)1-D pro fi le embedded in violet noise and its reconstruction;(b)time domain value of the extracted violet noise;(c)PSD of the extracted noise.
通過上述對SSA方法的描述可知,在SSA中惟一的參數為嵌入窗口長度L.對于上述仿真數據,我們選取窗口長度L=50,同時提取信號和噪聲.重建信號通過選取所需的分解模式并累加而獲得,如圖3(a)青綠線所示.同時,由于仿真信號不含有振蕩,因此可以將剩下的分解模式進一步累加,以對噪聲進行提取.獲得的噪聲如圖3(b)所示,其標準差為0.29μm,噪聲功率譜如圖3(c)所示,由于其功率譜密度斜率為f2,為紫色噪聲,表明SSA方法可以很好地提取噪聲的類型.此外,如果將提取后的信號和噪聲相加,即可獲得原始信號,這也說明了SSA方法對信號的分解是可逆的,不會造成原始信號的信息丟失.
使用2.1節中描述的雙光梳測距系統測定約0.5 m處的目標靶的距離.目標靶裝載于線性平移臺,以約2μm步長進行移動.圖4(a)中黑線所示為采集到的原始信號,測量更新速率為200 Hz,測量標準差為1.9968μm.
針對上述采集到的距離信號,選取L=40的嵌入窗口長度.發現40個分解成分中的第一個可以作為重建信號,如圖4(a)青綠線所示,其測量標準差為0.1522μm.通過對比原始數據的標準差,可以看出SSA方法對雙光梳測距系統精度具有顯著的提升作用,測距精度提升了13倍左右.
剩余的39個分解成分可以歸類為隨機噪聲,累加后重建為提取噪聲,時域波形如圖4(b)所示,其標準差為1.88μm,功率譜如圖4(c)所示,由于其功率譜密度為水平曲線,顯示出測距系統精度主要受到白噪聲的影響.
進而,還可以分析SSA方法在更遠的測距量程下的適用性.考慮遠距離雙梳測距,空氣折射率的估計誤差和飛秒光源的重復頻率的測量誤差分別會以10?7×D和10?9×D引入測量不確定度[9].因此,在大氣層內的遠距離測距不確定度主要由空氣折射率估計誤差主導;而在真空中,比如星間高精度基線測量,測量不確定度則完全由飛秒光源的重復頻率的測量誤差引入.以上兩種誤差因素均具有隨機噪聲的屬性,可以通過SSA方法分離誤差與信號,達到減小測量不確定度的目的.

圖4 雙光梳測距系統的直接測距結果與SSA分析 (a)實驗測得的激光器和目標之間距離的原始數據及SSA提取結果;(b)提取的噪聲時域波形;(c)噪聲功率譜密度Fig.4.Signal obtained from a dual-comb laser ranging system and data processing results by using SSA method:(a)The measured distance of a moving target and its reconstructed signal;(b)time domain value of the extracted noise;(c)PSD of the extracted noise.
本文應用SSA這種無模型方法處理含噪飛行時間法距離信號.通過SSA,可以在沒有任何運動模型和隨機過程假設的前提下同時提取信號和噪聲.這表明SSA對于重建測距信號和提取系統噪聲是一種有效的分析工具.在仿真實驗中,高精度地提取了含有紫色(功率譜密度斜率為f2)噪聲的一維階梯信號.實驗中,采用200 Hz更新速率的簡化的雙光梳測距系統追蹤移動合作靶的位置.通過SSA對信號的后處理,將合作靶的距離信息有效地從量子噪聲中提取出來,并獲得了0.1522μm的測距精度.值得注意的是,SSA方法本身對信號維度沒有任何要求,因此很容易擴展到高維情況.這意味著SSA輔助下的飛秒激光測距系統可以用于物體表面形貌測量,尤其對含有非連續臺階結構的大規模集成電路印刷版的缺陷檢測等具有潛在的應用價值.
[1]BobroffN 1993Meas.Sci.Technol.4 907
[2]Smullin L D,Fiocco G 1962Nature194 1267
[3]Minoshima K,Matsumoto H 2000Appl.Optics39 5512
[4]Coddington I,Swann W C,Nenadovic L,Newbury N R 2009Nat.Photonics3 351
[5]Lee J,Kim Y J,Lee K,Lee S,Kim S W 2010Nat.Photonics4 716
[6]Qin P,Chen W,Song Y J,Hu M L,Chai L,Wang C Y 2012Acta Phys.Sin.61 240601(in Chinese)[秦鵬,陳偉,宋有建,胡明列,柴路,王清月2012物理學報61 240601]
[7]Xing S J,Zhang F M,Cao S Y,Wang G W,Qu X H 2013Acta Phys.Sin.62 170603(in Chinese)[邢書劍,張福民,曹士英,王高文,曲興華2013物理學報62 170603]
[8]Zhang X S,Yi W M,Hu M H,Yang Z H,Wu G H 2016Acta Phys.Sin.65 080602(in Chinese)[張曉聲,易旺民,胡明皓,楊再華,吳冠豪2016物理學報65 080602]
[9]Shi H,Song Y,Liang F,Xu L,Hu M,Wang C 2015Opt.Express23 14057
[10]Broomhead D S,King G P 1986Physica D20 217
[11]Vautard R,Ghil M 1989Physica D35 395
[12]Vautard R,Yiou P,Ghil M 1992Physica D58 95
[13]Hassani H,Zhigljavsky A 2009J.Syst.Sci.Complex.22 372
[14]Hassani H,Webster A,Silva E S,Heravi S 2015Tourism Manage.46 322
[15]Chen Q,van Dam T,Sneeuw N,Collilieux X,Weigelt M,Rebischung P 2013J.Geodyn.72 25
[16]Zabalza J,Ren J,Wang Z,Marshall S,Wang J 2014IEEE Geosci.Remote S.11 1886
[17]Zabalza J,Ren J,Wang Z,Zhao H,Wang J,Marshall S 2015IEEE J.Sel.Top.Appl.8 2845
[18]Zabalza J,Ren J,Zheng J,Han J,Zhao H,Li S,Marshall S 2015IEEE T.Geosci.Remote53 4418
[19]Rafert J B,Zabalza J,Marshall S,Ren J 2016Appl.Spectrosc.70 1582