-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathgetBackbone.m
116 lines (100 loc) · 2.76 KB
/
getBackbone.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
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
function [bbSubs, bbLen, bbImg]=getBackbone(skelImg,somabw)
% "getBackbone" is used to get the longest path passing the pollenPos, from a connected skeleton bw image.
% It's like "getLongestPath" but there is difference: the path must pass the pollenPos.
% "bbSubs" should start with the end point closest to pollen grain, e.g.,
% bbSubs is ordered.
global gImg;
if isempty(find(skelImg,1))
error('Error: The skelImg is all black!');
end
[D vertices]=getDistMat(skelImg);
% If D==0, which means there is only one point in skelImg.
if length(D)==1
bbSubs=vertices(1,2:3);
bbLen=1;
bbImg=skelImg;
return;
end
% Find the end point closest to pollenPos, except for shortEP.
shortDist=inf;
% gImg=skelImg;
epIdx=0;
for i=1:size(vertices,1)
% nbr=nbr8(vertices(i,2:3));
% if size(nbr,1)==1 % end point vertex.
if vertices(i,4) && ~vertices(i,5) % is end point and not shortEP.
dist=euDist(vertices(i,2:3),pollenPos);
if dist<shortDist
shortDist=dist;
epIdx=i;
end
end
end
if epIdx==0 % If the pollenPos can't specify an end point, the longest path is used.
% [Y I]=max(D(:));
% bbLen=Y;
% [row col]=ind2sub(size(D),I);
% sp=vertices(row,2:3);
% ep=vertices(col,2:3);
error('epIdx is zero!!!');
else % Choose the longest path passing pollenPos.
[Y I]=max(D(:,epIdx));
bbLen=Y;
sp=vertices(epIdx,2:3); % The sp is the end point closest to pollen grain.
ep=vertices(I,2:3);
end
%% Get bbImg.
gImg=skelImg;
for i=1:size(vertices,1)
if i==I || i==epIdx
continue;
end
nbrs=nbr8(vertices(i,2:3));
if nbrs(1)~=0
tempImg=gImg; % Backup gImg.
% img(vertices(i,2),vertices(i,3))=0;
epTrace=traceToEJ(vertices(i,2:3),0);
gImg(epTrace(1),epTrace(2))=1; % Need to keep the ep.
if size(nbrs,1)==1 % If it's end point, no connectness check is needed.
continue;
end
% If backbone's two end points is not connected, reverse the
% erasing.
L=bwlabel(gImg,8);
if L(sp(1),sp(2))~=L(ep(1),ep(2))
gImg=tempImg; % restore gImg.
end
end
end
bbImg=gImg;
% Get the backbone indices sequence.
% Now img is the backbone img.
% Since sp is the closest end point to pollen grain, the bbSubs will start
% at the end where pollen grain resides.
bbSubs=getPathSubs(sp);
end
%%%%%
function subs=getPathSubs(sp)
% Get the path pixel subcripts from start point.
% subs: subs [row col] for backbone pixels in connection order, which is good for tracing.
% Now gImg is the backbone img.
global gImg;
subs=zeros(2,2);
gImg(sp(1),sp(2))=0;
len=1;
subs(len,:)=sp;
nbr1=nbr8(sp);
while nbr1(1)~=0
if size(nbr1,1)~=1
% When Ren-shape joint is met, first trace to its 4-nbr, then 8-nbr.
% fprintf(1,'Error: There are %d nbrs at sp
% %d\t%d.\n',size(nbr1,1),sp(1),sp(2));
nbr1=nbr1(1,:);
end
gImg(nbr1(1),nbr1(2))=0;
len=len+1;
subs(len,:)=nbr1;
% sp=nbr1;
nbr1=nbr8(nbr1);
end
end