-// This library is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
-// Lesser General Public License for more details.
-//
-// You should have received a copy of the GNU Lesser General Public
-// License along with this library; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
-//
-// See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org
+// Copyright (C) 2007-2023 CEA/DEN, EDF R&D, OPEN CASCADE
+//
+// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
+// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
#include "VTKViewer_ConvexTool.h"
-
-#include <vtkUnstructuredGrid.h>
-#include <vtkTriangle.h>
-#include <vtkConvexPointSet.h>
-#include <vtkMath.h>
+#include "VTKViewer_GeometryFilter.h"
#include <set>
-#include <iterator>
+#include <map>
#include <algorithm>
-#include <math.h>
-
-static float FACE_ANGLE_TOLERANCE=1.5;
+#include <iterator>
-typedef std::set<vtkIdType> TUIDS; // unique ids
-typedef std::map<vtkIdType,TUIDS> TPTOIDS; // id points -> unique ids
+#include <vtkUnstructuredGrid.h>
+#include <vtkGeometryFilter.h>
+#include <vtkDelaunay3D.h>
+#include <vtkGenericCell.h>
+#include <vtkPolyData.h>
+#include <vtkCellData.h>
+#include <vtkPoints.h>
+#include <vtkIdList.h>
+#include <vtkCell.h>
+#include <vtkPlane.h>
+#include <vtkMath.h>
+#include <vtkCellArray.h>
+#include <vtkTriangle.h>
+#include <vtkOrderedTriangulator.h>
-namespace CONVEX_TOOL
-{
#ifdef _DEBUG_
- static int MYDEBUG = 0;
+static int DEBUG_TRIA_EXECUTE = 0;
#else
- static int MYDEBUG = 0;
+static int DEBUG_TRIA_EXECUTE = 0;
#endif
-static void GetCenter(vtkUnstructuredGrid* theGrid,TCell theptIds,float *center)
+namespace
{
- float *p;
- center[0] = center[1] = center[2] = 0.0;
-
- int numIds=theptIds.size();
-
- // get the center of the cell
- for (int i=0; i < numIds; i++)
- {
- p = theGrid->GetPoint(theptIds[i]);
- for (int j=0; j < 3; j++)
- {
- center[j] += p[j];
- }
- }
- for (int j=0; j<3; j++)
+ typedef std::vector<vtkIdType> TConnectivities;
+
+ struct TPolygon
+ {
+ TConnectivities myConnectivities;
+ double myOrigin[3];
+ double myNormal[3];
+ TPolygon(const TConnectivities& theConnectivities,
+ double theOrigin[3],
+ double theNormal[3]):
+ myConnectivities(theConnectivities)
{
- center[j] /= numIds;
+ myOrigin[0] = theOrigin[0];
+ myOrigin[1] = theOrigin[1];
+ myOrigin[2] = theOrigin[2];
+
+ myNormal[0] = theNormal[0];
+ myNormal[1] = theNormal[1];
+ myNormal[2] = theNormal[2];
}
+ };
+
+ typedef std::vector<TPolygon> TPolygons;
+}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_Triangulator
+::VTKViewer_Triangulator():
+ myCellIds(vtkIdList::New()),
+ myPointIds(NULL),
+ myFaceIds(vtkIdList::New()),
+ myPoints(vtkPoints::New())
+{}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_Triangulator
+::~VTKViewer_Triangulator()
+{
+ myCellIds->Delete();
+ myFaceIds->Delete();
+ myPoints->Delete();
}
-static void ReverseIds(TCell &theIds)
+
+//----------------------------------------------------------------------------
+vtkPoints*
+VTKViewer_Triangulator
+::InitPoints(vtkUnstructuredGrid *theInput,
+ vtkIdType theCellId)
{
- int i;
- vtkIdType tmp;
- vtkIdType npts=theIds.size();
-
- for(i=0;i<(npts/2);i++){
- tmp = theIds[i];
- theIds[i] = theIds[npts-i-1];
- theIds[npts-i-1] = tmp;
+ myPoints->Reset();
+ myPoints->Modified(); // the VTK bug
+
+ vtkIdType aNumPts;
+ vtkIdType const *tmp;
+ theInput->GetCellPoints(theCellId, aNumPts, tmp);
+ myPointIds = const_cast<vtkIdType*>( tmp );
+ if ( aNumPts > 0 ) {
+ double anAbsoluteCoord[3];
+ myPoints->SetNumberOfPoints(aNumPts);
+ vtkPoints *anInputPoints = theInput->GetPoints();
+ for (auto aPntId = 0; aPntId < aNumPts; aPntId++) {
+ anInputPoints->GetPoint(myPointIds[aPntId], anAbsoluteCoord);
+ myPoints->SetPoint(aPntId, anAbsoluteCoord);
+ }
}
+
+ return myPoints;
}
-// caclulation of connected faces (faceId -> (faceId1,faceId2, ...))
-void GetFriends(const TPTOIDS p2faces,const TCellArray f2points,TPTOIDS& face2face_output)
+
+//----------------------------------------------------------------------------
+vtkIdType
+VTKViewer_Triangulator
+::GetNbOfPoints()
{
- TCellArray::const_iterator f2pIter = f2points.begin();
+ return myPoints->GetNumberOfPoints();
+}
- for( ; f2pIter!=f2points.end() ; f2pIter++ ){
- vtkIdType faceId = f2pIter->first;
- TCell face_points = f2pIter->second;
- int nb_face_points = face_points.size();
-
- vtkIdType id1;
- vtkIdType id2;
- TPTOIDS::const_iterator faces1;
- TPTOIDS::const_iterator faces2;
-
- id1 = face_points[0];
- faces1 = p2faces.find(id1);
-
- TUIDS output_faces;
-
- for(int i=1 ; i<nb_face_points ; i++ ){
- id2 = face_points[i];
+//----------------------------------------------------------------------------
+vtkIdType
+VTKViewer_Triangulator
+::GetPointId(vtkIdType thePointId)
+{
+ return thePointId;
+}
- faces2 = p2faces.find(id2);
-
- std::set_intersection(faces1->second.begin(), faces1->second.end(), faces2->second.begin(), faces2->second.end(),
- std::inserter(output_faces,output_faces.begin()));
-
- id1 = id2;
- faces1 = faces2;
- }
- id1 = face_points[0];
- faces1 = p2faces.find(id1);
- std::set_intersection(faces1->second.begin(), faces1->second.end(), faces2->second.begin(), faces2->second.end(),
- std::inserter(output_faces,output_faces.begin()));
-
- output_faces.erase(faceId); // erase the face id for which we found friends
- if(MYDEBUG){
- cout << "fId[" << faceId <<"]: ";
- std::copy(output_faces.begin(), output_faces.end(), std::ostream_iterator<vtkIdType>(cout, " "));
- cout << endl;
- }
-
- face2face_output[faceId] = output_faces;
- }
+//----------------------------------------------------------------------------
+double
+VTKViewer_Triangulator
+::GetCellLength()
+{
+ double aBounds[6];
+ myPoints->GetBounds(aBounds);
+
+ double aCoordDiff[3];
+ aCoordDiff[0] = (aBounds[1] - aBounds[0]);
+ aCoordDiff[1] = (aBounds[3] - aBounds[2]);
+ aCoordDiff[2] = (aBounds[5] - aBounds[4]);
+
+ return sqrt(aCoordDiff[0]*aCoordDiff[0] +
+ aCoordDiff[1]*aCoordDiff[1] +
+ aCoordDiff[2]*aCoordDiff[2]);
}
-bool IsConnectedFacesOnOnePlane( vtkUnstructuredGrid* theGrid,
- vtkIdType theFId1, vtkIdType theFId2,
- TUIDS FpIds1, TUIDS FpIds2 )
+
+//----------------------------------------------------------------------------
+void
+VTKViewer_Triangulator
+::GetCellNeighbors(vtkUnstructuredGrid *theInput,
+ vtkIdType theCellId,
+ vtkCell* theFace,
+ vtkIdList* theCellIds)
{
- bool status = false;
- TUIDS common_ids;
- std::set_intersection(FpIds1.begin(), FpIds1.end(), FpIds2.begin(), FpIds2.end(),
- std::inserter(common_ids,common_ids.begin()));
-
- /* Number of common ids = 2 (A1,A2)
-
-
- _ _ _ _ _ _ _ _ _ vectors:
- | \ / | v1 {A2,A1}
- \ / v2 {A1,B1}
- | | A2 | v3 {A1,C1}
- |
- | | |
- |
- | | A1 |
- / \
- |_ _ _ _ _/ \_ _ _ _ _|
- B2 B1 C1 C2
-
- */
- TUIDS::iterator common_iter = common_ids.begin();
- if(common_ids.size() == 2){
- TUIDS::iterator loc_id1_0 = FpIds1.find(*(common_iter));
- common_iter++;
- TUIDS::iterator loc_id1_1 = FpIds1.find(*(common_iter));
-
- TUIDS::iterator loc_id2_0 = FpIds1.begin();
- TUIDS::iterator loc_id2_1 = FpIds2.begin();
-
- vtkIdType A1 = *loc_id1_0;
- vtkIdType A2 = *loc_id1_1;
- vtkIdType B1;
- vtkIdType C1;
-
- for(;loc_id2_0!=FpIds1.end();loc_id2_0++)
- if(*loc_id2_0 != A1 && *loc_id2_0!= A2){
- B1 = *loc_id2_0;
- break;
- }
- for(;loc_id2_1!=FpIds2.end();loc_id2_1++)
- if(*loc_id2_1 != A1 && *loc_id2_1!= A2){
- C1 = *loc_id2_1;
- break;
- }
- if(MYDEBUG) cout <<endl;
- if(MYDEBUG) cout << "FId_1="<<theFId1<<" FId_2="<<theFId2<<endl;
- if(MYDEBUG) cout << " A1="<<A1<<" A2="<<A2<<" B1="<<B1<<" C1="<<C1<<" ->";
- float *p[4];
- float v1[3],v2[3],v3[3];
- p[0] = theGrid->GetPoint(A1);
- p[1] = theGrid->GetPoint(A2);
- p[2] = theGrid->GetPoint(B1);
- p[3] = theGrid->GetPoint(C1);
-
- for(int i=0;i<3;i++){
- v1[i] = p[1][i] - p[0][i];
- v2[i] = p[2][i] - p[0][i];
- v3[i] = p[3][i] - p[0][i];
- }
-
-
- float vec_b1[3];
- vtkMath::Cross(v2,v1,vec_b1);
- float vec_b2[3];
- vtkMath::Cross(v1,v3,vec_b2);
+ myFaceIds->Reset();
+ vtkIdList *anIdList = theFace->PointIds;
+ myFaceIds->InsertNextId(myPointIds[anIdList->GetId(0)]);
+ myFaceIds->InsertNextId(myPointIds[anIdList->GetId(1)]);
+ myFaceIds->InsertNextId(myPointIds[anIdList->GetId(2)]);
- float b1 = vtkMath::Norm(vec_b1);
+ theInput->GetCellNeighbors(theCellId, myFaceIds, theCellIds);
+}
- float b2 = vtkMath::Norm(vec_b2);
- float aCos = vtkMath::Dot(vec_b1,vec_b2)/(b1*b2);
-
- float angle=90.0;
- angle = aCos>=1.0 ? 0.0 : 180*acosf(aCos)/vtkMath::Pi();
-
- if( angle <= FACE_ANGLE_TOLERANCE)
- status = true;
- if (MYDEBUG){
- for(int k=0;k<4;k++){
- cout << " (";
- for(int j=0;j<2;j++){
- cout << p[k][j] << ",";
- }
- cout << p[k][2] << ") ";
- }
- cout << "angle="<<angle<<" status="<<status<<endl;
- }
-
- } else if (common_ids.size() >2){
- // not implemented yet
- if(MYDEBUG) cout << "Warning! VTKViewer_ConvexTool::IsConnectedFacesOnOnePlane()";
- } else {
- // one or no connection ... continue
- }
-
- return status;
+
+//----------------------------------------------------------------------------
+vtkIdType
+VTKViewer_Triangulator
+::GetConnectivity(vtkIdType thePntId)
+{
+ return myPointIds[thePntId];
}
-void GetAllFacesOnOnePlane( TPTOIDS theFaces, vtkIdType faceId,
- TUIDS &new_faces, TCell &new_faces_v2 )
+
+//----------------------------------------------------------------------------
+bool
+VTKViewer_Triangulator
+::Execute(vtkUnstructuredGrid *theInput,
+ vtkCellData* thInputCD,
+ vtkIdType theCellId,
+ int theShowInside,
+ int theAllVisible,
+ int theAppendCoincident3D,
+ const char* theCellsVisibility,
+ vtkPolyData *theOutput,
+ vtkCellData* theOutputCD,
+ int theStoreMapping,
+ std::vector<vtkIdType>& theVTK2ObjIds,
+ std::vector< std::vector<vtkIdType> >& theDimension2VTK2ObjIds,
+ bool theIsCheckConvex)
{
- if (new_faces.find(faceId) != new_faces.end()) return;
-
- new_faces.insert(new_faces.begin(),faceId);
- new_faces_v2.push_back(faceId);
-
- TPTOIDS::const_iterator aIter1 = theFaces.find(faceId);
- if(aIter1!=theFaces.end()){
- TUIDS::const_iterator aIter2 = (aIter1->second).begin();
- for(;aIter2!=(aIter1->second).end();aIter2++){
- if (new_faces.find(*aIter2) != new_faces.end()) continue;
- GetAllFacesOnOnePlane(theFaces,*aIter2,
- new_faces,new_faces_v2); // recurvise
+ vtkPoints *aPoints = InitPoints(theInput, theCellId);
+ vtkIdType aNumPts = GetNbOfPoints();
+ if(DEBUG_TRIA_EXECUTE) cout<<"Triangulator - aNumPts = "<<aNumPts<<"\n";
+
+ if(aNumPts == 0)
+ return true;
+
+ // To calculate the bary center of the cell
+ double aCellCenter[3] = {0.0, 0.0, 0.0};
+ {
+ double aPntCoord[3];
+ for (auto aPntId = 0; aPntId < aNumPts; aPntId++) {
+ aPoints->GetPoint(GetPointId(aPntId),aPntCoord);
+ if(DEBUG_TRIA_EXECUTE) cout<<"\taPntId = "<<GetPointId(aPntId)<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}\n";
+ aCellCenter[0] += aPntCoord[0];
+ aCellCenter[1] += aPntCoord[1];
+ aCellCenter[2] += aPntCoord[2];
}
+ aCellCenter[0] /= aNumPts;
+ aCellCenter[1] /= aNumPts;
+ aCellCenter[2] /= aNumPts;
}
- return;
-}
-void GetSumm(TCell v1,TCell v2,TCell &output)
-{
- output.clear();
-
- if(MYDEBUG) cout << "========================================="<<endl;
- if(MYDEBUG) cout << "v1:";
- if(MYDEBUG) std::copy(v1.begin(), v1.end(), std::ostream_iterator<vtkIdType>(cout, " "));
- if(MYDEBUG) cout << "\tv2:";
- if(MYDEBUG) std::copy(v2.begin(), v2.end(), std::ostream_iterator<vtkIdType>(cout, " "));
- if(MYDEBUG) cout << endl;
-
- TUIDS v1_set;
- std::copy(v1.begin(), v1.end(), std::inserter(v1_set,v1_set.begin()));
- TUIDS v2_set;
- std::copy(v2.begin(), v2.end(), std::inserter(v2_set,v2_set.begin()));
- TUIDS tmpIntersection;
- std::set_intersection(v1_set.begin(),v1_set.end(),v2_set.begin(),v2_set.end(), std::inserter(tmpIntersection,tmpIntersection.begin()));
- if(MYDEBUG) std::copy(tmpIntersection.begin(),tmpIntersection.end(), std::ostream_iterator<vtkIdType>(cout, " "));
- if(MYDEBUG) cout << endl;
-
- if(tmpIntersection.size() < 2)
- if(MYDEBUG) cout << __FILE__ << "[" << __LINE__ << "]: Warning ! Wrong ids" << endl;
+ double aCellLength = GetCellLength();
+ vtkIdType aNumFaces = GetNumFaces();
+
+ static double EPS = 1.0E-2;
+ double aDistEps = aCellLength/3.0 * EPS;
+ if(DEBUG_TRIA_EXECUTE) cout<<"\taNumFaces = "<<aNumFaces<<"; aCellLength = "<<aCellLength<<"; aDistEps = "<<aDistEps<<"\n";
+
+ // To initialize set of points that belong to the cell
+ typedef std::set<vtkIdType> TPointIds;
+ TPointIds anInitialPointIds;
+ for(vtkIdType aPntId = 0; aPntId < aNumPts; aPntId++)
+ anInitialPointIds.insert(GetPointId(aPntId));
- TCell::iterator v1_iter = v1.begin();
+ // To initialize set of points by face that belong to the cell and backward
+ typedef std::set<vtkIdType> TFace2Visibility;
+ TFace2Visibility aFace2Visibility;
- for(;v1_iter!=v1.end();v1_iter++){
-
- vtkIdType curr_id = *v1_iter;
-
- output.push_back(curr_id);
+ typedef std::set<TPointIds> TFace2PointIds;
+ TFace2PointIds aFace2PointIds;
+
+ for (auto aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
+ vtkCell* aFace = GetFace(aFaceId);
- if(tmpIntersection.find(curr_id) != tmpIntersection.end()){
- TCell::iterator v1_iter_tmp;
- v1_iter_tmp = v1_iter;
- v1_iter++;
-
- if(v1_iter==v1.end()) v1_iter=v1.begin();
-
- curr_id = *v1_iter;
-
- if(tmpIntersection.find(curr_id) != tmpIntersection.end()){
- TCell::iterator v2_iter = v2.begin();
- for(;v2_iter!=v2.end();v2_iter++){
- vtkIdType v2_id = *v2_iter;
- if(tmpIntersection.find(v2_id) == tmpIntersection.end())
- output.push_back(v2_id);
- }
- }
-
- v1_iter = v1_iter_tmp;
- curr_id = *v1_iter;
+ GetCellNeighbors(theInput, theCellId, aFace, myCellIds);
+ bool process = myCellIds->GetNumberOfIds() <= 0 ? true : theAppendCoincident3D;
+ if((!theAllVisible && !theCellsVisibility[myCellIds->GetId(0)]) ||
+ myCellIds->GetNumberOfIds() <= 0 || theShowInside || process)
+ {
+ TPointIds aPointIds;
+ vtkIdList *anIdList = aFace->PointIds;
+ aPointIds.insert(anIdList->GetId(0));
+ aPointIds.insert(anIdList->GetId(1));
+ aPointIds.insert(anIdList->GetId(2));
+ aFace2PointIds.insert(aPointIds);
+ aFace2Visibility.insert(aFaceId);
}
}
- if(MYDEBUG) cout << "Result: " ;
- if(MYDEBUG) std::copy(output.begin(),output.end(),std::ostream_iterator<vtkIdType>(cout, " "));
- if(MYDEBUG) cout << endl;
-}
-void GetPolygonalFaces(vtkUnstructuredGrid* theGrid,int cellId,TCellArray &outputCellArray)
-{
- if (theGrid->GetCellType(cellId) == VTK_CONVEX_POINT_SET){
- // get vtkCell type = VTK_CONVEX_POINT_SET
- if(vtkConvexPointSet* convex = static_cast<vtkConvexPointSet*>(theGrid->GetCell(cellId))){
- TCellArray f2points;
- float convex_center[3]; // convex center point coorinat
- int aNbFaces = convex->GetNumberOfFaces();
- int numPts = convex->GetNumberOfPoints();
- TCell convex_ids;
- TPTOIDS p2faces; // key=pointId , value facesIds set
+ ::TPolygons aPolygons;
+
+ for (auto aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
+ if(aFace2Visibility.find(aFaceId) == aFace2Visibility.end())
+ continue;
+
+ vtkCell* aFace = GetFace(aFaceId);
+
+ vtkIdList *anIdList = aFace->PointIds;
+ vtkIdType aNewPts[3] = {anIdList->GetId(0), anIdList->GetId(1), anIdList->GetId(2)};
+
+ // To initialize set of points for the plane where the trinangle face belong to
+ TPointIds aPointIds;
+ aPointIds.insert(aNewPts[0]);
+ aPointIds.insert(aNewPts[1]);
+ aPointIds.insert(aNewPts[2]);
+
+ // To get know, if the points of the trinagle were already observed
+ bool anIsObserved = aFace2PointIds.find(aPointIds) == aFace2PointIds.end();
+ if(DEBUG_TRIA_EXECUTE) {
+ cout<<"\taFaceId = "<<aFaceId<<"; anIsObserved = "<<anIsObserved;
+ cout<<"; aNewPts = {"<<aNewPts[0]<<", "<<aNewPts[1]<<", "<<aNewPts[2]<<"}\n";
+ }
+
+ if(!anIsObserved){
+ // To get coordinates of the points of the traingle face
+ double aCoord[3][3];
+ aPoints->GetPoint(aNewPts[0],aCoord[0]);
+ aPoints->GetPoint(aNewPts[1],aCoord[1]);
+ aPoints->GetPoint(aNewPts[2],aCoord[2]);
+
+ /* To calculate plane normal for face (aFace)
+
+
+ ^ aNormal
+ |
+ | ^ aVector01
+ | /
+ /_________> aVector02
+
+
+ */
+ double aVector01[3] = { aCoord[1][0] - aCoord[0][0],
+ aCoord[1][1] - aCoord[0][1],
+ aCoord[1][2] - aCoord[0][2] };
- for(int i=0;i<numPts;i++)
- convex_ids.push_back(convex->GetPointId(i));
+ double aVector02[3] = { aCoord[2][0] - aCoord[0][0],
+ aCoord[2][1] - aCoord[0][1],
+ aCoord[2][2] - aCoord[0][2] };
- GetCenter(theGrid,convex_ids,convex_center);
-
- for (vtkIdType faceId=0; faceId < aNbFaces; faceId++){
- vtkCell *aFace = convex->GetFace(faceId);
- int numFacePts = aFace->GetNumberOfPoints();
- TCell aIds;
-
- int i = 0;
- for(i=0;i<numFacePts;i++)
- aIds.push_back(aFace->GetPointId(i));
-
- float v_a[3],v_b[3],v_convex2face[3]; // vectors
- float *id_0,*id_1,*id_n;
- /*=============================================
- ,+- - - - _
- _ / id_n | v_b {id_0,id_n}
- v_b / _
- / | v_a {id_0,id_1}
- /
- / |
- + id_0
- \
- _ \ |
- v_a \
- \ id_1 |
- "+- - - -
-
- ============================================*/
- id_0 = theGrid->GetPoint(aIds[0]);
- id_1 = theGrid->GetPoint(aIds[1]);
- id_n = theGrid->GetPoint(aIds[numFacePts-1]);
-
- for(i=0;i<3;i++){
- v_a[i] = id_1[i] - id_0[i];
- v_b[i] = id_n[i] - id_0[i];
- v_convex2face[i] = id_0[i] - convex_center[i];
- }
-
- if (vtkMath::Determinant3x3(v_a,v_b,v_convex2face) < 0){
- ReverseIds(aIds);
- }
-
- for(i=0;i<(int)aIds.size();i++){
- TUIDS &acell = p2faces[aIds[i]];
- acell.insert(faceId);
- }
-
- f2points[faceId] = aIds;
+ vtkMath::Normalize(aVector01);
+ vtkMath::Normalize(aVector02);
- }
+ // To calculate the normal for the triangle
+ double aNormal[3];
+ vtkMath::Cross(aVector02,aVector01,aNormal);
- TPTOIDS face2face;
- GetFriends(p2faces,f2points,face2face);
+ vtkMath::Normalize(aNormal);
- TPTOIDS face2points;
+ // To calculate what points belong to the plane
+ // To calculate bounds of the point set
+ double aCenter[3] = {0.0, 0.0, 0.0};
+ {
+ TPointIds::const_iterator anIter = anInitialPointIds.begin();
+ TPointIds::const_iterator anEndIter = anInitialPointIds.end();
+ for(; anIter != anEndIter; anIter++){
+ double aPntCoord[3];
+ vtkIdType aPntId = *anIter;
+ aPoints->GetPoint(aPntId,aPntCoord);
+
+ double aVector0Pnt[3] = { aPntCoord[0] - aCoord[0][0],
+ aPntCoord[1] - aCoord[0][1],
+ aPntCoord[2] - aCoord[0][2] };
+
+
+ vtkMath::Normalize(aVector0Pnt);
+
+ double aNormalPnt[3];
+ // calculate aNormalPnt
+ {
+ double aCosPnt01 = vtkMath::Dot(aVector0Pnt,aVector01);
+ double aCosPnt02 = vtkMath::Dot(aVector0Pnt,aVector02);
+ if(aCosPnt01<-1)
+ aCosPnt01 = -1;
+ if(aCosPnt01>1)
+ aCosPnt01 = 1;
+ if(aCosPnt02<-1)
+ aCosPnt02 = -1;
+ if(aCosPnt02>1)
+ aCosPnt02 = 1;
+
+ double aDist01,aDist02;// deflection from Pi/3 angle (equilateral triangle)
+ double aAngPnt01 = fabs(acos(aCosPnt01));
+ double aAngPnt02 = fabs(acos(aCosPnt02));
+
+ /* check that triangle similar to equilateral triangle
+ AOC or COB ?
+ aVector0Pnt = (OC)
+ aVector01 = (OB)
+ aVector02 = (OA)
+
+ B
+ ^ aVector01 C
+ | ^ aVector0Pnt
+ | _____/
+ | ___/
+ |/________> aVector02
+ O A
+ */
+ aDist01 = fabs(aAngPnt01-(vtkMath::Pi())/3.0);
+ aDist02 = fabs(aAngPnt02-(vtkMath::Pi())/3.0);
+
+ // caculate a normal for best triangle
+ if(aDist01 <= aDist02)
+ vtkMath::Cross(aVector0Pnt,aVector01,aNormalPnt);
+ else
+ vtkMath::Cross(aVector0Pnt,aVector02,aNormalPnt);
+
+ }
+
+ vtkMath::Normalize(aNormalPnt);
+
+ if(DEBUG_TRIA_EXECUTE)
+ cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"};";
+
+ double aDist = vtkPlane::DistanceToPlane(aPntCoord,aNormal,aCoord[0]);
+ if(DEBUG_TRIA_EXECUTE) cout<<": aDist = "<<aDist;
+ if(fabs(aDist) < aDistEps){
+ aPointIds.insert(aPntId);
+ aCenter[0] += aPntCoord[0];
+ aCenter[1] += aPntCoord[1];
+ aCenter[2] += aPntCoord[2];
+ if(DEBUG_TRIA_EXECUTE) cout << "; Added = TRUE" << endl;
+ } else {
+ if(DEBUG_TRIA_EXECUTE) cout << "; Added = FALSE" << endl;
+ }
+ }
+ size_t aNbPoints = aPointIds.size();
+ aCenter[0] /= aNbPoints;
+ aCenter[1] /= aNbPoints;
+ aCenter[2] /= aNbPoints;
+ }
- // copy TCellArray::f2points to TPTOIDS::face2points
- for(TCellArray::iterator f2points_iter=f2points.begin();
- f2points_iter!=f2points.end();
- f2points_iter++){
-
- TUIDS tmp;
- for(TCell::iterator points_iter=(f2points_iter->second).begin();
- points_iter!=(f2points_iter->second).end();
- points_iter++)
- tmp.insert(*points_iter);
-
- face2points[f2points_iter->first] = tmp;
- } // end copy
-
+ //To sinchronize orientation of the cell and its face
+ double aVectorC[3] = { aCenter[0] - aCellCenter[0],
+ aCenter[1] - aCellCenter[1],
+ aCenter[2] - aCellCenter[2] };
+ vtkMath::Normalize(aVectorC);
- TPTOIDS new_face2faces; // which connected and in one plane
+ double aDot = vtkMath::Dot(aNormal,aVectorC);
+ if(DEBUG_TRIA_EXECUTE) {
+ cout<<"\t\taNormal = {"<<aNormal[0]<<", "<<aNormal[1]<<", "<<aNormal[2]<<"}";
+ cout<<"; aVectorC = {"<<aVectorC[0]<<", "<<aVectorC[1]<<", "<<aVectorC[2]<<"}\n";
+ cout<<"\t\taDot = "<<aDot<<"\n";
+ }
+ if(aDot > 0){
+ aNormal[0] = -aNormal[0];
+ aNormal[1] = -aNormal[1];
+ aNormal[2] = -aNormal[2];
+ }
- TPTOIDS::const_iterator aF2FIter = face2face.begin();
- for(;aF2FIter!=face2face.end();aF2FIter++){
- vtkIdType f_key = aF2FIter->first;
- TUIDS &faces = new_face2faces[f_key];
- //faces.insert(f_key);
- TUIDS f_friends = aF2FIter->second;
- TUIDS::const_iterator a_friends_iter = f_friends.begin();
- for(;a_friends_iter!=f_friends.end();a_friends_iter++){
- vtkIdType friend_id = *a_friends_iter;
- if( IsConnectedFacesOnOnePlane(theGrid,f_key,friend_id,
- (face2points.find(f_key))->second,
- (face2points.find(friend_id))->second)){
- faces.insert(friend_id);
- } // end if
-
- } // end a_friends_iter
- } // end aF2FIter
+ // To calculate the primary direction for point set
+ double aVector0[3] = { aCoord[0][0] - aCenter[0],
+ aCoord[0][1] - aCenter[1],
+ aCoord[0][2] - aCenter[2] };
+ vtkMath::Normalize(aVector0);
- if(MYDEBUG)
- {
- TPTOIDS::const_iterator new_face2face_iter = new_face2faces.begin();
- cout << "Connected faces and on plane:" << endl;
- for(;new_face2face_iter!=new_face2faces.end();new_face2face_iter++){
- cout << "Group ["<<new_face2face_iter->first<<"] :";
- TUIDS::const_iterator new_faces_iter = (new_face2face_iter->second).begin();
- for(;new_faces_iter!=(new_face2face_iter->second).end();new_faces_iter++)
- cout << " " << *new_faces_iter ;
- cout << endl;
- }
+ if(DEBUG_TRIA_EXECUTE) {
+ cout<<"\t\taCenter = {"<<aCenter[0]<<", "<<aCenter[1]<<", "<<aCenter[2]<<"}";
+ cout<<"; aVector0 = {"<<aVector0[0]<<", "<<aVector0[1]<<", "<<aVector0[2]<<"}\n";
}
- TPTOIDS output_newid2face;
- TCellArray output_newid2face_v2;
+ // To calculate the set of points by face those that belong to the plane
+ TFace2PointIds aRemoveFace2PointIds;
{
- TUIDS already_in;
- TUIDS already_in_tmp;
- int k=0;
- TPTOIDS::const_iterator new_face2face_iter = new_face2faces.begin();
- for(;new_face2face_iter!=new_face2faces.end();new_face2face_iter++){
- if(already_in.find(new_face2face_iter->first) != already_in.end())
- continue;
- if(new_face2face_iter->second.size() > 1)
- continue;
-
- TCell &tmp_v2 = output_newid2face_v2[k];
- tmp_v2.push_back(new_face2face_iter->first);
- already_in.insert(new_face2face_iter->first);
-
- TUIDS::const_iterator new_faces_iter = (new_face2face_iter->second).begin();
- for(;new_faces_iter!=(new_face2face_iter->second).end();new_faces_iter++){
- if(already_in.find(*new_faces_iter) != already_in.end()) continue;
- already_in.insert(*new_faces_iter);
-
- already_in_tmp.clear();
- already_in_tmp.insert(new_face2face_iter->first);
-
- TUIDS &tmp = output_newid2face[k];
- GetAllFacesOnOnePlane(new_face2faces,*new_faces_iter,
- already_in_tmp,tmp_v2);
-
- for(TUIDS::const_iterator aIter=already_in_tmp.begin();
- aIter!=already_in_tmp.end();
- aIter++)
- {
- already_in.insert(*aIter);
- tmp.insert(*aIter);
- }
- }
- k++;
- }
+ TFace2PointIds::const_iterator anIter = aFace2PointIds.begin();
+ TFace2PointIds::const_iterator anEndIter = aFace2PointIds.end();
+ for(; anIter != anEndIter; anIter++){
+ const TPointIds& anIds = *anIter;
+ TPointIds anIntersection;
+ std::set_intersection(aPointIds.begin(),aPointIds.end(),
+ anIds.begin(),anIds.end(),
+ std::inserter(anIntersection,anIntersection.begin()));
+
+
+ if(DEBUG_TRIA_EXECUTE) {
+ cout << "anIntersection:";
+ TPointIds::iterator aII = anIntersection.begin();
+ for(;aII!=anIntersection.end();aII++)
+ cout << *aII << ",";
+ cout << endl;
+ cout << "anIds :";
+ TPointIds::const_iterator aIIds = anIds.begin();
+ for(;aIIds!=anIds.end();aIIds++)
+ cout << *aIIds << ",";
+ cout << endl;
+ }
+ if(anIntersection == anIds){
+ aRemoveFace2PointIds.insert(anIds);
+ }
+ }
}
- if(MYDEBUG) {
- cout << "++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl;
- {
- TPTOIDS::const_iterator new_face2face_iter = output_newid2face.begin();
- for(;new_face2face_iter!=output_newid2face.end();new_face2face_iter++){
- cout << "Group ["<<new_face2face_iter->first<<"] :";
- TUIDS::const_iterator new_faces_iter = (new_face2face_iter->second).begin();
- for(;new_faces_iter!=(new_face2face_iter->second).end();new_faces_iter++)
- cout << " " << *new_faces_iter ;
- cout << endl;
- }
- }
- cout << "++++++++++++++++++++++++++ +++++++++++++++++++++++++++++++"<<endl;
- cout << "++++++++++++++++++++++++++ +++++++++++++++++++++++++++++++"<<endl;
- cout << "+++++++++++++++++++++++ ++ ++ ++++++++++++++++++++++++++++"<<endl;
- cout << "+++++++++++++++++++++++++ ++++++++++++++++++++++++++++++"<<endl;
- cout << "++++++++++++++++++++++++++ +++++++++++++++++++++++++++++++"<<endl;
- {
- TCellArray::const_iterator new_face2face_iter = output_newid2face_v2.begin();
- for(;new_face2face_iter!=output_newid2face_v2.end();new_face2face_iter++){
- cout << "Group ["<<new_face2face_iter->first<<"] :";
- TCell::const_iterator new_faces_iter = (new_face2face_iter->second).begin();
- for(;new_faces_iter!=(new_face2face_iter->second).end();new_faces_iter++)
- cout << " " << *new_faces_iter ;
- cout << endl;
- }
- }
+ // To remove from the set of points by face those that belong to the plane
+ {
+ TFace2PointIds::const_iterator anIter = aRemoveFace2PointIds.begin();
+ TFace2PointIds::const_iterator anEndIter = aRemoveFace2PointIds.end();
+ for(; anIter != anEndIter; anIter++){
+ const TPointIds& anIds = *anIter;
+ aFace2PointIds.erase(anIds);
+ }
}
- TCellArray output_new_face2ids;
-// {
-// TPTOIDS::const_iterator new_face2face_iter = output_newid2face.begin();
-// for(;new_face2face_iter!=output_newid2face.end();new_face2face_iter++){
-
-// vtkIdType new_faceId = new_face2face_iter->first;
-// TUIDS::const_iterator new_faces_iter = (new_face2face_iter->second).begin();
-// vtkIdType fId0 = *new_faces_iter;
-// TCellArray::const_iterator pIds0_iter = f2points.find(fId0);
-// TCell pIds0 = pIds0_iter->second;
-// TCell &output = output_new_face2ids[new_faceId];
-// new_faces_iter++;
-// if(new_face2face_iter->second.size() > 2 ){}
-// for(;new_faces_iter!=(new_face2face_iter->second).end();new_faces_iter++){
-
-// vtkIdType faceId = *new_faces_iter;
-// // find how much nodes in face (f2points)
-// TCellArray::const_iterator pIds_iter = f2points.find(faceId);
-// TCell pIds = pIds_iter->second;
-
-// GetSumm(pIds0,pIds,output);
-// pIds0 = output;
-
-// } // end new_faces_iter
-
-// } // new_face2face_iter
-// }
+ // To sort the planar set of the points accrding to the angle
{
- TCellArray::const_iterator new_face2face_iter = output_newid2face_v2.begin();
- for(;new_face2face_iter!=output_newid2face_v2.end();new_face2face_iter++){
-
- vtkIdType new_faceId = new_face2face_iter->first;
- TCell::const_iterator new_faces_iter = (new_face2face_iter->second).begin();
- vtkIdType fId0 = *new_faces_iter;
- TCellArray::const_iterator pIds0_iter = f2points.find(fId0);
- TCell pIds0 = pIds0_iter->second;
- TCell &output = output_new_face2ids[new_faceId];
- new_faces_iter++;
- if(new_face2face_iter->second.size() == 1 ){
- TCellArray::const_iterator pIds_iter = f2points.find(fId0);
- TCell pIds = pIds_iter->second;
- output = pIds;
- continue;
- }
- for(;new_faces_iter!=(new_face2face_iter->second).end();new_faces_iter++){
-
- vtkIdType faceId = *new_faces_iter;
- // find how much nodes in face (f2points)
- TCellArray::const_iterator pIds_iter = f2points.find(faceId);
- TCell pIds = pIds_iter->second;
-
- GetSumm(pIds0,pIds,output);
- pIds0 = output;
-
- } // end new_faces_iter
-
- } // new_face2face_iter
+ typedef std::map<double,vtkIdType> TSortedPointIds;
+ TSortedPointIds aSortedPointIds;
+
+ TPointIds::const_iterator anIter = aPointIds.begin();
+ TPointIds::const_iterator anEndIter = aPointIds.end();
+ for(; anIter != anEndIter; anIter++){
+ double aPntCoord[3];
+ vtkIdType aPntId = *anIter;
+ aPoints->GetPoint(aPntId,aPntCoord);
+ double aVector[3] = { aPntCoord[0] - aCenter[0],
+ aPntCoord[1] - aCenter[1],
+ aPntCoord[2] - aCenter[2] };
+ vtkMath::Normalize(aVector);
+
+ double aCross[3];
+ vtkMath::Cross(aVector,aVector0,aCross);
+ double aCr = vtkMath::Dot(aCross,aNormal);
+ bool aGreaterThanPi = aCr < 0;
+ double aCosinus = vtkMath::Dot(aVector,aVector0);
+ double anAngle = 0.0;
+ if(aCosinus >= 1.0){
+ aCosinus = 1.0;
+ } else if (aCosinus <= -1.0){
+ aCosinus = -1.0;
+ anAngle = vtkMath::Pi();
+ } else {
+ anAngle = acos(aCosinus);
+ if(aGreaterThanPi)
+ anAngle = 2*vtkMath::Pi() - anAngle;
+ }
+
+ if(DEBUG_TRIA_EXECUTE) {
+ cout << "\t\t\t vtkMath::Dot(aCross,aNormal)="<<aCr<<endl;
+ cout<<"\t\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}";
+ cout<<"; aGreaterThanPi = "<<aGreaterThanPi<<"; aCosinus = "<<aCosinus<<"; anAngle = "<<anAngle<<"\n";
+ }
+ aSortedPointIds[anAngle] = aPntId;
+ }
+
+ if(!aSortedPointIds.empty()){
+ size_t aNumFacePts = aSortedPointIds.size();
+ ::TConnectivities aConnectivities(aNumFacePts);
+ TSortedPointIds::const_iterator anIter = aSortedPointIds.begin();
+ TSortedPointIds::const_iterator anEndIter = aSortedPointIds.end();
+ if(DEBUG_TRIA_EXECUTE) cout << "Polygon:";
+ for(vtkIdType anId = 0; anIter != anEndIter; anIter++, anId++){
+ vtkIdType aPntId = anIter->second;
+ aConnectivities[anId] = GetConnectivity(aPntId);
+ if(DEBUG_TRIA_EXECUTE) cout << aPntId << ",";
+ }
+ if(DEBUG_TRIA_EXECUTE) cout << endl;
+ aPolygons.push_back(::TPolygon(aConnectivities,aCenter,aNormal));
+ }
}
-
- if(MYDEBUG) {
- cout << "++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl;
- cout << "++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl;
- cout << "++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl;
+ }
+ }
+ if(aPolygons.empty())
+ return true;
+
+ // To check, whether the polygons give a convex polyhedron or not
+ if(theIsCheckConvex){
+ size_t aNbPolygons = aPolygons.size();
+ for (size_t aPolygonId = 0; aPolygonId < aNbPolygons; aPolygonId++) {
+ ::TPolygon& aPolygon = aPolygons[aPolygonId];
+ double* aNormal = aPolygon.myNormal;
+ double* anOrigin = aPolygon.myOrigin;
+ if(DEBUG_TRIA_EXECUTE) {
+ cout<<"\taPolygonId = "<<aPolygonId<<"\n";
+ cout<<"\t\taNormal = {"<<aNormal[0]<<", "<<aNormal[1]<<", "<<aNormal[2]<<"}";
+ cout<<"; anOrigin = {"<<anOrigin[0]<<", "<<anOrigin[1]<<", "<<anOrigin[2]<<"}\n";
}
- outputCellArray = output_new_face2ids;//f2points;
+ for(vtkIdType aPntId = 0; aPntId < aNumPts; aPntId++){
+ double aPntCoord[3];
+ vtkIdType anId = GetPointId(aPntId);
+ aPoints->GetPoint(anId,aPntCoord);
+ double aDist = vtkPlane::Evaluate(aNormal,anOrigin,aPntCoord);
+ if(DEBUG_TRIA_EXECUTE) cout<<"\t\taPntId = "<<anId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}; aDist = "<<aDist<<"\n";
+ if(aDist < -aDistEps)
+ return false;
+ }
+ }
+ }
+
+
+ // To pass resulting set of the polygons to the output
+ {
+ size_t aNbPolygons = aPolygons.size();
+ for (size_t aPolygonId = 0; aPolygonId < aNbPolygons; aPolygonId++) {
+ ::TPolygon& aPolygon = aPolygons[aPolygonId];
+ if(DEBUG_TRIA_EXECUTE) cout << "PoilygonId="<<aPolygonId<<" | ";
+ TConnectivities& aConnectivities = aPolygon.myConnectivities;
+ if(DEBUG_TRIA_EXECUTE) {
+ for(size_t i=0;i<aConnectivities.size();i++)
+ cout << aConnectivities[i] << ",";
+ cout << endl;
+ }
+ int aNbPoints = (int)aConnectivities.size();
+ vtkIdType aNewCellId = theOutput->InsertNextCell(VTK_POLYGON,aNbPoints,&aConnectivities[0]);
+ if(theStoreMapping)
+ VTKViewer_GeometryFilter::InsertId( theCellId, VTK_POLYGON, theVTK2ObjIds, theDimension2VTK2ObjIds );
+ theOutputCD->CopyData(thInputCD,theCellId,aNewCellId);
+ }
+ }
+
+ if(DEBUG_TRIA_EXECUTE) cout<<"\tTriangulator - Ok\n";
+
+ return true;
+}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_OrderedTriangulator
+::VTKViewer_OrderedTriangulator():
+ myTriangulator(vtkOrderedTriangulator::New()),
+ myBoundaryTris(vtkCellArray::New()),
+ myTriangle(vtkTriangle::New())
+{
+ myBoundaryTris->Allocate(VTK_CELL_SIZE);
+ myTriangulator->PreSortedOff();
+}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_OrderedTriangulator
+::~VTKViewer_OrderedTriangulator()
+{
+ myTriangle->Delete();
+ myBoundaryTris->Delete();
+ myTriangulator->Delete();
+}
+
+
+//----------------------------------------------------------------------------
+vtkPoints*
+VTKViewer_OrderedTriangulator
+::InitPoints(vtkUnstructuredGrid *theInput,
+ vtkIdType theCellId)
+{
+ myBoundaryTris->Reset();
+
+ vtkPoints* aPoints = VTKViewer_Triangulator::InitPoints(theInput, theCellId);
+ vtkIdType aNumPts = myPoints->GetNumberOfPoints();
+ if ( aNumPts > 0 ) {
+ myTriangulator->InitTriangulation(0.0, 1.0, 0.0, 1.0, 0.0, 1.0, aNumPts);
+
+ double aBounds[6];
+ myPoints->GetBounds(aBounds);
+ double xSize, ySize, zSize;
+ xSize = aBounds[1] - aBounds[0];
+ ySize = aBounds[3] - aBounds[2];
+ zSize = aBounds[5] - aBounds[4];
+ double anAbsoluteCoord[3];
+ double aParamentrucCoord[3];
+ for (auto aPntId = 0; aPntId < aNumPts; aPntId++) {
+ myPoints->GetPoint(aPntId, anAbsoluteCoord);
+ aParamentrucCoord[0] = xSize==0. ? 0. : ((anAbsoluteCoord[0] - aBounds[0]) / xSize);
+ aParamentrucCoord[1] = ySize==0. ? 0. : ((anAbsoluteCoord[1] - aBounds[2]) / ySize);
+ aParamentrucCoord[2] = zSize==0. ? 0. : ((anAbsoluteCoord[2] - aBounds[4]) / zSize);
+ myTriangulator->InsertPoint(aPntId, anAbsoluteCoord, aParamentrucCoord, 0);
}
- } else {
- // not implemented
+
+ myTriangulator->Triangulate();
+ myTriangulator->AddTriangles(myBoundaryTris);
}
+
+ return aPoints;
+}
+
+
+//----------------------------------------------------------------------------
+vtkIdType
+VTKViewer_OrderedTriangulator
+::GetNumFaces()
+{
+ return myBoundaryTris->GetNumberOfCells();
+}
+
+
+//----------------------------------------------------------------------------
+vtkCell*
+VTKViewer_OrderedTriangulator
+::GetFace(vtkIdType theFaceId)
+{
+ vtkIdType aNumCells = myBoundaryTris->GetNumberOfCells();
+ if ( theFaceId < 0 || theFaceId >= aNumCells )
+ return NULL;
+
+ vtkIdType *aCells = myBoundaryTris->GetData()->GetPointer(0);
+
+ // Each triangle has three points plus number of points
+ vtkIdType *aCellPtr = aCells + 4*theFaceId;
+
+ myTriangle->PointIds->SetId(0, aCellPtr[1]);
+ myTriangle->Points->SetPoint(0, myPoints->GetPoint(aCellPtr[1]));
+
+ myTriangle->PointIds->SetId(1, aCellPtr[2]);
+ myTriangle->Points->SetPoint(1, myPoints->GetPoint(aCellPtr[2]));
+
+ myTriangle->PointIds->SetId(2, aCellPtr[3]);
+ myTriangle->Points->SetPoint(2, myPoints->GetPoint(aCellPtr[3]));
+
+ return myTriangle;
}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_DelaunayTriangulator
+::VTKViewer_DelaunayTriangulator():
+ myUnstructuredGrid(vtkUnstructuredGrid::New()),
+ myGeometryFilter(vtkGeometryFilter::New()),
+ myDelaunay3D(vtkDelaunay3D::New()),
+ myPolyData(NULL)
+{
+ myUnstructuredGrid->Initialize();
+ myUnstructuredGrid->Allocate();
+ myUnstructuredGrid->SetPoints(myPoints);
+
+ myDelaunay3D->SetInputData(myUnstructuredGrid);
+ myGeometryFilter->SetInputConnection(myDelaunay3D->GetOutputPort());
+ myPolyData = myGeometryFilter->GetOutput();
+}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_DelaunayTriangulator
+::~VTKViewer_DelaunayTriangulator()
+{
+ myUnstructuredGrid->Delete();
+ myGeometryFilter->Delete();
+ myDelaunay3D->Delete();
+}
+
+
+//----------------------------------------------------------------------------
+vtkPoints*
+VTKViewer_DelaunayTriangulator
+::InitPoints(vtkUnstructuredGrid *theInput,
+ vtkIdType theCellId)
+{
+ vtkPoints* aPoints = VTKViewer_Triangulator::InitPoints(theInput, theCellId);
+
+ myPoints->Modified();
+ myUnstructuredGrid->Modified();
+ myGeometryFilter->Update();
+
+ return aPoints;
+}
+
+
+//----------------------------------------------------------------------------
+vtkIdType
+VTKViewer_DelaunayTriangulator
+::GetNumFaces()
+{
+ return myPolyData->GetNumberOfCells();
+}
+
+
+//----------------------------------------------------------------------------
+vtkCell*
+VTKViewer_DelaunayTriangulator
+::GetFace(vtkIdType theFaceId)
+{
+ return myPolyData->GetCell(theFaceId);
}
