符艷軍,孫開(kāi)鋒
(1 西安外事學(xué)院,西安 710077; 2 西安精密機(jī)械研究所,西安 710075)
在基于視覺(jué)的飛行器導(dǎo)航與制導(dǎo)、遙感衛(wèi)星災(zāi)害監(jiān)控與環(huán)境監(jiān)測(cè)、醫(yī)學(xué)圖像分析等應(yīng)用中, 常常需要對(duì)不同成像傳感器獲取的異源圖像進(jìn)行匹配。
基于互信息測(cè)度的匹配方法是一種完全基于圖像灰度統(tǒng)計(jì)概率的匹配方法,不需要對(duì)原圖像間的灰度關(guān)系作任何假設(shè),非常適合于異源圖像的匹配,但互信息測(cè)度最大的缺陷是計(jì)算量大, 匹配耗時(shí)較長(zhǎng)[1]。為此,研究者們就圍繞減少互信息測(cè)度計(jì)算量以及改進(jìn)搜索策略等方面提出了相應(yīng)的加速方法。文獻(xiàn)[2-5]分別采用模擬退火、粒子群等優(yōu)化搜索算法提高匹配速度;文獻(xiàn)[6]采用灰度壓縮和粒子群優(yōu)化算法相結(jié)合加快匹配速度;文獻(xiàn)[7]通過(guò)提取多層次特征點(diǎn)減少了互信息測(cè)度的計(jì)算量;文獻(xiàn)[8-10]采用由粗到精的多層匹配方法對(duì)算法進(jìn)行加速;文獻(xiàn)[11-15]將分層匹配與各種優(yōu)化算法相結(jié)合提高匹配速度。
文中利用匹配過(guò)程中相鄰基準(zhǔn)子圖間灰度統(tǒng)計(jì)特征之間的相關(guān)性,通過(guò)差量法減少每一個(gè)匹配位置互信息的計(jì)算量來(lái)加快匹配速度。與已有的以灰度壓縮、特征提取或采取優(yōu)化搜索策略等快速互信息匹配方法相比,該方法既沒(méi)有減少參與匹配的像素?cái)?shù),也無(wú)需對(duì)圖像灰度等級(jí)進(jìn)行處理,不會(huì)對(duì)匹配精度造成影響。
圖像熵描述了圖像信源的平均信息量,反映了圖像灰度的統(tǒng)計(jì)信息,相似的兩幅圖像其圖像熵也相近。利用互信息進(jìn)行圖像匹配的實(shí)質(zhì)是:當(dāng)兩幅圖像在空間位置配準(zhǔn)時(shí),其重疊部分所對(duì)應(yīng)像素對(duì)的灰度互信息達(dá)到最大值。 兩幅圖像U和V的互信息I(U,V)可以用下式表示:
I(U,V)=H(U)+H(V)-H(U,V)
(1)
其中:H(U)、H(V) 分別表示圖像U及V的熵;H(U,V)表示兩幅圖像的聯(lián)合熵。熵通常由變量的概率密度來(lái)表示。針對(duì)灰度圖像,可以分別用兩幅圖像的灰度直方圖和聯(lián)合直方圖進(jìn)行估計(jì)。 計(jì)算公式分別為:
(2)
(3)
(4)
式(2)~式(4)中,hU(i)表示圖像U中灰度值為i的像素出現(xiàn)的次數(shù);hV(i)表示圖像V中灰度值為i的像素出現(xiàn)的次數(shù);hUV(i,j)表示同一像素位置上圖像U中灰度值為i、圖像V中灰度值為j的像素對(duì)出現(xiàn)的次數(shù);T為圖像的總像素?cái)?shù)。
研究發(fā)現(xiàn)[16],互信息本身的大小與待配準(zhǔn)兩圖像間的重疊度有一定的關(guān)聯(lián)性,為了消除這種關(guān)聯(lián)關(guān)系,文獻(xiàn)[17]提出了利用歸一化互信息(normalized mutual information, NMI)作為相似性測(cè)度。NMI能夠減少對(duì)圖像重疊部分的敏感性,實(shí)驗(yàn)表明它比標(biāo)準(zhǔn)的互信息方法更具魯棒性。NMI表達(dá)式為:
NMI(U,V)=(H(U)+H(V))/H(U,V)
(5)
在飛行器導(dǎo)航中,匹配的目的是確定飛行器的當(dāng)前位置,通過(guò)計(jì)算實(shí)測(cè)圖與基準(zhǔn)圖中每一個(gè)匹配位置對(duì)應(yīng)基準(zhǔn)子圖的互信息并找出最大互信息值對(duì)應(yīng)位置,即可確定基準(zhǔn)圖與實(shí)測(cè)圖的相對(duì)位置, 實(shí)現(xiàn)對(duì)飛行器的定位。
與基于特征的圖像匹配算法相比,基于最大歸一化互信息的圖像匹配最大缺點(diǎn)是計(jì)算速度慢。 每計(jì)算一次互信息測(cè)度,需要遍歷圖像中的所有像素,這使得計(jì)算量大大增加。由式(5)可以看出,歸一化互信息的計(jì)算量主要集中在兩幅圖像各自的熵及其聯(lián)合熵的計(jì)算。從式(2)~式(4)可以看出,圖像熵及聯(lián)合熵均為求和運(yùn)算,其計(jì)算量主要集中于每一個(gè)求和項(xiàng)的計(jì)算,因此,如果能減少每次匹配過(guò)程中需要計(jì)算的求和項(xiàng)的項(xiàng)數(shù), 則可以加快匹配速度。因?yàn)閳D像分布具有塊狀結(jié)構(gòu),在遍歷式搜索匹配過(guò)程中,各相鄰基準(zhǔn)子圖的灰度統(tǒng)計(jì)量具有很強(qiáng)的相關(guān)性,基準(zhǔn)圖像(i,j)位置處對(duì)應(yīng)子圖的灰度直方圖與其四鄰域位置處對(duì)應(yīng)子圖的灰度直方圖非常相近。
圖1 相鄰基準(zhǔn)子圖及其灰度直方圖
如圖1所示,圖1(a)為原始基準(zhǔn)圖像,圖1(b)為圖1(a)的一個(gè)子圖,圖1(c)為圖1(b)在原基準(zhǔn)圖上右移一個(gè)像素對(duì)應(yīng)的子圖,圖1(d)和圖1(e)分別為圖1(b)和圖1(c)的灰度直方圖。從圖中可以看出,相鄰兩個(gè)子圖的直方圖非常相似,也即兩相鄰子圖存在大量的灰度值其出現(xiàn)的頻率沒(méi)有發(fā)生變化,相應(yīng)的,在式(2)~式(4)的求和項(xiàng)中分別有很多項(xiàng)數(shù)沒(méi)有發(fā)生變化,因此, 當(dāng)前匹配位置處基準(zhǔn)子圖的熵以及與實(shí)測(cè)圖的聯(lián)合熵的計(jì)算可以前一匹配位置處基準(zhǔn)子圖的熵及與實(shí)測(cè)圖的聯(lián)合熵的計(jì)算為基準(zhǔn),通過(guò)比較相鄰兩個(gè)匹配位置處基準(zhǔn)子圖對(duì)應(yīng)的直方圖及聯(lián)合直方圖的變化情況,用差量法進(jìn)行計(jì)算。即在整個(gè)匹配過(guò)程中,除了要完整計(jì)算基準(zhǔn)圖上第一匹配位置(1,1)處對(duì)應(yīng)基準(zhǔn)子圖的熵以及與實(shí)測(cè)圖的聯(lián)合熵外,其它所有匹配位置的基準(zhǔn)子圖的熵及與實(shí)測(cè)圖的聯(lián)合熵的計(jì)算主要集中在與前一匹配位置相比灰度直方圖矩陣和聯(lián)合直方圖矩陣中發(fā)生變化的那些元素對(duì)應(yīng)的求和項(xiàng)上。據(jù)此, 文中互信息匹配算法的流程圖如圖2所示。具體步驟如下:
1)獲取實(shí)測(cè)圖像A和基準(zhǔn)圖像R,將兩幅圖像的灰度值調(diào)整到同一灰度區(qū)間。
圖2 文中匹配算法流程圖
2)按照式(2)或式(3)計(jì)算A的信息熵H(A)并存儲(chǔ)。
3)從基準(zhǔn)圖像R的左上角(1,1)點(diǎn)截取與實(shí)測(cè)圖像A大小相等的第1幅基準(zhǔn)子圖S1,統(tǒng)計(jì)S1各灰度值出現(xiàn)次數(shù)并存入初始灰度直方圖矩陣hS1,按照式(2)或式(3)計(jì)算S1的信息熵H(S1)并存儲(chǔ),最后將hS1和H(S1)分別賦值給hS0和H(S0)。
4)統(tǒng)計(jì)實(shí)測(cè)圖A與第1幅基準(zhǔn)子圖S1對(duì)應(yīng)像素位置的灰度對(duì)出現(xiàn)次數(shù)并存入初始聯(lián)合直方圖矩陣hAS1,按照式(4)計(jì)算初始聯(lián)合熵H(AS1)并保存,將hAS1和H(AS1)分別賦值給hAS0和H(AS0);再按照式(5)計(jì)算實(shí)測(cè)圖與基準(zhǔn)子圖的歸一化互信息值并存儲(chǔ)。
5)從基準(zhǔn)圖上(1,2)點(diǎn)位置開(kāi)始,逐行逐像素按照“Z”字形進(jìn)行迭代匹配。匹配過(guò)程中, 當(dāng)前匹配位置記為(i,j),當(dāng)j>1時(shí),當(dāng)前匹配位置的前一匹配位置應(yīng)為(i,j-1);當(dāng)j=1時(shí),則當(dāng)前匹配位置的前一匹配位置應(yīng)為(i-1,j),則當(dāng)前位置(i,j)處的歸一化互信息可按如下步驟計(jì)算:
①通過(guò)比較當(dāng)前匹配位置(i,j)處基準(zhǔn)子圖和前一匹配位置(i,j-1)或(i-1,j)處基準(zhǔn)子圖對(duì)應(yīng)列或行像素變化情況,按照差量法將前一基準(zhǔn)子圖的直方圖矩陣hS0更新為當(dāng)前子圖的直方圖矩陣hS。具體做法為: 在行方向匹配時(shí),當(dāng)前匹配位置的前一匹配位置為(i,j-1),則以hS0為基準(zhǔn),減去前一基準(zhǔn)子圖第一列對(duì)應(yīng)像素灰度值出現(xiàn)次數(shù),加上當(dāng)前基準(zhǔn)子圖最后一列對(duì)應(yīng)像素灰度值出現(xiàn)次數(shù),即得到hS;在列方向匹配時(shí),當(dāng)前匹配位置的前一匹配位置為(i-1,j),則以hS0為基準(zhǔn),減去前一基準(zhǔn)子圖第一行對(duì)應(yīng)像素灰度值出現(xiàn)次數(shù),加上當(dāng)前基準(zhǔn)子圖最后一行對(duì)應(yīng)像素灰度值出現(xiàn)次數(shù),即得到hS。
②找出hS與hS0的不同元素,分別組成矩陣hΔS和hΔS0,按照式(6)、式(7)分別計(jì)算hΔS和hΔS0對(duì)應(yīng)的信息熵ΔS及ΔS0;其中,normS0(S)為與hΔS0中元素對(duì)應(yīng)的前一基準(zhǔn)子圖熵值計(jì)算過(guò)程中的求和項(xiàng),因?yàn)檫@些求和項(xiàng)在前一匹配位置已經(jīng)計(jì)算過(guò),此次無(wú)需再重復(fù)計(jì)算,而只需計(jì)算這些項(xiàng)的和即可。
(6)
(7)
③以前一基準(zhǔn)子圖的信息熵H(S0)為基準(zhǔn),按照式(8)通過(guò)差量法將H(S0)更新為當(dāng)前基準(zhǔn)子圖的信息熵H(S)。
H(S)=H(S0)-ΔS0+ΔS
(8)
④計(jì)算實(shí)測(cè)圖與當(dāng)前基準(zhǔn)子圖的聯(lián)合直方圖矩陣hAS,找出hAS和前一基準(zhǔn)子圖聯(lián)合直方圖矩陣hAS0的不同元素,分別組成矩陣hΔAS和hΔAS0,按照式(8)、式(9)分別計(jì)算hΔAS和hΔAS0對(duì)應(yīng)的聯(lián)合熵ΔAS及ΔAS0;其中,normAS0為與hΔAS0中元素對(duì)應(yīng)的前一基準(zhǔn)子圖聯(lián)合熵計(jì)算過(guò)程中的求和項(xiàng),因?yàn)檫@些求和項(xiàng)在前一匹配位置已經(jīng)計(jì)算過(guò),此次無(wú)需再重復(fù)計(jì)算,而只需計(jì)算這些項(xiàng)的和即可。
(9)
(10)
⑤以前一基準(zhǔn)子圖的聯(lián)合熵H(AS0)為基準(zhǔn),按照式(11)通過(guò)差量法將H(AS0)更新為當(dāng)前基準(zhǔn)子圖的聯(lián)合熵H(AS)。
H(AS)=H(AS0)-ΔAS0+ΔAS
(11)
⑥按照式(12)計(jì)算實(shí)測(cè)圖與當(dāng)前基準(zhǔn)子圖的歸一化互信息。
NMI(A,S)=(H(A)+H(S))/H(A,S)
(12)
6) 遍歷整個(gè)基準(zhǔn)圖,選取歸一化互信息最大值對(duì)應(yīng)的位置作為最終匹配點(diǎn)。
文中算法與傳統(tǒng)方法相比,計(jì)算量的差別主要體現(xiàn)在每一個(gè)匹配位置處基準(zhǔn)子圖各統(tǒng)計(jì)量的計(jì)算上,為此只分析基準(zhǔn)子圖各統(tǒng)計(jì)量計(jì)算的時(shí)間復(fù)雜度。
設(shè)基準(zhǔn)圖尺寸為N×N,實(shí)測(cè)圖尺寸為n×n,則基準(zhǔn)子圖大小與實(shí)測(cè)圖大小一致均為n×n。
1)基準(zhǔn)子圖直方圖計(jì)算
每一個(gè)匹配位置處,傳統(tǒng)算法需要遍歷整個(gè)基準(zhǔn)子圖來(lái)計(jì)算該子圖的直方圖,而文中方法只須遍歷基準(zhǔn)子圖的第一行(或列)和最后一行(或列)即可獲得該子圖直方圖。對(duì)于整個(gè)匹配過(guò)程來(lái)說(shuō),傳統(tǒng)方法其時(shí)間復(fù)雜度為O((N-n)2n2),文中方法時(shí)間復(fù)雜度為O((N-n)2n)。
2)基準(zhǔn)子圖熵及聯(lián)合熵計(jì)算
由于相鄰基準(zhǔn)子圖其直方圖發(fā)生變化情況跟圖像內(nèi)容有關(guān),很難從理論上推導(dǎo)兩子圖間有多少個(gè)灰度值其出現(xiàn)頻率會(huì)發(fā)生變化,也很難推導(dǎo)出其聯(lián)合直方圖中有多少個(gè)灰度值對(duì)的頻率會(huì)發(fā)生變化。以直方圖為例,最壞情況時(shí)有2n(n為基準(zhǔn)子圖各邊長(zhǎng)所含像素個(gè)數(shù))個(gè)灰度值頻率發(fā)生變化,且2n不超過(guò)基準(zhǔn)圖灰度等級(jí)數(shù)。
為此,文中在標(biāo)準(zhǔn)圖像集上進(jìn)行了大量的仿真實(shí)驗(yàn),通過(guò)對(duì)仿真結(jié)果進(jìn)行統(tǒng)計(jì)得出:相鄰子圖間灰度直方圖中發(fā)生變化的元素平均為50%左右,聯(lián)合灰度直方圖中發(fā)生變化的元素平均為75%左右,這就意味著文中方法相比傳統(tǒng)方法,其基準(zhǔn)子圖熵值的計(jì)算量減少了近50%,聯(lián)合熵計(jì)算量減少了約25%,時(shí)間復(fù)雜度或時(shí)間頻度有所降低。部分實(shí)驗(yàn)結(jié)果列于表1~表4。
另一方面,因文中方法需要存儲(chǔ)前一匹配位置基準(zhǔn)子圖的多個(gè)統(tǒng)計(jì)量,算法的空間復(fù)雜度有所增加,但由于硬件技術(shù)的發(fā)展大大提高了計(jì)算機(jī)的存儲(chǔ)容量,使得存儲(chǔ)容量的局限性對(duì)于算法的影響大大降低。
為了驗(yàn)證算法的有效性,將文中方法應(yīng)用于SAR與可見(jiàn)光圖像匹配,以Matlab R2015為仿真平臺(tái),進(jìn)行兩類(lèi)實(shí)驗(yàn):一類(lèi)是通過(guò)比較每個(gè)匹配位置處文中方法與傳統(tǒng)遍歷式互信息匹配方法在計(jì)算量方面的差別以及總耗時(shí)來(lái)說(shuō)明文中方法的有效性;另一類(lèi)是將文中方法與多分辨率分層匹配相結(jié)合,與單純的多分辨率方法進(jìn)行匹配耗時(shí)及匹配精度的比較,其中的多分辨率分層匹配采用d3小波基,小波分解級(jí)數(shù)為1。限于篇幅,文中只列出部分結(jié)果。
匹配過(guò)程中設(shè): 基準(zhǔn)子圖直方圖計(jì)算中需要統(tǒng)計(jì)的像素灰度值的次數(shù)為M,基準(zhǔn)子圖熵值計(jì)算中需要計(jì)算的求和項(xiàng)的項(xiàng)數(shù)為M1,基準(zhǔn)子圖與實(shí)測(cè)圖聯(lián)合熵計(jì)算中需要計(jì)算的求和項(xiàng)項(xiàng)數(shù)為M2。
選取可見(jiàn)光與SAR圖像進(jìn)行匹配,如圖3所示,其中圖3(a)為可見(jiàn)光基準(zhǔn)圖,大小為422×358;圖3(b)為SAR實(shí)測(cè)圖像1,大小為150×150,圖3(c)為SAR實(shí)測(cè)圖像2,大小為100×90,圖3(d)為基于文中差量法匹配時(shí)實(shí)測(cè)圖在基準(zhǔn)圖上的匹配定位結(jié)果;表1為差量法與傳統(tǒng)標(biāo)準(zhǔn)算法列方向匹配時(shí)部分匹配位置計(jì)算量統(tǒng)計(jì)結(jié)果的比較,表2為差量法與傳統(tǒng)標(biāo)準(zhǔn)算法行方向匹配時(shí)部分匹配位置計(jì)算量統(tǒng)計(jì)結(jié)果的比較,表3為文中差量法與傳統(tǒng)標(biāo)準(zhǔn)方法匹配結(jié)果和匹配耗時(shí)的比較,表4為差量法結(jié)合多分辨率分層匹配法與單純多分辨率分層匹配法的匹配結(jié)果和匹配耗時(shí)的比較。
圖3 SAR與可見(jiàn)光圖像匹配
由圖3可以看出,采用歸一化互信息作為匹配測(cè)度,可以實(shí)現(xiàn)異源圖像之間的匹配;由表1~表3可以看出,相比傳統(tǒng)方法,文中差量法能夠有效減少除第一個(gè)匹配位置外其它所有匹配位置互信息測(cè)度的計(jì)算量,總匹配耗時(shí)平均減少20%以上;由表4可以看出,文中差量法可以與多分辨率分層匹配技術(shù)相結(jié)合,其匹配精度與單純的多分辨率分層匹配算法精度相同(在有些情況下匹配精度有所下降),但其匹配耗時(shí)平均下降20%左右。
與已有的以灰度壓縮、特征提取或采取優(yōu)化搜索策略等快速互信息匹配方法不同,文中基于圖像的塊狀結(jié)構(gòu)特性以及匹配過(guò)程中相鄰子圖統(tǒng)計(jì)特性的相關(guān)性, 通過(guò)差量法減少互信息測(cè)度的計(jì)算量來(lái)加快匹配速度, 其互信息測(cè)度計(jì)算精度及匹配精度與原始的基于標(biāo)準(zhǔn)互信息計(jì)算方法的計(jì)算精度及匹配精度相同, 但計(jì)算量大為減少。 另外,如果參與匹配的兩幅圖像為高分辨率清晰圖像,則可以把文中方法與灰度壓縮、多分辨率分層匹配等方法相結(jié)合,即先對(duì)待匹配的圖像進(jìn)行灰度壓縮或多分辨率分解,然后采用遍歷式方法進(jìn)行搜索,而在每一個(gè)搜索點(diǎn)處采用文中的差量法計(jì)算互信息測(cè)度,這樣可以進(jìn)一步減少匹配耗時(shí)。
表1 SAR與可見(jiàn)光圖像在列方向匹配時(shí)部分位置計(jì)算量統(tǒng)計(jì)結(jié)果
表2 SAR與可見(jiàn)光圖像在行方向匹配時(shí)部分位置計(jì)算量統(tǒng)計(jì)結(jié)果
表3 差量法與傳統(tǒng)標(biāo)準(zhǔn)法在SAR與可見(jiàn)光圖像上的匹配結(jié)果
表4 差量法結(jié)合多分辨率分層算法與單純多分辨率分層算法在SAR與可見(jiàn)光圖像上的匹配結(jié)果