aimsalgo 6.0.0
Neuroimaging image processing
meshmorphomat_d.h
Go to the documentation of this file.
1/* This software and supporting documentation are distributed by
2 * Institut Federatif de Recherche 49
3 * CEA/NeuroSpin, Batiment 145,
4 * 91191 Gif-sur-Yvette cedex
5 * France
6 *
7 * This software is governed by the CeCILL-B license under
8 * French law and abiding by the rules of distribution of free software.
9 * You can use, modify and/or redistribute the software under the
10 * terms of the CeCILL-B license as circulated by CEA, CNRS
11 * and INRIA at the following URL "http://www.cecill.info".
12 *
13 * As a counterpart to the access to the source code and rights to copy,
14 * modify and redistribute granted by the license, users are provided only
15 * with a limited warranty and the software's author, the holder of the
16 * economic rights, and the successive licensors have only limited
17 * liability.
18 *
19 * In this respect, the user's attention is drawn to the risks associated
20 * with loading, using, modifying and/or developing or reproducing the
21 * software by the user in light of its specific status of free software,
22 * that may mean that it is complicated to manipulate, and that also
23 * therefore means that it is reserved for developers and experienced
24 * professionals having in-depth computer knowledge. Users are therefore
25 * encouraged to load and test the software's suitability as regards their
26 * requirements in conditions enabling the security of their systems and/or
27 * data to be ensured and, more generally, to use and operate it in the
28 * same conditions as regards security.
29 *
30 * The fact that you are presently reading this means that you have had
31 * knowledge of the CeCILL-B license and that you accept its terms.
32 */
33
34
35#ifndef AIMS_DISTANCEMAP_MESHMORPHOMAT_D_H
36#define AIMS_DISTANCEMAP_MESHMORPHOMAT_D_H
37
38#include <aims/distancemap/meshmorphomat.h>
39#include <aims/distancemap/stlsort.h>
40#include <aims/distancemap/meshdistance.h>
41#include <aims/connectivity/meshcc_d.h>
42#include <stack>
43#include <float.h>
44
45namespace aims
46{
47 namespace meshdistance
48 {
49
50template<class T>
51Texture<T> LineariseSulci( const AimsSurface<3,Void> & mesh,
52 const Texture<T> & inittex,
53 const T & label,
54 const T & Back, const T & For,
55 unsigned & min, unsigned & max, unsigned nbDil)
56{
57
58 Texture<float> dist;
59 Texture<T> tex = inittex;
60 bool allowUnreached = false;
61 const std::vector<Point3df> & vert = mesh.vertex();
62 unsigned n = vert.size(),imax,i,j,jmax,nbpt=0;
63 TimeTexture<float> distance(11,n);
64 float distMax = 0;
65 std::vector<unsigned> nodeMax;
66 const std::vector< AimsVector<uint,3> > & poly = mesh.polygon();
67 ASSERT (n == tex.nItem() );
68
69
70 //init object map
71 for (i=0;i<n;++i)
72 {
73 if (inittex.item(i) == label)
74 {
75 tex.item(i) = Back;
76 ++nbpt;
77 }
78 else
79 tex.item(i) = For;
80 }
81
82 if (nbpt == 0 || nbpt == 1)
83 {
84 std::cout << "The label " << label << " has only one node ! \n";
85 return(inittex);
86 }
87 //ncc = AimsMeshLabelNbConnectedComponent(mesh,tex,Back);
88 //ASSERT ( ncc == 1 );
89
90 // neighbours map
91
92 std::map<unsigned, std::set<unsigned> > neighbours;
93 unsigned v1, v2, v3;
94
95 for( i=0; i<poly.size(); ++i )
96 {
97 v1 = poly[i][0];
98 v2 = poly[i][1];
99 v3 = poly[i][2];
100 if(tex.item(v1)!=For
101 && tex.item(v2)!=For)
102 {
103 neighbours[v1].insert( v2 );
104 neighbours[v2].insert( v1 );
105 }
106 if(tex.item(v1)!=For
107 && tex.item(v3)!=For)
108 {
109 neighbours[v1].insert( v3 );
110 neighbours[v3].insert( v1 );
111 }
112 if(tex.item(v2)!=For
113 && tex.item(v3)!=For)
114 {
115 neighbours[v2].insert( v3 );
116 neighbours[v3].insert( v2 );
117 }
118 }
119
120 //init distance map
121 i=0;
122 while (tex.item(i) == For)
123 ++ i;
124 tex.item(i) = label;
125
126 dist = MeshDistance<T>(mesh, tex, allowUnreached);
127 distance[0] = dist;
128
129 for (i=0;i<n;++i)
130 if (distance[0].item(i) >= distMax)
131 {
132 distMax = distance[0].item(i);
133 imax = i;
134 }
135
136 nodeMax.push_back(imax);
137
138 //Distance map
139 for (j=1;j<5;++j)
140 {
141 distMax = 0;
142 for (i=0;i<n;++i)
143 {
144 if (inittex.item(i) > 0)
145 tex.item(i) = Back;
146 else
147 tex.item(i) = For;
148 }
149 tex.item(imax) = 1;
150 dist = MeshDistance<T>(mesh, tex, allowUnreached);
151 distance[j] = dist;
152 for (i=0;i<n;++i)
153 if (dist.item(i) > distMax )
154 {
155 distMax = dist.item(i);
156 imax = i;
157 }
158 nodeMax.push_back(imax);
159 }
160
161 Texture<T> sulci(n);
162 std::set<unsigned>::iterator in,en;
163 float dmin;
164 min = imax = nodeMax[4];
165 max = jmax = nodeMax[3];
166 std::stack<unsigned> line;
167 i=imax;
168 while (i != jmax)
169 {
170 line.push(i);
171 dmin = FLT_MAX;
172 for (in = neighbours[i].begin(), en = neighbours[i].end(); in !=en ; ++in)
173 if (dist.item(*in) < dmin)
174 {
175 i = *in;
176 dmin = dist.item(*in);
177 }
178 sulci.item(i) = label;
179 }
180
181 unsigned inc = 0, l = line.size();
182
183 std::cout << "Size of the line: " << l << " nbDil: " << nbDil << std::endl;
184 if (nbDil != 0)
185 while (!line.empty())
186 {
187 ++inc;
188 i = line.top();
189 line.pop();
190 if (inc <= nbDil || inc >= l - nbDil)
191 sulci.item(i) = Back;
192 }
193
194 for (i=0;i<n;++i)
195 {
196 if (inittex.item(i) != label)
197 sulci.item(i) = inittex.item(i);
198 }
199
200 return (sulci);
201}
202
203
204
205
206
207
208
209
210template<class T>
211Texture<T> SplitSulci( const AimsSurface<3,Void> & mesh,
212 const Texture<T> & inittex,
213 const T & Back, const T & For ,
214 const T label,
215 const unsigned min, const unsigned max)
216{
217
218
219 const std::vector<Point3df> & vert = mesh.vertex();
220 unsigned ncc, n = vert.size(),i;
221 Texture<T> dilat(n),line(n),split(n);
222 const std::vector< AimsVector<uint,3> > & poly = mesh.polygon();
223
224
225 ASSERT (n = inittex.nItem() );
226
227 // Forbidden neighbours map
228
229 std::set<unsigned> neighbours;
230 unsigned v1, v2, v3;
231 std::set<unsigned>::iterator is,es;
232
233 i = min;
234 v1 = poly[i][0];
235 v2 = poly[i][1];
236 v3 = poly[i][2];
237 neighbours.insert(v1);
238 neighbours.insert(v2);
239 neighbours.insert(v3);
240 i = max;
241 v1 = poly[i][0];
242 v2 = poly[i][1];
243 v3 = poly[i][2];
244 neighbours.insert(v1);
245 neighbours.insert(v2);
246 neighbours.insert(v3);
247 neighbours.insert(min);
248 neighbours.insert(max);
249
250 for (i=0;i<n;++i)
251 {
252 if (inittex.item(i) == label)
253 dilat.item(i) = label;
254 else
255 dilat.item(i) = Back;
256 }
257
258
259
260 for (is=neighbours.begin(), es=neighbours.end();is != es; ++is)
261 //if (dilat.item(*is) == Back)
262 dilat.item(*is) = For;
263
264
265 dilat = MeshDilation<T>(mesh,dilat,Back,For,1,true);
266
267 for (i=0;i<n;++i)
268 split.item(i) = dilat.item(i);
269
270
271 for (i=0;i<n;++i)
272 if (inittex.item(i) == label || dilat.item(i) == For )
273 split.item(i) = Back;
274
275
276 ncc = AimsMeshLabelNbConnectedComponent(mesh,split,label);
277 std::cout << "ncc: " << ncc << std::endl;
278
279
280 return(split);
281
282}
283
284
285
286//closing and thinning of the sulci
287template<class T>
288Texture<T> CloseSulci( const AimsSurface<3,Void> & mesh,
289 const Texture<T> & sulctex,
290 float nbDilMax,bool connexity,
291 const T & Back, const T & For,
292 const std::vector<std::list<unsigned> > & neigho)
293{
294
295 const std::vector<Point3df> & vert = mesh.vertex();
296 unsigned inc=0,nbLabel,i,/*nbDil,*/n = vert.size();
297 Texture<T> sulctemp(n),ccsulcus(n);
298 Texture<T> tex(n);
299 short ncc;
300 T label;
301 //std::set<short> labsulci;
302 std::set<T> cc_sulci_labels;//labels of the connex components
303 typename std::set<T>::iterator ilab,elab;
304 bool hole=true;
305
306 // init
307 for (i=0; i<n; ++i)
308 {
309 cc_sulci_labels.insert(sulctex.item(i));
310 tex.item(i) = sulctex.item(i);
311 }
312
313 cc_sulci_labels.erase(Back);
314 cc_sulci_labels.erase(For);
315 nbLabel = cc_sulci_labels.size();
316 std::cout << "\nNb of labels: " << nbLabel << std::endl;
317
318
319 for ( ilab = cc_sulci_labels.begin(), elab=cc_sulci_labels.end(); ilab!=elab; ++ilab)
320 {
321 ++inc;
322 //nbDil=0;
323 label = *ilab;
324 std::cout << inc << "/" << nbLabel << std::endl;
325 for (i=0; i<n; ++i)
326 {
327 if (sulctex.item(i) == label)
328 ccsulcus.item(i) = label;
329 else
330 ccsulcus.item(i) = Back;
331 }
332 ncc = AimsMeshLabelNbConnectedComponent(mesh,ccsulcus,label);
333 std::cout << "Nb initial of cc: " << ncc << std::endl;
334
335 if (ncc == 1)
336 {
337 hole = HasHole(ccsulcus,mesh,label);
338
339 if (hole)
340 std::cout << "With hole(s) \n";
341 else
342 std::cout << "Without hole \n";
343 }
344 float nbd=1;
345 unsigned NB_PT_MAX = 10,nbpt;
346
347 //homotopic closing
348 while ( (ncc != 1 || hole) && nbd <= nbDilMax)
349 {
350 sulctemp = MeshDilation<T>( mesh, ccsulcus,Back,For,nbd,connexity);
351 ncc = AimsMeshLabelNbConnectedComponent(mesh,sulctemp,label);
352 std::cout << "After Dilation (" << nbd << ") : cc: " << ncc << std::endl ;
353 if (ncc == 1)
354 {
355 hole = HasHole(sulctemp,mesh,label);
356 if (hole)
357 std::cout << "With hole(s) \n";
358 else
359 std::cout << "Without hole \n";
360 }
361 nbd += 1;
362 }
363
364 if (ncc != 1)
365 std::cerr <<"Not able to close the projected connnected component, but I continue.\n";
366 if (nbd == 1)
367 sulctemp = ccsulcus;
368 ccsulcus = MeshDilation<T>( mesh, ccsulcus,Back,For,1,connexity);
369 std::cout << "Skeletonization" << std::endl;
370 ccsulcus = MeshSkeletization( mesh,sulctemp,label,Back,neigho);
371 nbpt = NbOfPoint(sulctemp,label);
372 if ( nbpt > NB_PT_MAX)
373 ccsulcus = MeshSkeletization( mesh,ccsulcus,label,Back,neigho);//2 steps for the ebarbulage
374 else
375 std::cout << "Not enough point for the ebarbulage ( " << nbpt << " < " << NB_PT_MAX << ")" << std::endl;
376 ncc = AimsMeshLabelNbConnectedComponent(mesh,ccsulcus,label);
377 std::cout << "ncc: " << ncc << std::endl ;
378 hole = HasHole(ccsulcus,mesh,label);
379 if (hole)
380 std::cout << "With hole(s) \n";
381 else
382 std::cout << "Without hole \n";
383
384 //old thinning
385 /*
386 if (ncc != 1 && nbDil > 0)
387 {
388 std::cout << "Impossible to close the sulci of label " << label << std::endl;
389 ccsulcus = AimsMeshFilterConnectedComponent( mesh, ccsulcus,label );
390 }
391 if (AimsMeshLabelNbConnectedComponent(mesh,ccsulcus,label) == 1)
392 {
393 std::cout << "Boundaries thinning \n";
394 ccsulcus = LineariseSulci<short>( mesh, ccsulcus ,label,Back,For,min,max,nbDil);
395 //std::cout << "Split sulci \n";
396 //ccsulcus = SplitSulci(mesh,ccsulcus,Back,For,label,min,max);
397
398 }
399 */
400
401 //check if the node is free
402 for (i=0;i<n;++i)
403 if (ccsulcus.item(i) == label)
404 {
405 if (sulctex.item(i)==label || sulctex.item(i) == Back )
406 tex.item(i)=label;
407 else
408 std::cout << "The point " << i << " cannot change its label \n";
409 }
410
411
412 for (i=0;i<n;++i)
413 if (tex.item(i) == label)
414 {
415 if (ccsulcus.item(i) == label)
416 tex.item(i) = label;
417 else
418 tex.item(i) = Back ;
419 }
420 }
421
422 return(tex);
423}
424
425
426
427//id but no thinning for label in the set label_forbidden
428template<class T>
429Texture<T> CloseSulci( const AimsSurface<3,Void> & mesh,
430 const Texture<T> & sulctex,
431 float nbDilMax,bool connexity,
432 const T & Back, const T & For,
433 const std::vector<std::list<unsigned> > & neigho,
434 const std::set<T> & label_forbidden)
435{
436
437 const std::vector<Point3df> & vert = mesh.vertex();
438 unsigned inc=0,nbLabel,i,/*nbDil,*/n = vert.size();
439 Texture<T> sulctemp(n),ccsulcus(n);
440 Texture<T> tex(n);
441 short ncc;
442 T label;
443 //std::set<short> labsulci;
444 std::set<T> cc_sulci_labels;//labels of the connex components
445 typename std::set<T>::iterator ilab,elab;
446 bool hole=true;
447
448 // init
449 for (i=0; i<n; ++i)
450 {
451 cc_sulci_labels.insert(sulctex.item(i));
452 tex.item(i) = sulctex.item(i);
453 }
454
455 cc_sulci_labels.erase(Back);
456 cc_sulci_labels.erase(For);
457 nbLabel = cc_sulci_labels.size();
458 std::cout << "\nNb of labels: " << nbLabel << std::endl;
459
460 typename std::set<T>::iterator lfe = label_forbidden.end();
461
462 for ( ilab = cc_sulci_labels.begin(), elab=cc_sulci_labels.end(); ilab!=elab; ++ilab)
463 {
464 ++inc;
465 //nbDil=0;
466 label = *ilab;
467 std::cout << inc << "/" << nbLabel << std::endl;
468 for (i=0; i<n; ++i)
469 {
470 if (sulctex.item(i) == label)
471 ccsulcus.item(i) = label;
472 else
473 ccsulcus.item(i) = Back;
474 }
475 ncc = AimsMeshLabelNbConnectedComponent(mesh,ccsulcus,label);
476 std::cout << "Nb initial of cc: " << ncc << std::endl;
477
478 if (ncc == 1)
479 {
480 hole = HasHole(ccsulcus,mesh,label);
481
482 if (hole)
483 std::cout << "With hole(s) \n";
484 else
485 std::cout << "Without hole \n";
486 }
487 float nbd=1;
488 unsigned NB_PT_MAX = 10,nbpt;
489
490 //homotopic closing
491 while ( (ncc != 1 || hole) && nbd <= nbDilMax)
492 {
493 sulctemp = MeshDilation<T>( mesh, ccsulcus,Back,For,nbd,connexity);
494 ncc = AimsMeshLabelNbConnectedComponent(mesh,sulctemp,label);
495 std::cout << "After Dilation (" << nbd << ") : cc: " << ncc << std::endl ;
496 if (ncc == 1)
497 {
498 hole = HasHole(sulctemp,mesh,label);
499 if (hole)
500 std::cout << "With hole(s) \n";
501 else
502 std::cout << "Without hole \n";
503 }
504 nbd += 1;
505 }
506
507 if (ncc != 1)
508 std::cerr <<"Not able to close the projected connnected component, but I continue.\n";
509 if (nbd == 1)
510 sulctemp = ccsulcus;
511 ccsulcus = MeshDilation<T>( mesh, ccsulcus,Back,For,1,connexity);
512 std::cout << "Skeletonization" << std::endl;
513
514 if (label_forbidden.find(*ilab) == lfe)
515 {
516 ccsulcus = MeshSkeletization( mesh,sulctemp,label,Back,neigho);
517 nbpt = NbOfPoint(sulctemp,label);
518 if ( nbpt > NB_PT_MAX)
519 ccsulcus = MeshSkeletization( mesh,ccsulcus,label,Back,neigho);//2 steps for the ebarbulage
520 else
521 std::cout << "Not enough point for the ebarbulage ( " << nbpt << " < " << NB_PT_MAX << ")" << std::endl;
522 ncc = AimsMeshLabelNbConnectedComponent(mesh,ccsulcus,label);
523 std::cout << "ncc: " << ncc << std::endl ;
524 hole = HasHole(ccsulcus,mesh,label);
525 if (hole)
526 std::cout << "With hole(s) \n";
527 else
528 std::cout << "Without hole \n";
529 }
530 else
531 std::cout << "No thinning for the INSULA\n";
532
533 //old thinning
534 /*
535 if (ncc != 1 && nbDil > 0)
536 {
537 std::cout << "Impossible to close the sulci of label " << label << std::endl;
538 ccsulcus = AimsMeshFilterConnectedComponent( mesh, ccsulcus,label );
539 }
540 if (AimsMeshLabelNbConnectedComponent(mesh,ccsulcus,label) == 1)
541 {
542 std::cout << "Boundaries thinning \n";
543 ccsulcus = LineariseSulci<short>( mesh, ccsulcus ,label,Back,For,min,max,nbDil);
544 //std::cout << "Split sulci \n";
545 //ccsulcus = SplitSulci(mesh,ccsulcus,Back,For,label,min,max);
546
547 }
548 */
549
550 //check if the node is free
551 for (i=0;i<n;++i)
552 if (ccsulcus.item(i) == label)
553 {
554 if (sulctex.item(i)==label || sulctex.item(i) == Back )
555 tex.item(i)=label;
556 else
557 std::cout << "The point " << i << " cannot change its label \n";
558 }
559
560
561 for (i=0;i<n;++i)
562 if (tex.item(i) == label)
563 {
564 if (ccsulcus.item(i) == label)
565 tex.item(i) = label;
566 else
567 tex.item(i) = Back ;
568 }
569 }
570
571
572 return(tex);
573}
574
575
576namespace
577{
578
579 template<class T> inline
580 bool SimplePoint( const Texture<T> & tex,
581 const std::list<unsigned> & neigho)
582 {
583
584 std::list<unsigned>::const_iterator il=neigho.begin(), el=neigho.end();
585 short l=tex.item(*il);
586 int nb=0;//nb chgmt of value
587
588 while(il != el)
589 {
590 if (tex.item(*il) != l)
591 {
592 ++nb;
593 l = tex.item(*il);
594 }
595 ++il;
596 }
597 return(nb == 1 || nb == 2);
598 }
599
600}
601
602template<class T> inline
603bool HasHole( const Texture<T> &tex,
604 const AimsSurface<3,Void> & mesh, const T label )
605{
606
607 const std::vector<Point3df> & vert = mesh.vertex();
608 const std::vector< AimsVector<uint,3> > & poly = mesh.polygon();
609 unsigned i, n = vert.size();
610 unsigned S=0,A=0,F=0;
611
612 ASSERT( tex.nItem() == n );
613
614 // neighbours map
615
616 std::set<std::set<unsigned> > neighbours;
617 std::set<unsigned> temp;
618 unsigned v1, v2, v3;
619
620 for( i=0; i<poly.size(); ++i )
621 {
622 v1 = poly[i][0];
623 v2 = poly[i][1];
624 v3 = poly[i][2];
625
626 if(tex.item(v1)==label && tex.item(v2)==label)
627 {
628 temp.clear();
629 temp.insert(v1);
630 temp.insert(v2);
631 neighbours.insert(temp);
632 }
633
634 if(tex.item(v1)==label && tex.item(v3)==label)
635 {
636 temp.clear();
637 temp.insert(v1);
638 temp.insert(v3);
639 neighbours.insert(temp);
640 }
641
642 if(tex.item(v3)==label && tex.item(v2)==label)
643 {
644 temp.clear();
645 temp.insert(v2);
646 temp.insert(v3);
647 neighbours.insert(temp);
648 }
649
650 if(tex.item(v3)==label && tex.item(v2)==label && tex.item(v1)==label)
651 ++F;
652 }
653
654 //Nb de sommets S
655 for (i=0;i<n;++i)
656 if (tex.item(i) == label)
657 ++S;
658
659 //Nb d'aretes A
660 A = neighbours.size();
661
662
663 return( (S - A + F != 1) );
664
665}
666
667
668namespace
669{
670
671 template<class T> inline
672 bool ImmortalPoint( const Texture<T> & tex,
673 const std::list<unsigned> & neigho )
674 {
675
676 std::list<unsigned>::const_iterator il=neigho.begin(), el=neigho.end();
677 short l=tex.item(*il);
678 int nb=0;//nb chgmt of value
679
680 while(il != el)
681 {
682 if (tex.item(*il) != l)
683 {
684 ++nb;
685 l = tex.item(*il);
686 }
687 ++il;
688 }
689 il = neigho.begin();
690 if (tex.item(*il) != l)
691 ++nb;
692
693 return(nb > 2 || nb ==1 );
694 }
695
696}
697
698
699
700template<class T> inline
701unsigned NbOfPoint( const Texture<T> tex, const T label )
702{
703 unsigned nb = 0;
704 for (unsigned i=0; i < tex.nItem(); ++i)
705 if (tex.item(i) == label)
706 ++nb;
707
708 return(nb);
709
710}
711
712
713template<class T>
714Texture<T> MeshSkeletization( const AimsSurface<3,Void> & mesh,
715 const Texture<T> & inittex,
716 T label, T Back,
717 const std::vector<std::list<unsigned> > & neigho )
718{
719
720 unsigned i,n = mesh.vertex().size();
721 Texture<short> skelet(n),temp(n);
722 Texture<short> immortal(n); //0 : still alive , 1: immortal
723 Texture<T> tex(n);
724 std::list<unsigned>::const_iterator il,el;
725
726 for (i=0; i<n; ++i)
727 if (inittex.item(i) == label)
728 skelet.item(i) = 1;
729
730 //init
731
732
733 //Def initial immortal point map
734 for (i=0; i<n; ++i)
735 if ( skelet.item(i) == 1 && ImmortalPoint(skelet,neigho[i]) )
736 immortal.item(i) = 1;
737
738 for (i=0; i<n; ++i)
739 if (skelet.item(i) == 1)
740 temp.item(i) = 0;
741 else
742 temp.item(i) = 1;
743
744 Texture<float> dist;
745 std::multimap<float,unsigned> dfront;
746 std::multimap<float,unsigned>::iterator id,ed;
747 dist = MeshDistance(mesh,temp,true);
748
749
750 // node of the front are ordered according to their distance from the background
751 for (i=0; i<n; ++i)
752 if (skelet.item(i) == 1 && immortal.item(i) == 0)
753 for (il=neigho[i].begin(), el=neigho[i].end();il!=el;++il)
754 if (skelet.item(*il) == 0 )
755 dfront.insert( std::pair<float,unsigned>(dist.item(i),i) );
756
757
758 while (dfront.size()>0)
759 {
760
761 for (id=dfront.begin(), ed=dfront.end();id != ed; ++id)
762 if ( immortal.item(id->second) == 0 )
763 {
764 skelet.item(id->second) = 0;
765 //upgrade immortal map
766 for (il=neigho[id->second].begin(), el=neigho[id->second].end();il!=el;++il)
767 if ( skelet.item(*il) == 1 && ImmortalPoint(skelet,neigho[*il]) )
768 immortal.item(*il) = 1;
769 }
770
771 //upgrade front
772
773 dfront.clear();
774 for (i=0; i<n; ++i)
775 if (skelet.item(i) == 1 && immortal.item(i) == 0)
776 for (il=neigho[i].begin(), el=neigho[i].end();il!=el;++il)
777 if (skelet.item(*il) == 0 )
778 dfront.insert( std::pair<float,unsigned>(dist.item(i),i) );
779 }
780
781 for (i=0; i<n; ++i)
782 if ( skelet.item(i) == 1)
783 tex.item(i) =label;
784 else
785 tex.item(i) = Back;
786
787 return(tex);
788
789}
790
791template <class T>
792Texture<T> MeshDilation( const AimsSurface<3,Void> & mesh,
793 const Texture<T> & inittex,
794 const T & Back, const T & For,
795 const float dist,bool connexity )
796{
797
798 Texture<T> tex;
799 tex = MeshVoronoiT<T>( mesh, inittex, Back, For, dist, connexity,true );
800 return(tex);
801
802}
803
804template <class T>
805Texture<T> MeshErosion( const AimsSurface<3,Void> & mesh,
806 const Texture<T> & inittex,
807 const T & Back, const T & For,
808 const float dist,bool connexity )
809{
810
811 Texture<T> tex;
812 tex = MeshVoronoiT<T>( mesh, inittex, Back,For,dist, connexity,false );
813 return(tex);
814}
815
816 }
817}
818
819#endif
820
ASSERT
#define ASSERT(EX)
AimsSurface::polygon
const std::vector< AimsVector< uint, D > > & polygon() const
AimsSurface::vertex
const std::vector< Point3df > & vertex() const
Texture::nItem
size_t nItem() const
Texture::item
const T & item(int n) const
TimeTexture::item
const T & item(int n) const
AimsMeshLabelNbConnectedComponent
unsigned AimsMeshLabelNbConnectedComponent(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, const T lab)
Definition meshcc_d.h:167
meshcc_d.h
aims::meshdistance::SplitSulci
Texture< T > SplitSulci(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, const T &Back, const T &For, const T label, const unsigned min, const unsigned max)
Definition meshmorphomat_d.h:211
aims::meshdistance::MeshErosion
Texture< T > MeshErosion(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, const T &Back, const T &For, const float dist, bool connexity)
Definition meshmorphomat_d.h:805
aims::meshdistance::MeshDilation
Texture< T > MeshDilation(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, const T &Back, const T &For, const float dist, bool connexity)
Definition meshmorphomat_d.h:792
aims::meshdistance::min
float min(float x, float y)
Definition thickness.h:106
aims::meshdistance::CloseSulci
Texture< T > CloseSulci(const AimsSurface< 3, Void > &mesh, const Texture< T > &sulctex, float nbDilMax, bool connexity, const T &Back, const T &For, const std::vector< std::list< unsigned > > &neigho, const std::set< T > &label_forbidden)
Close the sulcus map (dilation and skeletization).
Definition meshmorphomat_d.h:429
aims::meshdistance::max
float max(float x, float y)
Definition thickness.h:98
aims::meshdistance::HasHole
bool HasHole(const Texture< T > &tex, const AimsSurface< 3, Void > &mesh, const T label)
Definition meshmorphomat_d.h:603
aims::meshdistance::LineariseSulci
Texture< T > LineariseSulci(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, const T &label, const T &Back, const T &For, unsigned &min, unsigned &max, unsigned nbDil)
Definition meshmorphomat_d.h:51
aims::meshdistance::MeshSkeletization
Texture< T > MeshSkeletization(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, T label, T Back, const std::vector< std::list< unsigned > > &neigho)
Skeletization of the components of the innitex with label label.
Definition meshmorphomat_d.h:714
aims::meshdistance::MeshVoronoiT
Texture< T > MeshVoronoiT(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, const T &Back, const T &For, float dist, bool connexity, bool object)
Compute a geodesic voronoi diagram (dist = MAX_FLOAT,object=true ) of objects defined in inittex.
Definition meshvoronoi_d.h:59
aims::meshdistance::NbOfPoint
unsigned NbOfPoint(const Texture< T > tex, const T label)
Definition meshmorphomat_d.h:701
aims::meshdistance::MeshDistance
Texture< float > MeshDistance(const AimsSurface< 3, Void > &mesh, const Texture< T > &inittex, bool allowUnreached, float max_dist=FLT_MAX)
Computes a distance texture over a mesh.
Definition meshdistance_d.h:42