-
Notifications
You must be signed in to change notification settings - Fork 12
/
Test_Skel2Graph3D.m
77 lines (59 loc) · 1.87 KB
/
Test_Skel2Graph3D.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
clear all;
close all;
% load example binary skeleton image
load skel
w = size(skel,1);
l = size(skel,2);
h = size(skel,3);
% initial step: condense, convert to voxels and back, detect cells
[~,node,link] = Skel2Graph3D(skel,0);
% total length of network
wl = sum(cellfun('length',{node.links}));
skel2 = Graph2Skel3D(node,link,w,l,h);
[~,node2,link2] = Skel2Graph3D(skel2,0);
% calculate new total length of network
wl_new = sum(cellfun('length',{node2.links}));
% iterate the same steps until network length changed by less than 0.5%
while(wl_new~=wl)
wl = wl_new;
skel2 = Graph2Skel3D(node2,link2,w,l,h);
[A2,node2,link2] = Skel2Graph3D(skel2,0);
wl_new = sum(cellfun('length',{node2.links}));
end;
% display result
figure();
hold on;
for i=1:length(node2)
x1 = node2(i).comx;
y1 = node2(i).comy;
z1 = node2(i).comz;
if(node2(i).ep==1)
ncol = 'c';
else
ncol = 'y';
end;
for j=1:length(node2(i).links) % draw all connections of each node
if(node2(node2(i).conn(j)).ep==1)
col='k'; % branches are black
else
col='r'; % links are red
end;
if(node2(i).ep==1)
col='k';
end;
% draw edges as lines using voxel positions
for k=1:length(link2(node2(i).links(j)).point)-1
[x3,y3,z3]=ind2sub([w,l,h],link2(node2(i).links(j)).point(k));
[x2,y2,z2]=ind2sub([w,l,h],link2(node2(i).links(j)).point(k+1));
line([y3 y2],[x3 x2],[z3 z2],'Color',col,'LineWidth',2);
end;
end;
% draw all nodes as yellow circles
plot3(y1,x1,z1,'o','Markersize',9,...
'MarkerFaceColor',ncol,...
'Color','k');
end;
axis image;axis off;
set(gcf,'Color','white');
drawnow;
view(-17,46);