結合RANSAC算法與傳統相似變換模型的穩定性分析

劉忠賀，李宗春，郭迎鋼，何 華，趙文斌

(1.信息工程大學，河南 鄭州 450001；2.中國科學院上海高等研究院，上海 201210)

加速器裝置運營期間，磁鐵等關鍵元件會受到外界諸多因素影響發生位移，如不能及時發現并做出調整，則裝置運轉性能會受到影響，嚴重時會導致粒子運行軌道發生偏轉[1-2]。元件是否發生變形只能依靠控制網對其進行監測和調整，所以定期對控制網進行復測并對網點進行穩定性分析是一項非常重要的工作[3]。

德國測量學者Pelzer提出的平均間隙法[4]利用統計檢驗的方法對控制網幾何圖形作整體檢驗，判斷控制網在兩期觀測間是否發生顯著性變化，依次尋找變形點并將其剔除，直到檢驗通過[5-6]。基于平差后的兩期坐標之間存在一定的平移、旋轉、尺度信息，可采用傳統相似變換(Traditional Similarity Transformation, TST)法計算兩期控制網之間的變換參數，進而求出變形點的位移信息，但易受到位移較大點的影響，計算的模型參數不準確，穩定性分析結果不可靠。基于此，陳永奇[7]提出了迭代加權相似變換法(Iterative Weighted Similarity Transformation, IWST)，以所有點位移的一次范數最小為目標函數[8]，采用變權迭代的方式來計算點的位移量，在TST法基礎上，增加對變形點的抗差性，提高穩定點判別的正確率。Duchnowski Robert[9]將基于中值的穩健估計理論應用到控制網的穩定性分析中，在一定程度上抵抗了粗差的影響。Krzysztof Nowel[10]提出觀測值差異分析變形的穩健估計方法，是通過消除觀測數據中存在的粗差來提高穩定性分析的可靠性。A.R.Amiri-Simkooei等人[11]同時將兩期觀測數據進行平差解算，將每一個點當作變形點，其變形量設為待求參數進行求解，依次剔除位移最大的點，然后進行整體檢驗，直到檢驗通過。該方法計算繁瑣，類似于單點檢驗方法，并且有時會出現最大變形點判斷錯誤的情況。郭迎鋼[12]在加速器隧道控制網穩定性分析過程中加入了網點的空間靈敏度橢球，并分析了控制網的可監測性，結合整體穩定性檢驗和平均間隙法對隧道控制網進行了穩定性分析，但在求解網點變形量時未顧及較大變形點對變換參數的影響。

TST模型在求解相似變換參數時易受到較大變形點的影響，不具有抗差性，若能通過一定的方法將變形較大的點剔除，利用剩余較為穩定的點求解變換參數，則可以得到更加準確的網點變形信息。隨機抽樣一致性(Random Sample Consensus，RANSAC)算法[13]可用于數據的提純，減少或避免異常點對模型的干擾[14-17]。文中結合RANSAC算法與TST模型，使該模型在應用時具備抗差性，并通過試驗進行驗證。

1 函數模型

1.1 TST模型

兩期控制網之間可通過公共點求解轉換參數，七參數坐標轉換模型可表示為：

(1)

式中:(X，Y，Z)與(x，y，z)為同一點在兩期觀測下的坐標；兩期觀測之間繞X，Y，Z軸的旋轉角分別為ω，φ，κ；(Δx，Δy，Δz)為平移量；m為尺度。兩期觀測間旋轉矩陣R為：

(2)

將式(1)線性化，得：

(3)

式中：(vX,vY,vZ)為殘差；(X0，Y0，Z0)為第二期坐標系下的坐標(x，y，z)根據兩期之間轉換參數初值轉換至第一期坐標系下坐標。

TST模型可表示為：

(4)

式中：H表示系數矩陣；t=[dmdωdφdκdΔxdΔydΔz]T；d為位移向量。

1.2 IWST模型

由于TST模型在求解變換參數時容易受到變形點的影響，所以在其基礎上，陳永奇提出了IWST模型[7]。在求解變換參數時，以網點位移向量d的一次范數最小為準則，根據位移的大小不斷調整權值，通過迭代求解最終的位移。目標函數可表示為：

∑|di|=min(i=1,2,…,3n).

(5)

初始計算時，取權陣P(k=1)=I為單位陣，通過TST模型計算出位移向量d(k=1)。

對于第k+1次迭代，第i點權陣為：

(6)

為了避免權值無限大，分母上加了一個微小量ε。迭代過程如以下形式：

(7)

其中

S(k)=I-H(HTP(k)H)HTP(k).

(8)

式中：ΔX為兩期坐標差；Q11，Q22分別為第一期和第二期平差后的協因數陣。

當前后兩次迭代后的位移量之差小于某一閾值μ時，停止迭代，即:

|d(k+1)-d(k)|<μ.

(9)

最后，根據最終位移量，利用T檢驗法對網點的穩定性進行判斷。

1.3 文中方法

1.3.1RANSAC算法思想

RANSAC算法是一種通過迭代計算過程來尋求樣本中滿足某一正確模型對應的最大內點集，并利用該集合中樣本重新估算模型的過程[18]。基本過程：

1)給定一個數據集，從中隨機抽取t個樣本并計算待求參數模型初值；

2)設置閾值ε，利用求出的參數模型初值驗證集合中剩余樣本的模型誤差是否超限。閾值范圍內的樣本歸為內點；否則，歸為外點；

3)統計內點集中樣本數量，再次隨機抽取t個樣本，重復上述步驟；

4)經過r次迭代次數后，選取迭代過程中的最大內點集，判斷該集合樣本數量是否達到k值，若滿足條件，利用該集合計算參數模型并篩選內點集，計算最終的模型；否則，算法失敗。

上述過程涉及到的參數根據實際情況設定，代表含義如下：

t：計算模型參數所需最少樣本數量；

ε：決定一致集樣本數量多少；

k：表征正確模型對應一致集的最少樣本數；

r：與外點所占樣本集比例有關，決定模型的可靠程度。

1.3.2 穩定性分析流程

最小二乘估計是一種無偏估計，在無系統誤差、粗差的情況下，可得到最優化參數估計結果，且計算過程簡便。然而，該估計方法對粗差非常敏感，當觀測數據中存在較大變形點時，變形量極易發生偏移、轉移或分配到其他點的坐標信息中，對結果造成影響。為了得到準確的相似變換參數，必須將控制網中的變形點予以剔除，使其不參與變換參數的計算，通過穩定的點來估計正確的轉換模型，在此基礎上分析網點的穩定性，求解其絕對位移量。前面已經介紹了利用RANSAC算法提高坐標系之間轉換參數精度的方法，將其應用到控制網的穩定性分析中來，對TST方法進行改進。基本過程如下：

1)設隧道控制網兩期觀測誤差方程為：

(10)

采用重心基準，分別求解兩期觀測坐標X1,X2；單位權中誤差σ1,σ2；協因數陣Q11,Q22。

2)經過平差后的兩期坐標之間存在一定的平移、旋轉及尺度縮放關系，通過公共點轉換的方法對各參數進行求解，計算過程中采用RANSAC算法篩選相對穩定點，具體流程見圖1。

篩選流程主要分為兩步，首先通過迭代篩選出最大內點集，然后利用該內點集重新估算參數模型，迭代篩選最終的穩定點集。

3)利用篩選出的公共點求解變換參數，將所有點代入該參數模型，計算網點的位移量。

4)通過T檢驗法對每一個點進行穩定性判斷，構造統計量。

(11)

(12)

式中：f1,f2為兩期觀測自由度。

選擇一定顯著性水平α，對網點穩定性進行判斷。

Ti≤Fα(wi,f1+f2).

(13)

滿足式(3)的點即為穩定點，反之，則為不穩定點。

TST的改進之處在于第2步求解轉換參數過程中引入了RANSAC算法，使得TST模型在應用時具有穩健性，通過篩選并剔除某些變形點，進而提高轉換參數的精度，使得兩期坐標間的轉換關系更加準確可靠。

2 試驗與分析

2.1 模擬試驗

為了驗證文中方法的可行性，以粒子加速器隧道控制網為例，設計了模擬試驗。

2.1.1 數據

在4 m×2.5 m×4 m范圍內模擬了一組數據，

圖1 相對穩定點篩選流程

包含S1～S6共6個測站以及A1～A16共16個控制點，為了使測量條件更加貼近實際，點位分布均勻且測站完全在控制點所包含的空間內部，網型及試驗場景如圖2所示。

圖2 試驗場景

理論數據見表1，并在各點的坐標分量上加入方向、大小不同的位移量，作為第二期理論數據。

根據表1，對A1,A4,A7,A9,A11,A13,A15共7個點的某一坐標分量加入了位移，考慮到5 m范圍內，點位測量精度不超過0.05 mm，所以加入的位移量最小為0.05 mm，最大為3 mm，變形點占比43.75%。兩期數據的理論坐標設計完成后，按照測角精度±(15+6×10-6·D)μm，測距精度為±0.5·Dμm/m的標稱精度向數據添加隨機誤差，作為兩期對比試驗數據。

2.1.2 試驗與分析

分別采用TST模型、IWST模型及文中方法對兩期模擬數據進行處理，分析控制網的變形及穩定情況。計算過程中，取顯著性水平α=0.05，在對各點進行T檢驗時，F0.05(3,420)=2.626 1，3種方法分析結果見表2。

表1 網點及測站設計坐標 mm

由表2可知，在16個點組成的控制網中，TST模型受到變形點的影響，最終穩定性分析效果較差，準確率達到50%。結合RANSAC算法與TST模型，最終結果接近實際情況，分析變形點結果準確率接近90%，誤判了A5點。對于IWST模型，誤判了兩個點，其分析結果劣于文中方法。為了檢驗3種方法求解變換參數的準確性，對比其效果優劣，將3種方法計算的變形點位移同理論值進行對比，見表3。

表2 不同模型網點穩定性分析結果對比

表3 變形點位移量與理論值對比 mm

從表3可以發現，從整體效果來看，文中方法所計算的變形點位移量最接近真實變形情況。對于所有變形點，文中方法計算的A1點位移與理論值相差最大，達到0.025 mm，分析其原因，可能是在模擬數據時，添加的隨機誤差最大不超過0.05 mm，而在該點上隨機誤差大小接近0.025 mm，文中方法雖然是一種穩健估計方法，但作用對象是平差后的點坐標，對于平差后坐標本身所帶有的誤差并沒有消除作用。然而，剩余變形點變形量均與理論變形值相差甚小。

2.2 實例分析

2.2.1 數據來源

文中采用2017年與2019年的SSRF復測數據，對兩期控制網同名點進行穩定性分析。采用重心基準，分別對兩期觀測數據進行平差，結果見表4。

表4 平差結果

表4中自由度表示多余觀測個數。經查詢，兩期測量數據中同名點數量共65個，平差后的同名點分布如圖3所示。

圖3 同名點分布圖

從圖3可以看出，同名點為儲存環部分網點，以這些點為基礎，采用不同方法分析其穩定性。

2.2.2 試驗與分析

利用TST模型、IWST模型及文中方法對65個同名點進行處理。TST模型和IWST模型直接對所有點進行處理，文中方法處理的基本過程為：①篩選同名點中部分相對穩定點，將距離閾值設為0.3 mm，最終得到35個相對穩定點；②利用這些點計算變換參數；③計算剩余點在此變換參數下的位移量。計算時，取顯著性水平α=0.05，3種方法分析結果見表5。

表5 不同方法穩定性分析結果

從表5可以發現，3種方法分析結果略有不同。3種方法均得出BD512點位變形點，且位移量基本相等，在3 mm左右；相比IWST模型，TST模型及文中方法分析結果中，均多出一個變形點，且不是同一點，其位移量較小。3種方案處理后的65個同名點變形量如圖4所示。

圖4 不同模型所得變形量

從圖4直觀地發現，3種方案處理結果均發現第42號點(BD512)發生了變形，且十分明顯，而其余點變形相對較小。除BD512點外，1號點(BD303)變形最大，在2.5 mm左右，向后依次減小，趨于平穩，3種方法所得結果一致，說明控制網各網點之間相對穩定。為了進一步驗證3種方法處理結果的正確性，利用SA軟件對兩期數據進行處理，將兩期數據中的同名點進行最佳擬合變換，得到的變形量見圖5。

圖5 最佳擬合變換

從圖5中可以看出，第42號點發生了較大變形，與以上3種方法處理結果一致。綜合幾種方法，該點可能發生變形，需要到實地勘測變形情況，確定其是否真正發生變形，對該點坐標信息進行更新并分析變形原因。

模擬試驗雖然初步驗證文中方法的有效性，但由于條件受限，加上TST模型及文中方法分析結果中所多出的變形點，對于實例分析部分，并未對相應變形點進行實地考察，不能確定點位是否真正發生了變形，只能從數據分析結果來判斷，是文中的不足之處。

3 結束語

以粒子加速器隧道控制網為例，通過模擬試驗與實例分析，驗證文中方法在隧道控制網穩定性分析中有效性(亦可推廣至其他場景)，并得出以下3點結論。

1)結合RANSAC算法與TST模型，使得該模型在計算相似變換參數時具備穩健性，能夠避免較大變形點對變換參數的影響；

2)同TST模型及IWST模型相比，文中方法穩定點判別準確率最高，且計算得到的位移量與實際情況更加相符；

3)文中方法在TST模型基礎上提高了坐標轉換參數的精度和可靠性，在判斷控制網穩定性方面更具可行性。