]> Creatis software - creaVtk.git/blob - bbtk_creaVtk_PKG/src/bbcreaVtkCleanMeshWithPatch.cxx
#3513 CleanMeshWithPatch
[creaVtk.git] / bbtk_creaVtk_PKG / src / bbcreaVtkCleanMeshWithPatch.cxx
1 //===== 
2 // Before editing this file, make sure it's a file of your own (i.e.: it wasn't generated from xml description; if so : your modifications will be lost)
3 //===== 
4 #include "bbcreaVtkCleanMeshWithPatch.h"
5 #include "bbcreaVtkPackage.h"
6
7 #include <vtkAppendPolyData.h>
8 #include <vtkStaticPointLocator.h>
9 #include <vtkDijkstraGraphGeodesicPath.h>
10 #include <vtkIdList.h>
11 #include <vtkCharArray.h>
12 #include <vtkPolyDataConnectivityFilter.h>
13 #include <vtkPointData.h>
14 #include <vtkCellData.h>
15 #include <vtkFillHolesFilter.h>
16 #include <vtkPolyDataNormals.h>
17 #include <vtkIdFilter.h>
18 #include <vtkFeatureEdges.h>
19 #include <vtkCleanPolyData.h>
20 #include <vtkTriangleFilter.h>
21 #include <vtkStripper.h>
22 #include <vtkTriangleStrip.h>
23 #include <vtkRuledSurfaceFilter.h>
24
25
26
27 // #include <vtkImprintFilter.h>
28
29
30 namespace bbcreaVtk
31 {
32
33 BBTK_ADD_BLACK_BOX_TO_PACKAGE(creaVtk,CleanMeshWithPatch)
34 BBTK_BLACK_BOX_IMPLEMENTATION(CleanMeshWithPatch,bbtk::AtomicBlackBox);
35 //===== 
36 // Before editing this file, make sure it's a file of your own (i.e.: it wasn't generated from xml description; if so : your modifications will be lost)
37 //===== 
38 void CleanMeshWithPatch::Process()
39 {
40 // THE MAIN PROCESSING METHOD BODY
41 //   Here we simply set the input 'In' value to the output 'Out'
42 //   And print out the output value
43 // INPUT/OUTPUT ACCESSORS ARE OF THE FORM :
44 //    void bbSet{Input|Output}NAME(const TYPE&)
45 //    const TYPE& bbGet{Input|Output}NAME() const 
46 //    Where :
47 //    * NAME is the name of the input/output
48 //      (the one provided in the attribute 'name' of the tag 'input')
49 //    * TYPE is the C++ type of the input/output
50 //      (the one provided in the attribute 'type' of the tag 'input')
51 //    bbSetOutputOut( bbGetInputIn() );
52 //    std::cout << "Output value = " <<bbGetOutputOut() << std::endl;
53   
54     std::vector<int>    lstIndexs   = bbGetInputLstIndexs();
55     if (( lstIndexs.size()>=3 ) && (bbGetInputMesh()!=NULL) ) {
56         
57         
58         // Step 1. ----Get perimeter CtrlPoints-----
59         std::vector<double> lstX        = bbGetInputLstX();
60         std::vector<double> lstY        = bbGetInputLstY();
61         std::vector<double> lstZ        = bbGetInputLstZ();
62         std::vector<double> perimeterX;
63         std::vector<double> perimeterY;
64         std::vector<double> perimeterZ;
65         std::vector<double> tmpX;
66         std::vector<double> tmpY;
67         std::vector<double> tmpZ;
68         int iGroup,sizeGroups = lstIndexs.size();
69         int i,iGeneral=0,size;
70         for (iGroup=0 ;  iGroup<sizeGroups ; iGroup++)
71         {
72             if (iGroup==0)
73             {
74                 for (i=iGeneral ;  i<iGeneral+lstIndexs[iGroup] ; i++)
75                 {
76                     perimeterX.push_back( lstX[i] );
77                     perimeterY.push_back( lstY[i] );
78                     perimeterZ.push_back( lstZ[i] );
79                 } // for first group
80             } // if iGroup == 0
81             if ((iGroup!=0) && (iGroup!=sizeGroups-1))
82             {
83                 i=iGeneral+(lstIndexs[iGroup]-1);
84                 perimeterX.push_back( lstX[i] );
85                 perimeterY.push_back( lstY[i] );
86                 perimeterZ.push_back( lstZ[i] );
87                 tmpX.insert( tmpX.begin() , lstX[iGeneral] );
88                 tmpY.insert( tmpY.begin() , lstY[iGeneral] );
89                 tmpZ.insert( tmpZ.begin() , lstZ[iGeneral] );
90             }
91             if (iGroup==sizeGroups-1)
92             {
93                 for (i=iGeneral+lstIndexs[iGroup]-1 ;  i>=iGeneral ; i--)
94                 {
95                     perimeterX.push_back( lstX[i] );
96                     perimeterY.push_back( lstY[i] );
97                     perimeterZ.push_back( lstZ[i] );
98                 } // for i  last group
99             } // if iGroup == 0
100             iGeneral=iGeneral+lstIndexs[iGroup];
101         } // for
102         size=tmpX.size();
103         for (i=0; i<size ; i++)
104         {
105             perimeterX.push_back( tmpX[i] );
106             perimeterY.push_back( tmpY[i] );
107             perimeterZ.push_back( tmpZ[i] );
108         } // for i tmpX
109         // Step 2. ----Get Ids in Mesh of perimeter CtrlPoints-----
110         std::vector<double>    spc;
111         if (bbGetInputSpacing().size()<3)
112         {
113             spc.push_back(1);
114             spc.push_back(1);
115             spc.push_back(1);
116         } else {
117             spc = bbGetInputSpacing();
118         }//
119         vtkStaticPointLocator   *pointLocator;
120         pointLocator   = vtkStaticPointLocator::New();
121         pointLocator->SetDataSet( bbGetInputMesh() );
122         pointLocator->BuildLocator();
123         long int id;
124         size=perimeterX.size();
125         std::vector<long int> lstIdsCtrlPoints;
126         for (i=0;i<size;i++)
127         {
128             id = pointLocator->FindClosestPoint( perimeterX[i]*spc[0] ,perimeterY[i]*spc[1] , perimeterZ[i]*spc[2] ) ;
129             lstIdsCtrlPoints.push_back( id );
130         } // for i  perimeter
131         // Step 3 Get geodesic path from segments
132         long int id1,id2;
133         std::vector<long int> lstIdsGeodesicPerimter;
134         vtkAppendPolyData* append = vtkAppendPolyData::New();
135         for (i=0;i<size;i++)
136         {
137             id1=lstIdsCtrlPoints[i];
138             if (i<size-1)
139             {
140                 id2 = lstIdsCtrlPoints[i+1];
141             } else {
142                 id2 = lstIdsCtrlPoints[0];
143             }// if
144             vtkDijkstraGraphGeodesicPath *pathFilter    = vtkDijkstraGraphGeodesicPath::New();
145             pathFilter->SetInputData( bbGetInputMesh() );
146             pathFilter->SetStartVertex( id1 );
147             pathFilter->SetEndVertex( id2 );
148             pathFilter->Update();
149             append->AddInputData( pathFilter->GetOutput() );
150             append->Update();
151             
152             // Step 3.1 Puts ids of perimeter Subsegment in vector  lstIdsGeodesicPerimter
153             vtkIdList   *IdsList            = pathFilter->GetIdList();
154             long int    iIdLst,sizeIdLst    = IdsList->GetNumberOfIds();
155             for (iIdLst=0; iIdLst<sizeIdLst; iIdLst++)
156             {
157                 lstIdsGeodesicPerimter.push_back( IdsList->GetId(iIdLst) );
158             }// for iIdLst
159             
160             pathFilter->Delete();
161         } // for i  perimeter
162         append->Update();
163 //        bbSetOutputOut( append->GetOutput() );
164
165         // Step 4. -- Define scalars ---
166         vtkCharArray    *scalarsArray   = vtkCharArray::New();
167         vtkPoints       *points         = bbGetInputMesh()->GetPoints();
168         long int        iMPs,sizeMPs    = points->GetNumberOfPoints();
169         scalarsArray->SetNumberOfValues(sizeMPs);
170         for (iMPs=0; iMPs<sizeMPs; iMPs++)
171         {
172             scalarsArray->SetValue(iMPs,15);
173         } // for iMPs
174         long int    iIdLst,sizeIdLst    = lstIdsGeodesicPerimter.size();
175         for (iIdLst=0; iIdLst<sizeIdLst; iIdLst++)
176         {
177             scalarsArray->SetValue( lstIdsGeodesicPerimter[iIdLst]  ,100 );
178         }// for iIdLst
179         scalarsArray->SetName("scalarsPerimeterPatch");
180         bbGetInputMesh()->GetPointData()->AddArray( scalarsArray );
181         bbGetInputMesh()->GetPointData()->SetActiveScalars( "scalarsPerimeterPatch" );
182         bbGetInputMesh()->GetCellData()->SetActiveScalars( "scalarsPerimeterPatch" );
183
184         // Step 5.  -- Connectivity Filter
185         vtkPolyDataConnectivityFilter *connectivity = vtkPolyDataConnectivityFilter::New();
186         connectivity->SetInputData( bbGetInputMesh() );
187         connectivity->SetScalarConnectivity(true);
188
189         // just perimeter mesh
190 //        connectivity->SetExtractionModeToPointSeededRegions();
191 //        connectivity->InitializeSeedList();
192 //        connectivity->AddSeed( lstIdsGeodesicPerimter[0] );
193 //        connectivity->SetScalarRange(90,110);
194 //        connectivity->SetFullScalarConnectivity(false);
195
196         //connectivity->SetExtractionModeToCellSeededRegions();
197         //connectivity->InitializeSeedList();
198         //connectivity->AddSeed(0);
199         //connectivity->SetScalarRange(10,20);
200         //connectivity->SetFullScalarConnectivity(true);
201
202         //connectivity->SetExtractionModeToPointSeededRegions();
203         //connectivity->InitializeSeedList();
204         //connectivity->AddSeed(0);
205         //connectivity->SetScalarRange(10,20);
206         //connectivity->SetFullScalarConnectivity(true);
207
208         // OK
209
210         connectivity->SetExtractionModeToLargestRegion();
211         connectivity->SetScalarRange(10,20);
212         connectivity->SetFullScalarConnectivity(true);
213
214         
215         connectivity->Update();
216 //        bbSetOutputOut( connectivity->GetOutput() );
217
218         printf("EED CleanMeshWithPatch::Process   connectivity GetNumberOfPoints  %ld \n   ", connectivity->GetOutput()->GetNumberOfPoints() );
219
220 /*
221         vtkAppendPolyData* append3 = vtkAppendPolyData::New();
222         append3->AddInputData( connectivity->GetOutput()  );
223         append3->AddInputData( bbGetInputPatch()  );
224         append3->Update();
225
226         vtkStripper *stripper = vtkStripper::New();
227         stripper->SetInputData( append3->GetOutput() );
228         stripper->Update();
229         bbSetOutputOut( stripper->GetOutput() );
230 */
231         
232 /*
233         vtkFeatureEdges *edges1 = vtkFeatureEdges::New();
234         edges1->SetInputData( connectivity->GetOutput() );
235         edges1->BoundaryEdgesOn();
236         edges1->ManifoldEdgesOff();
237         edges1->NonManifoldEdgesOff();
238         edges1->FeatureEdgesOff();
239         edges1->Update();
240 //        bbSetOutputOut( edges1->GetOutput() );
241
242         vtkFeatureEdges *edges2 = vtkFeatureEdges::New();
243         edges2->SetInputData( bbGetInputPatch() );
244         edges2->BoundaryEdgesOn();
245         edges2->ManifoldEdgesOff();
246         edges2->NonManifoldEdgesOff();
247         edges2->FeatureEdgesOff();
248         edges2->Update();
249 //        bbSetOutputOut( edges1->GetOutput() );
250
251         vtkAppendPolyData* append3 = vtkAppendPolyData::New();
252         append3->AddInputData( edges1->GetOutput()  );
253         append3->AddInputData( edges2->GetOutput()  );
254         append3->Update();
255 //        bbSetOutputOut( append3->GetOutput() );
256         vtkRuledSurfaceFilter *ruledSurfaceFilter = vtkRuledSurfaceFilter::New();
257         ruledSurfaceFilter->SetInputData( append3->GetOutput() );
258         ruledSurfaceFilter->Update();
259         bbSetOutputOut( ruledSurfaceFilter->GetOutput() );
260
261  */
262         
263         
264         // ok
265 //        bbSetOutputOut( connectivity->GetOutput() );
266
267         
268         vtkCleanPolyData *connectivityClean = vtkCleanPolyData::New();
269         connectivityClean->SetInputData( connectivity->GetOutput() );
270         connectivityClean->Update();
271         
272         //step 6.1.  -- connect both meshes     Patch->ConnectivityClean    --
273         Stretch( connectivityClean->GetOutput() , bbGetInputPatch() );
274         //step 6.2.  -- connect both meshes     ConnectivityClean->Patch    --
275         Stretch(  bbGetInputPatch(), connectivityClean->GetOutput()  );
276
277      /*
278         //step 6.1.  -- connect both meshes     Patch->ConnectivityClean    --
279         vtkIdFilter *idFilter = vtkIdFilter::New();
280         idFilter->SetInputData( bbGetInputPatch() );
281         //idFilter->SetIdsArrayName("ids");
282         //idFilter->SetCellIdsArrayName("ids");
283         idFilter->SetPointIdsArrayName("ids");
284         idFilter->SetPointIds(true);
285         idFilter->SetCellIds(false);
286 //# Available for vtk>=8.3:
287 //           #idFilter.SetPointIdsArrayName(arrayName)
288 //           #idFilter.SetCellIdsArrayName(arrayName)
289         idFilter->Update();
290         vtkFeatureEdges *edges = vtkFeatureEdges::New();
291         edges->SetInputData( idFilter->GetOutput() );
292         edges->BoundaryEdgesOn();
293         edges->ManifoldEdgesOff();
294         edges->NonManifoldEdgesOff();
295         edges->FeatureEdgesOff();
296         edges->Update();
297         vtkIdTypeArray* arrayIds = vtkIdTypeArray::SafeDownCast(edges->GetOutput()->GetPointData()->GetArray("ids"));
298         long int iIds,sizeArrayIds    = edges->GetOutput()->GetNumberOfPoints();
299 printf("\n EED CleanMeshWithPatch::Process  sizeArrayIds=%ld \n", sizeArrayIds);
300         
301         
302         vtkPolyData             *newMesh        = connectivityClean->GetOutput();
303         vtkStaticPointLocator   *pointLocator2  = vtkStaticPointLocator::New();
304         pointLocator2->SetDataSet( newMesh );
305         pointLocator2->BuildLocator();
306         vtkPoints *pointsPatch    = bbGetInputPatch()->GetPoints();
307         vtkPoints *pointsNewMesh  = newMesh->GetPoints();
308         double pP[3], pPP[3];
309         long int idNewMesh;
310         for (iIds=0; iIds<sizeArrayIds ; iIds++ )
311         {
312             pointsPatch->GetPoint( arrayIds->GetValue(iIds) , pP );
313             idNewMesh = pointLocator2->FindClosestPoint(pP);
314             pointsNewMesh->GetPoint( idNewMesh , pPP );
315             pointsPatch->SetPoint( arrayIds->GetValue(iIds) , pPP );
316  printf(" iIds=%ld   pP=%f %f %f     idNewMesh=%ld   pPP=%f %f %f \n", arrayIds->GetValue(iIds) , pP[0], pP[1], pP[2], idNewMesh, pPP[0], pPP[1], pPP[2]  );
317         } // for iIds
318         bbGetInputPatch()->Modified();
319
320         //step 6.2.  -- connect both meshes     ConnectivityClean->Patch    --
321         vtkIdFilter *idFilter2 = vtkIdFilter::New();
322         idFilter2->SetInputData( newMesh );
323         //idFilter->SetIdsArrayName("ids");
324         //idFilter->SetCellIdsArrayName("ids");
325         idFilter2->SetPointIdsArrayName("ids");
326         idFilter2->SetPointIds(true);
327         idFilter2->SetCellIds(false);
328 //# Available for vtk>=8.3:
329 //           #idFilter.SetPointIdsArrayName(arrayName)
330 //           #idFilter.SetCellIdsArrayName(arrayName)
331         idFilter2->Update();
332         vtkFeatureEdges *edges2 = vtkFeatureEdges::New();
333         edges2->SetInputData( idFilter2->GetOutput() );
334         edges2->BoundaryEdgesOn();
335         edges2->ManifoldEdgesOff();
336         edges2->NonManifoldEdgesOff();
337         edges2->FeatureEdgesOff();
338         edges2->Update();
339         vtkIdTypeArray* arrayIds2 = vtkIdTypeArray::SafeDownCast(edges2->GetOutput()->GetPointData()->GetArray("ids"));
340         long int iIds2,sizeArrayIds2    = edges2->GetOutput()->GetNumberOfPoints();
341 printf("\n EED CleanMeshWithPatch::Process  sizeArrayIds2=%ld \n", sizeArrayIds2);
342 //        vtkPolyData             *newMesh        = connectivityClean->GetOutput();
343         vtkStaticPointLocator   *pointLocator3  = vtkStaticPointLocator::New();
344         pointLocator3->SetDataSet( bbGetInputPatch() );
345         pointLocator3->BuildLocator();
346         vtkPoints *pointsNewMesh2   = newMesh->GetPoints();
347         vtkPoints *pointsPatch2     = bbGetInputPatch()->GetPoints();
348         double pP2[3], pPP2[3];
349         long int idPatch;
350         for (iIds2=0; iIds2<sizeArrayIds2 ; iIds2++ )
351         {
352             pointsNewMesh2->GetPoint( arrayIds2->GetValue(iIds2) , pP2 );
353             idPatch = pointLocator3->FindClosestPoint(pP2);
354             pointsPatch2->GetPoint( idPatch , pPP2 );
355             pointsNewMesh2->SetPoint( arrayIds2->GetValue(iIds) , pPP2 );
356  printf(" iIds=%ld   pP=%f %f %f     idNewMesh=%ld   pPP=%f %f %f \n", arrayIds2->GetValue(iIds) , pP2[0], pP2[1], pP2[2], idPatch, pPP2[0], pPP2[1], pPP2[2]  );
357         } // for iIds
358         newMesh->Modified();
359 */
360         
361         
362         // step 7. -- Append meshes--
363         vtkAppendPolyData* append2 = vtkAppendPolyData::New();
364         append2->AddInputData( connectivityClean->GetOutput()  );
365         append2->AddInputData( bbGetInputPatch()  );
366         append2->Update();
367
368         vtkCleanPolyData *append2Clean = vtkCleanPolyData::New();
369         append2Clean->SetInputData( append2->GetOutput() );
370         append2Clean->Update();
371         
372 //        vtkTriangleFilter *triangles = vtkTriangleFilter::New();
373 //        triangles->SetInputData( append2Clean->GetOutput() );
374 //        triangles->Update();
375 //        bbSetOutputOut( triangles->GetOutput() );
376
377         
378         //step 8.  -- Clean  --
379
380         vtkPolyDataNormals *normals =vtkPolyDataNormals::New();
381   //      normals->SetInputConnection( triangles->GetOutputPort() );
382        normals->SetInputConnection( append2Clean->GetOutputPort() );
383         normals->ConsistencyOn();
384         normals->SplittingOff();
385         normals->Update();
386        bbSetOutputOut( normals->GetOutput() );
387
388 //        vtkFillHolesFilter *fillHoles = vtkFillHolesFilter::New();
389 //        fillHoles->SetInputData( normals->GetOutput() );
390 //        fillHoles->SetHoleSize(1000);
391 //        fillHoles->Update();
392 //        bbSetOutputOut( fillHoles->GetOutput() );
393         
394         
395         /*
396         // Step 4. -- vtkImprintFilter ---
397         vtkImprintFilter *imp = vtkImprintFilter::New();
398 //        imp->SetTargetData( append->GetOutput() );
399 //        imp->SetImprintData( bbGetInputMesh() );
400         
401 //        imp->SetTargetData( bbGetInputMesh() );
402 //        imp->SetImprintData( append->GetOutput() );
403
404 //        imp->SetTargetData( bbGetInputMesh() );
405 //        imp->SetImprintData( bbGetInputPatch() );
406         
407           imp->SetTargetData(  bbGetInputMesh() );
408           imp->SetImprintData( bbGetInputPatch()  );
409
410    //     imp->SetOutputTypeToProjectedImprint();
411         imp->SetOutputTypeToTargetCells();
412     //    imp->SetOutputTypeToImprintedCells();
413     //    imp->SetOutputTypeToImprintedRegion();
414     //    imp->SetOutputTypeToMergedImprint();
415         
416         imp->SetTolerance(500);
417         imp->Update();
418          bbSetOutputOut( imp->GetOutput() );
419 */
420         
421     } else {
422         bbSetOutputOut( NULL );
423     } // if LstIndexs size >=3
424         
425 }
426
427
428
429 void CleanMeshWithPatch::Stretch(vtkPolyData *newMesh, vtkPolyData *patchMesh)
430 {
431     
432     vtkIdFilter *idFilter = vtkIdFilter::New();
433     idFilter->SetInputData( patchMesh );
434     //idFilter->SetIdsArrayName("ids");
435     //idFilter->SetCellIdsArrayName("ids");
436     idFilter->SetPointIdsArrayName("ids");
437     idFilter->SetPointIds(true);
438     idFilter->SetCellIds(false);
439 //# Available for vtk>=8.3:
440 //           #idFilter.SetPointIdsArrayName(arrayName)
441 //           #idFilter.SetCellIdsArrayName(arrayName)
442     idFilter->Update();
443     vtkFeatureEdges *edges = vtkFeatureEdges::New();
444     edges->SetInputData( idFilter->GetOutput() );
445     edges->BoundaryEdgesOn();
446     edges->ManifoldEdgesOff();
447     edges->NonManifoldEdgesOff();
448     edges->FeatureEdgesOff();
449     edges->Update();
450     vtkIdTypeArray* arrayIds = vtkIdTypeArray::SafeDownCast(edges->GetOutput()->GetPointData()->GetArray("ids"));
451     long int iIds,sizeArrayIds    = edges->GetOutput()->GetNumberOfPoints();
452     
453     
454 //    vtkPolyData             *newMesh        = connectivityClean->GetOutput();
455     vtkStaticPointLocator   *pointLocator2  = vtkStaticPointLocator::New();
456     pointLocator2->SetDataSet( newMesh );
457     pointLocator2->BuildLocator();
458     vtkPoints *pointsPatch    = patchMesh->GetPoints();
459     vtkPoints *pointsNewMesh  = newMesh->GetPoints();
460     double pP[3], pPP[3];
461     long int idNewMesh;
462     for (iIds=0; iIds<sizeArrayIds ; iIds++ )
463     {
464         pointsPatch->GetPoint( arrayIds->GetValue(iIds) , pP );
465         idNewMesh = pointLocator2->FindClosestPoint(pP);
466         pointsNewMesh->GetPoint( idNewMesh , pPP );
467         pointsPatch->SetPoint( arrayIds->GetValue(iIds) , pPP );
468     } // for iIds
469     patchMesh->Modified();
470 }
471
472
473
474 //===== 
475 // Before editing this file, make sure it's a file of your own (i.e.: it wasn't generated from xml description; if so : your modifications will be lost)
476 //===== 
477 void CleanMeshWithPatch::bbUserSetDefaultValues()
478 {
479 //  SET HERE THE DEFAULT INPUT/OUTPUT VALUES OF YOUR BOX 
480 //    Here we initialize the input 'In' to 0
481    bbSetInputMesh(NULL);
482    bbSetInputPatch(NULL);
483    bbSetOutputOut(NULL);
484 }
485
486 //===== 
487 // Before editing this file, make sure it's a file of your own (i.e.: it wasn't generated from xml description; if so : your modifications will be lost)
488 //===== 
489 void CleanMeshWithPatch::bbUserInitializeProcessing()
490 {
491 //  THE INITIALIZATION METHOD BODY :
492 //    Here does nothing 
493 //    but this is where you should allocate the internal/output pointers 
494 //    if any
495 }
496
497 //===== 
498 // Before editing this file, make sure it's a file of your own (i.e.: it wasn't generated from xml description; if so : your modifications will be lost)
499 //===== 
500 void CleanMeshWithPatch::bbUserFinalizeProcessing()
501 {
502 //  THE FINALIZATION METHOD BODY :
503 //    Here does nothing 
504 //    but this is where you should desallocate the internal/output pointers 
505 //    if any
506 }
507
508 } // EO namespace bbcreaVtk
509
510