1 // SALOME OBJECT : kernel of SALOME component
3 // Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org
24 // File : VTKViewer_GeometryFilter.cxx
25 // Author : Michael ZORIN
29 #include "VTKViewer_GeometryFilter.h"
31 #include <vtkSmartPointer.h>
32 #include <vtkCellArray.h>
33 #include <vtkCellData.h>
34 #include <vtkGenericCell.h>
35 #include <vtkHexahedron.h>
36 #include <vtkMergePoints.h>
37 #include <vtkObjectFactory.h>
38 #include <vtkPointData.h>
39 #include <vtkPolyData.h>
40 #include <vtkPyramid.h>
41 #include <vtkStructuredGrid.h>
43 #include <vtkUnsignedCharArray.h>
44 #include <vtkUnstructuredGrid.h>
50 #include <vtkDelaunay3D.h>
51 #include <vtkGeometryFilter.h>
65 #define USE_ROBUST_TRIANGULATION
67 //----------------------------------------------------------------------------
68 vtkCxxRevisionMacro(VTKViewer_GeometryFilter, "$Revision$");
69 vtkStandardNewMacro(VTKViewer_GeometryFilter);
71 VTKViewer_GeometryFilter
72 ::VTKViewer_GeometryFilter():
79 VTKViewer_GeometryFilter
80 ::~VTKViewer_GeometryFilter()
84 //----------------------------------------------------------------------------
86 VTKViewer_GeometryFilter
89 vtkDataSet *input= this->GetInput();
90 vtkIdType numCells=input->GetNumberOfCells();
97 if (input->GetDataObjectType() == VTK_UNSTRUCTURED_GRID){
98 this->UnstructuredGridExecute();
101 vtkGeometryFilter::Execute();
105 //----------------------------------------------------------------------------
107 VTKViewer_GeometryFilter
108 ::UnstructuredGridExecute()
110 vtkUnstructuredGrid *input= (vtkUnstructuredGrid *)this->GetInput();
111 vtkCellArray *Connectivity = input->GetCells();
113 if ( Connectivity == NULL )
115 vtkDebugMacro(<<"Nothing to extract");
124 vtkPoints *p = input->GetPoints();
125 vtkIdType numCells=input->GetNumberOfCells();
126 vtkPointData *pd = input->GetPointData();
127 vtkCellData *cd = input->GetCellData();
128 vtkPolyData *output = this->GetOutput();
129 vtkPointData *outputPD = output->GetPointData();
131 #ifdef USE_ROBUST_TRIANGULATION
132 vtkUnstructuredGrid* anUnstructuredGrid = vtkUnstructuredGrid::New();
133 vtkPoints* aDelaunayPoints = vtkPoints::New();
135 vtkDelaunay3D* aDelaunay3D = vtkDelaunay3D::New();
136 aDelaunay3D->SetInput(anUnstructuredGrid);
138 vtkGeometryFilter* aGeometryFilter = vtkGeometryFilter::New();
139 aGeometryFilter->SetInput(aDelaunay3D->GetOutput());
142 vtkCellData *outputCD = output->GetCellData();
143 vtkGenericCell *cell = vtkGenericCell::New();
146 vtkIdList *cellIds = vtkIdList::New();
147 vtkIdList *faceIds = vtkIdList::New();
151 int faceId, *faceVerts, numFacePts;
153 int PixelConvert[4], aNewPts[VTK_CELL_SIZE];
155 unsigned char updateLevel = (unsigned char)(output->GetUpdateGhostLevel());
156 unsigned char *cellGhostLevels = 0;
163 vtkDebugMacro(<<"Executing geometry filter for unstructured grid input");
165 vtkDataArray* temp = 0;
168 temp = cd->GetArray("vtkGhostLevels");
170 if ( (!temp) || (temp->GetDataType() != VTK_UNSIGNED_CHAR)
171 || (temp->GetNumberOfComponents() != 1))
173 vtkDebugMacro("No appropriate ghost levels field available.");
177 cellGhostLevels = ((vtkUnsignedCharArray*)temp)->GetPointer(0);
180 // Determine nature of what we have to do
181 if ( (!this->CellClipping) && (!this->PointClipping) &&
182 (!this->ExtentClipping) )
190 cellVis = new char[numCells];
193 // Just pass points through, never merge
194 output->SetPoints(input->GetPoints());
195 outputPD->PassData(pd);
197 outputCD->CopyAllocate(cd,numCells,numCells/2);
199 output->Allocate(numCells/4+1,numCells);
201 // Loop over the cells determining what's visible
204 for (cellId=0, Connectivity->InitTraversal();
205 Connectivity->GetNextCell(npts,pts);
209 if ( this->CellClipping && cellId < this->CellMinimum ||
210 cellId > this->CellMaximum )
216 for (i=0; i < npts; i++)
218 x = p->GetPoint(pts[i]);
219 if ( (this->PointClipping && (pts[i] < this->PointMinimum ||
220 pts[i] > this->PointMaximum) ) ||
221 (this->ExtentClipping &&
222 (x[0] < this->Extent[0] || x[0] > this->Extent[1] ||
223 x[1] < this->Extent[2] || x[1] > this->Extent[3] ||
224 x[2] < this->Extent[4] || x[2] > this->Extent[5] )) )
228 }//point/extent clipping
230 }//if point clipping needs checking
232 }//if not all visible
234 // Loop over all cells now that visibility is known
235 // (Have to compute visibility first for 3D cell boundarys)
236 int progressInterval = numCells/20 + 1;
238 myVTK2ObjIds.clear();
239 myVTK2ObjIds.reserve(numCells);
241 for (cellId=0, Connectivity->InitTraversal();
242 Connectivity->GetNextCell(npts,pts);
245 //Progress and abort method support
246 if ( !(cellId % progressInterval) )
248 vtkDebugMacro(<<"Process cell #" << cellId);
249 this->UpdateProgress ((float)cellId/numCells);
252 // Handle ghost cells here. Another option was used cellVis array.
253 if (cellGhostLevels && cellGhostLevels[cellId] > updateLevel)
254 { // Do not create surfaces in outer ghost cells.
258 if (allVisible || cellVis[cellId]) //now if visible extract geometry
260 //special code for nonlinear cells - rarely occurs, so right now it
262 vtkIdType aCellType = input->GetCellType(cellId);
269 case VTK_POLY_VERTEX:
270 newCellId = output->InsertNextCell(aCellType,npts,pts);
272 myVTK2ObjIds.push_back(cellId); //apo
274 outputCD->CopyData(cd,cellId,newCellId);
279 newCellId = output->InsertNextCell(aCellType,npts,pts);
281 myVTK2ObjIds.push_back(cellId);
282 outputCD->CopyData(cd,cellId,newCellId);
288 newCellId = output->InsertNextCell(aCellType,npts,pts);
290 myVTK2ObjIds.push_back(cellId);
291 outputCD->CopyData(cd,cellId,newCellId);
294 case VTK_TRIANGLE_STRIP:
295 newCellId = output->InsertNextCell(aCellType,npts,pts);
297 myVTK2ObjIds.push_back(cellId);
298 outputCD->CopyData(cd,cellId,newCellId);
302 newCellId = output->InsertNextCell(aCellType,npts,pts);
304 myVTK2ObjIds.push_back(cellId);
305 outputCD->CopyData(cd,cellId,newCellId);
308 case VTK_CONVEX_POINT_SET: {
309 //cout<<"cellId = "<<cellId<<"\n";
313 #ifdef USE_ROBUST_TRIANGULATION
314 aPoints = aDelaunayPoints;
315 anUnstructuredGrid->Initialize();
316 anUnstructuredGrid->Allocate();
317 anUnstructuredGrid->SetPoints(aDelaunayPoints);
320 input->GetCellPoints(cellId,aNumPts,aPts);
323 aDelaunayPoints->SetNumberOfPoints(aNumPts);
324 vtkPoints *anInputPoints = input->GetPoints();
325 for (int aPntId = 0; aPntId < aNumPts; aPntId++) {
326 anInputPoints->GetPoint(aPts[aPntId],aPntCoord);
327 aDelaunayPoints->SetPoint(aPntId,aPntCoord);
331 input->GetCell(cellId,cell);
332 aPoints = input->GetPoints();
333 aNumPts = cell->GetNumberOfPoints();
335 // To calculate the bary center of the cell
336 float aCellCenter[3] = {0.0, 0.0, 0.0};
339 for (int aPntId = 0; aPntId < aNumPts; aPntId++) {
340 #ifdef USE_ROBUST_TRIANGULATION
341 aPoints->GetPoint(aPntId,aPntCoord);
343 aPoints->GetPoint(cell->GetPointId(aPntId),aPntCoord);
345 //cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}\n";
346 aCellCenter[0] += aPntCoord[0];
347 aCellCenter[1] += aPntCoord[1];
348 aCellCenter[2] += aPntCoord[2];
350 aCellCenter[0] /= aNumPts;
351 aCellCenter[1] /= aNumPts;
352 aCellCenter[2] /= aNumPts;
355 #ifdef USE_ROBUST_TRIANGULATION
356 aGeometryFilter->Update();
357 vtkPolyData* aPolyData = aGeometryFilter->GetOutput();
359 float aCellLength = aPolyData->GetLength();
360 int aNumFaces = aPolyData->GetNumberOfCells();
362 float aCellLength = sqrt(cell->GetLength2());
363 int aNumFaces = cell->GetNumberOfFaces();
366 static float EPS = 1.0E-5;
367 float aDistEps = aCellLength * EPS;
369 // To initialize set of points that belong to the cell
370 typedef std::set<vtkIdType> TPointIds;
371 TPointIds anInitialPointIds;
372 for(vtkIdType aPntId = 0; aPntId < aNumPts; aPntId++){
373 #ifdef USE_ROBUST_TRIANGULATION
374 anInitialPointIds.insert(aPntId);
376 anInitialPointIds.insert(cell->GetPointId(aPntId));
380 // To initialize set of points by face that belong to the cell and backward
381 typedef std::set<vtkIdType> TFace2Visibility;
382 TFace2Visibility aFace2Visibility;
384 typedef std::set<TPointIds> TFace2PointIds;
385 TFace2PointIds aFace2PointIds;
387 for (int aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
388 #ifdef USE_ROBUST_TRIANGULATION
389 vtkCell* aFace = aPolyData->GetCell(aFaceId);
391 vtkCell* aFace = cell->GetFace(aFaceId);
393 vtkIdList *anIdList = aFace->PointIds;
394 aNewPts[0] = anIdList->GetId(0);
395 aNewPts[1] = anIdList->GetId(1);
396 aNewPts[2] = anIdList->GetId(2);
398 #ifdef USE_ROBUST_TRIANGULATION
400 faceIds->InsertNextId(aPts[aNewPts[0]]);
401 faceIds->InsertNextId(aPts[aNewPts[1]]);
402 faceIds->InsertNextId(aPts[aNewPts[2]]);
403 input->GetCellNeighbors(cellId, faceIds, cellIds);
405 input->GetCellNeighbors(cellId, anIdList, cellIds);
407 if((!allVisible && !cellVis[cellIds->GetId(0)]) ||
408 cellIds->GetNumberOfIds() <= 0 ||
412 aPointIds.insert(aNewPts[0]);
413 aPointIds.insert(aNewPts[1]);
414 aPointIds.insert(aNewPts[2]);
416 aFace2PointIds.insert(aPointIds);
417 aFace2Visibility.insert(aFaceId);
421 for (int aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
422 if(aFace2Visibility.find(aFaceId) == aFace2Visibility.end())
425 #ifdef USE_ROBUST_TRIANGULATION
426 vtkCell* aFace = aPolyData->GetCell(aFaceId);
428 vtkCell* aFace = cell->GetFace(aFaceId);
430 vtkIdList *anIdList = aFace->PointIds;
431 aNewPts[0] = anIdList->GetId(0);
432 aNewPts[1] = anIdList->GetId(1);
433 aNewPts[2] = anIdList->GetId(2);
435 // To initialize set of points for the plane where the trinangle face belong to
437 aPointIds.insert(aNewPts[0]);
438 aPointIds.insert(aNewPts[1]);
439 aPointIds.insert(aNewPts[2]);
441 // To get know, if the points of the trinagle were already observed
442 bool anIsObserved = aFace2PointIds.find(aPointIds) == aFace2PointIds.end();
443 //cout<<"\taFaceId = "<<aFaceId<<"; anIsObserved = "<<anIsObserved;
444 //cout<<"; aNewPts = {"<<aNewPts[0]<<", "<<aNewPts[1]<<", "<<aNewPts[2]<<"}\n";
447 // To get coordinates of the points of the traingle face
449 aPoints->GetPoint(aNewPts[0],aCoord[0]);
450 aPoints->GetPoint(aNewPts[1],aCoord[1]);
451 aPoints->GetPoint(aNewPts[2],aCoord[2]);
453 // To calculate plane normal
454 float aVector01[3] = { aCoord[1][0] - aCoord[0][0],
455 aCoord[1][1] - aCoord[0][1],
456 aCoord[1][2] - aCoord[0][2] };
458 float aVector02[3] = { aCoord[2][0] - aCoord[0][0],
459 aCoord[2][1] - aCoord[0][1],
460 aCoord[2][2] - aCoord[0][2] };
463 vtkMath::Cross(aVector02,aVector01,aCross21);
465 vtkMath::Normalize(aCross21);
467 // To calculate what points belong to the plane
468 // To calculate bounds of the point set
469 float aCenter[3] = {0.0, 0.0, 0.0};
471 TPointIds::const_iterator anIter = anInitialPointIds.begin();
472 TPointIds::const_iterator anEndIter = anInitialPointIds.end();
473 for(; anIter != anEndIter; anIter++){
475 vtkIdType aPntId = *anIter;
476 aPoints->GetPoint(aPntId,aPntCoord);
477 float aDist = vtkPlane::DistanceToPlane(aPntCoord,aCross21,aCoord[0]);
478 //cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}; aDist = "<<aDist<<"\n";
479 if(fabs(aDist) < aDistEps){
480 aPointIds.insert(aPntId);
481 aCenter[0] += aPntCoord[0];
482 aCenter[1] += aPntCoord[1];
483 aCenter[2] += aPntCoord[2];
486 int aNbPoints = aPointIds.size();
487 aCenter[0] /= aNbPoints;
488 aCenter[1] /= aNbPoints;
489 aCenter[2] /= aNbPoints;
492 // To calculate the primary direction for point set
493 float aVector0[3] = { aCoord[0][0] - aCenter[0],
494 aCoord[0][1] - aCenter[1],
495 aCoord[0][2] - aCenter[2] };
497 //To sinchronize orientation of the cell and its face
498 float aVectorC[3] = { aCenter[0] - aCellCenter[0],
499 aCenter[1] - aCellCenter[1],
500 aCenter[2] - aCellCenter[2] };
501 vtkMath::Normalize(aVectorC);
503 float aDot = vtkMath::Dot(aCross21,aVectorC);
504 //cout<<"\t\taCross21 = {"<<aCross21[0]<<", "<<aCross21[1]<<", "<<aCross21[2]<<"}";
505 //cout<<"; aVectorC = {"<<aVectorC[0]<<", "<<aVectorC[1]<<", "<<aVectorC[2]<<"}\n";
506 //cout<<"\t\taDot = "<<aDot<<"\n";
508 aCross21[0] = -aCross21[0];
509 aCross21[1] = -aCross21[1];
510 aCross21[2] = -aCross21[2];
513 vtkMath::Normalize(aVector0);
515 //cout<<"\t\taCenter = {"<<aCenter[0]<<", "<<aCenter[1]<<", "<<aCenter[2]<<"}";
516 //cout<<"; aVector0 = {"<<aVector0[0]<<", "<<aVector0[1]<<", "<<aVector0[2]<<"}\n";
518 // To calculate the set of points by face those that belong to the plane
519 TFace2PointIds aRemoveFace2PointIds;
521 TFace2PointIds::const_iterator anIter = aFace2PointIds.begin();
522 TFace2PointIds::const_iterator anEndIter = aFace2PointIds.end();
523 for(; anIter != anEndIter; anIter++){
524 const TPointIds& anIds = *anIter;
525 TPointIds anIntersection;
526 std::set_intersection(aPointIds.begin(),aPointIds.end(),
527 anIds.begin(),anIds.end(),
528 std::inserter(anIntersection,anIntersection.begin()));
530 if(anIntersection == anIds){
531 aRemoveFace2PointIds.insert(anIds);
536 // To remove from the set of points by face those that belong to the plane
538 TFace2PointIds::const_iterator anIter = aRemoveFace2PointIds.begin();
539 TFace2PointIds::const_iterator anEndIter = aRemoveFace2PointIds.end();
540 for(; anIter != anEndIter; anIter++){
541 const TPointIds& anIds = *anIter;
542 aFace2PointIds.erase(anIds);
546 // To sort the planar set of the points accrding to the angle
548 typedef std::map<float,vtkIdType> TSortedPointIds;
549 TSortedPointIds aSortedPointIds;
551 TPointIds::const_iterator anIter = aPointIds.begin();
552 TPointIds::const_iterator anEndIter = aPointIds.end();
553 for(; anIter != anEndIter; anIter++){
555 vtkIdType aPntId = *anIter;
556 aPoints->GetPoint(aPntId,aPntCoord);
557 float aVector[3] = { aPntCoord[0] - aCenter[0],
558 aPntCoord[1] - aCenter[1],
559 aPntCoord[2] - aCenter[2] };
560 vtkMath::Normalize(aVector);
563 vtkMath::Cross(aVector,aVector0,aCross);
564 bool aGreaterThanPi = vtkMath::Dot(aCross,aCross21) < 0;
565 float aCosinus = vtkMath::Dot(aVector,aVector0);
570 static float a2Pi = 2.0 * vtkMath::Pi();
571 float anAngle = acos(aCosinus);
572 //cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}";
573 //cout<<"; aGreaterThanPi = "<<aGreaterThanPi<<"; aCosinus = "<<aCosinus<<"; anAngle = "<<anAngle<<"\n";
575 anAngle = a2Pi - anAngle;
576 aSortedPointIds[anAngle] = aPntId;
577 //cout<<"\t\t\tanAngle = "<<anAngle<<"\n";
579 if(!aSortedPointIds.empty()){
580 aCellType = VTK_POLYGON;
581 int numFacePts = aSortedPointIds.size();
582 std::vector<vtkIdType> aConnectivities(numFacePts);
583 TSortedPointIds::const_iterator anIter = aSortedPointIds.begin();
584 TSortedPointIds::const_iterator anEndIter = aSortedPointIds.end();
585 for(vtkIdType anId = 0; anIter != anEndIter; anIter++, anId++){
586 vtkIdType aPntId = anIter->second;
587 #ifdef USE_ROBUST_TRIANGULATION
588 aConnectivities[anId] = aPts[aPntId];
590 aConnectivities[anId] = aPntId;
593 newCellId = output->InsertNextCell(aCellType,numFacePts,&aConnectivities[0]);
595 myVTK2ObjIds.push_back(cellId);
596 outputCD->CopyData(cd,cellId,newCellId);
605 for (faceId = 0; faceId < 4; faceId++)
608 faceVerts = vtkTetra::GetFaceArray(faceId);
609 faceIds->InsertNextId(pts[faceVerts[0]]);
610 faceIds->InsertNextId(pts[faceVerts[1]]);
611 faceIds->InsertNextId(pts[faceVerts[2]]);
612 aCellType = VTK_TRIANGLE;
614 input->GetCellNeighbors(cellId, faceIds, cellIds);
615 if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ||
616 (!allVisible && !cellVis[cellIds->GetId(0)]) )
618 for ( i=0; i < numFacePts; i++)
619 aNewPts[i] = pts[faceVerts[i]];
620 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
622 myVTK2ObjIds.push_back(cellId);
623 outputCD->CopyData(cd,cellId,newCellId);
629 for (faceId = 0; faceId < 6; faceId++)
632 faceVerts = vtkVoxel::GetFaceArray(faceId);
633 faceIds->InsertNextId(pts[faceVerts[0]]);
634 faceIds->InsertNextId(pts[faceVerts[1]]);
635 faceIds->InsertNextId(pts[faceVerts[2]]);
636 faceIds->InsertNextId(pts[faceVerts[3]]);
637 aCellType = VTK_QUAD;
639 input->GetCellNeighbors(cellId, faceIds, cellIds);
640 if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ||
641 (!allVisible && !cellVis[cellIds->GetId(0)]) )
643 for ( i=0; i < numFacePts; i++)
644 aNewPts[i] = pts[faceVerts[PixelConvert[i]]];
645 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
647 myVTK2ObjIds.push_back(cellId);
648 outputCD->CopyData(cd,cellId,newCellId);
653 case VTK_HEXAHEDRON: {
654 for (faceId = 0; faceId < 6; faceId++)
657 faceVerts = vtkHexahedron::GetFaceArray(faceId);
658 faceIds->InsertNextId(pts[faceVerts[0]]);
659 faceIds->InsertNextId(pts[faceVerts[1]]);
660 faceIds->InsertNextId(pts[faceVerts[2]]);
661 faceIds->InsertNextId(pts[faceVerts[3]]);
662 aCellType = VTK_QUAD;
664 input->GetCellNeighbors(cellId, faceIds, cellIds);
665 if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ||
666 (!allVisible && !cellVis[cellIds->GetId(0)]) )
668 for ( i=0; i < numFacePts; i++)
669 aNewPts[i] = pts[faceVerts[i]];
670 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
672 myVTK2ObjIds.push_back(cellId);
673 outputCD->CopyData(cd,cellId,newCellId);
679 for (faceId = 0; faceId < 5; faceId++)
682 faceVerts = vtkWedge::GetFaceArray(faceId);
683 faceIds->InsertNextId(pts[faceVerts[0]]);
684 faceIds->InsertNextId(pts[faceVerts[1]]);
685 faceIds->InsertNextId(pts[faceVerts[2]]);
686 aCellType = VTK_TRIANGLE;
688 if (faceVerts[3] >= 0)
690 faceIds->InsertNextId(pts[faceVerts[3]]);
691 aCellType = VTK_QUAD;
694 input->GetCellNeighbors(cellId, faceIds, cellIds);
695 if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ||
696 (!allVisible && !cellVis[cellIds->GetId(0)]) )
698 for ( i=0; i < numFacePts; i++)
699 aNewPts[i] = pts[faceVerts[i]];
700 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
702 myVTK2ObjIds.push_back(cellId);
703 outputCD->CopyData(cd,cellId,newCellId);
709 for (faceId = 0; faceId < 5; faceId++)
712 faceVerts = vtkPyramid::GetFaceArray(faceId);
713 faceIds->InsertNextId(pts[faceVerts[0]]);
714 faceIds->InsertNextId(pts[faceVerts[1]]);
715 faceIds->InsertNextId(pts[faceVerts[2]]);
716 aCellType = VTK_TRIANGLE;
718 if (faceVerts[3] >= 0)
720 faceIds->InsertNextId(pts[faceVerts[3]]);
721 aCellType = VTK_QUAD;
724 input->GetCellNeighbors(cellId, faceIds, cellIds);
725 if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ||
726 (!allVisible && !cellVis[cellIds->GetId(0)]) )
728 for ( i=0; i < numFacePts; i++)
729 aNewPts[i] = pts[faceVerts[i]];
730 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
732 myVTK2ObjIds.push_back(cellId);
733 outputCD->CopyData(cd,cellId,newCellId);
739 case VTK_QUADRATIC_EDGE:
740 case VTK_QUADRATIC_TRIANGLE:
741 case VTK_QUADRATIC_QUAD:
742 case VTK_QUADRATIC_TETRA:
743 case VTK_QUADRATIC_HEXAHEDRON:
744 if(!myIsWireframeMode){
745 input->GetCell(cellId,cell);
746 vtkIdList *pts = vtkIdList::New();
747 vtkPoints *coords = vtkPoints::New();
748 vtkIdList *cellIds = vtkIdList::New();
751 if ( cell->GetCellDimension() == 1 ) {
752 aCellType = VTK_LINE;
754 cell->Triangulate(0,pts,coords);
755 for (i=0; i < pts->GetNumberOfIds(); i+=2) {
756 aNewPts[0] = pts->GetId(i);
757 aNewPts[1] = pts->GetId(i+1);
758 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
760 myVTK2ObjIds.push_back(cellId);
761 outputCD->CopyData(cd,cellId,newCellId);
764 else if ( cell->GetCellDimension() == 2 ) {
765 aCellType = VTK_TRIANGLE;
767 cell->Triangulate(0,pts,coords);
768 for (i=0; i < pts->GetNumberOfIds(); i+=3) {
769 aNewPts[0] = pts->GetId(i);
770 aNewPts[1] = pts->GetId(i+1);
771 aNewPts[2] = pts->GetId(i+2);
772 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
774 myVTK2ObjIds.push_back(cellId);
775 outputCD->CopyData(cd,cellId,newCellId);
778 else //3D nonlinear cell
780 aCellType = VTK_TRIANGLE;
782 for (int j=0; j < cell->GetNumberOfFaces(); j++){
783 vtkCell *face = cell->GetFace(j);
784 input->GetCellNeighbors(cellId, face->PointIds, cellIds);
785 if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ) {
786 face->Triangulate(0,pts,coords);
787 for (i=0; i < pts->GetNumberOfIds(); i+=3) {
788 aNewPts[0] = pts->GetId(i);
789 aNewPts[1] = pts->GetId(i+1);
790 aNewPts[2] = pts->GetId(i+2);
791 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
793 myVTK2ObjIds.push_back(cellId);
794 outputCD->CopyData(cd,cellId,newCellId);
805 case VTK_QUADRATIC_EDGE: {
806 aCellType = VTK_POLY_LINE;
813 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
815 myVTK2ObjIds.push_back(cellId);
817 outputCD->CopyData(cd,cellId,newCellId);
820 case VTK_QUADRATIC_TRIANGLE: {
821 aCellType = VTK_POLYGON;
831 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
833 myVTK2ObjIds.push_back(cellId);
835 outputCD->CopyData(cd,cellId,newCellId);
838 case VTK_QUADRATIC_QUAD: {
839 aCellType = VTK_POLYGON;
851 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
853 myVTK2ObjIds.push_back(cellId);
855 outputCD->CopyData(cd,cellId,newCellId);
858 case VTK_QUADRATIC_TETRA: {
859 aCellType = VTK_POLYGON;
862 //---------------------------------------------------------------
870 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
872 myVTK2ObjIds.push_back(cellId);
874 outputCD->CopyData(cd,cellId,newCellId);
876 //---------------------------------------------------------------
884 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
886 myVTK2ObjIds.push_back(cellId);
888 outputCD->CopyData(cd,cellId,newCellId);
890 //---------------------------------------------------------------
898 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
900 myVTK2ObjIds.push_back(cellId);
902 outputCD->CopyData(cd,cellId,newCellId);
904 //---------------------------------------------------------------
912 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
914 myVTK2ObjIds.push_back(cellId);
916 outputCD->CopyData(cd,cellId,newCellId);
920 case VTK_QUADRATIC_HEXAHEDRON: {
921 aCellType = VTK_POLYGON;
924 //---------------------------------------------------------------
928 aNewPts[3] = pts[17];
930 aNewPts[5] = pts[12];
932 aNewPts[7] = pts[16];
934 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
936 myVTK2ObjIds.push_back(cellId);
938 outputCD->CopyData(cd,cellId,newCellId);
940 //---------------------------------------------------------------
944 aNewPts[3] = pts[18];
946 aNewPts[5] = pts[13];
948 aNewPts[7] = pts[17];
950 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
952 myVTK2ObjIds.push_back(cellId);
954 outputCD->CopyData(cd,cellId,newCellId);
956 //---------------------------------------------------------------
958 aNewPts[1] = pts[10];
960 aNewPts[3] = pts[19];
962 aNewPts[5] = pts[14];
964 aNewPts[7] = pts[18];
966 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
968 myVTK2ObjIds.push_back(cellId);
970 outputCD->CopyData(cd,cellId,newCellId);
972 //---------------------------------------------------------------
974 aNewPts[1] = pts[11];
976 aNewPts[3] = pts[16];
978 aNewPts[5] = pts[15];
980 aNewPts[7] = pts[19];
982 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
984 myVTK2ObjIds.push_back(cellId);
986 outputCD->CopyData(cd,cellId,newCellId);
988 //---------------------------------------------------------------
994 aNewPts[5] = pts[10];
996 aNewPts[7] = pts[11];
998 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1000 myVTK2ObjIds.push_back(cellId);
1002 outputCD->CopyData(cd,cellId,newCellId);
1004 //---------------------------------------------------------------
1005 aNewPts[0] = pts[4];
1006 aNewPts[1] = pts[12];
1007 aNewPts[2] = pts[5];
1008 aNewPts[3] = pts[13];
1009 aNewPts[4] = pts[6];
1010 aNewPts[5] = pts[14];
1011 aNewPts[6] = pts[7];
1012 aNewPts[7] = pts[15];
1014 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1016 myVTK2ObjIds.push_back(cellId);
1018 outputCD->CopyData(cd,cellId,newCellId);
1029 vtkDebugMacro(<<"Extracted " << input->GetNumberOfPoints() << " points,"
1030 << output->GetNumberOfCells() << " cells.");
1032 #ifdef USE_ROBUST_TRIANGULATION
1033 anUnstructuredGrid->Delete();
1034 aDelaunayPoints->Delete();
1036 aDelaunay3D->Delete();
1037 aGeometryFilter->Delete();
1052 //----------------------------------------------------------------------------
1054 VTKViewer_GeometryFilter
1055 ::SetInside(int theShowInside)
1057 if(myShowInside == theShowInside)
1060 myShowInside = theShowInside;
1065 VTKViewer_GeometryFilter
1068 return myShowInside;
1072 //----------------------------------------------------------------------------
1074 VTKViewer_GeometryFilter
1075 ::SetWireframeMode(int theIsWireframeMode)
1077 if(myIsWireframeMode == theIsWireframeMode)
1080 myIsWireframeMode = theIsWireframeMode;
1085 VTKViewer_GeometryFilter
1086 ::GetWireframeMode()
1088 return myIsWireframeMode;
1092 //----------------------------------------------------------------------------
1094 VTKViewer_GeometryFilter
1095 ::SetStoreMapping(int theStoreMapping)
1097 if(myStoreMapping == theStoreMapping)
1100 myStoreMapping = theStoreMapping;
1105 VTKViewer_GeometryFilter
1108 return myStoreMapping;
1112 //----------------------------------------------------------------------------
1113 vtkIdType VTKViewer_GeometryFilter::GetElemObjId(int theVtkID){
1114 if(myVTK2ObjIds.empty() || theVtkID > myVTK2ObjIds.size()) return -1;
1115 #if defined __GNUC_2__
1116 return myVTK2ObjIds[theVtkID];
1118 return myVTK2ObjIds.at(theVtkID);