張曉林++趙志丹++王曉玲

摘 要：光鑷系統(tǒng)中光阱剛度系數(shù)的快速、準確、高精度測量是光鑷進行粒子定量操控的關鍵因素。本文給出了一種利用數(shù)字圖像相關匹配的方法實現(xiàn)亞像素的位移測量。為提高測量精度，此方法在捕獲過程中采集多幅圖像并以相鄰前一幅圖像作為模板去匹配下一幅圖像，以此類推，最后將每個增量位移場疊加作為最終的位移，同時在亞像素計算時采用較高精度的梯度法；為了提高計算速度，在整像素點匹配時，采用了動態(tài)閾值序貫相似法。最后，通過模擬仿真實驗，對此算法的有效性進行了驗證，結(jié)果表明本算法比普通亞像素算法快1～2倍。

關鍵詞：亞像素位移測量 數(shù)字圖像相關法 序貫相似法 動態(tài)閾值 梯度法

中圖分類號：TP2 文獻標識碼：A 文章編號：1672-3791（2014）01（b）-0001-04

光鑷（Optical tweezers）又稱為單光束梯度力光阱，是一種利用高度聚焦的激光束形成的三維梯度勢阱來捕獲、操控微小粒子的技術[1]。光鑷自1986年由Arthur Ashkin[2]發(fā)明以來，以其非接觸、地損傷等優(yōu)點，已被廣泛應用于物理學中的激光冷卻、膠體化學、生物醫(yī)學尤其是分子生物學等領域[3]，成為一項重要的研究工具。

光鑷的一個重要功能為微小力的測量，對光阱的剛度進行標定是光鑷測力的重要環(huán)節(jié)。標定光阱剛度有許多方法，常用的有流體力學法、熱運動分析法、功率譜法和外加周期驅(qū)動力法等[4～6]，文獻[7]詳細分析比較這四種方法的優(yōu)缺點，其中流體力學法和熱運動分析法均需要CCD跟蹤微粒運行軌跡，在微粒運行一段期間內(nèi)拍攝大量的圖像，然后再對些圖像進行后期處理。由此可知在硬件條件一定的情況下圖像亞像素分析對測量精度有著重大的意義。

亞像素位移測量的算法主要有如下幾種：亞像素灰度插值法[8]、曲面擬合法[9]、相關系數(shù)插值法，牛頓-拉普森[10]（Newton-Rapshon，簡稱N-R）、基于梯度的方法；頻率相關法，后驗概率算法，神經(jīng)網(wǎng)絡方法和基于迭代的最小二乘法[11]。這些算法測量精度所稱精度能到0.005～0.1pixel。常用的算法就三種，下面就簡單的總結(jié)這三種算法的優(yōu)缺性：亞像素灰度插值法，計算量大，精度較低，一般較少直接使用；牛頓-拉普森是基于最優(yōu)化的思想，建立合理的位移和變形模型后然后進多次迭代后然后求出其中的參數(shù)，由于迭代過程中要用到灰度插值，以及灰度的梯度插值，因此目前這個算法是精度最高的，但耗時也最長。相關系數(shù)曲面擬合法，有著很強的抗噪性能，但計算精度比相關系數(shù)插值法略低。梯度法（又稱微區(qū)統(tǒng)計特性梯度法），其基本思想是微小物體的近似剛性位移后微小變形前后點對應點的灰度值保持不變。潘兵等人，指出梯度法與曲面擬合法具有相同的效率，且精度優(yōu)于相關系數(shù)曲面擬合[12]。

本文用圖像相關法[12]并利用位移場的連續(xù)性，設計了一套分步計算整像素位移、壓像素位移，最后將兩者疊加作為最終的位移，如圖1。在具體的計算過程中將集標準化協(xié)方差相關法，動態(tài)閾值序貫相似法，梯度法眾多優(yōu)點集于一體。最后通過仿真實驗驗證算法的有效性和實驗。

1 本文算法的計算流程

1.1 整像素位移測量

為了提高計算的效率和速度這里采用了動態(tài)閾值序貫相似法（SSDA）[13]。下面給出動態(tài)閾值SSDA算法大致流程：

（1）定義絕對誤差：

其中，，。

（2）在相鄰圖像中前一幅中取選定一塊大小和位置合適的圖像作為模板中心坐標為；

（3）確定在后面一幅圖像中的搜索范圍（即子圖的遍歷范圍）；

（4）在后面一幅圖像中，計算模板圖像與初始位置子圖中所有像素點的的累加值，并將其作為閾值的初始值；

（5）計算模板和一個位置子圖中對應點的并累加記作；

（6）在計算并累加過程中比較與的大小，若在計算完每一行或一列后就立刻比較與，若，則停止計算，并將圖像子移動到下一個位置，重復（4）計算，加快匹配速度；

（7）若再遍歷模板圖像與該位置子圖的所有像素點后，有，則用T更新，并記錄此子圖中心點的坐標。

下面給出算法的具體流程如圖2所示。

1.2 亞像素位移的測量

為了獲取更高的精度，需要在正像素結(jié)果的基礎上進一步進行亞像素位移的求解，設變形前的圖像為，變形后的圖像為，分別為對應于原圖像中所求位移點在變形的圖像中對應點的整像素位移，為對應于整像素位移結(jié)果的亞像素位移。

當選物體作微小位移時，且物體表面上任意一點在周圍的領域內(nèi)的元面積足夠小，則小面元可以看成近似剛體運動，亦元面內(nèi)所有的點的均勻相同的位移量。根據(jù)數(shù)字圖像基本假設，在微區(qū)內(nèi)，和有下面的關系：

同時，考慮到相關搜索對應變的不敏感性，設真實位移為：

定義是微區(qū)內(nèi)的刻畫與的相似程度函數(shù)：

在微區(qū)內(nèi)，和相似程度最大，應該滿足（4）式，此時應該取駐值。

將（4）式待入中并進一步表示的函數(shù)如下：

對應真實的微小變形應有：

本文中選取了公式（9）作為相似程度函數(shù)：

式中，為模板在點點處的灰度值，是變形子區(qū)在點處灰度值；分別是模板區(qū)域與變形子區(qū)的中像素灰度值的平均值。詳細的推導過程，見文獻[14]。

1.3 位移場的疊加

位移場疊加時考慮相鄰圖像之間的位置傳遞。例如，根據(jù)圖像1和圖像2，可以算出圖像2相對圖像1的位移增量為 ，其中是以圖像1為參考系的坐標。同理由圖像2、圖像3，得出，其中是以圖像2為坐標系，簡單推導后得到圖像3相對圖像1的位移場為：

據(jù)此類推，可以得到任意一個圖像相對于第一副圖像的位移。

2 數(shù)值模擬

采用斯坦福大學Peng Zhou等人提出的算法[15]，生成標準散斑圖如圖3，每一幅圖大小為512×512，散斑的光強程高斯分布，散斑尺寸大小為4個像素，散斑數(shù)為1200，在水平方向移動0～0.1pixel像素內(nèi)以0.01步長生成9幅圖像，在0.1～1像素位移范圍內(nèi)以0.1為步長，生成9幅圖像。同時在每幅圖像上取上5個不同位置的采用41×41模板，然后經(jīng)行統(tǒng)計分析，如圖4所示。endprint

由圖4可以看出，在理想的條件下，本文的亞像素有0.005pixel精度，單從精度要求上，滿足光鑷光鑷中對亞像素速位移精度的要求。同時發(fā)現(xiàn)在0～0.5pixel時，計算位移值大于設定位移值，而在0.5～1pixel時計算位移小于設定的位移值。設定位移在0.5pixel時出現(xiàn)反轉(zhuǎn)現(xiàn)象。

第二組實驗，生成10幅散斑圖像，每幅圖像沿y軸上移動移動3.25 pixel，在x軸上外加一個隨機0.1*rand的小位移量作為步進電機的運作時的擾動，同時兼顧電磁噪聲和其他噪聲，在后續(xù)的圖像中加均值為0，方差為0.01+0.01*rand高斯噪聲和噪聲密度0.01+rand*0.01的椒鹽噪聲如圖5所示。表1，采用相鄰的模板匹配得到的不同幀數(shù)上的匹配點的坐標。圖5，微粒的不同時刻的相對位置顯示。

由第一幅到最后一幅在y軸理論位移量為256+9×3.25=285.250；本文提出的算法最后結(jié)果為：285.195，誤差值為0.055pixel；而直接用第一幅和最后一幅相關算出的位移為：285.097；與理論值相差為0.1530pixel。顯然本文提出的算法有較強的抗干擾性。

第三組實驗，是比較在整像素點搜索采用動態(tài)閾值的SSDA的和普通的搜索的時間，而亞像素點的計算采用相同算法。對比，實驗結(jié)果如下表2。表格中的運行時間是在CPU為Intel（R） Core（TM）2 T5870，主頻為2.00 GHz處理器，內(nèi)存大小為2 G，計算5次匹配所畫的時間。實際上程序運行的時間取決于程序效率，編程語言、以及計算機硬件設備。本文算法基本比普通亞像素算法快1～2倍。倍數(shù)相差不明顯的原因是：亞像素計算的時間占整個計算時間很大的部分。

3 結(jié)論

由上面的分析可以得出，本文給出一種了利用圖像相關分析法通過整數(shù)像素位移場計算、亞像素位移計算、和位移場的疊加，來實現(xiàn)大位移場的高精度的亞像素位移測量的方法。在整像素的計算時采的動態(tài)的閾值的序貫相似法能使整個匹配過程所花時間節(jié)省1～2倍；第二組實驗中看出，粒子運動時提高采集圖像的頻率能的改善測量精度，本文直接用起始位置和終點位置的圖像計算出的位移誤差是采集多張圖像和計算相鄰的位移最后按1.3所述疊加后位移出差的近似9倍。此外，雖然本文盡力考慮了各種的干擾影響，步進電機的震動，電磁脈沖的等，但是實際實驗中的影響遠遠不如此，如細胞各自的布朗運動、焦平面的變化、光照的不均勻等。下一步工作應在光鑷捕獲微粒實驗中檢驗本算法的性能。

參考文獻

[1] 余娜，蔡志崗，梁業(yè)旺，等.光鑷系統(tǒng)的組建及光阱效應的觀察[J].大學物理，2010，29（3）：59-62.

[2] Ashkin A，Dziedzic JM， Bjorkholm JE，Chu S.Observation of a single-beam gradient force optical trap for dielectric particles[J].Opt Lett，1986，11（5）：288-290.

[3] 任洪亮，莊禮輝，李銀妹.雙光鑷測量膠體微粒間相互作用勢[J].中國激光，2008，35（1）：151-155.

[4] Florin E.-L，Pralle A，Stelzer E.H.K，et al.Photonic force microscope calibration by theraml noise analysis[J].J.Appl.Phys.A.1998，66：S75-S78.

[5] John Bechhoefer，Scott WilsonFaster cheaper，safer optical tweezers for the undergraduate laboratory[J].Am.J.Phys.2002，70（4）：393-400.

[6] 龔契，陳洪濤，李銀妹.四種光阱剛度測量法的實驗研究與比較[J].中國科大學報，2005，35：601.

[7] Bamea D I Silveman H E A class of algorithms for digital image reg istration[J].IEEE，1972，C-21（2）：179-186.

[8] 于起峰.基于圖像的精密測量與運動測量[M].北京：科學出版社，2002.

[9] 王懷文，亢一郎，謝和平.數(shù)字散斑相關方法與應用研究進展[J].力學進展，2005，32（2）：195-203.

[10] Bruck HA McNeil SR ，Sutton MA，et，al.Bigital Image correlation using Newton-Rapshon method of partial differential correction[J].Exprimental Mechanics，1989，29（3）：261-267.

[11] 潘兵，吳大方，國寶橋.數(shù)字體圖像相關方法中基于迭代最小二乘法的物體內(nèi)部變形測量[J].實驗力學，2011，26（6）：666-673.

[12] 潘兵，謝惠民，戴福隆.數(shù)字圖像相關中亞像素位移測量算法的研究[J].力學學報，2007，36（2）：245-252.

[13] 杜德生，葉建平，等.一種新的自適應閾值SSDA算法[J].現(xiàn)在電子技術，2010，6：135-139.

[14] 張軍，金光昌，馬少鵬，等.基于微區(qū)域統(tǒng)計特性的數(shù)字散斑相關測量亞像素位移梯度算法[J].光學技術，2003，29（4）：467-472.

[15] Zhou P，Goodson KE.Subpixel displacement and deformation gradient measurement using digital image speckle correlation（DISC）[J].Opt.Eng.2001，40（8）：1613-1620（August 2001）.endprint

由圖4可以看出，在理想的條件下，本文的亞像素有0.005pixel精度，單從精度要求上，滿足光鑷光鑷中對亞像素速位移精度的要求。同時發(fā)現(xiàn)在0～0.5pixel時，計算位移值大于設定位移值，而在0.5～1pixel時計算位移小于設定的位移值。設定位移在0.5pixel時出現(xiàn)反轉(zhuǎn)現(xiàn)象。

第二組實驗，生成10幅散斑圖像，每幅圖像沿y軸上移動移動3.25 pixel，在x軸上外加一個隨機0.1*rand的小位移量作為步進電機的運作時的擾動，同時兼顧電磁噪聲和其他噪聲，在后續(xù)的圖像中加均值為0，方差為0.01+0.01*rand高斯噪聲和噪聲密度0.01+rand*0.01的椒鹽噪聲如圖5所示。表1，采用相鄰的模板匹配得到的不同幀數(shù)上的匹配點的坐標。圖5，微粒的不同時刻的相對位置顯示。

由第一幅到最后一幅在y軸理論位移量為256+9×3.25=285.250；本文提出的算法最后結(jié)果為：285.195，誤差值為0.055pixel；而直接用第一幅和最后一幅相關算出的位移為：285.097；與理論值相差為0.1530pixel。顯然本文提出的算法有較強的抗干擾性。

第三組實驗，是比較在整像素點搜索采用動態(tài)閾值的SSDA的和普通的搜索的時間，而亞像素點的計算采用相同算法。對比，實驗結(jié)果如下表2。表格中的運行時間是在CPU為Intel（R） Core（TM）2 T5870，主頻為2.00 GHz處理器，內(nèi)存大小為2 G，計算5次匹配所畫的時間。實際上程序運行的時間取決于程序效率，編程語言、以及計算機硬件設備。本文算法基本比普通亞像素算法快1～2倍。倍數(shù)相差不明顯的原因是：亞像素計算的時間占整個計算時間很大的部分。

3 結(jié)論

由上面的分析可以得出，本文給出一種了利用圖像相關分析法通過整數(shù)像素位移場計算、亞像素位移計算、和位移場的疊加，來實現(xiàn)大位移場的高精度的亞像素位移測量的方法。在整像素的計算時采的動態(tài)的閾值的序貫相似法能使整個匹配過程所花時間節(jié)省1～2倍；第二組實驗中看出，粒子運動時提高采集圖像的頻率能的改善測量精度，本文直接用起始位置和終點位置的圖像計算出的位移誤差是采集多張圖像和計算相鄰的位移最后按1.3所述疊加后位移出差的近似9倍。此外，雖然本文盡力考慮了各種的干擾影響，步進電機的震動，電磁脈沖的等，但是實際實驗中的影響遠遠不如此，如細胞各自的布朗運動、焦平面的變化、光照的不均勻等。下一步工作應在光鑷捕獲微粒實驗中檢驗本算法的性能。

參考文獻

[1] 余娜，蔡志崗，梁業(yè)旺，等.光鑷系統(tǒng)的組建及光阱效應的觀察[J].大學物理，2010，29（3）：59-62.

[2] Ashkin A，Dziedzic JM， Bjorkholm JE，Chu S.Observation of a single-beam gradient force optical trap for dielectric particles[J].Opt Lett，1986，11（5）：288-290.

[3] 任洪亮，莊禮輝，李銀妹.雙光鑷測量膠體微粒間相互作用勢[J].中國激光，2008，35（1）：151-155.

[4] Florin E.-L，Pralle A，Stelzer E.H.K，et al.Photonic force microscope calibration by theraml noise analysis[J].J.Appl.Phys.A.1998，66：S75-S78.

[5] John Bechhoefer，Scott WilsonFaster cheaper，safer optical tweezers for the undergraduate laboratory[J].Am.J.Phys.2002，70（4）：393-400.

[6] 龔契，陳洪濤，李銀妹.四種光阱剛度測量法的實驗研究與比較[J].中國科大學報，2005，35：601.

[7] Bamea D I Silveman H E A class of algorithms for digital image reg istration[J].IEEE，1972，C-21（2）：179-186.

[8] 于起峰.基于圖像的精密測量與運動測量[M].北京：科學出版社，2002.

[9] 王懷文，亢一郎，謝和平.數(shù)字散斑相關方法與應用研究進展[J].力學進展，2005，32（2）：195-203.

[10] Bruck HA McNeil SR ，Sutton MA，et，al.Bigital Image correlation using Newton-Rapshon method of partial differential correction[J].Exprimental Mechanics，1989，29（3）：261-267.

[11] 潘兵，吳大方，國寶橋.數(shù)字體圖像相關方法中基于迭代最小二乘法的物體內(nèi)部變形測量[J].實驗力學，2011，26（6）：666-673.

[12] 潘兵，謝惠民，戴福隆.數(shù)字圖像相關中亞像素位移測量算法的研究[J].力學學報，2007，36（2）：245-252.

[13] 杜德生，葉建平，等.一種新的自適應閾值SSDA算法[J].現(xiàn)在電子技術，2010，6：135-139.

[14] 張軍，金光昌，馬少鵬，等.基于微區(qū)域統(tǒng)計特性的數(shù)字散斑相關測量亞像素位移梯度算法[J].光學技術，2003，29（4）：467-472.

[15] Zhou P，Goodson KE.Subpixel displacement and deformation gradient measurement using digital image speckle correlation（DISC）[J].Opt.Eng.2001，40（8）：1613-1620（August 2001）.endprint

由圖4可以看出，在理想的條件下，本文的亞像素有0.005pixel精度，單從精度要求上，滿足光鑷光鑷中對亞像素速位移精度的要求。同時發(fā)現(xiàn)在0～0.5pixel時，計算位移值大于設定位移值，而在0.5～1pixel時計算位移小于設定的位移值。設定位移在0.5pixel時出現(xiàn)反轉(zhuǎn)現(xiàn)象。

第二組實驗，生成10幅散斑圖像，每幅圖像沿y軸上移動移動3.25 pixel，在x軸上外加一個隨機0.1*rand的小位移量作為步進電機的運作時的擾動，同時兼顧電磁噪聲和其他噪聲，在后續(xù)的圖像中加均值為0，方差為0.01+0.01*rand高斯噪聲和噪聲密度0.01+rand*0.01的椒鹽噪聲如圖5所示。表1，采用相鄰的模板匹配得到的不同幀數(shù)上的匹配點的坐標。圖5，微粒的不同時刻的相對位置顯示。

由第一幅到最后一幅在y軸理論位移量為256+9×3.25=285.250；本文提出的算法最后結(jié)果為：285.195，誤差值為0.055pixel；而直接用第一幅和最后一幅相關算出的位移為：285.097；與理論值相差為0.1530pixel。顯然本文提出的算法有較強的抗干擾性。

第三組實驗，是比較在整像素點搜索采用動態(tài)閾值的SSDA的和普通的搜索的時間，而亞像素點的計算采用相同算法。對比，實驗結(jié)果如下表2。表格中的運行時間是在CPU為Intel（R） Core（TM）2 T5870，主頻為2.00 GHz處理器，內(nèi)存大小為2 G，計算5次匹配所畫的時間。實際上程序運行的時間取決于程序效率，編程語言、以及計算機硬件設備。本文算法基本比普通亞像素算法快1～2倍。倍數(shù)相差不明顯的原因是：亞像素計算的時間占整個計算時間很大的部分。

3 結(jié)論

由上面的分析可以得出，本文給出一種了利用圖像相關分析法通過整數(shù)像素位移場計算、亞像素位移計算、和位移場的疊加，來實現(xiàn)大位移場的高精度的亞像素位移測量的方法。在整像素的計算時采的動態(tài)的閾值的序貫相似法能使整個匹配過程所花時間節(jié)省1～2倍；第二組實驗中看出，粒子運動時提高采集圖像的頻率能的改善測量精度，本文直接用起始位置和終點位置的圖像計算出的位移誤差是采集多張圖像和計算相鄰的位移最后按1.3所述疊加后位移出差的近似9倍。此外，雖然本文盡力考慮了各種的干擾影響，步進電機的震動，電磁脈沖的等，但是實際實驗中的影響遠遠不如此，如細胞各自的布朗運動、焦平面的變化、光照的不均勻等。下一步工作應在光鑷捕獲微粒實驗中檢驗本算法的性能。

參考文獻

[1] 余娜，蔡志崗，梁業(yè)旺，等.光鑷系統(tǒng)的組建及光阱效應的觀察[J].大學物理，2010，29（3）：59-62.

[2] Ashkin A，Dziedzic JM， Bjorkholm JE，Chu S.Observation of a single-beam gradient force optical trap for dielectric particles[J].Opt Lett，1986，11（5）：288-290.

[3] 任洪亮，莊禮輝，李銀妹.雙光鑷測量膠體微粒間相互作用勢[J].中國激光，2008，35（1）：151-155.

[4] Florin E.-L，Pralle A，Stelzer E.H.K，et al.Photonic force microscope calibration by theraml noise analysis[J].J.Appl.Phys.A.1998，66：S75-S78.

[5] John Bechhoefer，Scott WilsonFaster cheaper，safer optical tweezers for the undergraduate laboratory[J].Am.J.Phys.2002，70（4）：393-400.

[6] 龔契，陳洪濤，李銀妹.四種光阱剛度測量法的實驗研究與比較[J].中國科大學報，2005，35：601.

[7] Bamea D I Silveman H E A class of algorithms for digital image reg istration[J].IEEE，1972，C-21（2）：179-186.

[8] 于起峰.基于圖像的精密測量與運動測量[M].北京：科學出版社，2002.

[9] 王懷文，亢一郎，謝和平.數(shù)字散斑相關方法與應用研究進展[J].力學進展，2005，32（2）：195-203.

[10] Bruck HA McNeil SR ，Sutton MA，et，al.Bigital Image correlation using Newton-Rapshon method of partial differential correction[J].Exprimental Mechanics，1989，29（3）：261-267.

[11] 潘兵，吳大方，國寶橋.數(shù)字體圖像相關方法中基于迭代最小二乘法的物體內(nèi)部變形測量[J].實驗力學，2011，26（6）：666-673.

[12] 潘兵，謝惠民，戴福隆.數(shù)字圖像相關中亞像素位移測量算法的研究[J].力學學報，2007，36（2）：245-252.

[13] 杜德生，葉建平，等.一種新的自適應閾值SSDA算法[J].現(xiàn)在電子技術，2010，6：135-139.

[14] 張軍，金光昌，馬少鵬，等.基于微區(qū)域統(tǒng)計特性的數(shù)字散斑相關測量亞像素位移梯度算法[J].光學技術，2003，29（4）：467-472.

[15] Zhou P，Goodson KE.Subpixel displacement and deformation gradient measurement using digital image speckle correlation（DISC）[J].Opt.Eng.2001，40（8）：1613-1620（August 2001）.endprint