上官晉太，連瑋

（長治學(xué)院 計(jì)算機(jī)系，山西 長治 046011）

0 引言

點(diǎn)匹配方法是非剛體圖像配準(zhǔn)中的一個(gè)最基本的方法，廣泛地應(yīng)用于計(jì)算機(jī)視覺、圖像分析和模式識別的過程中［1-6］。比如，圖像融合［7］就是典型的點(diǎn)匹配問題。與其他基于特征的匹配方法相比，點(diǎn)匹配方法更加簡單，且點(diǎn)匹配方法也是其他基于特征的匹配方法的基礎(chǔ)，因此，其在非剛體圖像配準(zhǔn)中受到了越來越多的重視。非剛體點(diǎn)匹配涉及兩個(gè)基本的過程，第一個(gè)是兩個(gè)點(diǎn)集中點(diǎn)與點(diǎn)的對應(yīng)關(guān)系的確定，第二個(gè)是點(diǎn)集之間的非剛體變換。這兩個(gè)問題不是互相獨(dú)立的，它們之間相互影響，一旦其中之一被確定了，對另一個(gè)的求解就是非常簡單的事情。對于點(diǎn)匹配問題的解決，一般有兩種途徑，一種是將對應(yīng)關(guān)系和空間變換分別處理，只解決其中之一，另一種是將兩個(gè)問題聯(lián)合處理，我們可以將點(diǎn)之間的對應(yīng)關(guān)系和空間變換看成是兩個(gè)變量，用迭代的方法先保持一個(gè)變量不變，估計(jì)另一個(gè)變量，然后在得到估計(jì)量的基礎(chǔ)上，保持這個(gè)估計(jì)量不變，對先前不變的那個(gè)變量進(jìn)行新的估計(jì)，在交替的過程中兩個(gè)變量相互改進(jìn)，直到得到最終的優(yōu)化解。因此，我們可以將第一種途徑稱之為獨(dú)立估計(jì)法，第二種稱之為聯(lián)合估計(jì)法。到目前為止，在已發(fā)表的文獻(xiàn)中大部分只專注于解決剛體變換問題，比如，已提出的慣性矩法［8］，此法先求出數(shù)據(jù)的“質(zhì)心”和“主軸”，“質(zhì)心”可以用來解決空間變換中的位移問題，而“主軸”提供數(shù)據(jù)集的整體的指向，由此可以解決數(shù)據(jù)集的旋轉(zhuǎn)角度問題。同樣，只用于解決剛體變換問題的方法還有Hough 變換［9］，Hausdauff 距離［10］等，這一類的方法屬于獨(dú)立估計(jì)法?，F(xiàn)實(shí)生活中，非剛體配準(zhǔn)的應(yīng)用場合更廣泛一些，比如醫(yī)學(xué)圖像配準(zhǔn)就是典型的非剛體圖像配準(zhǔn)［11-12］。非剛體的配準(zhǔn)，由于最優(yōu)解的搜索空間維度過大，求解過程過于復(fù)雜，所以一般利用聯(lián)合估計(jì)法。在聯(lián)合估計(jì)法中，兩個(gè)點(diǎn)集之間對應(yīng)關(guān)系的確定是解決問題的關(guān)鍵所在。常用的方法有形狀上下文方法和基于迭代最近點(diǎn)的方法。

形狀上下文算法是通過一個(gè)基于輪廓的形狀描述子，來表達(dá)點(diǎn)集中一個(gè)點(diǎn)和其余點(diǎn)的相對位置關(guān)系，點(diǎn)與點(diǎn)的對應(yīng)關(guān)系通過最小化形狀上下文距離之和來實(shí)現(xiàn)。在利用形狀上下文進(jìn)行非剛體點(diǎn)匹配時(shí)，賈耕云等［13］提出了一種具有對稱不變性的改進(jìn)形狀上下文特征提取與匹配的算法，該算法能夠有效地在互相對稱的相似圖像間建立匹配，從而使匹配更加穩(wěn)定。師碩等［14］針對人臉匹配在光照、姿態(tài)、表情等背景因素影響下匹配準(zhǔn)確率低的問題，提出了一種基于SURF 和形狀上下文的人臉匹配算法，此方法的應(yīng)用有效地增加了匹配點(diǎn)對數(shù)目，提高了人臉圖像匹配的準(zhǔn)確率?；烘I等［15］提出的改進(jìn)算法，引入雙局部二進(jìn)制模式紋理特征權(quán)重值，利用梯度信息使得該算法對噪聲、光照變化具有較高的魯棒性。這些方法均對形狀上下文匹配算法進(jìn)行了一定程度的改進(jìn)，但都未解決傳統(tǒng)算法容易產(chǎn)生點(diǎn)對之間的誤匹配問題，特別是當(dāng)采樣點(diǎn)在柵格交界處，即在相鄰兩個(gè)柵格邊界位置時(shí)，在有輕微形變的情況下，就會引起階躍變化，從而導(dǎo)致誤匹配。

另一種常用的點(diǎn)匹配的算法是迭代最近點(diǎn)算法（Iterative Closest Point，ICP）［16］，此算法通過在兩個(gè)點(diǎn)集間進(jìn)行雙向的最近點(diǎn)確認(rèn)，得到對應(yīng)矩陣，由對應(yīng)矩陣得到模板點(diǎn)集到目標(biāo)點(diǎn)集的映射位置，從而建立點(diǎn)對應(yīng)關(guān)系。石愛軍等［17］提出了遺傳算法結(jié)合自適應(yīng)閾值約束的ICP 算法，此算法可以提高配準(zhǔn)精度，同時(shí)減小配準(zhǔn)時(shí)間，能在一定程度上滿足實(shí)時(shí)性的要求。張春雷等［18］提出了一種將三點(diǎn)法與迭代最近點(diǎn)算法結(jié)合的配準(zhǔn)策略，其配準(zhǔn)精度和穩(wěn)定性均優(yōu)于傳統(tǒng)算法，且具備高效、易于操作的特點(diǎn)。以上這些改進(jìn)的迭代最近點(diǎn)算法在性能上都有一些提高，但是均無法解決迭代最近點(diǎn)算法易于收斂于局部最小值點(diǎn)的缺點(diǎn)。本文針對這一個(gè)問題，提出了在對應(yīng)矩陣上加適當(dāng)隨機(jī)擾動的方法，來產(chǎn)生更加合理的點(diǎn)對應(yīng)關(guān)系，從而可以避免非剛體點(diǎn)匹配過早地陷入局部最小值點(diǎn)。

1 兩個(gè)點(diǎn)集中點(diǎn)對應(yīng)關(guān)系的確定

在聯(lián)合估計(jì)法中，有一種方法最為常見，那就是基于ICP 的聯(lián)合估計(jì)法。它的原理就是先在兩個(gè)點(diǎn)集之間初步確定對應(yīng)關(guān)系，在對應(yīng)關(guān)系確定后進(jìn)行空間變換，空間變換后可得到新的估計(jì)點(diǎn)集，用此點(diǎn)集與目標(biāo)點(diǎn)集再一次進(jìn)行對應(yīng)關(guān)系的確定，這樣周而復(fù)始地進(jìn)行對應(yīng)關(guān)系和空間變換的輪流迭代，最后收斂于一個(gè)全局最優(yōu)點(diǎn)。圖1 是要配準(zhǔn)的兩個(gè)點(diǎn)集，分別為模板點(diǎn)集Q 和目標(biāo)點(diǎn)集P［3］。配準(zhǔn)的任務(wù)是將模板點(diǎn)集經(jīng)過空間變換匹配到目標(biāo)點(diǎn)集上去。

圖1 （a）模板點(diǎn)集Q 和（b）目標(biāo)點(diǎn)集PFig. 1 （a）Template point set Q and（b）Target point set P

首先，對模板點(diǎn)集中的每一個(gè)點(diǎn)，在目標(biāo)點(diǎn)集中找到離自己最近的點(diǎn)，這樣可以形成一個(gè)模板點(diǎn)集到目標(biāo)點(diǎn)集的對應(yīng)矩陣M1，矩陣中的行號為模板點(diǎn)集中對應(yīng)點(diǎn)的標(biāo)號，矩陣中的列號 為 目 標(biāo) 點(diǎn) 集 中 對 應(yīng) 點(diǎn) 的 標(biāo) 號 ，若M1(i，j)=1，說明模板點(diǎn)集中的第i點(diǎn)在目標(biāo)點(diǎn)集中找到的最近點(diǎn)為第j個(gè)點(diǎn)。再對目標(biāo)點(diǎn)集中的每個(gè)點(diǎn)，尋找在模板點(diǎn)集中的最近點(diǎn)，形成目標(biāo)點(diǎn)集到模板點(diǎn)集的對應(yīng)矩陣M2，M2中的行號和列號規(guī)定和M1中是相同的。將兩者進(jìn)行平均得到兩個(gè)點(diǎn)集之間的雙向?qū)?yīng)矩陣M。 有M=0.5(M1+M2)。 為 了 觀 察 和 敘 述方便，對圖1 中的兩個(gè)點(diǎn)集分別進(jìn)行等間隔抽樣，得到兩個(gè)均由8 個(gè)點(diǎn)形成的點(diǎn)集來說明雙向?qū)?yīng)矩陣M的形成過程。抽樣后形成的點(diǎn)集如圖2 所示。用ICP 法可以得到模板點(diǎn)集到目 標(biāo) 點(diǎn) 集 的 對 應(yīng) 矩 陣M1，其 中，除M1(1，5)、M1(2，4)、M1(3，3)、M1(4，7)、M1(5，8)、M1(6，8)、M1(7，8)和M1(8，8)為1 外 其 余 都 為0，說明模板點(diǎn)集中的第1 個(gè)點(diǎn)在目標(biāo)點(diǎn)集中找到的與自己最近的點(diǎn)為第5 點(diǎn)，以此類推，模板點(diǎn)集中的第8 點(diǎn)在目標(biāo)點(diǎn)集中找到的與自己最近的點(diǎn)也為第8 點(diǎn)。在這里可以看到，兩個(gè)點(diǎn)集之間并非一一對應(yīng)關(guān)系，模板點(diǎn)集中的5、6、7、8 這4 個(gè)點(diǎn)在目標(biāo)點(diǎn)集中的最近點(diǎn)都是第8點(diǎn)。這樣的多對一，或一對多的關(guān)系是非剛體點(diǎn)匹配中比較難解決的問題之一，這也會造成匹配困難。同樣的方法，求目標(biāo)點(diǎn)集到模板點(diǎn)集 的 對 應(yīng) 矩 陣M2，可 知，M2矩 陣 中 除M2(1，1)、M2(1，5)、M2(2，2)、M2(2，3)、M2(2，4)、M2(5，6)、M2(5，7) 和M2(5，8) 為1 外 其 余 都 為0，說明目標(biāo)點(diǎn)集中的第1、5 點(diǎn)在模板點(diǎn)集中找到的與自己最近的點(diǎn)都為第1 點(diǎn)。目標(biāo)點(diǎn)集中的2、3、4 點(diǎn)在模板點(diǎn)集中找到的最近的點(diǎn)都是第2 點(diǎn)。目標(biāo)點(diǎn)集中的6、7、8 點(diǎn)在模板點(diǎn)集中找到的最近點(diǎn)都是第5 點(diǎn)。從圖2（b）中可以看到兩個(gè)點(diǎn)集中的對應(yīng)關(guān)系。由此形成的雙向?qū)?yīng)矩陣M，再對M按行進(jìn)行歸一化處理可得M′，結(jié)果如下所示：

圖2 抽樣后兩個(gè)點(diǎn)集形成的對應(yīng)關(guān)系（a）抽樣后的兩個(gè)點(diǎn)集；（b）兩個(gè)點(diǎn)集之間的點(diǎn)對應(yīng)關(guān)系Fig. 2 Correspondence between two sampled point sets（a）Two sampled point sets；（b）The correspondence between two point sets

由歸一化后的雙向?qū)?yīng)矩陣M′，可以求出變換后的新的點(diǎn)集Q′。變換公式如下：

如果直接利用M′矩陣進(jìn)行計(jì)算，所得到的點(diǎn)將會出現(xiàn)歸并現(xiàn)象，點(diǎn)數(shù)將會減少，這樣就無法在兩個(gè)點(diǎn)集之間建立起一一對應(yīng)的關(guān)系，匹配的結(jié)果很可能會陷入局部極小值，無法得到全局最優(yōu)解。通過觀察矩陣M′可以找出產(chǎn)生點(diǎn)歸并現(xiàn)象的原因，在此矩陣中M′(6，8) 、M′(7，8)和M′(8，8)都為1，這使得通過公式（1）產(chǎn)生的映射會出現(xiàn)3 點(diǎn)并1 點(diǎn)的結(jié)果。為了避免這樣的情況發(fā)生，可以在矩陣M上加上一個(gè)隨機(jī)擾動。圖3 為不加擾動與加擾動的對比情況。

從圖3 的（a）中可以看到，如果在對應(yīng)矩陣M上不加擾動，利用歸一化后的矩陣M′進(jìn)行變換，得到的新點(diǎn)集Q′只有6 個(gè)點(diǎn)，用這6 個(gè)點(diǎn)與目標(biāo)點(diǎn)集的8個(gè)點(diǎn)進(jìn)行匹配是不合理的。在圖3 的（b）中加上了隨機(jī)擾動，這樣可以恢復(fù)到8 個(gè)對應(yīng)點(diǎn)的狀態(tài)。這樣做有利于匹配結(jié)果收斂于全局最優(yōu)點(diǎn)。

圖3 對應(yīng)矩陣M 加隨機(jī)擾動與不加擾動情況對比（a）沒有加隨機(jī)擾動的情況；（b）加隨機(jī)擾動的情況Fig. 3 Comparison between the situation with and without random perturbations on the corresponding matrix M（a）The case without perturbation；（b）The case with perturbation

2 空間變換

非剛體匹配中，最常用的空間變換是薄板樣條插值法。此方法源自數(shù)據(jù)插值和函數(shù)逼近理論［19］，它尋找一個(gè)通過所有的控制點(diǎn)的彎曲最小的光滑曲面。設(shè)：

其中rij=|pi－pj|代表兩點(diǎn)之間的距離。再定義如下矩陣：

求得w和a后，已知點(diǎn)坐標(biāo)(x，y)可求得變換后的坐標(biāo)(x′，y′)，二者之間關(guān)系為：

3 配準(zhǔn)實(shí)驗(yàn)及結(jié)果

用圖1 中的兩個(gè)點(diǎn)集進(jìn)行匹配實(shí)驗(yàn)，對應(yīng)關(guān)系矩陣的求解和空間變換交替進(jìn)行，每一步的迭代都是向最終解的逼近，經(jīng)過多次的運(yùn)算，最后得到一個(gè)全局最優(yōu)解。圖4 是匹配過程及結(jié)果，在圖（a）中可以看到兩個(gè)點(diǎn)集之間的對應(yīng)關(guān)系，在目標(biāo)點(diǎn)集曲率大的地方，出現(xiàn)了多對一或者一對多的情況。圖（b）是最終匹配結(jié)果，除曲線的頂點(diǎn)處和末尾處外，其余部分都得到了較好的配準(zhǔn)。圖（c）是空間變換情況，顯示了原點(diǎn)集的網(wǎng)格圖（虛線部分）和變換后的網(wǎng)格圖（實(shí)線部分）之間的對應(yīng)關(guān)系。

圖4 匹配實(shí)驗(yàn)與結(jié)果（a）點(diǎn)集之間的對應(yīng)關(guān)系；（b）配準(zhǔn)結(jié)果；（c）空間變換Fig.4 Registration test and result（a）Correspondence；（b）Registration result；（c）Transformation

不同的隨機(jī)擾動會對最后匹配結(jié)果產(chǎn)生不同的影響。圖5 是在雙向?qū)?yīng)矩陣M上加不同程度隨機(jī)擾動后的匹配結(jié)果?？梢钥吹剑簣D（a）為沒有加隨機(jī)擾動的情況，模板點(diǎn)集的點(diǎn)普遍下移，造成兩個(gè)點(diǎn)集在頂點(diǎn)部分沒有得到很好的配準(zhǔn)，圖（b）是在雙向?qū)?yīng)矩陣M上加0.26 倍的均值為0、方差為1的正態(tài)分布的隨機(jī)擾動時(shí)的情況，在整個(gè)范圍內(nèi)，模板點(diǎn)集和目標(biāo)點(diǎn)集得到了比較好的匹配。當(dāng)所加的隨機(jī)擾動進(jìn)一步增大到0.6 倍時(shí)，匹配結(jié)果如圖（c）所示，此時(shí)的匹配結(jié)果，甚至比不加隨機(jī)擾動時(shí)還要差，這是因?yàn)檫^量的隨機(jī)擾動，造成了對應(yīng)關(guān)系的模糊，由這種模糊引起的對應(yīng)的不確定性會引起匹配效果的惡變。

圖5 不同程度隨機(jī)擾動的匹配結(jié)果對比（a）無隨機(jī)擾動；（b）隨機(jī)擾動系數(shù)為0.26；（c）隨機(jī)擾動系數(shù)為0.6Fig.5 Comparison of registration results with different degree of perturbation（a）Without any perturbation；（b）The perturbation coefficient is 0.26；（c）The perturbation coefficient is 0.6

為了定量地研究在對應(yīng)矩陣上加不同程度的隨機(jī)擾動對匹配結(jié)果的影響，本文中先定義配準(zhǔn)誤差如下：

其中，dxi為變形后的模板點(diǎn)集Q′中的每個(gè)點(diǎn)與其在目標(biāo)點(diǎn)集P中找到的最近點(diǎn)之間的距離，dyj為目標(biāo)點(diǎn)集P中的每個(gè)點(diǎn)與其在模板點(diǎn)集Q′中找到的最近點(diǎn)之間的距離。R和S分別為模板點(diǎn)集和目標(biāo)點(diǎn)集的點(diǎn)數(shù)。我們以步長為0.1，范圍從0 到1 取11 個(gè)值，作為加在對應(yīng)矩陣M上的均值為0、方差為1 的正態(tài)分布隨機(jī)擾動的系數(shù)，每一種情況下運(yùn)行30次，最后計(jì)算配準(zhǔn)誤差的均值和方差，結(jié)果如表1所示。

由表中數(shù)據(jù)可知，當(dāng)沒有在對應(yīng)矩陣M上加隨機(jī)擾動時(shí)（此時(shí)擾動系數(shù)為0），配準(zhǔn)的結(jié)果是固定不變的，每次配準(zhǔn)的誤差都為0.019 3，此時(shí)方差自然為0。對于隨機(jī)擾動系數(shù)不為0 的情況，其配準(zhǔn)結(jié)果是一個(gè)概率分布問題，從表1 中可以看出，配準(zhǔn)誤差的均值最小值點(diǎn)在0.2 到0.3 之間，當(dāng)隨機(jī)擾動系數(shù)大于0.5 以后，配準(zhǔn)誤差急劇增加，這是因?yàn)檫^大的隨機(jī)擾動掩蓋了點(diǎn)之間對應(yīng)關(guān)系的本質(zhì)特性。為了進(jìn)一步找到最小值點(diǎn)的位置，再取步長為0.01，范圍從0.21 到0.3 取10 個(gè)值，作為隨機(jī)擾動的系數(shù)，每一種情況下同樣運(yùn)行30 次，最后計(jì)算配準(zhǔn)誤差的均值和方差，結(jié)果如表2 所示。從表中可以看到，當(dāng)所加正態(tài)分布的隨機(jī)擾動系數(shù)為0.26時(shí)，配準(zhǔn)誤差的均值最小為0.014 4，說明此時(shí)的配準(zhǔn)效果最好，其結(jié)果也好于沒有加隨機(jī)擾動的情況。

表1 擾動系數(shù)從0到1，步長為0.1時(shí)配準(zhǔn)誤差的情況Table 1 The registration error when the perturbation coefficient ranges from 0 to 1 and the step size is 0.1

表2 擾動系數(shù)從0.21到0.3，步長為0.01時(shí)配準(zhǔn)誤差的情況Table 2 Registration error when the perturbation coefficient ranges from 0.21 to 0.3 and the step size is 0.01

為了討論此方法對不同數(shù)據(jù)圖集的適用性，本文選擇不同的點(diǎn)集進(jìn)行配準(zhǔn)實(shí)驗(yàn)，圖6 是本文用于配準(zhǔn)實(shí)驗(yàn)的點(diǎn)集，分別為點(diǎn)集1 到點(diǎn)集8。圖中用“○”和“+”分別代表模板點(diǎn)集和目標(biāo)點(diǎn)集。

圖6 用于配準(zhǔn)實(shí)驗(yàn)的點(diǎn)集Fig. 6 Point sets used for registration experiments

和前面一樣，采取由粗到細(xì)的方法，先以步長為0.1，從0 到1 選取隨機(jī)擾動系數(shù)，每一種情況下運(yùn)行30 次，最后觀察其配準(zhǔn)誤差的均值和方差，在可能最小值點(diǎn)附近，再以步長為0.01 取10 個(gè)值，作為隨機(jī)擾動系數(shù)進(jìn)行進(jìn)一步配準(zhǔn)實(shí)驗(yàn)，以配準(zhǔn)誤差最小為原則，兼顧考慮方差影響，找到對于此點(diǎn)集的最優(yōu)隨機(jī)擾動系數(shù)。對于不同的點(diǎn)集其結(jié)果如表3 所示。從表中可以看到，對于不同的點(diǎn)集最“合適”的擾動系數(shù)并不相同，這與點(diǎn)集組成的復(fù)雜程度、模板點(diǎn)集和目標(biāo)點(diǎn)集之間的相對位置及二者之間的形變有關(guān)。

表3 不同點(diǎn)集的配準(zhǔn)實(shí)驗(yàn)結(jié)果Table 3 Registration experiment results of different sets

為了和沒有加隨機(jī)擾動時(shí)的情況進(jìn)行比較，表4 中列出了圖6 中的8 個(gè)點(diǎn)集，在隨機(jī)擾動系數(shù)為0時(shí)的配準(zhǔn)誤差。對比表3 和表4 可以看出，同序號的點(diǎn)集在其對應(yīng)矩陣上加上合理的隨機(jī)擾動，其配準(zhǔn)誤差小于無隨機(jī)擾動時(shí)的配準(zhǔn)誤差。表中同時(shí)列出了二者之比，這些比值的均值為0.551 6，即加上適當(dāng)?shù)碾S機(jī)擾動后，配準(zhǔn)誤差可以減小到原來的一半左右。這說明一般情況下，此算法對不同點(diǎn)集的配準(zhǔn)精度會有一定程度的提升效果。

表4 擾動系數(shù)為0時(shí)不同點(diǎn)集的配準(zhǔn)實(shí)驗(yàn)結(jié)果Table 4 Registration experiment results of different sets when the perturbation coefficient is 0

4 結(jié)論

在聯(lián)合估計(jì)法的非剛體點(diǎn)匹配過程中，我們應(yīng)用了對應(yīng)關(guān)系和空間變換交替迭代的方法來完成整個(gè)配準(zhǔn)過程。對應(yīng)關(guān)系和空間變換在這個(gè)迭代過程中相互改進(jìn)，目的是使最終的配準(zhǔn)結(jié)果能收斂于一個(gè)全局最優(yōu)點(diǎn)。在這個(gè)過程中，我們首先確定對應(yīng)關(guān)系，對應(yīng)關(guān)系的確定是下一步空間變換的基礎(chǔ)。由于兩個(gè)點(diǎn)集的點(diǎn)與點(diǎn)之間關(guān)系的復(fù)雜性會導(dǎo)致對應(yīng)關(guān)系中多對一或者一對多的情況出現(xiàn)，這就為下一步空間變換帶來了一定的難度。如果不進(jìn)行適當(dāng)?shù)母深A(yù)，最終的配準(zhǔn)會導(dǎo)致得到的是一個(gè)局部最小值點(diǎn)。我們的辦法是在得到的對應(yīng)矩陣M上加上一個(gè)適當(dāng)?shù)碾S機(jī)擾動，然后再對矩陣M按行進(jìn)行歸一化處理，得到矩陣M′，這樣得到的變換后的點(diǎn)集與目標(biāo)點(diǎn)集進(jìn)行配準(zhǔn)可以得到更好的結(jié)果。本文對多組已知點(diǎn)集加不同程度的隨機(jī)擾動，通過多輪仿真實(shí)驗(yàn)，結(jié)果表明此方法有效。