rsvs3D  0.0.0
Codes for the c++ implementation of the 3D RSVS
RSVSalgorithm.cpp
1 #include "RSVSalgorithm.hpp"
2 
3 #include <cstdlib>
4 #include <iostream>
5 #include <vector>
6 
7 #include "snake.hpp"
8 #include "snakeengine.hpp"
9 #include "warning.hpp"
10 
11 using namespace std;
12 
13 std::vector<int> FindSpawnVerts(const mesh &meshin, vector<int> &vertList, vector<int> &voluOutList, int outerBorder)
14 {
15  // Function which identifies spawn points
16  // Spawn points are:
17  // - Any point part of a cell which touches the void
18  // Put that is not part of a surface that is on the void
19  // itself.
20  // - Points That are on the border of two cells one with some VOS
21  // One without.
22 
23  int ii, ni, jj, nj, surfSurb, vert_i;
24  vector<int> offBorderVert, internalSurf, voluIndIntern;
25 
26  if (outerBorder == 1)
27  { // spawn at outer borders
28  meshin.GetOffBorderVert(offBorderVert, voluOutList, 1);
29  meshin.ElmOnParentBound(internalSurf, voluIndIntern, false, true);
30  }
31  else if (outerBorder == 0)
32  {
33  meshin.GetOffBorderVert(offBorderVert, voluOutList, 0);
34  meshin.ElmOnParentBound(internalSurf, voluIndIntern, false, false);
35  }
36  else
37  {
38  RSVS3D_ERROR_ARGUMENT("outerBorder has an unknown value");
39  }
40  ni = meshin.verts.size();
41  vertList.reserve(ni);
42  vertList.clear();
43  ni = internalSurf.size();
44 
45  cout << "Vertices : " << ni << " (interal), ";
46  for (ii = 0; ii < ni; ++ii)
47  {
48  surfSurb = meshin.surfs.find(internalSurf[ii]);
49  nj = meshin.surfs(surfSurb)->edgeind.size();
50  for (jj = 0; jj < nj; jj++)
51  {
52  vert_i = meshin.edges.isearch(meshin.surfs(surfSurb)->edgeind[jj])->vertind[0];
53  if (!meshin.verts.isearch(vert_i)->isBorder)
54  {
55  vertList.push_back(meshin.verts.isearch(vert_i)->index);
56  }
57  vert_i = meshin.edges.isearch(meshin.surfs(surfSurb)->edgeind[jj])->vertind[1];
58  if (!meshin.verts.isearch(vert_i)->isBorder)
59  {
60  vertList.push_back(meshin.verts.isearch(vert_i)->index);
61  }
62  }
63  }
64 
65  sort(offBorderVert);
66  unique(offBorderVert);
67  ni = offBorderVert.size();
68  cout << ni << " (border); " << endl;
69  for (ii = 0; ii < ni; ++ii)
70  {
71  vertList.push_back(offBorderVert[ii]);
72  }
73 
74  sort(vertList);
75  unique(vertList);
76 
77  ni = voluIndIntern.size();
78  for (ii = 0; ii < ni; ++ii)
79  {
80  voluOutList.push_back(voluIndIntern[ii]);
81  }
82 
83  sort(voluOutList);
84  unique(voluOutList);
85  return (offBorderVert);
86 }
87 
88 void SpawnSnakeAndMove(snake &snakein, std::vector<int> vertSpawn)
89 {
90  vector<int> isImpact;
91  vector<double> dt;
92  int nVe, nE, nS, nVo;
93  snakein.snakeconn.size(nVe, nE, nS, nVo);
94  int nNew = vertSpawn.size();
95  // cout << "vertices to output " << ni << endl;
96  snakein.reserve(nVe + nNew * 15, nE + nNew * 15, nS + nNew * 15, nVo + nNew * 1);
97  for (int ii = 0; ii < nNew; ++ii)
98  {
99  SpawnAtVertex(snakein, vertSpawn[ii]);
100  // cout << ii << " " ;
101  }
102  // cout << " vertices Dones" << endl;
103  // Move to half distances
104  int nTotSnax = snakein.snaxs.size();
105  for (int ii = 0; ii < nVe; ++ii)
106  {
107  snakein.snaxs[ii].v = 0.0;
108  snakein.snaxs[ii].d = 0.5;
109  }
110  for (int ii = nVe; ii < nTotSnax; ++ii)
111  {
112  snakein.snaxs[ii].v = 1.0;
113  }
114  snakein.CalculateTimeStep(dt, 0.51);
115  snakein.UpdateDistance(dt);
116  snakein.PrepareForUse();
117  snakein.UpdateCoord();
118  snakein.PrepareForUse();
119  snakein.SnaxImpactDetection(isImpact);
120  MergeAllContactVertices(snakein, isImpact);
121  snakein.PrepareForUse();
122  CleanupSnakeConnec(snakein);
123  snakein.PrepareForUse();
124 }
125 
126 void SpawnRSVS(snake &snakein, int outerBorder)
127 {
128  // Function which handles
129  // - spawning
130  // - growing the snake
131  // - Identifying snaxels to remove
132  // - update connectivity:
133  // + find snaxsurfs in invalid snake volus(surfs)
134  // + find snaxedges in these surfs
135  // + find snaxels
136  // * invalid snakevolus=[border cell, empty cell]
137  int ii, ni;
138  vector<int> vertSpawn, borderVertSpawn;
139  vector<int> voluSnaxDelete;
140  vector<int> isImpact;
141  vector<double> dt;
142  snakein.snakemesh()->SetBorders();
143  // snakein.snakemesh()->disp();
144  borderVertSpawn = FindSpawnVerts(*(snakein.snakemesh()), vertSpawn, voluSnaxDelete, outerBorder);
145  ni = vertSpawn.size();
146  // cout << "vertices to output " << ni << endl;
147  snakein.reserve(ni * 15, ni * 15, ni * 15, ni * 15);
148  for (ii = 0; ii < ni; ++ii)
149  {
150  SpawnAtVertex(snakein, vertSpawn[ii]);
151  // cout << ii << " " ;
152  }
153  // cout << " vertices Dones" << endl;
154  // Move to half distances
155  ni = snakein.snaxs.size();
156  for (ii = 0; ii < ni; ++ii)
157  {
158  snakein.snaxs[ii].v = 1.0;
159  }
160  snakein.CalculateTimeStep(dt, 0.6);
161  snakein.UpdateDistance(dt);
162  snakein.PrepareForUse();
163  snakein.UpdateCoord();
164  snakein.PrepareForUse();
165  // snakein.displight();
166  // detect and remove impacts
167  snakein.SnaxImpactDetection(isImpact);
168  MergeAllContactVertices(snakein, isImpact);
169  snakein.PrepareForUse();
170  CleanupSnakeConnec(snakein);
171  snakein.PrepareForUse();
172  // snakein.displight();
173  // Remove one of the 'snakes'
174  if (snakein.Check3D())
175  {
176  RemoveSnakeInVolu(snakein, voluSnaxDelete, outerBorder);
177  }
178  else
179  {
180  RemoveSnakeInSurf(snakein, voluSnaxDelete, outerBorder);
181  }
182  snakein.PrepareForUse();
183  // snakein.displight();
184  // Second spawn phase
185  if (outerBorder == 1)
186  {
187  SpawnSnakeAndMove(snakein, borderVertSpawn);
188  }
189 
190  snakein.OrientFaces();
191  cout << "Initialisation DONE!" << endl;
192 }
193 
194 void RemoveSnakeInVolu(snake &snakein, vector<int> &voluInd, int outerBorder)
195 {
196  int ii, ni, jj, nj, nBlocks;
197  vector<int> delSurf, delEdge, delSnax, tempSurf, tempEdge, tempSnax, subSurf, subEdge, vertBlocks;
198  vector<bool> isBlockDel;
199 
200  delSurf.reserve(snakein.snaxsurfs.size());
201  delEdge.reserve(snakein.snaxedges.size());
202  delSnax.reserve(snakein.snaxs.size());
203 
204  ni = voluInd.size();
205 
206  for (ii = 0; ii < ni; ++ii)
207  {
208  tempSurf.clear();
209  snakein.snaxsurfs.findsiblings(voluInd[ii], tempSurf);
210  tempEdge = ConcatenateVectorField(snakein.snakeconn.surfs, &surf::edgeind, tempSurf);
211  subEdge = snakein.snakeconn.edges.find_list(tempEdge);
212 
213  tempSnax = ConcatenateVectorField(snakein.snakeconn.edges, &edge::vertind, subEdge);
214  nj = tempSnax.size();
215  for (jj = 0; jj < nj; ++jj)
216  {
217  delSnax.push_back(tempSnax[jj]);
218  }
219  }
220  // cout << "nSnax " << delSnax.size() << endl;
221  // cout << "Find snax to del" << endl;
222  vertBlocks.clear();
223  nBlocks = snakein.snakeconn.ConnectedVertex(vertBlocks);
224  bool delVertInVolu = outerBorder > 0;
225  isBlockDel.assign(nBlocks, !delVertInVolu);
226  ni = delSnax.size();
227  for (ii = 0; ii < ni; ++ii)
228  {
229  isBlockDel[vertBlocks[snakein.snakeconn.verts.find(delSnax[ii])] - 1] = delVertInVolu;
230  }
231  // cout << "Find All snax to del" << endl;
232  delSnax.clear();
233  delSurf.clear();
234  delEdge.clear();
235  ni = vertBlocks.size();
236 
237  for (ii = 0; ii < ni; ++ii)
238  {
239  if (isBlockDel[vertBlocks[ii] - 1])
240  {
241  delSnax.push_back(snakein.snaxs(ii)->index);
242 
243  subEdge = snakein.snakeconn.edges.find_list(snakein.snakeconn.verts(ii)->edgeind);
244  nj = snakein.snakeconn.verts(ii)->edgeind.size();
245  for (jj = 0; jj < nj; ++jj)
246  {
247  delEdge.push_back(snakein.snakeconn.verts(ii)->edgeind[jj]);
248  }
249 
250  tempSurf = ConcatenateVectorField(snakein.snakeconn.edges, &edge::surfind, subEdge);
251  nj = tempSurf.size();
252  for (jj = 0; jj < nj; ++jj)
253  {
254  delSurf.push_back(tempSurf[jj]);
255  }
256  }
257  }
258 
259  ni = delSurf.size();
260  for (ii = 0; ii < ni; ++ii)
261  {
262  snakein.snakeconn.RemoveIndex(3, delSurf[ii]);
263  }
264 
265  sort(delSnax);
266  unique(delSnax);
267  sort(delEdge);
268  unique(delEdge);
269  sort(delSurf);
270  unique(delSurf);
271 
272  // ni=delEdge.size();
273  // for(ii=0; ii<ni; ++ii){snakein.snakeconn.RemoveIndex(2,delEdge[ii]);}
274 
275  // snakein.displight();
276  snakein.snaxs.remove(delSnax);
277  snakein.snaxedges.remove(delEdge);
278  snakein.snaxsurfs.remove(delSurf);
279  snakein.snakeconn.verts.remove(delSnax);
280  snakein.snakeconn.edges.remove(delEdge);
281  snakein.snakeconn.surfs.remove(delSurf);
282 
283  snakein.snakeconn.TightenConnectivity();
284  snakein.HashArray();
285  snakein.snakeconn.TestConnectivityBiDir(__PRETTY_FUNCTION__);
286  snakein.ForceCloseContainers();
287  snakein.PrepareForUse();
288  // snakein.displight();
289  if (outerBorder > 0)
290  {
291  snakein.Flip();
292  }
293  // cout << "Before Assignement of internal verts" << endl;
294  snakein.AssignInternalVerts();
295  // cout << "After Assignement of internal verts" << endl;
296 }
297 
298 void RemoveSnakeInSurf(snake &snakein, vector<int> &voluInd, int outerBorder)
299 {
300  int ii, ni, jj, nj, nBlocks;
301  vector<int> delSurf, delEdge, delSnax, tempSurf, tempEdge, tempSnax, subSurf, subEdge, vertBlocks;
302  vector<bool> isBlockDel;
303 
304  delSurf.reserve(snakein.snaxsurfs.size());
305  delEdge.reserve(snakein.snaxedges.size());
306  delSnax.reserve(snakein.snaxs.size());
307 
308  ni = voluInd.size();
309 
310  for (ii = 0; ii < ni; ++ii)
311  {
312  tempSurf.clear();
313  snakein.snaxedges.findsiblings(voluInd[ii], subEdge);
314  // tempEdge=ConcatenateVectorField(snakein.snakeconn.surfs,
315  // &surf::edgeind,tempSurf);
316  // subEdge=snakein.snakeconn.edges.find_list(tempEdge);
317 
318  tempSnax = ConcatenateVectorField(snakein.snakeconn.edges, &edge::vertind, subEdge);
319  nj = tempSnax.size();
320  for (jj = 0; jj < nj; ++jj)
321  {
322  delSnax.push_back(tempSnax[jj]);
323  }
324  }
325  // cout << "nSnax " << delSnax.size() << endl;
326  // cout << "Find snax to del" << endl;
327  vertBlocks.clear();
328  nBlocks = snakein.snakeconn.ConnectedVertex(vertBlocks);
329  bool delVertInVolu = outerBorder > 0;
330  isBlockDel.assign(nBlocks, !delVertInVolu);
331  ni = delSnax.size();
332  for (ii = 0; ii < ni; ++ii)
333  {
334  isBlockDel[vertBlocks[snakein.snakeconn.verts.find(delSnax[ii])] - 1] = delVertInVolu;
335  }
336  // cout << "Find All snax to del" << endl;
337  delSnax.clear();
338  delSurf.clear();
339  delEdge.clear();
340  ni = vertBlocks.size();
341 
342  for (ii = 0; ii < ni; ++ii)
343  {
344  if (isBlockDel[vertBlocks[ii] - 1])
345  {
346  delSnax.push_back(snakein.snaxs(ii)->index);
347 
348  subEdge = snakein.snakeconn.edges.find_list(snakein.snakeconn.verts(ii)->edgeind);
349  nj = snakein.snakeconn.verts(ii)->edgeind.size();
350  for (jj = 0; jj < nj; ++jj)
351  {
352  delEdge.push_back(snakein.snakeconn.verts(ii)->edgeind[jj]);
353  }
354 
355  // tempSurf=ConcatenateVectorField(snakein.snakeconn.edges,
356  // &edge::surfind,subEdge); nj=tempSurf.size(); for(jj=0; jj<nj; ++jj){
357  // delSurf.push_back(tempSurf[jj]);
358  // }
359  }
360  }
361 
362  ni = delEdge.size();
363  for (ii = 0; ii < ni; ++ii)
364  {
365  snakein.snakeconn.RemoveIndex(2, delEdge[ii]);
366  }
367 
368  sort(delSnax);
369  unique(delSnax);
370  sort(delEdge);
371  unique(delEdge);
372  // sort(delSurf);
373  // unique(delSurf);
374 
375  // ni=delEdge.size();
376  // for(ii=0; ii<ni; ++ii){snakein.snakeconn.RemoveIndex(2,delEdge[ii]);}
377 
378  // snakein.displight();
379  snakein.snaxs.remove(delSnax);
380  snakein.snaxedges.remove(delEdge);
381  // snakein.snaxsurfs.remove(delSurf);
382  snakein.snakeconn.verts.remove(delSnax);
383  snakein.snakeconn.edges.remove(delEdge);
384  // snakein.snakeconn.surfs.remove(delSurf);
385 
386  snakein.snakeconn.TightenConnectivity();
387  snakein.HashArray();
388  snakein.snakeconn.TestConnectivityBiDir(__PRETTY_FUNCTION__);
389  snakein.ForceCloseContainers();
390  snakein.PrepareForUse();
391  // snakein.displight();
392  if (outerBorder > 0)
393  {
394  snakein.Flip();
395  }
396  // cout << "Before Assignement of internal verts" << endl;
397  snakein.AssignInternalVerts();
398  // cout << "After Assignement of internal verts" << endl;
399 }
Functions which are part of the RSVS algorithm but not core to the snaking process.
Class for mesh handling.
Definition: mesh.hpp:592
int ConnectedVertex(std::vector< int > &vertBlock) const
Return in a vector for each vertex a block number which it is part of.
Definition: mesh.cpp:4633
Definition: snake.hpp:83
Provides the core restricted surface snake container.
Functions needed to evolve the r-surface snake.
Provides the error and warning system used by the RSVS3D project.
#define RSVS3D_ERROR_ARGUMENT(M)
Throw a invalid_argument.
Definition: warning.hpp:148