Matlab DIP(瓦)ch11表示與描述練習

     這一章主要是練習一些圖像的表示和描述，內容主要包括圖像的一些常用表示方法，比如說鏈碼，MPP，標記，邊界片段，骨骼等。另一方面是圖像的描述算子，比如說邊界描述算子，區域描述算子等。有些課本練習的代碼暫時還沒有搞定。

     其練習過程和結果如下所示：

%%  數組單元概念
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1102(a)(noisy_circular_stroke).tif');
imshow(f);
[a b c d]=image_stats(f);
a=image_stats1(f);
b=image_stats2(f);

%% Freeman鏈碼及其某些變體
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1102(a)(noisy_circular_stroke).tif');
subplot(231),imshow(f);
title('freeman原圖');

h=fspecial('average',9);%9*9平滑算子
g=imfilter(f,h,'replicate');%濾波
subplot(232),imshow(g);
title('freeman中值濾波后圖像');

g=im2bw(g,0.5);%將原圖像轉換為二值圖像
subplot(233),imshow(g);
title('freeman圖像閾值處理后');

%函數boundaries(BW,CONN,DIR)表示的是跟蹤目標的邊界，返回值為一個p*1的數組單元，p為目標的個數
%每一個單元又是一個Q*2的矩陣，即Q個點的x,y坐標
B=boundaries(g);

%函數cellfun(FUN,C)是指對一個單元數組中的每一個單元應用函數FUN
d=cellfun('length',B);%求B中每一個目標邊界的長度,所以返回值d是一個向量

[max_d,k]=max(d);%返回向量d中最大的值，存在max_d中，k為其索引。
v=B{k(1)};%如果最大邊界不止一條，則取出其中的一條即可。v是一個坐標數組
[M,N]=size(g);

%函數bound2im(b,M,N,x0,y0)是生成一副二值圖像,大小為M*N，x0和y0是b中最小的x和y軸坐標
g=bound2im(v,M,N,min(v(:,1)),min(v(:,2)));

subplot(234),imshow(g);
title('取出最大邊界后圖像');

% [S, SU] = bsubsamp(B, GRIDSEP)指的是對邊界B進行子采樣，每個GRIDSEP個點采樣一次
%輸出S是采樣后的值，輸出SU是歸一化的邊界，什么意思？
[s,su]=bsubsamp(v,50);
g2=bound2im(s,M,N,min(s(:,1)),min(s(:,2)));
subplot(235),imshow(g2);
title('蔣采樣后的邊界圖像');

%connectpoly(X,Y)函數是按照坐標(X,Y)順時針或者逆時針連接成多邊形
cn=connectpoly(s(:,1),s(:,2));
g2=bound2im(cn,M,N,min(cn(:,1)),min(cn(:,2)));
subplot(236),imshow(g2);
title('連接成多邊形圖像');
%freeman鏈碼過程如下：

%函數fchcode為計算邊界的freeman鏈碼
c=fchcode(su);
%鏈碼結果：

%% 獲得一個區域的邊界的細胞牆(木有搞定)
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1107(a)(mapleleaf).tif');
subplot(231),imshow(f);
title('原圖像');

B=bwperim(f,8);%返回一個包含二進制圖像的邊界圖像，采用4近鄰
subplot(232),imshow(B);
title('二進制邊界圖像');

Q=qtdecomp(B,0,2);%對B圖進行四叉樹分解，最小尺寸為2*2，用0做閾值
subplot(233),imshow(Q);
title('四叉樹分解后圖像')

%% 使用函數minperpoly
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1107(a)(mapleleaf).tif');
subplot(231),imshow(f);
title('原圖像');

b=boundaries(f,4,'cw');%找出B的邊界，順時針4鄰接
b=b{1};%因為這幅圖比較特殊，它只有一個邊界物體
[M,N]=size(f);
xmin=min(b(:,1));
ymin=min(b(:,2));
bim=bound2im(b,M,N,xmin,ymin);%將邊界物體轉換成邊界圖像
subplot(232),imshow(bim);
title('找出的邊界圖像');

% [X, Y] = minperpoly(F, CELLSIZE)為計算二值圖像的最小周長多邊形
[x,y]=minperpoly(f,2)%2為用於形成邊界方形的大小，x和y是MPP上頂點的坐標

b2=connectpoly(x,y);%將頂點依次連接成多邊形
b2=bound2im(b2,M,N,xmin,ymin);
subplot(233),imshow(b2);
title('使用大小為2*2的方形邊界單元獲得的MPP');

[x,y]=minperpoly(f,3)%2為用於形成邊界方形的大小，x和y是MPP上頂點的坐標
b3=connectpoly(x,y);%將頂點依次連接成多邊形
b3=bound2im(b3,M,N,xmin,ymin);
subplot(234),imshow(b3);
title('使用大小為3*3的方形邊界單元獲得的MPP');

[x,y]=minperpoly(f,4)%2為用於形成邊界方形的大小，x和y是MPP上頂點的坐標
b4=connectpoly(x,y);%將頂點依次連接成多邊形
b4=bound2im(b4,M,N,xmin,ymin);
subplot(235),imshow(b4);
title('使用大小為4*4的方形邊界單元獲得的MPP');

[x,y]=minperpoly(f,8)%2為用於形成邊界方形的大小，x和y是MPP上頂點的坐標
b8=connectpoly(x,y);%將頂點依次連接成多邊形
b8=bound2im(b8,M,N,xmin,ymin);
subplot(236),imshow(b8);
title('使用大小為8*8的方形邊界單元獲得的MPP');
%使用minperpoly過程如下：

%% 標記(木有搞定)
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1111(a)(boundary_sq.tif');
imshow(f);
title('標記原圖');

[st,angle,x0,y0]=signature(f);
plot(angle,st)

%% 計算一個區域的骨骼
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1113(a)(chromo_original).tif');
subplot(231),imshow(f);
whos f;%其實是個二值圖像
title('人體染色體原始圖');

f=im2double(f);%不再是二值圖像了
h=fspecial('gaussian',25,25);
g=imfilter(f,h,'replicate');
subplot(232),imshow(g);
title('高斯平滑后的圖像');%為什么平滑后那些雜點就沒了呢，至少平滑到那點時像素值它自己的權值最大的啊？

g=im2bw(g,1.5*graythresh(g));%其中的graythresh()函數為自動閾值，處理后g為二值圖像了
subplot(233),imshow(g);
title('經閾值處理后');

s=bwmorph(g,'skel',Inf);%提取骨骼,bwmorph為二進制形態學處理函數
subplot(234),imshow(s);
title('骨骼化后圖像');

s1=bwmorph(s,'spur',8);%刪除毛刺
subplot(235),imshow(s1);
title('8次去毛刺后');

s2=bwmorph(s1,'spur',7);%刪除毛刺
subplot(236),imshow(s2);
title('再7次去毛刺后');
%提取染色體骨骼圖過程如下：

%% 傅里葉描述子
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1113(a)(chromo_original).tif');
subplot(331),imshow(f);
whos f;%其實是個二值圖像
title('人體染色體原始圖');

f=im2double(f);%不再是二值圖像了
h=fspecial('gaussian',15,9);
g=imfilter(f,h,'replicate');
g=im2bw(g,0.7);%其中的graythresh()函數為自動閾值，處理后g為二值圖像了
subplot(332),imshow(g);
title('經閾值處理后');

b=boundaries(g);
b=b{1};
bim=bound2im(b,344,270);%344,270為圖片的尺寸
subplot(333),imshow(bim);
title('提取的邊界圖');

z=frdescp(b);%對邊界坐標b進行傅里葉變換系數(有多少個點就有多少個系數)，將b的坐標點看成是復平面中的某個復數
z14=ifrdescp(z,546);%用50%的描述子進行逆變換
z546im=bound2im(z14,344,270);
subplot(334),imshow(z546im);
title('546個描述子恢復后');

z110=ifrdescp(z,110);%用50%的描述子進行逆變換
z110im=bound2im(z110,344,270);
subplot(335),imshow(z110im);
title('110個描述子恢復后');

z56=ifrdescp(z,56);%用50%的描述子進行逆變換
z56im=bound2im(z56,344,270);
subplot(336),imshow(z56im);
title('56個描述子恢復后');

z28=ifrdescp(z,28);%用50%的描述子進行逆變換
z28im=bound2im(z28,344,270);
subplot(337),imshow(z28im);
title('28個描述子恢復后');

z14=ifrdescp(z,14);%用50%的描述子進行逆變換
z14im=bound2im(z14,344,270);
subplot(338),imshow(z14im);
title('14個描述子恢復后');

z8=ifrdescp(z,8);%用50%的描述子進行逆變換
z8im=bound2im(z8,344,270);
subplot(339),imshow(z8im);
title('8個描述子恢復后');
%傅里葉描述子圖過程：

%% 函數regionprops的使用
clc
clear
B =  [1     1     1     0     0     0     0     0
      1     1     1     0     1     1     0     0
      1     1     1     0     1     1     0     0
      1     1     1     0     0     0     0     0
      1     1     1     0     0     0     1     0
      1     1     1     0     0     0     1     0
      1     1     1     0     0     1     1     0
      1     1     1     0     0     0     0     0];
  B=bwlabel(B);%將矩陣B變成標記矩陣，默認為8連通，將連通區域用不同的整數給標記出來
  
  D=regionprops(B,'area','boundingbox')
  w=[D.Area]
  NR=length(w)
  
  V=cat(1,D.BoundingBox)%cat函數表示沿着第一個方向級聯存放數組
%該函數的使用過程如下：

%% 統計紋理度量(木有搞定)
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1119(a)(superconductor-with-window).tif');
sA=imhist(f);
statxture(sA)%計算圖像的6個統計度量量

%% 計算頻譜紋理
clc
clear

f1=imread('..\images\dipum_images_ch11\Fig1121(a)(random_matches).tif');
subplot(241),imshow(f1);
title('隨機排列火柴');

f2=imread('..\images\dipum_images_ch11\Fig1121(b)(ordered_matches).tif');
subplot(245),imshow(f2);
title('規則排列火柴');

[srad1,sang1,s1]=specxture(f1);
subplot(242),imshow(s1);
title('隨機排列火柴頻譜圖');

subplot(243),plot(srad1);
title('隨機排列火柴幅度圖');

subplot(244),plot(sang1);
title('隨機排列火柴角度圖');

[srad2,sang2,s2]=specxture(f2);
subplot(246),imshow(s2);
title('規則排列火柴頻譜圖');

subplot(247),plot(srad2);
title('規則排列火柴幅度圖');

subplot(248),plot(sang2);
title('規則排列火柴角度圖');
%計算頻譜紋理圖過程如下：

%% 不變矩
clc
clear
f=imread('..\images\dipum_images_ch11\Fig1123(a)(Original_Padded_to_568_by_568).tif');
imshow(f);
title('不變矩原圖');
%不變矩原始圖像如下：

%padarray函數為矩陣0填充，圖像原始大小為400*400，但是需要旋轉45度，所以最大變為568*568
fp=padarray(f,[84 84],'both');%因此四個方向每個方向都要擴充84個像素

fhs=f(1:2:end,1:2:end);%表示從第一個到最后一個，每2個采樣一次
figure,subplot(221),imshow(fhs);
title('半大圖像');
fhsp=padarray(fhs,[184 184],'both');

fm=fliplr(f);
subplot(222),imshow(fm);
title('鏡面圖像');
fmp=padarray(fm,[84 84],'both');

fr2=imrotate(f,2,'bilinear');%旋轉2度，雙線性插值
subplot(223),imshow(fr2);
title('旋轉2度后圖像');

fr45=imrotate(f,45,'bilinear');
subplot(224),imshow(fr45);
title('旋轉45度后圖像');
%不變矩變換圖像如下：

phiorig=abs(log(invmoments(f)))
phihalf=abs(log(invmoments(fhs)))
phimirror=abs(log(invmoments(fm)))
phirot2=abs(log(invmoments(fr2)))
phirot45=abs(log(invmoments(fr45)))
%這5幅圖像的7個不變矩圖像如下：

%圖像主分量
clc
clear
f1=imread('..\images\dipum_images_ch11\Fig1125(a)(WashingtonDC_Band1_512).tif');
subplot(231),imshow(f1);
title('可見藍色波段圖像');

f2=imread('..\images\dipum_images_ch11\Fig1125(b)(WashingtonDC_Band2_512).tif');
subplot(232),imshow(f2);
title('可見綠色波段圖像');

f3=imread('..\images\dipum_images_ch11\Fig1125(c)(WashingtonDC_Band3_512).tif');
subplot(233),imshow(f3);
title('可見紅色波段圖像');

f4=imread('..\images\dipum_images_ch11\Fig1125(d)(WashingtonDC_Band4_512).tif');
subplot(234),imshow(f4);
title('近紅外光圖像');

f5=imread('..\images\dipum_images_ch11\Fig1125(e)(WashingtonDC_Band5_512).tif');
subplot(235),imshow(f5);
title('中紅外光圖像');

f6=imread('..\images\dipum_images_ch11\Fig1125(f)(WashingtonDC_Band6_512).tif');
subplot(236),imshow(f6);
title('熱紅外光圖像');
%6幅多譜圖像顯示結果如下：

S=cat(3,f1,f2,f3,f4,f5,f6);%將6幅圖像在第三維中連接起來
[X,R]=imstack2vectors(S);%將S組織到數組X中
p=princomp(X,6);%princomp用來獲得6幅圖像的主分量

g1=p.Y(:,1);%計算第1主維分量
g1=reshape(g1,512,512);
figure,subplot(231),imshow(g1,[]);

g2=p.Y(:,2);%計算第2主維分量
g2=reshape(g2,512,512);
subplot(232),imshow(g2,[]);

g3=p.Y(:,3);%計算第3主維分量
g3=reshape(g3,512,512);
subplot(233),imshow(g3,[]);

g4=p.Y(:,4);%計算第4主維分量
g4=reshape(g4,512,512);
subplot(234),imshow(g4,[]);

g5=p.Y(:,5);%計算第5主維分量
g5=reshape(g5,512,512);
subplot(235),imshow(g5,[]);

g6=p.Y(:,6);%計算第6主維分量
g6=reshape(g6,512,512);
subplot(236),imshow(g6,[]);
%主分量分解過程如下：

 

     圖像的表示和描述對后面做計算機視覺處理起來很大的作用。當然針對不同的研究對象其表示和描述的特點也不同。