轉自：http://blog.sina.com.cn/s/blog_4b9b714a0100d5k5.html

?圖像配準的方法

1? 基于特征的圖像配準

基于特征的圖像配準首先提取圖像信息的特征，然后以這些特征為模型進行配準。特征提取的結果是一含有特征的表和對圖像的描述，每個特征由一組屬性表示，對屬性的進一步描述包括邊緣的定向和弧度、區域的大小等。局部特征之間存在著相互關系，如幾何關系、輻射度量關系、拓撲關系等。可以用這些局部特征之間的關系描述全局特征。通常基于局部特征配準大多都是基于點、線或邊緣的，而全局特征的配準則是利用局部特征之間的關系進行配準的方法。

由于圖像的特征點比圖像的像素點要少很多，因此大大減少了配準過程的計算量，但特征提取方法的計算代價通常較大，不便于實時應用。特征點的配準度量值對位置的變化比較敏感，可以大大提高配準的精確程度。對于紋理較少的圖像區域提取的特征的密度通常比較稀少，局部特征的提取就比較困難。特征點的提取過程可以減少噪聲的影響，對灰度變化、圖像形變和遮擋等都有較好的適應能力。因此，在圖像配準領域得到了廣泛應用。基于特征的圖像配準方法有兩個重要環節：特征提取和特征配準。

2? 基于互信息的圖像配準

醫學圖像配準技術從基于特征的配準方法發展到基于統計的配準方法有其突破性的意義。與基于特征的配準方法相比，基于統計的配準方法的突出優點為魯棒性好、配準精度高、人工干預少。基于統計的配準方法通常是指最大互信息的圖像配準方法。

基于互信息的圖像配準是用兩幅圖像的聯合概率分布與完全獨立時的概率分布的廣義距離來估計互信息，并作為多模態醫學圖像配準的測度。當兩幅基于共同的解剖結構的圖像達到最佳配準時，它們的對應像素的灰度互信息應為最大。由于基于互信息的配準對噪聲比較敏感，首先，通過濾波和分割等方法對圖像進行預處理。然后進行采樣、變換、插值、優化從而達到配準的目的。

基于互信息的配準技術屬于基于像素相似性的方法。它基于圖像中所有的像素進行配準，基于互信息的圖像配準引入了信息論中的概念，如熵、邊緣熵、聯合熵和互信息等，可使配準精度達到亞像素級的高精度。????

基于互信息只依賴于圖像本身的信息，不需要對圖像進行特征點提取和組織分類等預處理，是一種自動而有效的配準算法。該算法可靠，對圖像中的幾何失真、灰度不均勻和數據的缺失等不敏感。不依賴于任何成像設備，可應用于多模態醫學圖像配準。基于互信息的圖像配準也有其缺點，它運算量大，對噪聲敏感，要求待配準圖像間聯合概率分布函數必須是嚴格正性的。

【例1】計算兩幅圖像的互信息。

function mi = MI(a,b)

%Caculate MI of a and b in the region of the overlap part

%計算重疊部分

[Ma,Na] = size(a);

[Mb,Nb] = size(b);

M=min(Ma,Mb);

N=min(Na,Nb);

%初始化直方圖數組

hab = zeros(256,256);

ha = zeros(1,256);

hb = zeros(1,256);

%歸一化

if max(max(a))~=min(min(a))

a = (a-min(min(a)))/(max(max(a))-min(min(a)));

else

a = zeros(M,N);

end

if max(max(b))-min(min(b))

b = (b-min(min(b)))/(max(max(b))-min(min(b)));

else

b = zeros(M,N);

end

a = double(int16(a*255))+1;

b = double(int16(b*255))+1;

%統計直方圖

for i=1:M

for j=1:N

indexx = a(i,j);

indexy = b(i,j) ;

hab(indexx,indexy) = hab(indexx,indexy)+1;%聯合直方圖

ha(indexx) = ha(indexx)+1;%a圖直方圖

hb(indexy) = hb(indexy)+1;%b圖直方圖

end

end

%計算聯合信息熵

hsum = sum(sum(hab));

index = find(hab~=0);

p = hab/hsum;

Hab = sum(sum(-p(index).*log(p(index))));

%計算a圖信息熵

hsum = sum(sum(ha));

index = find(ha~=0);

p = ha/hsum;

Ha = sum(sum(-p(index).*log(p(index))));

%計算b圖信息熵

hsum = sum(sum(hb));

index = find(hb~=0);

p = hb/hsum;

Hb = sum(sum(-p(index).*log(p(index))));

%計算a和b的互信息

mi = Ha+Hb-Hab;

%計算a和b的歸一化互信息

%mi = hab/(Ha+Hb);

3? 醫學圖像配準

目前，圖像配準技術在醫學領域應用得相當廣泛。醫學圖像配準具有很重要的臨床應用價值。對各種使用相同或者不同的成像手段所獲得的醫學圖像進行配準不僅可以應用于醫療診斷，也可以應用于手術計劃的制定、治療計劃的制定、病理情況的跟蹤和治療效果的評價等各個方面。

下面給出一個具體的圖像配準的例子來形象說明配準的應用。

?

(1) 讀取圖像數據。因為源圖像以矩陣形式存在一個二進制的文件里，用fread可將其讀取到變量矩陣中。將讀取文件編制成一個子函數(RTIread.m)。源代碼如下：

function imMatrix=RTIread(FILENAME,SIZE)

%RTIread Read the image matrix from binary "Registration Test Image" file.

% imMatrix=RTIread(FILENAME,SIZE) opens the file FILENAME, and reads the

% number of elements specified by SIZE.

%

% FILENAME is a string containing the name of the file to be opened.

% Valid entries for SIZE are:

% N read N elements into a column vector.

% inf read to the end of the file.

% [M,N] read elements to fill an M-by-N matrix, in column order.

% N can be inf, but M can't.

%

% It returns the image matrix.

fid=fopen(FILENAME,'r');

imMatrix=fread(fid,SIZE,'uint8=>uint8');

fclose(fid);

這里選取了兩張256×256的圖片，文件名為sag1.tif和sag2.tif。

?

運行以下代碼讀取圖像矩陣：

base=RTIread('sag1.tif',[256,256]);

input=RTIread(sag2.tif,[256,256]);

(2) 選取匹配點。根據預定的配準方法，選定足夠的匹配點對。運行下列代碼：

cpselect(input,base); %please select 15 points for test.

?

注意：

這里允許選擇15個點進行測試。

操作很簡單，只需注意選點要均勻布開，以增加其代表性。選定完畢，再選擇File→Save Points to Workspace命令將數據保存到工作區中。Workspace立刻多出兩個N×2的數組(其中N為選定的匹配點對數)，分別為input_points和base_points。

(3) 利用十字相關法調整選定了的匹配點。如果對上面一個步驟的選點工作不滿意，可以進行這一步。如果感覺尚可，則跳過這一步，以減小運算量。

運行代碼：

input_points_corr = cpcorr(input_points,base_points,input,base);

%優化選擇點的作用

input_points_corr為優化后在輸入圖片的對應匹配點。

(4) 計算變換公式的參數。利用cp2tform選定變換類型(配準方法)，計算變換參數。以下只需選定一種即可。

% (1) not Fine-tune points

Tlinear = cp2tform(input_points,base_points,'linear conformal');

Taffine = cp2tform(input_points,base_points,'affine');

Tprojective = cp2tform(input_points,base_points,'projective');

Tpolynomial2 = cp2tform(input_points,base_points,'polynomial',2);

Tpolynomial3 = cp2tform(input_points,base_points,'polynomial',3);

Tpolynomial4 = cp2tform(input_points,base_points,'polynomial',4);

Tpiecewise = cp2tform(input_points,base_points,'piecewise linear');

Tlwm = cp2tform(input_points,base_points,'lwm');

% (2)Fine-tune points

fTlinear = cp2tform(input_points_corr,base_points,'linear conformal');

fTaffine = cp2tform(input_points_corr,base_points,'affine');

fTprojective = cp2tform(input_points_corr,base_points,'projective');

fTpolynomial2 = cp2tform(input_points_corr,base_points,'polynomial',2);

fTpolynomial3 = cp2tform(input_points_corr,base_points,'polynomial',3);

fTpolynomial4 = cp2tform(input_points_corr,base_points,'polynomial',4);

fTpiecewise = cp2tform(input_points_corr,base_points,'piecewise linear');

fTlwm = cp2tform(input_points_corr,base_points,'lwm');

(5) 變換圖像。

% Transform the unregistered image to bring it into alignment.

title('image registration polynomial method');

subplot(2,2,1);

imshow(base);

title('Base image');

subplot(2,2,2);

imshow(input);

title('Input image');

subplot(2,2,3);

imshow(imtransform(input,Tpolynomial2));

title('registered image');

subplot(2,2,4);

imshow(imtransform(input,fTpolynomial2));

title('registered image(fine-tune points)');

?

小結：

(1) 選擇適當的方法來建立轉換參數，并非算法越復雜越好，應參考成像因素(退化因素)。

(2) 利用cpselect選擇匹配點，cpselect可以返回一個GUI句柄。在使用中，應當盡量選擇恰當的匹配區域，這樣能夠得到較好的效果。

上面的這個例子典型地說明了圖像配準的步驟：提取圖像特征，根據特征對進行幾何變換，執行整個變換。