泛函回歸代理及條件期望配準(zhǔn)的機械擺動測量

鄭思凡，陳平平，蘇凱雄，吳永春

（1. 福州大學(xué)物理與信息工程學(xué)院，福建福州350116；2. 泉州黎明職業(yè)大學(xué)智能制造工程學(xué)院，福建泉州362000；3. 數(shù)字電視智能化技術(shù)國家地方聯(lián)合工程研究中心，福建福州350116）

1 引 言

機械擺動測量在精密制造與智能制造的精密閉環(huán)控制具有廣泛的應(yīng)用，如五軸數(shù)控機床擺動式回轉(zhuǎn)軸角位置測量［1］，固體發(fā)動機柔性擺動噴管擺角測量［2］，起重機吊鉤擺動角測量［3］，偏擺儀指針偏擺檢測等。相應(yīng)地，工程上也出現(xiàn)了各種基于光電技術(shù)的高精度測角儀［4］，但一般對測量相對位置與光照環(huán)境要求較高，測量對象單一，無法滿足廣泛應(yīng)用需求，因而如何運用圖像運動感知技術(shù)對上述機械擺件進行逐幀提取分割，對擺角進行視覺測量則具有很大的工業(yè)實際意義。

隨著圖像采集及視頻壓縮傳輸軟硬件技術(shù)的發(fā)展，高幀率，高分辨率自帶光學(xué)防抖的CMOS 與CCD 攝像機的出現(xiàn)為高速的機械運動非接觸式視覺測量提供了硬件實現(xiàn)條件，如三星6 400 萬像素ISOCELL Bright GW1 CMOS 感光元，其單個像素面積僅0.8 μm，支持Tetracell 像素合成技術(shù)和Remosaic 還原，可直接輸出480 FPS 的1 080 pixel 圖像［5］。當(dāng)前國內(nèi)也有關(guān)于三維運動圖像測量分析系統(tǒng)的相關(guān)軟件產(chǎn)品，如北京現(xiàn)代富博科技有限公司的MIAS3.0［6］，但是產(chǎn)品需要手工標(biāo)定跟蹤像素，在自動化方面存在提升空間。

當(dāng)前 SLAM［7］（Simultaneous Localization and Mapping）地圖重建技術(shù)及其衍生品視覺里程計［8］（Visual Odometry）技術(shù)的發(fā)展為視覺測速提供了許多成熟的算法框架，二者的核心框架均在于通過相鄰幀對應(yīng)目標(biāo)點云配準(zhǔn)完成對相機的位姿運動參數(shù)的估計。對于SLAM 技術(shù)的視覺測速，當(dāng)前也有比較成熟的研究成果，如圖優(yōu)化工具g2o（General Graph Optimization）［9］所集成的基于指數(shù)映射的李代數(shù)擾動模型的集束調(diào)整（Bundle Adjustment）［10］算法。該算法將有約束的旋轉(zhuǎn)矩陣?yán)钊耗繕?biāo)函數(shù)，轉(zhuǎn)化為無約束的李代數(shù)矢量優(yōu)化問題，從而可以用普通牛頓迭代法求解。但上述研究僅局限于相機本身運動估計，如何將點云配準(zhǔn)結(jié)合運動分割，用于工業(yè)上機械剛體點云的運動參數(shù)估計與測量則有較大研究空間。

為實現(xiàn)對弱結(jié)構(gòu)張量特征點光流的運動分割，T. Brox，J. Malik 等提出通過對變分稠密光流軌跡時空相似度拉氏譜圖譜聚類的運動分割算法［11-13］。此類算法是采用ε-NN 或K-NN 近鄰矩陣來構(gòu)建，其超參數(shù)ε 或K的取值由數(shù)據(jù)集訓(xùn)練后統(tǒng)一設(shè)定，無法根據(jù)實際光流稠密度不同自適應(yīng)調(diào)整。為此Roberto Tron，Ren′e Vidal 等根據(jù)SFM（Structure From Motion）［14］三維重建理論為基礎(chǔ)，將運動分割轉(zhuǎn)化為剛體光流軌跡，聚類到2～4 維的低維超平面子空間里來完成運動分割。代表性算法有廣義主成分分析、局部子空間相似度聚類、低秩分解、低秩子空間聚類、多段學(xué)習(xí)等［15］。上述算法充分挖掘了剛體光流軌跡蘊含的子空間結(jié)構(gòu)信息，克服了光流空間分散及密度不均造成的分割尺度誤差。但是各子空間的維度及個數(shù)參數(shù)仍然需要事先指定。有鑒于此，von 等于2013 年提出了稀疏子空間聚類［16］，將傳統(tǒng)的基于ε-NN 或K-NN 近鄰矩陣替換為L1范數(shù)圖。其中各光流軌跡為圖的節(jié)點，邊的權(quán)重則是以其他軌跡作為字典基進行稀疏學(xué)習(xí)后的自線性表達系數(shù)組成。這樣利用高維數(shù)據(jù)的稀疏性使得同一流形內(nèi)部取得不為零的表達系數(shù)，而不同流形的節(jié)點其代表的軌跡互表達系數(shù)為零。從而以特征選擇的方式自適應(yīng)表達了光流軌跡密度不均的流形結(jié)構(gòu)，解決了ε-NN 或K 自適應(yīng)調(diào)整問題。論文［17-19］等通過擴展修正不同的字典基表達系數(shù)的正則約束，改進SSC 聚類的類內(nèi)一致性及稀疏性。但是，上述算法均是針對塊狀移動物體進行分割提取，在機械擺動的視覺運動分割中，因背景復(fù)雜性及運動對象的材質(zhì)（如金屬件），存在大量弱結(jié)構(gòu)張量特征點及被中斷為多段零散的光流軌跡，利用傳統(tǒng)的時空相似度譜聚類易形成碎片化的塊狀軌跡群，難以一次性完成所有軌跡的子空間聚類，而是需要分段多次進行，不僅計算量大且在光流遮擋處密度過于稀疏影響配準(zhǔn)精度，因此SSC在機械擺動的運動分割方面存在較大的局限性。

為此本文提出了一種以高斯隨機過程泛函回歸代理的弱結(jié)構(gòu)張量機械擺動測量算法。首先獲取機械擺動件LDOF 變分光流；其次通過構(gòu)建新的鄰接矩陣與相似度度量獲取過分割狀態(tài)的弧狀軌跡群；再者以時間為參變量分別對軌跡群x，y坐標(biāo)進行高斯回歸，學(xué)習(xí)出該軌跡群的平均軌跡，此平均軌跡代理該軌跡群作為稀疏子空間聚類的種子樣本，一次性完成聚類；這樣即克服了碎片化導(dǎo)致的分段多次聚類，也通過SSC 的模型選擇得到正確的運動物體數(shù)目從而克服了過分割的局限性，最后對每幀運動分割后的機械擺件特征點云匹配關(guān)系建模為各匹配對應(yīng)點為中心的GMM 分布，將該分布的數(shù)據(jù)最大似然參數(shù)估計作為運動位姿矩陣估計值，通過提取旋轉(zhuǎn)矩陣歐拉角完成擺角計算及測量。為證明所提算法的有效性，本文以理想雙搖桿模型的6 種不同照度下的車輛刮水器擺桿為對象，結(jié)合基于機器視覺的車輛安全自動化日檢工程項目于安檢現(xiàn)場提取了30 組數(shù)據(jù)實驗。分析比較了所提算法與當(dāng)前算法對擺桿的光流軌跡聚類角位移誤差。實驗結(jié)果表明：本算法能完整學(xué)習(xí)出等長軌跡并在滿足一定的環(huán)境亮度情況下，角位移測算值與人為標(biāo)定的回歸值均方差小于10%。同時解決了當(dāng)前算法存在的過分割，稠密度不均導(dǎo)致的配準(zhǔn)誤差及多次分段與迭代的運算量問題。

本文首先介紹測量對象四連桿機構(gòu)動力學(xué)原理及其理想角位移曲線；其次介紹稠密變分光流種類及被測刮刷總成運動光流軌跡的獲??；第三部分介紹傳統(tǒng)塊狀譜聚類運動分割存在的過分割及碎片化兩個局限性現(xiàn)象及原因并通過修改傳統(tǒng)塊狀軌跡群的相似度度量與近鄰矩陣構(gòu)建方法引出本文的弧狀超流素軌跡群概念及實現(xiàn)；第四部分介紹超流素泛函回歸代理參數(shù)的訓(xùn)練及推斷；第五部分介紹算法的聚類性能比較；第六部分為算法角位移測量誤差效果比較；最后為總結(jié)與展望。

2 四連桿機構(gòu)動力學(xué)原理

當(dāng)前客運車輛較常用的刮水器總成結(jié)構(gòu)由兩部分四連桿組成，為研究方便，可在solidworks建模如圖1 所示。圖中L1為左搖桿、L2為機架、L3為右搖桿、L4為耦合連桿，L1～4構(gòu)成雙搖桿機構(gòu)，L5為雨刮擺桿，L6為輸入連桿，L7為電機驅(qū)動曲柄，三者構(gòu)成主刮曲柄擺桿系統(tǒng)；通過設(shè)定左右擺桿長度起到控制輸入動桿與輸出桿速度和相位的作用。

為了構(gòu)建仿真算例，觀察左右搖桿角位移變化曲線，取L1=L2=75 cm，L3=L4=35 cm 并設(shè)仿真時間設(shè)置為5 s，電機轉(zhuǎn)速r=30 rmps，動畫幀率為30 FPS。算例求解完成后，得左搖桿的角速度變化曲線分別如圖2 所示。

圖1 刮水器連桿結(jié)構(gòu)Fig.1 Wiper link structure

圖2 左搖桿角位移變化曲線Fig.2 Angular velocity curve of left Rocker

由圖可見四連桿刮刷運動過程呈周期性變化，且在啟動與結(jié)束處速度存在收斂的過渡過程。若按照圖像逐幀分割進行還原，則會呈現(xiàn)如圖3 圖樣。

圖3 刮桿逐幀運動分割圖樣Fig.3 Fame-by-frame motion segmentation of Wiper

根據(jù)矢量運算原理，同樣可以得到右搖桿的角位移曲線，因為不同擺桿其在SFM 模型對應(yīng)的運動矩陣平移分量不同，所以左右擺桿處于不同仿射空間，在理想情況下，可以直接使用子空間聚類對滿幅等長軌跡一次性完成聚類分割，相應(yīng)地，各幀分割形成的點云經(jīng)過配準(zhǔn)并提取位姿旋轉(zhuǎn)矩陣的歐拉角則可一次性求得角位移曲線，但在復(fù)雜背景下將呈現(xiàn)為多段碎片化分割區(qū)域，無法一次完成分割。

3 稠密變分光流定義及計算

為了解決傳統(tǒng)光流亮度一致性方程梯度逆陣病態(tài)帶來的孔徑問題，Horn-Schunck 等提出了數(shù)據(jù)項與正則項聯(lián)合變分的稠密光流［11］。其中數(shù)據(jù)項主要完成亮度一致性約束，具體形式如式（1）所示：

其中：T(x，y)為參考圖像，I(x，y)為當(dāng)前圖像。u(x，y)，v(x，y)是圖像上每一點的偏移量，φ(x)為代價函數(shù)，可取絕對值，二次函數(shù)等一些非負(fù)對稱且正半軸單調(diào)遞增的函數(shù)。

為了對純色區(qū)域的弱結(jié)構(gòu)張量像素實現(xiàn)幀間稠密匹配，在數(shù)據(jù)項加入一個平滑項，使得相鄰兩個特征點的偏移量相差不能太大：

為研究方便，本文采用華睿A3200CU000 面陣高速工業(yè)相機進行采集。在幀率為120 FPS，分辨率720 pixel 的設(shè)定下，在安檢臺采集客車雨刮器上升動作視頻180 幀分析。當(dāng)設(shè)定空間采樣率為8、像素灰度結(jié)構(gòu)張量閾值為0.4 的情況下，可以采用前述OpenCV3.4 提供的BroxOptical-Flow 類完成實時LDOF 光流計算并生成相應(yīng)的flo 文件，為了保證軌跡的魯棒性，在各相鄰幀的光流連線形成軌跡前，必須對各光流進行前向與后向方向一致性檢查［14］，這里設(shè)定一致性系數(shù)0.06，一致性偏移閾值為4，可得雨刮的軌跡分布如圖4 所示。

圖4 刮水器總的光流軌跡分布Fig.4 Total optical flow trajectory distribution of wiper

4 當(dāng)前運動分割局限性及原因分析

這里先考慮以軌跡時空相似度譜聚類前述光流軌跡的運動分割算法［11］，其中軌跡的相似度以共有幀的位移及速度范數(shù)距離的RBF 核來構(gòu)建，此相似度進一步作為ε-NN 近鄰矩陣的邊權(quán)重A(tri，trj)即：

其中：wv(tri，trj)為軌跡tri與軌跡trj在幀t的速度相似度，定義如下：

其中：Lrw為隨機矩陣，其元素pij=wij di可以看成是節(jié)點i通過隨機行走到達節(jié)點j的概率［21］，其中wij為節(jié)點i與節(jié)點j之間相似度權(quán)值，di為連接到節(jié)點i的邊總數(shù)。這樣通過特征分解可得特征矢量如式（8）：

最后取特征值小于0.15 對應(yīng)的特征矢量矩陣，構(gòu)成行矢量空間進行k-mean 聚類。圖5 展示了特征矢量對應(yīng)的特征值從小到大排列。

圖5 SC 特征值排列Fig.5 Eigenvalue index of Spectral Clustering

由圖可見，前5 個特征矢量特征值小于0.15，因此直接對前述光流軌跡譜聚類得到如下的聚類結(jié)果，如圖6 所示。

圖6 譜聚類運動分割結(jié)果Fig.6 Spectral clustering motion segmentation results

圖7 拉氏特征映射低維嵌入分布圖Fig.7 Low-dimensional embedding of Laplace Eigenmaps

相應(yīng)地，圖7 展示了軌跡相似度在1～5 各維度的低維嵌入分布。由上述可見，直接對軌跡相似度矩陣譜聚類能夠得到較好的類內(nèi)一致性。相似度矩陣也顯示出塊狀結(jié)構(gòu)，各低維嵌入的同類數(shù)據(jù)也保持相鄰的拓?fù)潢P(guān)系。

但是因其構(gòu)建相似度矩陣所采用統(tǒng)一固定的鄰域ε 參數(shù)，與尺度參數(shù)η無法反映同一擺桿不同段光流線速度存在的差異，導(dǎo)致其受到軌跡密度與速度分布不均的影響，將兩根擺桿過分割為5 根。因此，算法沒有利用同一剛體的軌跡處于同一子空間的線性表達關(guān)系的內(nèi)在信息，得到準(zhǔn)確的聚類超參數(shù)K，出現(xiàn)了過分割現(xiàn)象。

下面介紹分析此算法的另一個局限性即碎片化現(xiàn)象，為此在前述得到軌跡相似度拉氏矩陣的情況下，直接取K=24 得到超流素并以顏色區(qū)分，如圖8 所示（彩圖見期刊電子版）。由圖可以看出，傳統(tǒng)譜聚類運動分割在K 取較大值時，因為遮擋及軌跡的中斷將使得分割結(jié)果碎片化，無法理想將線速度相近的軌跡聚類為同一軌跡群。

圖8 SC 軌跡聚類碎片化現(xiàn)象Fig.8 fragmentation phenomenon caused by SC

下面從相似度的距離度量定義及鄰接矩陣的構(gòu)建方式兩個方面來說明上述聚類形成碎片化的原因并提出弧狀軌跡群的分割方法。

第一，在構(gòu)建鄰接矩陣方面，以圖9 為例，參考文獻［11］提出以軌跡公共時間幀里兩條軌跡相應(yīng)X，Y位移差最大值作為ε，則圖9 中軌跡A與B，C與D兩兩之間滿足條件構(gòu)成鄰接關(guān)系，而A與C，B與D因為中間出現(xiàn)了軌跡中斷，沒有公共時間幀，于是鄰接關(guān)系為零，不構(gòu)成鄰接關(guān)系，在后續(xù)的聚類結(jié)果里被聚類為不同類別，E，F(xiàn)，G，H同樣如此，從而出現(xiàn)了碎片化。

圖9 SC 軌跡聚類碎片化現(xiàn)象Fig.9 Fragmentation phenomenon caused by SC

于是本文在原鄰接關(guān)系的基礎(chǔ)上，將共有時間幀推廣到在軌跡終點處且時間差Δt在10 幀以內(nèi)且軌跡起點在ε位移半徑以內(nèi)的另一條軌跡建立軌跡鄰接關(guān)系，這樣就使得A與C，B與D分別構(gòu)成鄰接關(guān)系為后面將他們聚類為同一類別創(chuàng)造了基礎(chǔ)。

其次，從公式（4）定義的時空相似度來看，圖10 中A與B的相似度大于A與C和B與D，同理E與F相似度大于E與G和F與H（塊狀軌跡群相似度＞弧狀軌跡群相似度），因此在K值較小時A與C，B與D還有可能通過某條中間橋接軌跡分別聚為同一類，但是在K取較大時A與B，C與D各自聚為同一類形成塊狀軌跡群，而AC，BD則聚為不同軌跡群，從而造成K值較大時出現(xiàn)碎片化。

圖10 稀疏子空間運動分割結(jié)果Fig.10 Results of Sparse subspace clustering

為了使A與C聚為同一類形成弧狀軌跡群，而A與B形成不同軌跡群，考慮到同一半徑處的像素運動軌跡不管被中斷為幾段，各段的平均曲率或線速度絕對值相似，因此這里將軌跡相似度距離度量修改為各軌跡的平均曲率，使得AC的相似度大于AB的相似度，同理，EG相似度大于EF相似度，這樣在K 取較大值時準(zhǔn)確構(gòu)建弧狀軌跡群，為后面的GP回歸及Lasso 回歸提供了正確的種子樣本。

綜合上述，這里ε取8，K取48，將各段軌跡經(jīng)立方樣條spline 擬合后對函數(shù)結(jié)構(gòu)求一二階導(dǎo)數(shù)所得平均曲率作為近鄰矩陣的邊權(quán)重A(tri，trj)后再譜聚類可得48 個弧狀軌跡群如圖11。

由圖可見，上述分割方式解決了直接譜聚類在弧向方向的碎片化問題。各弧狀軌跡群在徑向方向則形成高維的線性表達關(guān)系，這一特征可以用子空間聚類方式再進一步聚類，從而可以以無監(jiān)督的機器學(xué)習(xí)方式經(jīng)過兩層回歸對原始零散的軌跡群進行完整聚類。

圖11 弧狀軌跡群聚類對碎片化的改進結(jié)果Fig.11 Improved results by hierarchical clustering

5 超流素的定義及泛函回歸參數(shù)的學(xué)習(xí)與樣本推斷

5.1 超流素定義

超流素是在軌跡元（Tracklet）［20-22］基礎(chǔ)上經(jīng)過一層聚類的軌跡群，代表一段視頻中光流軌跡時空相似度最小分割單位，將其作為視頻分割的中層表達形式由文獻［23］提出。為了與相關(guān)文獻保持名稱一致，下面以超流素的概念作為本文軌跡群的定義。

5.2 高斯回歸定義及新值統(tǒng)計推斷

高斯回歸代理作為一種泛函回歸［24］，是基于假定采樣點之間的距離或方向可存在變化的空間相關(guān)性基礎(chǔ)上，將此相關(guān)性以統(tǒng)計模型進行學(xué)習(xí)擬合，利用統(tǒng)計推斷對新測試值做出最佳預(yù)測。其本質(zhì)可以看為一種貝葉斯推斷，因此如果數(shù)據(jù)X似然部分采用高斯似然作為核函數(shù)，則可以利用共軛先驗特性獲取新數(shù)據(jù)協(xié)方差∑*及平均值μ*的閉式解。由此通過貝葉斯推斷可得新測試樣本X*的泛函分布f*如下：

其中：X為樣本數(shù)據(jù)，y為樣本對應(yīng)的訓(xùn)練值，K為數(shù)據(jù)核函數(shù)矩陣，此矩陣可以分為樣本數(shù)據(jù)與測試數(shù)據(jù)兩部分，分別以腳標(biāo)*與y標(biāo)出，定義為：K=κ(X，X)；K*=κ(X，X*)；K**=κ(X*，X*)。

將前述超流素的坐標(biāo)位移容限15 看為高斯回歸的數(shù)據(jù)噪聲協(xié)方差，則可得核函數(shù)矩陣元素對應(yīng)核函數(shù)如下：

5.3 高斯回歸超參數(shù)學(xué)習(xí)

設(shè)θ=(?，σ2f)則可得上述MAP 估計對超參數(shù)的梯度如下：

令此梯度為零可得超參數(shù)閉式解，而上述高斯過程對新數(shù)據(jù)的MAP 估計等價于無限維的RKHS 空間二次范數(shù)正則化泛函的線性回歸，根據(jù)里斯表示定理（Riesz Representation Theorem）可以由數(shù)據(jù)樣本張成的有限維空間表示如下：

其中：α?y為各已知樣本與測試點距離的加權(quán)系數(shù)是測試點與各樣本的核函數(shù)矢量，其元素為κ(xi，x*)，即測試數(shù)據(jù)與樣本映射到RKHS 希爾伯特空間的內(nèi)積距離。

為反映各軌跡時間方面的同步與次序信息，這里對超流素軌跡的回歸以時間為中間參變量，分別對軌跡群的橫縱坐標(biāo)進行插值。如前述取超參數(shù)=15，=20，?=480 為初值，利用公式（14）對每個超流素的超參數(shù)重新學(xué)習(xí)后，由公式（9）對所有時間段橫坐標(biāo)值重新推斷，則可得回歸的平均橫坐標(biāo)值如圖12 所示。

圖12 超流素橫坐標(biāo)隨時間的高斯回歸Fig.12 Gaussian regression of Superfloxel’s x coordinate

相應(yīng)地，y坐標(biāo)的回歸平均軌跡如圖13所示。

圖13 超流素縱坐標(biāo)隨時間的高斯回歸Fig.13 Gaussian regression of Superfloxel’s y coordinate

5.4 算法實現(xiàn)及時間復(fù)雜度

為了考量算法時間復(fù)雜度，這里假定樣本數(shù)為r。由公式（17）可見，為了算法穩(wěn)定性，不宜直接對Ky求逆，這里采用對Ky做Cholesky 分解為Ky=LLT則算法具體實現(xiàn)如下：

其中步驟1 的Cholesky 分解時間復(fù)雜度為O(r3)，步驟2 求α?xí)r間復(fù)雜度為O(r2)，步驟3的矩陣乘法需要O(r)，步驟4 和步驟5 運算復(fù)雜度為O(r2)。為減少運算量，在實踐中只需要求出步驟3 所表示的坐標(biāo)值為回歸軌跡坐標(biāo)，步驟4～6 所求的回歸置信度可以省去。為進一步降低運算量，這里采用子集回歸的稀疏高斯回歸（SoD）［25］算法對數(shù)據(jù)核Gram 矩陣K分塊如下：

其中：Km m由m個引導(dǎo)點（Inducing Point）構(gòu)成活動集I（Active Set）對應(yīng)的分塊矩陣，則經(jīng)過SVD分解后原來核Gram 矩陣K可用其前最大m個特征值及相應(yīng)特征向量構(gòu)成的m階降秩矩陣近似，其定義如下：

相應(yīng)地，公式（17）核函數(shù)由測試點與子集樣本的核函數(shù)線性組合近似表達為：

其中：j為引導(dǎo)點序號，cij為線性表達系數(shù)，其構(gòu)成矩陣C由式（19）確定：

其閉式解為：

Copt=Knm K-1mm，這樣新的子集回歸對應(yīng)的平均值與協(xié)方差近似值如下：

這樣總的時間復(fù)雜度降為O(m2r)，這里取m值為幀數(shù)，并分別對前述48 個超流素進行坐標(biāo)插值回歸，可得回歸后的軌跡形狀如圖14 所示。

圖14 超流素高斯回歸平均軌跡Fig.14 Gaussian regression of Superfloxel

6 聚類性能比較

將上述回歸的軌跡作為種子樣本，假設(shè)軌跡長度為F幀，軌跡數(shù)為N，將軌跡每個點坐標(biāo)按照先橫坐標(biāo)后縱坐標(biāo)次序排列為矩陣Y的每一列。

則SSC 通過如下目標(biāo)函數(shù)求解線性自表達系數(shù)矩陣X：

其中：E為字典重構(gòu)誤差，‖E‖為誤差E 的Frobenius 范數(shù)，當(dāng)假定重構(gòu)誤差為高斯分布時此范數(shù)反映誤差代價，λ為權(quán)衡參數(shù)，平衡軌跡誤差與解的稀疏性及泛函的凸性，因為本文是針對回歸后的軌跡聚類，故誤差相比理想情況下的軌跡大，故對λ取較大的值0.1 以克服回歸帶來的誤差。

通過引入拉格朗日乘子Δ ∈RN×N及懲罰系數(shù)μ可將上式轉(zhuǎn)化為無約束目標(biāo)函數(shù)如下：

文獻［16］通過ADMM（交替方向乘子法）算法輪流優(yōu)化變量E，X，Δ 及μ將上述目標(biāo)函數(shù)轉(zhuǎn)化為兩個子優(yōu)化問題如下：

其中：θ=‖Y‖?X(Xk，Ek，Δk，μk)為目標(biāo)函數(shù)L對X的偏導(dǎo)函數(shù)，k為迭代次數(shù)，則公式（25）變?yōu)闃?biāo)準(zhǔn)LASSO 問題，可以采用通用迭代收縮閾值算法（Generalized Iterated Shrinkage Algorithm）［26］求解。另外通過引入輔助變量，即假設(shè)B=YX(k)-Y-Δ(k)/μ(k)，則公式（27）可改寫為：

一樣轉(zhuǎn)化為針對每個數(shù)據(jù)εj的標(biāo)準(zhǔn)LASSO問題從而同上求解。由式（26），式（27）可見ADMM 算法每進行一次迭代都需要求N 階逆陣，其運算復(fù)雜度為O(N3)，在N取1 000，m取100，超流素K取50 時，則r為2 000，可見ADMM 的單步復(fù)雜度O(N3)大于稀疏高斯回歸O(m2r)，但是高斯回歸不需要迭代，因此本文算法可減少一階運算量。

最后將所求的矩陣X元素作為L1 范數(shù)圖節(jié)點間的權(quán)值并預(yù)先取K=2 進行譜聚類可得如圖15 所示結(jié)果。

圖15 代理軌跡的SSC 聚類結(jié)果Fig.15 SSC clustering result of the surrogate trajectories

由圖可見，回歸軌跡代理的SSC 聚類（RSSSC）準(zhǔn)確標(biāo)識了左右擺桿，但上述的聚類數(shù)目為事先取定，依據(jù)前述譜圖理論，聚類個數(shù)的值可由特征值為零的個數(shù)決定，相應(yīng)圖5，這里將特征矢量對應(yīng)的特征值從小到大排列如圖16所示。

相應(yīng)的48 條種子軌跡相似度拉氏圖矩陣及1～2 維的拉氏映射嵌入分別如圖17～圖18 所示。由圖可見由線性自表達系數(shù)構(gòu)成L1 范數(shù)圖構(gòu)成的軌跡相似度圖相比原來的相似度矩陣，在特征值分布和低維嵌入數(shù)據(jù)的類間距離，更接近超參數(shù)K真實值2。

圖16 回歸軌跡代理的SSC 特征值排列Fig.16 Eigenvalue index of the regression surrogate SSC

從而說明了經(jīng)過超流素回歸軌跡代理的SSC 聚類在K-mean 聚類前可以通過矩陣干擾理論或模型選擇等方法獲取最佳聚類數(shù)目，克服了前述譜聚類過分割缺陷，同時這一優(yōu)勢使得本算法可以從光流軌跡中發(fā)現(xiàn)真正擺桿數(shù)目，在機械故障診斷中（如檢測因疲勞而停止擺動的賈卡針數(shù)目）有一定實際意義。

圖17 回歸軌跡代理的SSC 的相似度矩陣Fig.17 Affinity matrix of the regression surrogate SSC

圖18 回歸軌跡代理的SSC 低維嵌入分布圖Fig.18 Low-dimensional embedding of RS-SSC

7 角位移曲線的測量及比較

除了類內(nèi)一致性性能，下面從聚類稠密度提高的優(yōu)勢說明經(jīng)過代理回歸后SSC 聚類相比直接稀疏子空間聚類可提高后續(xù)擺件分割點云配準(zhǔn)得到的角位移精度。

7.1 非種子軌跡的分割及點云提取

代理回歸SSC 得到各幀種子軌跡聚類標(biāo)簽后，將種子軌跡標(biāo)簽在超流素內(nèi)部擴散至回歸前的所有軌跡樣本，使得運動分割后各種分割點云稠密度達到光流采樣的密度. 在第125 幀將擴散后非種子軌跡的像素點云及所屬聚類標(biāo)簽用顏色標(biāo)出如圖19 所示（彩圖見期刊電子版）。

圖19 回歸軌跡代理的SSC 聚類運動分割結(jié)果Fig.19 Motion segmentation results of the RS-SSC

7.2 最大條件期望點云配準(zhǔn)及角位移測量

由前述可以看出，因為光流阻斷，各幀分割的運動擺桿點云個數(shù)不同，不能直接采用ICP 等配準(zhǔn)算法求位姿，故這里采用論文［27］基于條件期望的最大似然參數(shù)估計來求點云間的位姿變換矩陣，為減少匹配離群點個數(shù)，這里取個數(shù)最大的點云為模型數(shù)據(jù)，表示為X={Xi}1≤i≤n，而待匹配各幀點云為Y={Yj}1≤j≤m，模型經(jīng)過位姿參數(shù)配準(zhǔn)后的點云為：μ(Xi；Θ)=RXi+t，Θ：={R，t}。假定數(shù)據(jù)點云每個點以其匹配的相應(yīng)模型點云中的點為中心的高斯分布，則配準(zhǔn)問題轉(zhuǎn)化為以匹配關(guān)系矩陣Z為隱含變量并以觀察到的數(shù)據(jù)為條件的GMM 最大似然參數(shù)估計問題，被估計的參數(shù)為Ψ={Θ，Σ1，…，Σn}，其中Σn為各點分布協(xié)方差，Θ：={R，t}為位姿參數(shù)。

估計的目標(biāo)函數(shù)即條件數(shù)據(jù)最大似然通過隱含變量Z可分解為：

其中，對應(yīng)各數(shù)據(jù)的隱含變量矩陣為Z={Zj}，1≤j≤m，矩陣元素取值為該點所匹配對應(yīng)模型點的序號，可定義該隨機變量取值分布如下：

αji=P(Zj=i∣Yj)，引入克羅內(nèi)克離散取值定義后可將上述隨機變量取值視為條件期望：

從而可以采用E-M 架構(gòu)對目標(biāo)函數(shù)求解，其中E 步為求解該匹配關(guān)系隱含變量的后驗估計如式（31）：

其中?3D為失配離群點集，定義如下：

其次在M 步分別對另外兩個運動參數(shù)估計如下：

為保證全局最優(yōu)，求解前，先對擺桿點云進行RANSAC 擬合得到的轉(zhuǎn)角結(jié)合平移分量為零構(gòu)成的旋轉(zhuǎn)矩陣作為參數(shù)R初值。

7.3 角位移測量誤差比較

為方便與其他兩種算法比較，將視頻分割為4 段長度為35 的視頻共140 幀數(shù)據(jù)，分別進行SSC 運動分割以保證每段有一定數(shù)量的種子軌跡. 為簡單起見，點云所在平面為XOY平面，Z坐標(biāo)為0，采用右手系確定Z軸方向后，通過反正切公式提取旋轉(zhuǎn)矩陣R的Z軸歐拉角為擺桿擺角測量值，基準(zhǔn)值采用人工標(biāo)定，標(biāo)定方法為在matlab2017 里figure 對象的鼠標(biāo)點擊回調(diào)函數(shù)里實現(xiàn)了取點及對雨刮根部主擺杠RANSAC 擬合與反正切求角，最后可得三種算法各幀擺角測量結(jié)果如圖20 所示（彩圖見期刊電子版）。圖中不同顏色點為SSC，SC，RS-SSC 各算法的測量值，曲線為測量值相應(yīng)的高斯回歸平均值，圖中SSC，SC，RS-SSC 對標(biāo)定回歸曲線的均方差分別為16.96，11.71，8.79。因此，RS-SSC 最接近標(biāo)定值，這是因為SSC 只能選長度大于35 的光流，而RS-SSC 可以將回歸前超流素所包含的所有的被阻斷的光流軌跡經(jīng)過標(biāo)簽擴散后一同參與點云配準(zhǔn)，其像素的稠密度大于前者，而SC 在擺桿接近處軌跡會出現(xiàn)混疊導(dǎo)致部分異常點降低了一部分精度。

圖20 SC，SSC，RS-SSC 各算法的擺角曲線檢測結(jié)果Fig.20 Angle curve detection of SC，SSC，RS-SSC

為了在工程現(xiàn)場考量算法的測量精度，結(jié)合客運車輛的機器視覺自動安檢項目，在現(xiàn)場采集了6 種不同環(huán)境亮度同一車型客車雨刮擺桿的5組角位移曲線作為測量對象，將人為標(biāo)定回歸值作為基準(zhǔn)則測量誤差平均值如下：

由表1 可知，隨著環(huán)境照度的提高，角位移檢測精確率均有所提高。這是因為環(huán)境亮度的提高增加了擺桿前景對比度，這樣減少了軌跡中斷現(xiàn)象使得種子軌跡數(shù)目增加，同時因前景特征點結(jié)構(gòu)張量增加也增加了總體軌跡的數(shù)目與稠密度?？傮w上看，當(dāng)環(huán)境亮度大于700 lux 時，代理回歸的SSC 的擺角測量誤差不超過10%，工程上可以作為此安檢項目中刮桿卡扣及其驅(qū)動電機的速度是否正常無故障的判斷依據(jù)。

表1 各算法在不同照度下的擺角檢測誤差Tab.1 Angle detection error under different algorithm and illuminance

8 結(jié) 論

本文在利用圖像視覺對單色的機械擺件的運動參數(shù)測量中，針對當(dāng)前兩類運動分割的缺陷，提出了一種通過超流素回歸平均軌跡代理的SSC 聚類的運動分割算法。該算法因采用軌跡自表達系數(shù)構(gòu)成的L1 范數(shù)圖取代傳統(tǒng)的變分光流的時空相似度拉式譜圖進行聚類，從而克服了選擇超參數(shù)ε、K的困難并以無監(jiān)督機器學(xué)習(xí)方式通過子空間結(jié)構(gòu)的信息挖掘出運動物體數(shù)目K，為運動分割用于機械擺動故障檢測提供了算法依據(jù)。

另一方面，本算法利用高斯泛函回歸學(xué)習(xí)出因背景復(fù)雜而被阻斷的軌跡并代理超流素軌跡群一次性完成SSC 運動分割，相比不回歸直接分時間片段的SSC 算法降低了運算量，增加了子空間聚類種子樣本的稠密度，結(jié)合SLAM 點云配準(zhǔn)技術(shù)提高了逐幀角位移測量精度. 所提出的算法提供了一個將圖像學(xué)與圖形學(xué)有機結(jié)合的機械運動測量框架，具有廣闊的研究空間及工業(yè)機器視覺測量應(yīng)用價值。