2

I have a large 3D triangulation (~2000 faces, ~1000 vertices), and I want to "crop" it, i.e. create a new tirangulation contains only the edges connecting a given list of vertex indices.

I work in MATLAB, and I tried to export the points and connectivity matrix to variables, and crop them manually:

% tri is the triangulation to be cropped, and inds is the vector of vertex 
% indices to remain in the final triangulation.

%% Create list of final vertices
p = tri.Points(inds, :);

%% Create new connectivity matrix
% start by cropping the connectivity matrix to contain only the final vertices
conn = [];
for j = 1:height(tri.ConnectivityList)
    if all( ismember(tri.ConnectivityList(j,:),inds) )
        conn = [conn; tri.ConnectivityList(j,:)];
    end
end

% then fix the indices to the new ones, since the index of each point changed
for i = 1:height(conn)
    for j = 1:width(conn)
        conn(i,j) = find(inds==conn(i,j),1);
    end
end

%% Create final triangulation
tri = triangulation(conn, p);

This solution works, but it is rather slow. Specifically, the cropping part using the ismember() function is a strong bottleneck.

After searching for a bit I found out somtimes people use containers.Map() to do similar jobs, but I couldn't understand if and how is this relevant to my case.

Any help is appreciated!

EDIT: here is an example conn matrix to work with:

1   2   3
4   2   5
2   4   3
6   5   7
8   9   10
11  12  7
13  10  14
15  16  17
14  16  15
18  19  17
20  21  22
23  24  25
24  26  25
27  21  24
20  28  21
27  24  23
29  13  12
10  23  14
14  23  25
30  29  31
31  11  1
30  31  32
33  34  8
34  35  9
11  7   2
7   36  6
17  19  37
36  17  37
38  39  40
41  22  39
41  42  22
42  20  22
43  44  28
45  20  42
40  39  35
39  22  35
46  33  8
47  48  49
47  50  51
49  52  50
53  54  18
55  53  18
55  56  53
57  58  53
53  58  59
25  56  55
25  26  56
24  60  61
62  63  57
56  64  57
65  66  67
67  68  65
68  69  62
68  62  64
61  64  56
65  64  61
69  70  62
19  18  48
48  54  49
71  72  49
71  73  72
71  58  73
59  71  54
59  58  71
73  74  75
76  77  74
63  78  79
37  47  80
2   7   5
81  82  6
81  6   80
83  81  80
5   6   82
46  84  85
29  84  13
85  84  29
84  46  13
46  8   13
8   10  13
8   34  9
7   12  36
36  12  15
12  14  15
13  14  12
36  15  17
14  25  16
17  16  18
16  25  18
9   86  10
10  86  27
10  27  23
9   22  86
86  21  27
22  21  86
24  21  60
21  28  60
29  12  11
87  31  1
87  32  31
31  29  11
30  85  29
88  39  89
40  35  34
11  2   1
9   35  22
38  89  39
41  39  88
20  45  43
41  88  42
43  28  20
47  37  48
51  80  47
47  49  50
49  72  52
18  25  55
53  59  54
56  26  61
60  28  66
61  60  65
60  66  65
90  57  63
64  62  57
64  65  68
26  24  61
63  91  92
62  70  78
63  62  78
63  92  90
53  56  57
58  57  90
37  19  48
18  54  48
49  54  71
73  58  74
74  58  76
58  90  76
75  74  77
91  63  79
80  51  83
44  66  28
6   36  80
80  36  37

pt is a points array of 92 points in 3D (92 by 3 matrix).

DeadlosZ
  • 165
  • 5
  • Are you able to share the typical/approximate size of your `tri.ConnectivityList` as well as the typical length of `inds` compared to your original number of vertices? – horchler Jun 30 '23 at 15:29
  • 1
    I think that you should be able to replace your first `for` loop with a single vectorized call to `ismember` and logical indexing: `conn = tri.ConnectivityList(all(ismember(tri.ConnectivityList, inds), 2), :);`. If you can verify that works and provide any other details to help me assess your case vis-a-vis performance I'll try to write an answer with some thoughts on optimization. – horchler Jun 30 '23 at 16:08
  • This looks promising, I'll try that soon. Also, I added the example `conn` matrix you requested. – DeadlosZ Jun 30 '23 at 18:09
  • @horchler I tried it, and it worked perfectly! About 60 times faster. You may want to write an answer so I can accpept it? – DeadlosZ Jun 30 '23 at 18:27
  • to DeadlosZ: depending on the choice of nodes to remove and the initial network, for instance the surface of a ball, chop-all first and then rebuild all the open connections may have more than one solution, so it's a way of cropping netwotks that may add uncertainty. But if you interleave chopping and rebuilding, it may take more time, but you avoid ending up with nodes apart at same distance and then having to choose, or equivalently, having more that 1 solution. – John BG Jun 30 '23 at 19:50

1 Answers1

2

Your first for loop can be replaced by a vectorized call to ismember and use of logical indexing:

conn = tri.ConnectivityList(all(ismember(tri.ConnectivityList, inds), 2), :);

In this case, all is applied across the columns of the logical output from ismember. This results in a logical vector the same length as tri.ConnectivityList that can be used to extract just the valid rows.

There are likely other ways to accomplish this. Notably, more general code might check if the number of points being "cropped" is more or less than half the number of original points. In that case, it might be more efficient to use logic to "crop-in" or "crop-out" points.

horchler
  • 18,384
  • 4
  • 37
  • 73