Matlab DIP(瓦)ch11表示与描述练习

     这一章主要是练习一些图像的表示和描述,内容主要包括图像的一些常用表示方法,比如说链码,MPP,标记,边界片段,骨骼等。另一方面是图像的描述算子,比如说边界描述算子,区域描述算子等。有些课本练习的代码暂时还没有搞定。

     其练习过程和结果如下所示:

  1 %%  数组单元概念
2 clc
3 clear
4 f=imread('..\images\dipum_images_ch11\Fig1102(a)(noisy_circular_stroke).tif');
5 imshow(f);
6 [a b c d]=image_stats(f);
7 a=image_stats1(f);
8 b=image_stats2(f);
9
10 %% Freeman链码及其某些变体
11 clc
12 clear
13 f=imread('..\images\dipum_images_ch11\Fig1102(a)(noisy_circular_stroke).tif');
14 subplot(231),imshow(f);
15 title('freeman原图');
16
17 h=fspecial('average',9);%9*9平滑算子
18 g=imfilter(f,h,'replicate');%滤波
19 subplot(232),imshow(g);
20 title('freeman中值滤波后图像');
21
22 g=im2bw(g,0.5);%将原图像转换为二值图像
23 subplot(233),imshow(g);
24 title('freeman图像阈值处理后');
25
26 %函数boundaries(BW,CONN,DIR)表示的是跟踪目标的边界,返回值为一个p*1的数组单元,p为目标的个数
27 %每一个单元又是一个Q*2的矩阵,即Q个点的x,y坐标
28 B=boundaries(g);
29
30 %函数cellfun(FUN,C)是指对一个单元数组中的每一个单元应用函数FUN
31 d=cellfun('length',B);%求B中每一个目标边界的长度,所以返回值d是一个向量
32
33 [max_d,k]=max(d);%返回向量d中最大的值,存在max_d中,k为其索引。
34 v=B{k(1)};%如果最大边界不止一条,则取出其中的一条即可。v是一个坐标数组
35 [M,N]=size(g);
36
37 %函数bound2im(b,M,N,x0,y0)是生成一副二值图像,大小为M*N,x0和y0是b中最小的x和y轴坐标
38 g=bound2im(v,M,N,min(v(:,1)),min(v(:,2)));
39
40 subplot(234),imshow(g);
41 title('取出最大边界后图像');
42
43 % [S, SU] = bsubsamp(B, GRIDSEP)指的是对边界B进行子采样,每个GRIDSEP个点采样一次
44 %输出S是采样后的值,输出SU是归一化的边界,什么意思?
45 [s,su]=bsubsamp(v,50);
46 g2=bound2im(s,M,N,min(s(:,1)),min(s(:,2)));
47 subplot(235),imshow(g2);
48 title('蒋采样后的边界图像');
49
50 %connectpoly(X,Y)函数是按照坐标(X,Y)顺时针或者逆时针连接成多边形
51 cn=connectpoly(s(:,1),s(:,2));
52 g2=bound2im(cn,M,N,min(cn(:,1)),min(cn(:,2)));
53 subplot(236),imshow(g2);
54 title('连接成多边形图像');
55%freeman链码过程如下:


56
57 %函数fchcode为计算边界的freeman链码
58 c=fchcode(su);
59%链码结果:


60
61 %% 获得一个区域的边界的细胞墙(木有搞定)
62 clc
63 clear
64 f=imread('..\images\dipum_images_ch11\Fig1107(a)(mapleleaf).tif');
65 subplot(231),imshow(f);
66 title('原图像');
67
68 B=bwperim(f,8);%返回一个包含二进制图像的边界图像,采用4近邻
69 subplot(232),imshow(B);
70 title('二进制边界图像');
71
72 Q=qtdecomp(B,0,2);%对B图进行四叉树分解,最小尺寸为2*2,用0做阈值
73 subplot(233),imshow(Q);
74 title('四叉树分解后图像')
75
76 %% 使用函数minperpoly
77 clc
78 clear
79 f=imread('..\images\dipum_images_ch11\Fig1107(a)(mapleleaf).tif');
80 subplot(231),imshow(f);
81 title('原图像');
82
83 b=boundaries(f,4,'cw');%找出B的边界,顺时针4邻接
84 b=b{1};%因为这幅图比较特殊,它只有一个边界物体
85 [M,N]=size(f);
86 xmin=min(b(:,1));
87 ymin=min(b(:,2));
88 bim=bound2im(b,M,N,xmin,ymin);%将边界物体转换成边界图像
89 subplot(232),imshow(bim);
90 title('找出的边界图像');
91
92 % [X, Y] = minperpoly(F, CELLSIZE)为计算二值图像的最小周长多边形
93 [x,y]=minperpoly(f,2)%2为用于形成边界方形的大小,x和y是MPP上顶点的坐标
94
95 b2=connectpoly(x,y);%将顶点依次连接成多边形
96 b2=bound2im(b2,M,N,xmin,ymin);
97 subplot(233),imshow(b2);
98 title('使用大小为2*2的方形边界单元获得的MPP');
99
100 [x,y]=minperpoly(f,3)%2为用于形成边界方形的大小,x和y是MPP上顶点的坐标
101 b3=connectpoly(x,y);%将顶点依次连接成多边形
102 b3=bound2im(b3,M,N,xmin,ymin);
103 subplot(234),imshow(b3);
104 title('使用大小为3*3的方形边界单元获得的MPP');
105
106 [x,y]=minperpoly(f,4)%2为用于形成边界方形的大小,x和y是MPP上顶点的坐标
107 b4=connectpoly(x,y);%将顶点依次连接成多边形
108 b4=bound2im(b4,M,N,xmin,ymin);
109 subplot(235),imshow(b4);
110 title('使用大小为4*4的方形边界单元获得的MPP');
111
112 [x,y]=minperpoly(f,8)%2为用于形成边界方形的大小,x和y是MPP上顶点的坐标
113 b8=connectpoly(x,y);%将顶点依次连接成多边形
114 b8=bound2im(b8,M,N,xmin,ymin);
115 subplot(236),imshow(b8);
116 title('使用大小为8*8的方形边界单元获得的MPP');
117%使用minperpoly过程如下:


118
119 %% 标记(木有搞定)
120 clc
121 clear
122 f=imread('..\images\dipum_images_ch11\Fig1111(a)(boundary_sq.tif');
123 imshow(f);
124 title('标记原图');
125
126 [st,angle,x0,y0]=signature(f);
127 plot(angle,st)
128
129 %% 计算一个区域的骨骼
130 clc
131 clear
132 f=imread('..\images\dipum_images_ch11\Fig1113(a)(chromo_original).tif');
133 subplot(231),imshow(f);
134 whos f;%其实是个二值图像
135 title('人体染色体原始图');
136
137 f=im2double(f);%不再是二值图像了
138 h=fspecial('gaussian',25,25);
139 g=imfilter(f,h,'replicate');
140 subplot(232),imshow(g);
141 title('高斯平滑后的图像');%为什么平滑后那些杂点就没了呢,至少平滑到那点时像素值它自己的权值最大的啊?
142
143 g=im2bw(g,1.5*graythresh(g));%其中的graythresh()函数为自动阈值,处理后g为二值图像了
144 subplot(233),imshow(g);
145 title('经阈值处理后');
146
147 s=bwmorph(g,'skel',Inf);%提取骨骼,bwmorph为二进制形态学处理函数
148 subplot(234),imshow(s);
149 title('骨骼化后图像');
150
151 s1=bwmorph(s,'spur',8);%删除毛刺
152 subplot(235),imshow(s1);
153 title('8次去毛刺后');
154
155 s2=bwmorph(s1,'spur',7);%删除毛刺
156 subplot(236),imshow(s2);
157 title('再7次去毛刺后');
158%提取染色体骨骼图过程如下:


159
160 %% 傅里叶描述子
161 clc
162 clear
163 f=imread('..\images\dipum_images_ch11\Fig1113(a)(chromo_original).tif');
164 subplot(331),imshow(f);
165 whos f;%其实是个二值图像
166 title('人体染色体原始图');
167
168 f=im2double(f);%不再是二值图像了
169 h=fspecial('gaussian',15,9);
170 g=imfilter(f,h,'replicate');
171 g=im2bw(g,0.7);%其中的graythresh()函数为自动阈值,处理后g为二值图像了
172 subplot(332),imshow(g);
173 title('经阈值处理后');
174
175 b=boundaries(g);
176 b=b{1};
177 bim=bound2im(b,344,270);%344,270为图片的尺寸
178 subplot(333),imshow(bim);
179 title('提取的边界图');
180
181 z=frdescp(b);%对边界坐标b进行傅里叶变换系数(有多少个点就有多少个系数),将b的坐标点看成是复平面中的某个复数
182 z14=ifrdescp(z,546);%用50%的描述子进行逆变换
183 z546im=bound2im(z14,344,270);
184 subplot(334),imshow(z546im);
185 title('546个描述子恢复后');
186
187 z110=ifrdescp(z,110);%用50%的描述子进行逆变换
188 z110im=bound2im(z110,344,270);
189 subplot(335),imshow(z110im);
190 title('110个描述子恢复后');
191
192 z56=ifrdescp(z,56);%用50%的描述子进行逆变换
193 z56im=bound2im(z56,344,270);
194 subplot(336),imshow(z56im);
195 title('56个描述子恢复后');
196
197 z28=ifrdescp(z,28);%用50%的描述子进行逆变换
198 z28im=bound2im(z28,344,270);
199 subplot(337),imshow(z28im);
200 title('28个描述子恢复后');
201
202 z14=ifrdescp(z,14);%用50%的描述子进行逆变换
203 z14im=bound2im(z14,344,270);
204 subplot(338),imshow(z14im);
205 title('14个描述子恢复后');
206
207 z8=ifrdescp(z,8);%用50%的描述子进行逆变换
208 z8im=bound2im(z8,344,270);
209 subplot(339),imshow(z8im);
210 title('8个描述子恢复后');
211%傅里叶描述子图过程:


212
213 %% 函数regionprops的使用
214 clc
215 clear
216 B = [1 1 1 0 0 0 0 0
217 1 1 1 0 1 1 0 0
218 1 1 1 0 1 1 0 0
219 1 1 1 0 0 0 0 0
220 1 1 1 0 0 0 1 0
221 1 1 1 0 0 0 1 0
222 1 1 1 0 0 1 1 0
223 1 1 1 0 0 0 0 0];
224 B=bwlabel(B);%将矩阵B变成标记矩阵,默认为8连通,将连通区域用不同的整数给标记出来
225
226 D=regionprops(B,'area','boundingbox')
227 w=[D.Area]
228 NR=length(w)
229
230 V=cat(1,D.BoundingBox)%cat函数表示沿着第一个方向级联存放数组
231%该函数的使用过程如下:


232
233 %% 统计纹理度量(木有搞定)
234 clc
235 clear
236 f=imread('..\images\dipum_images_ch11\Fig1119(a)(superconductor-with-window).tif');
237 sA=imhist(f);
238 statxture(sA)%计算图像的6个统计度量量
239
240 %% 计算频谱纹理
241 clc
242 clear
243
244 f1=imread('..\images\dipum_images_ch11\Fig1121(a)(random_matches).tif');
245 subplot(241),imshow(f1);
246 title('随机排列火柴');
247
248 f2=imread('..\images\dipum_images_ch11\Fig1121(b)(ordered_matches).tif');
249 subplot(245),imshow(f2);
250 title('规则排列火柴');
251
252 [srad1,sang1,s1]=specxture(f1);
253 subplot(242),imshow(s1);
254 title('随机排列火柴频谱图');
255
256 subplot(243),plot(srad1);
257 title('随机排列火柴幅度图');
258
259 subplot(244),plot(sang1);
260 title('随机排列火柴角度图');
261
262 [srad2,sang2,s2]=specxture(f2);
263 subplot(246),imshow(s2);
264 title('规则排列火柴频谱图');
265
266 subplot(247),plot(srad2);
267 title('规则排列火柴幅度图');
268
269 subplot(248),plot(sang2);
270 title('规则排列火柴角度图');
271%计算频谱纹理图过程如下:


272
273 %% 不变矩
274 clc
275 clear
276 f=imread('..\images\dipum_images_ch11\Fig1123(a)(Original_Padded_to_568_by_568).tif');
277 imshow(f);
278 title('不变矩原图');
279%不变矩原始图像如下:


280
281 %padarray函数为矩阵0填充,图像原始大小为400*400,但是需要旋转45度,所以最大变为568*568
282 fp=padarray(f,[84 84],'both');%因此四个方向每个方向都要扩充84个像素
283
284 fhs=f(1:2:end,1:2:end);%表示从第一个到最后一个,每2个采样一次
285 figure,subplot(221),imshow(fhs);
286 title('半大图像');
287 fhsp=padarray(fhs,[184 184],'both');
288
289 fm=fliplr(f);
290 subplot(222),imshow(fm);
291 title('镜面图像');
292 fmp=padarray(fm,[84 84],'both');
293
294 fr2=imrotate(f,2,'bilinear');%旋转2度,双线性插值
295 subplot(223),imshow(fr2);
296 title('旋转2度后图像');
297
298 fr45=imrotate(f,45,'bilinear');
299 subplot(224),imshow(fr45);
300 title('旋转45度后图像');
301%不变矩变换图像如下:


302
303 phiorig=abs(log(invmoments(f)))
304 phihalf=abs(log(invmoments(fhs)))
305 phimirror=abs(log(invmoments(fm)))
306 phirot2=abs(log(invmoments(fr2)))
307 phirot45=abs(log(invmoments(fr45)))
308%这5幅图像的7个不变矩图像如下:


309
310 %图像主分量
311 clc
312 clear
313 f1=imread('..\images\dipum_images_ch11\Fig1125(a)(WashingtonDC_Band1_512).tif');
314 subplot(231),imshow(f1);
315 title('可见蓝色波段图像');
316
317 f2=imread('..\images\dipum_images_ch11\Fig1125(b)(WashingtonDC_Band2_512).tif');
318 subplot(232),imshow(f2);
319 title('可见绿色波段图像');
320
321 f3=imread('..\images\dipum_images_ch11\Fig1125(c)(WashingtonDC_Band3_512).tif');
322 subplot(233),imshow(f3);
323 title('可见红色波段图像');
324
325 f4=imread('..\images\dipum_images_ch11\Fig1125(d)(WashingtonDC_Band4_512).tif');
326 subplot(234),imshow(f4);
327 title('近红外光图像');
328
329 f5=imread('..\images\dipum_images_ch11\Fig1125(e)(WashingtonDC_Band5_512).tif');
330 subplot(235),imshow(f5);
331 title('中红外光图像');
332
333 f6=imread('..\images\dipum_images_ch11\Fig1125(f)(WashingtonDC_Band6_512).tif');
334 subplot(236),imshow(f6);
335 title('热红外光图像');
336%6幅多谱图像显示结果如下:


337
338 S=cat(3,f1,f2,f3,f4,f5,f6);%将6幅图像在第三维中连接起来
339 [X,R]=imstack2vectors(S);%将S组织到数组X中
340 p=princomp(X,6);%princomp用来获得6幅图像的主分量
341
342 g1=p.Y(:,1);%计算第1主维分量
343 g1=reshape(g1,512,512);
344 figure,subplot(231),imshow(g1,[]);
345
346 g2=p.Y(:,2);%计算第2主维分量
347 g2=reshape(g2,512,512);
348 subplot(232),imshow(g2,[]);
349
350 g3=p.Y(:,3);%计算第3主维分量
351 g3=reshape(g3,512,512);
352 subplot(233),imshow(g3,[]);
353
354 g4=p.Y(:,4);%计算第4主维分量
355 g4=reshape(g4,512,512);
356 subplot(234),imshow(g4,[]);
357
358 g5=p.Y(:,5);%计算第5主维分量
359 g5=reshape(g5,512,512);
360 subplot(235),imshow(g5,[]);
361
362 g6=p.Y(:,6);%计算第6主维分量
363 g6=reshape(g6,512,512);
364 subplot(236),imshow(g6,[]);
365%主分量分解过程如下:

 

     图像的表示和描述对后面做计算机视觉处理起来很大的作用。当然针对不同的研究对象其表示和描述的特点也不同。

 

 

 

作者:tornadomeet 出处:http://www.cnblogs.com/tornadomeet 欢迎转载或分享,但请务必声明文章出处。 (新浪微博:tornadomeet,欢迎交流!)
posted on 2012-04-01 18:58  tornadomeet  阅读(5960)  评论(7编辑  收藏  举报

阿萨德发斯蒂芬