圖像配准是對取自不同時間、不同傳感器或者不同視角的同一場景的兩幅圖像或者多幅圖像匹配的過程,它被廣泛地應用在遙感圖像、醫學影像、三維重構、機器人視覺等諸多領域中.而匹配就是在一幅大圖像中搜尋目標,已知該圖中有要找的目標,且該目標同模板有相同的尺寸、方向和圖像,通過一定的算法可以在圖中找到目標,確定其坐標位.利用模板匹配可以在一幅圖像中找到已知的物體.比如抓拍到了一張射門的照片,要在該照片中找到足球的位置.這時就可以采用模板匹配的方法.
目前已經開發了許多處理各種不同種類的數據和問題的技術,總的來講,按照圖像的不同可分為三類:
第一類:由於獲取圖像的不同,造成圖像之間的錯位.為了配准這些圖像,可以采用空域的方法來去除這些差別.
第二類:由於獲取圖像的不同而造成的,不過這主要是由於光照或者天氣條件不同而造成的差別.這樣的圖像在亮度上有所區別,也有的會在空域中產生差別,例如透視變形.
第三類:圖像的差別是由於目標的運動,生長或者其他景物的變化而引起的.
第二類和第三類圖像差別並不是直接通過配准消除的,但這給圖像配准帶來了很大的困難.對於圖像的配准技術,可以把它看作是以下4點選擇的組合.
·特征空間
·搜索空間
·搜索策略
·相似度測量
特征空間指的是從圖像中提取出來的用來匹配的信息;搜索空間則是指用來校准圖像的圖像變換集;搜索策略決定如何在這個空間中選擇下一個變換,如何測試並搜索出最優的變換;相似度測量則決定了每一個配准測試中的相關特性.
在圖像配准中,特征空間、相似度測量、搜索空間以及搜索策略等的選擇都會影響到最后配准的精確度.所有的圖像配准方法都可以認為是這些選擇的組合.這個框架對於我們理解現在存在的各種各樣的配准方法的特點以及其之間的關系是相當有用的.
因為主要要用到匹配技術,所以下面我們先介紹一下圖像的匹配,再介紹一些圖像配准中所用到的理論基礎和配准方法.
1.圖像的模板匹配
1.1 模板匹配的基本概念
圖像匹配技術是數字圖像處理領域的一項重要研究,並已在計算機視覺、虛擬現實場景生成、航空航天遙感測量、醫學影像分析、光學和雷達跟蹤、景物制導等領域得到了廣泛的應用.根據已知模式(模板圖),到另一幅圖中搜索相匹配的子圖像的過程,稱為模板匹配.
模板就是一幅已知的小圖像.模板匹配就是在一幅大圖像中搜尋目標,已知該圖中
有要找的目標,且該目標同模板有相同的尺寸、方向和圖像,通過一定的算法可以在圖中找到目標,確定其坐標位置.
以8 位圖像(其1個像素由1個字節描述)為例,模板T( 個像素)疊放在被搜索圖 ( 個像素)上平移,模板覆蓋被搜索圖的那塊區域叫子圖 . 為子圖左上角在被搜索圖 上的坐標.搜索范圍是:
通過比較 和 的相似性,完成模板匹配過程.
注意:圖像的數據是從下到上、從左到右排列的.
可以用下式衡量 和 相似性:
=
上式的第一項為子圖的能量,第三項為模板的能量,都與模板匹配無關.第二項是模板和子圖的互相關,隨 而改變.當模板和子圖匹配時,該項有極大值.將其歸一化,得模板匹配的相關系數:
當模板和子圖完全一樣時,相關系數 .在被搜索圖S 中完成全部搜索后,找出 的最大值 ,其對應的子圖 即為匹配目標.顯然,用這種公式做圖像匹配計算量大、速度較慢.
另一種算法是衡量T 和Sij 的誤差,其公式為:
為最小值處即為匹配目標.為提高計算速度,取一個誤差閾值 ,當 > 時就停止該點的計算,繼續下一點計算.
試驗結果如下
結果表明:被搜索圖越大,匹配速度越慢;模板越小,匹配速度越快.誤差法速度較快,閾值的大小對匹配速度影響大,和模板的尺寸有關.
1.2 改進模板匹配算法
在誤差算法的基礎上設計了二次匹配誤差算法:
第一次匹配是粗略匹配.取模板的隔行隔列數據,即四分之一的模板數據,在被搜索圖上進行隔行隔列掃描匹配,即在原圖的四分之一范圍內匹配.由於數據量大幅度減少,匹配速度顯著提高.
為了合理的給出一個誤差閾值 ,我們設計了一個確定誤差閾值 的准則:
式中: 為各點平均的最大誤差,一般取40~50 即可; 為模板的長和寬.
第二次匹配是精確匹配.在第一次誤差最小點 的鄰域內,即在對角點為 , 的矩形內,進行搜索匹配,得到最后結果.
下表是相關法、誤差法、二次匹配誤差法這三種模板匹配算法對兩幅圖像進行模板匹配的結果比較,二次匹配誤差法的速度比其它算法快了10倍左右.
使用二次匹配誤差法對256*256 像素的被搜索圖進行模板匹配的結果如下:
從上圖結果可看出,第一次匹配位置是偶數,因為是隔行隔列進行搜索的;第二次則是精確位置.
1.3 二次匹配誤差法的主要代碼
#define AvEthreshold 40 // 各點平均誤差
int Ethreshold; // 誤差閾值
//計算誤差閾值
Ethreshold=AvEthreshold*((lTemplateHeight+1)/2)*((lTemplateWidth+1)/2);
// 第一次粗略匹配,找出誤差最小位置(nMaxHeight,nMaxWidth)
// 僅使用模板中隔行隔列的數據,在被搜索圖中隔行隔列匹配
nMinError = 99999999;
for (i = 0; i < lHeight - lTemplateHeight +1; i=i+2)
{
for(j = 0; j < lWidth - lTemplateWidth + 1; j=j+2)
{
nError = 0;
for (m = 0; m < lTemplateHeight; m=m+2)
{
for(n = 0; n < lTemplateWidth; n=n+2)
{
// 指向被搜索圖像倒數第i+m行,第j+n個象素的指針
lpSrc = (char *)lpDIBBits + lLineBytes * (i+m) + (j+n);
// 指向模板圖像倒數第m行,第n個象素的指針
lpTemplateSrc = (char *)lpTemplateDIBBits + lTemplateLineBytes * m + n;
pixel = (unsigned char)*lpSrc;
templatepixel = (unsigned char)*lpTemplateSrc;
nDelta=(int)pixel-templatepixel;
if (nDelta<0) // 做絕對值運算
{nError = nError-nDelta;}
else{nError = nError+nDelta;}
}
if (nError>(Ethreshold)) break; // 誤差大於閾值,進入下一點計算
}
if (nError < nMinError) //與最小誤差比較
{
nMinError = nError; // 記錄最小誤差及其坐標
nMaxWidth = j;
nMaxHeight = i;
}
}
}
2 圖像的配准技術
2.1 圖像配准理論基礎
我們先介紹圖像配准中經常用到的圖像變換、相似性測度、插值算法以及最小二乘法.
2.1.1 圖像變換
將一幅圖像與另一幅圖像對准,常需對一幅圖像進行一系列的變換,這些變換可分為剛體變換、仿射變換、投影變換和非線形變換.
一.剛體變換
如果第一幅圖像中的兩點間的距離經變換到第二幅圖像中后仍保持不變,則這種變換稱為剛體變換.剛體變換可分解為平移、旋轉和反轉(鏡像).在二維空間中,點 經過剛體變換到點 的變換公式為:
其中 為旋轉角, 為平移向量.
二.仿射變換
經過變換后第一幅圖像上的直線映射到第二幅圖像仍為直線,並且保持平衡關系,這樣的變換稱為仿射變換.仿射變換可以分解為線形(矩陣)變換和平移變換.在2D空間中,變換公式為:
其中 為實矩陣.
三. 投影變換
經過變換后第一幅圖像上的直線映射到第二幅圖像仍為直線,並且保持平行關系,這樣的變換稱為投影變換.投影變換可用高維空間上的線形(矩陣)變換來表示.變換公式為:
四. 非線形變換
非線形變換可把直線變為曲線.在2D空間中,可以用以下公式表示:
其中,F表示把第一幅圖像映射到第二幅圖像上的任意一種函數形式.典型的非線形變換如多項式變換,在2D空間中,多項式函數可寫為如下形式:
非線形變換比較適用於那些具有全局性形變的圖像配准問題,以及整體近似剛體但局部有形變的配准情況.
2.1.2 相似性測度
相似性測度主要有互相關相似性測度、序貫相似性測度、相位相關測度以及一些變形等.
一 互相關相似性測度
最基本的相似性測度就是互相關相似性測度了.假設有模板T和搜索圖S, 表示模板覆蓋下的那塊搜索圖,其中i,j表示位置.則其互相關相似性測度為:
或者歸一化為:
二 基於傅立葉分析的測度
基於傅立葉分析的測度主要有兩種,一個是將空域中的互相關在頻域中進行計算,另一個是相位相關法.
由基於傅立葉分析中的相關定理可知,兩個函數在定義域中的卷積等於它在頻域中的乘積,而相關則是卷積的一種特定形式.因此,可以利用快速傅立葉變換來求解相似度.
對 和 分別求DFT,可得到相關的離散傅立葉變換 為:
其中 和 分別是 和 的傅立葉變換.其中*為共扼運算的符號.
對 求傅立葉反變換,就可以得出空間域中的相關函數 為:
相位相關法可以用力計算兩幅圖像的平移.圖像的傅立葉變換由實部和虛部組成,即:
其中 為虛部算子.這也可以表示為:
則
相位相關則是基於傅立葉變換的這個性質來建立的.兩幅圖像的平移可以認為是其傅立葉變換的角度差別 .這樣計算兩幅圖像的互功率譜就可以得到兩幅圖像的角度差別:
其中 和 分別表示 和 上的平移量.
三 序貫相似性檢測
序貫相似性檢測算法實際上是基於互相關法,是為了快速計算而設計的算法.其基本思路就是設計一個閥值,在計算相似性的時候,累計誤差大於閥值則停止計算,這樣可以減少在誤配准點上的計算量.
在序貫相似性檢測算法的基礎上,也有不少的變形,例如基於排序的序貫相似性檢測算法,分層的序貫相似性檢測算法等.
2.1.3 插值
在圖像配准中,還需要涉及到許多的插值算法.由於粗配后所得出的像素坐標位置可能不在整數像素上,因此需要用灰度插值的方法對像素值進行估計.灰度插值的發建立在信號抽樣理論基礎上.現實世界中的圖像可以看作是一個連續函數 ,數字圖像是對這個函數的采樣,恢復連續函數可以用一個理想低通濾波器對離散化的圖像進行濾波,在空域上就表現為與Sinc函數的卷積.但是Sinc函數計算比較復雜,一般用比較容易計算的函數來逼近,常用的方法包括最近鄰法、雙線形法,但是它們的精度對於細化比較豐富的圖像來說還不夠.在本問中,采用3次插值的方法用當前像素周圍的4×4個像素進行插值,用如下所示的核函數對圖像進行卷積操作:
下圖是 函數和Sinc函數的圖像.從圖中可以看出 函數是Sinc函數的3次良好逼近.
2.1.4 最小二乘法
最小二乘法(LSM)是一種有效的統計方法,它可以用多個配准點找出最優參數解.設對應配准點的坐標分別為 和 .
由仿射模型可得如下方程組:
配准需要解決的問題就是求出線形最小二乘解 和 ,其中 和 為仿射變換參數,以6參數的仿射變換為例
,
,
.
可以證明,只要 的秩大於6,則存在 和 的最小二乘解:
2.2 圖像配准中常用的技術
圖像配准涉及到的技術比較多,總的說來,配准的方法主要有互相關法、傅立葉變換法、點映射法和彈性模型法.前兩種方法前面的理論基礎已經進行了介紹,下面我們就介紹一下點映射法和基於彈性模型的配准方法,再回頭敘述一下配准着需要注意的問題,包括特征空間的選擇,相似性的選擇以及搜索空間和策略的選擇.
2.2.1 點映射
點映射或者是基准點映射技術是在不知道兩幅圖像的映射方式時最常采用的配准方法.例如,如果圖像是在同一個場景中從不同的角度按照平滑的景深變化拍攝到的,則這兩幅圖像會由於透視變形而有所區別.我們不能確定正確的透視變換,因為通常並不知道在場景中真實的景深.但我們能在兩幅圖像中找到基准點,並利用一個普通的變換來進行配准.然而,如果場景並不是由光滑的表面組成,而是由很大的景深變化,這樣兩幅圖像之間就會產生變形以及遮擋等問題,物體將出現在圖像的不同位置,而且一些變形也會是局部的.當更多的變形是局部的,這就給全局的點映射方法圖像配准帶來了很大的困難.在這個情況下,利用局部的變換,例如局部點映射方法將會比較合適.
點映射方法通常有3個步驟組成:
第一步:計算圖像中的特征.
第二步:在參考圖像上找到特征點,也就是經常說的控制點,並在待配准圖像中找到相對應的特征點.
第三步:利用這些配准的特征點計算出空間映射參數,這個空間通常是2D的多項式函數.通常,還需要在一幅圖像中進行重采樣,對另一幅圖像應用空間映射和插值.
點映射方法經常不能穩定的得到配准,因此,點映射方法還經常利用各個階段之間的反饋來找到最優的變換.
前面提到的控制點在圖像配准中起到了重要的作用.在點配准后,點映射的方法就只剩下插值或者逼近了.這樣,點配准的精度就確定了最后配准的精度.
很多特征都可以用來作為控制點,我們可以將它們分為內在的和外部的特征點.內在的控制點指的是圖像中不依賴於圖像數據本身的一些點,它們通常是為了配准的目的放入場景中的標記點,並且很容易進行識別.例如在醫學圖像中,就經常往患者的皮膚或者其他不會產生變形和移位的位置上放置一些特征點,以便進行配准.外在的控制點指的是那些從數據中得到點.這些點可以是手工得到的,也可以是自動獲取的.手工的控制點也就是利用人的交互得到的點,例如一些可鑒別的基准點或者解剖結構.這些點一般都是選取為剛性、穩定並且在數據中很容易點擊得到等.當然,這需要由在其領域的專門知識.還有許多應用自動定位控制點.用來自動定位的控制點典型的有角點、直線的交點、曲線中的局部最大曲率點、具有局部最大曲率的窗口的中心、閉合區域的重心等等.這些特征通常都是在配准的兩幅圖像中唯一的,並且對於局部的變形表現更魯棒一些.
得到控制點后,就可以對這些控制點進行配准了.另一幅圖像中配准的控制點可以用手工點擊得到,也可以利用互相關等方法進行自動獲取.
2.2.2 基於彈性模型的配准
在很多圖像配准中,沒有直接應用插值來計算控制點之間的映射,而是采用了基於彈性模型的方法來進行.這些方法對圖像中的變形利用了一個彈性的變形來模型化.從另一個角度來說,配准變換就是一個彈性的材料經過最小的彎曲和拉伸變換的結果.彎曲和拉伸的量由彈性材料的能量狀態來表示.然而,插值的方法也和此類似,因為能量的最小化需要滿足彈性模型的限制,而這可以用樣條來解決.
通常,這些方法是逼近圖像之間的配准,盡管它們有時也利用特征,但是它們沒有包括點配准這個步驟.圖像或者物體是模型化為一個彈性的整體,並且兩幅圖像之間點或者特征點的相似性是用整體的"拉伸"的外部的力來表示的.這些方法通常用硬度或者是光滑度的約束上給了用戶很多的靈活性.最后確定的最小能量的狀態將決定定義配准的變形變換.但其問題在於最小能量狀態的求取上通常包括迭代的計算過程.
2.2.3 特征空間的選擇
配准的過程中,第一步就是特征空間的選擇.有許多的特征可供選擇,這可能就是圖像本身的亮度,同時也有許多其他類似的選擇,這包括邊沿、曲線、表面.顯著的特征也可以用來配准,例如角點、直線交點,高曲率的點.也可以是統計特征,例如不變矩、重心等.高層的結果和語義描述也可以用來配准.
顯著的特征通常是指在圖像中可以容易辨別處理的有意義特征的一些特殊的像素.統計特征指的是一個區域(這個區域通常是通過一個分割的預處理得到)的測度,它表示了這個區域的估計.特征空間是圖像配准,同時也是幾乎所有的高層圖像處理或者計算機視覺的基礎.
對於圖像配准,特征空間的選擇將影響:
·傳感器和場景中的數據什么特征的敏感的(通常,特征是選擇那些減少傳感器噪聲和其他變形,例如亮度的變換等).
·圖像是什么特征將會配准(例如在配准結構比紋理特征更有利).
·計算代價,通過減少計算的相似性,或者從另一個角度來說,增加預計算的必要性的代價.
2.2.4 相似性測度的選擇
圖像配准的第二步就是設計或者選擇相似性測度.這個階段和配准特征的選擇有很大關系.對於內部結構,也就是說圖像的不變特征,通常是通過特征空間和通過相似性測度提取得到的.典型的相似性測度有:經過或者不經過預先濾波的互相關(例如,配准濾波或者統計相關)、差別的絕對值(這是很高效的方法)、傅立葉不變特征(例如相位相關).利用曲線和表面作為特征空間的方法則需要取在最近鄰域差別的平方和.結構或者語義的方法則需要高度依賴於其特征的測度.例如,"自由"圖形之間的熵的最小變化就是用來作為一個結構模式識別的相似性測度.
相似性測度的選擇是一個圖像配准中最重要的步驟之一,它將決定如何確定配准變換.而且,其配准的程度最后應轉化為配准或者不配准.
2.2.5 搜索空間和策略的選擇
由於很多的配准特征和相似度測量方法都需要大量的計算量,因此圖像配准中最后一步設計就是一個選擇搜索策略的問題.搜索空間通常是將要找到的配准最優變換的變換集.我們可以在特征空間上利用相似性測度計算每一個變換.然而,在很多情況下,例如利用相關作為相似性測度的方法.減少測量計算的數量是很重要的.誤配准位置越多,計算量越大,這個問題就越嚴重.在許多情況下,搜索空間就是所有變換的空間.通常,可以將這個集合從影響的大小和搜索空間的復雜度分為全局或者局部的變換.在有些情況下,可以利用一些可以得到的信息去掉不可能配准的搜索子空間,從而達到減少計算量的目的.
通常的搜索策略包括分層或多分辨技術、判決序列、松弛、廣義Hough變換、線性規划、樹和圖像配准、動態規划以及啟發式的搜索等等.
3.圖像配准及其主要VC++代碼
上面我們介紹了圖像配准的理論,下面我們利用VC++來實現圖像的配准.我們給出一個半自動的基於特征的圖像配准MFC,在程序中,需要手工選取特征點,程序將自動尋找到相匹配的特征點,然后自動計算 仿射變換參數,並將兩幅圖像拼接,顯示匹配圖像.
配准后:
下面我們介紹一個在圖像配准中用到的一個主要函數TemplateMatch()的VC++源代碼:
BOOL CDlgRecMatch::TemplateMatch(CDib* pDibSrc, CDib* pDibTemplate)
{
// 指向源圖像的指針
LPBYTE lpSrc,lpTemplateSrc;
// 指向緩存圖像的指針
LPBYTE lpDst;
//循環變量
long i;
long j;
long m;
long n;
//中間結果
double dSigmaST;
double dSigmaS;
double dSigmaT;
//相似性測度
double R;
//最大相似性測度
double dbMaxR;
//最大相似性出現位置
int nMaxWidth;
int nMaxHeight;
//像素值
unsigned char unchPixel;
unsigned char unchTemplatePixel;
// 獲得圖象數據存儲的高度和寬度
CSize sizeSaveImage;
sizeSaveImage = pDibSrc->GetDibSaveDim();
// 獲得模板圖象數據存儲的高度和寬度
CSize sizeSaveTemplate;
sizeSaveTemplate = pDibTemplate->GetDibSaveDim();
// 暫時分配內存,以保存新圖像
CDib* pDibNew;
pDibNew = new CDib;
// 如果分配內存失敗,則推出
if(!CopyDIB(pDibSrc,pDibNew)){
// 釋放已分配內存
pDibNew->Empty();
// 返回
return FALSE;
}
// 初始化新分配的內存
lpDst = (LPBYTE)pDibNew->m_lpImage;
// 圖象的高度
int nImageHeight ;
nImageHeight = pDibSrc->m_lpBMIH->biHeight;
// 圖象的寬度
int nImageWidth;
nImageWidth = pDibSrc->m_lpBMIH->biWidth;
// 模板圖象的高度
int nTemplateHeight;
nTemplateHeight = pDibTemplate->m_lpBMIH->biHeight;
// 模板圖象的寬度
int nTemplateWidth;
nTemplateWidth = pDibTemplate->m_lpBMIH->biWidth;
//計算dSigmaT
dSigmaT = 0;
for (n = 0;n < nTemplateHeight ;n++)
{
for(m = 0;m < nTemplateWidth ;m++)
{
// 指向模板圖像倒數第j行,第i個象素的指針
lpTemplateSrc = (LPBYTE)pDibTemplate->m_lpImage + sizeSaveTemplate.cx * n + m;
unchTemplatePixel = (unsigned char)*lpTemplateSrc;
dSigmaT += (double)unchTemplatePixel*unchTemplatePixel;
}
}
//找到圖像中最大相似性的出現位置
dbMaxR = 0.0;
for (j = 0;j < nImageHeight - nTemplateHeight +1 ;j++)
{
for(i = 0;i < nImageWidth - nTemplateWidth + 1;i++)
{
dSigmaST = 0;
dSigmaS = 0;
for (n = 0;n < nTemplateHeight ;n++)
{
for(m = 0;m < nTemplateWidth ;m++)
{
// 指向源圖像倒數第j+n行,第i+m個象素的指針
lpSrc = (LPBYTE)pDibSrc->m_lpImage + sizeSaveImage.cx * (j+n) + (i+m);
// 指向模板圖像倒數第n行,第m個象素的指針
lpTemplateSrc = (LPBYTE)pDibTemplate->m_lpImage + sizeSaveTemplate.cx * n + m;
unchPixel = (unsigned char)*lpSrc;
unchTemplatePixel = (unsigned char)*lpTemplateSrc;
dSigmaS += (double)unchPixel*unchPixel;
dSigmaST += (double)unchPixel*unchTemplatePixel;
}
}
//計算相似性
R = dSigmaST / ( sqrt(dSigmaS)*sqrt(dSigmaT));
//與最大相似性比較
if (R > dbMaxR)
{
dbMaxR = R;
nMaxWidth = i;
nMaxHeight = j;
}
}
}
// 對目標圖象的象素進行賦值
for(i=0; i<nImageHeight; i++)
for( j=0; j<nImageWidth; j++){
lpDst[i*sizeSaveImage.cx +j] /=2;
}
//將最大相似性出現區域部分復制到目標圖像
for (n = 0;n < nTemplateHeight ;n++)
{
for(m = 0;m < nTemplateWidth ;m++)
{
lpTemplateSrc = (LPBYTE)pDibTemplate->m_lpImage + sizeSaveTemplate.cx * n + m;
lpDst = (LPBYTE)pDibNew->m_lpImage + sizeSaveImage.cx * (n+nMaxHeight) + (m+nMaxWidth);
*lpDst = *lpTemplateSrc;
}
}
// 復制圖像
memcpy(pDibSrc->m_lpImage, pDibNew->m_lpImage, nImageWidth * nImageHeight);
// 釋放內存
pDibNew->Empty();
// 返回
return TRUE;
}
其中CDib*pDibSrc指向CDib類的指針,含有待匹配圖象信息,CDib*pDibTemplate指向CDib類的指針,含有模板圖象信息 該函數將對圖象進行模板匹配操作.需要注意的是,此程序只處理256灰度級的圖象
