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
|
function loops = extractloops(edges)
%
% loops=extractloops(edges)
%
% extract individual loop or polyline segment from a collection of edges
%
% author: Qianqian Fang, <q.fang at neu.edu>
% date: 2007/11/21
%
% input:
% edges: two column matrix recording the starting/ending
% points of all edge segments
%
% output:
% loops: output, a single vector separated by NaN, each segment
% is a 3D polyline or loop consisted of node IDs
%
% example:
% edges=[1 2;2 3;1 4;3 4;7 3;1 9;5 6;6 7;10 9; 8 10;1 8;9 3;11 11;11 12];
% loops=extractloops(edges)
%
% -- this function is part of iso2mesh toolbox (http://iso2mesh.sf.net)
%
loops = [];
edges(edges(:, 1) == edges(:, 2), :) = []; % remove degenerated edges
loops = [loops, edges(1, :)];
loophead = edges(1, 1);
loopend = edges(1, end);
edges(1, :) = [];
while (~isempty(edges))
idx = [find(edges(:, 1) == loopend)', find(edges(:, 2) == loopend)'];
if (length(idx) > 1) % when a node with multiple connection found
idx = idx(1); % take the first connection and continue
end
if (isempty(idx)) % when an open-line segment gets to one end
% when both open ends are found
if (isempty([find(edges(:, 1) == loophead)', find(edges(:, 2) == loophead)']))
loops = [loops, nan];
loops = [loops, edges(1, :)];
loophead = edges(1, 1);
loopend = edges(1, end);
edges(1, :) = [];
else % only the first open end is found, flip and trace the other
[loophead, loopend] = deal(loopend, loophead);
lp = fliplr(loops);
seg = find(isnan(lp), 1);
if (isempty(seg))
loops = lp;
else
loops = [loops(1:end - seg(1) + 1) lp(1:seg(1) - 1)];
end
end
continue
end
if (length(idx) == 1) % tracing along a single line thread
idx = idx(1);
ed = edges(idx, :);
ed(ed == loopend) = [];
newend = ed(1);
if (newend == loophead) % when a loop is found
loops = [loops loophead nan];
edges(idx, :) = [];
if (size(edges, 1) == 0)
break
end
loops = [loops, edges(1, :)];
loophead = edges(1, 1);
loopend = edges(1, end);
edges(1, :) = [];
continue
else
loops = [loops, newend];
end
loopend = newend;
edges(idx, :) = [];
end
end
|
