From: apo <apo@opencascade.com>
Date: Mon, 20 Mar 2006 09:18:21 +0000 (+0000)
Subject: Fix for
X-Git-Tag: OCC_3_2_0a2_TC~11
X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=bd1675247d4845dfcfbd108c651afc6b10bc37d4;p=modules%2Fgui.git

Fix for
Bug IPAL11800
 	Fatal SIGSEGV on create presentation on cells for external file from PAL11736
Bug IPAL11667
  	SIGSEGV after trying to display mesh from "polyedres.med" in Mesh and Post-Pro modules.
---

diff --git a/src/VTKViewer/VTKViewer_ConvexTool.cxx b/src/VTKViewer/VTKViewer_ConvexTool.cxx
index 679e48428..464b27ce0 100644
--- a/src/VTKViewer/VTKViewer_ConvexTool.cxx
+++ b/src/VTKViewer/VTKViewer_ConvexTool.cxx
@@ -29,862 +29,552 @@
 
 #include "VTKViewer_ConvexTool.h"
 
-#include <vtkUnstructuredGrid.h>
-#include <vtkTriangle.h>
-#include <vtkConvexPointSet.h>
-#include <vtkUnstructuredGridWriter.h>
-#include <vtkMath.h>
-#include <vtkSmartPointer.h>
-
 #include <set>
-#include <iterator>
-#include <algorithm>
-#include <math.h>
-
-typedef vtkUnstructuredGrid TInput;
+#include <map>
 
-typedef std::set<vtkIdType> TUIDS; // unique ids 
-typedef std::map<vtkIdType,TUIDS> TPTOIDS; // id points -> unique ids
-
-namespace CONVEX_TOOL
-{
-  // all pairs
-  typedef std::pair<vtkIdType,vtkIdType> TPair;
-  typedef std::set<TPair> TSet;
-  
-  void
-  WriteToFile(vtkUnstructuredGrid* theDataSet, const std::string& theFileName)
-  {
-    vtkUnstructuredGridWriter* aWriter = vtkUnstructuredGridWriter::New();
-    aWriter->SetFileName(theFileName.c_str());
-    aWriter->SetInput(theDataSet);
-    aWriter->Write();
-    aWriter->Delete();
-  }
+#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>
 
-static float FACE_ANGLE_TOLERANCE=1.5;
-#define EPS 1.0e-38
-#define EPS_T 1.0e-3
-
-#ifdef _DEBUG_
-  static int MYDEBUG = 0;
-  static int MYDEBUG_REMOVE = 0;
-#else
-  static int MYDEBUG = 0;
-  static int MYDEBUG_REMOVE = 0;
-#endif
-
-/*! \fn static void GetCenter(TInput* theGrid,TCell theptIds,float *center)
- * \brief Calculation of geometry center.
- * \param theGrid - TInput cell.
- * \param theptIds - point ids.
- * \retval center - output array[3] with coordinates of geometry center.
- */
-static void GetCenter(vtkPoints* thePoints,float center[3])
+//----------------------------------------------------------------------------
+namespace
 {
-  float p[3];
-  center[0] = center[1] = center[2] = 0.0;
+  typedef std::vector<vtkIdType> TConnectivities;
 
-  int numPts = thePoints->GetNumberOfPoints();
-  if (numPts == 0) return;
-
-  // get the center of the cell
-  for (int i = 0; i < numPts; i++)
+  struct TPolygon
   {
-    thePoints->GetPoint(i, p);
-    for (int j = 0; j < 3; j++)
+    TConnectivities myConnectivities;
+    float myOrigin[3];
+    float myNormal[3];
+    TPolygon(const TConnectivities& theConnectivities,
+	     float theOrigin[3],
+	     float theNormal[3]):
+      myConnectivities(theConnectivities)
     {
-      center[j] += p[j];
+      myOrigin[0] = theOrigin[0];
+      myOrigin[1] = theOrigin[1];
+      myOrigin[2] = theOrigin[2];
+
+      myNormal[0] = theNormal[0];
+      myNormal[1] = theNormal[1];
+      myNormal[2] = theNormal[2];
     }
-  }
-  for (int j = 0; j < 3; j++)
-  {
-    center[j] /= numPts;
-  }
+  };
+
+  typedef std::vector<TPolygon> TPolygons;
 }
 
-/*! \fn static void ReverseIds(TCell &theIds)
- * \brief Reverse ids.
- * \param theIds - points ids.
- * \retval theIds - example input:(1,2,3,4) -> output:(4,3,2,1)
- */
-static void ReverseIds(TCell &theIds)
+//----------------------------------------------------------------------------
+VTKViewer_Triangulator
+::VTKViewer_Triangulator():
+  myInput(NULL),
+  myCellId(-1),
+  myShowInside(-1),
+  myAllVisible(-1),
+  myCellsVisibility(NULL),
+  myCellIds(vtkIdList::New())
+{}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_Triangulator
+::~VTKViewer_Triangulator()
 {
-  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;
-  }
+  myCellIds->Delete();
 }
 
-/*! \fn void GetFriends(const TPTOIDS p2faces,const TCellArray f2points,TPTOIDS& face2face_output)
- * \brief Caclulation of connected faces (faceId -> (faceId1,faceId2, ...))
- * \param p2faces - point to faces ids map.
- * \param f2points - faces to points ids map.
- * \retval face2face_output - faces to faces ids map.
- */
-void GetFriends(const TPTOIDS p2faces,const TCellArray f2points,TPTOIDS& face2face_output)
-{
-  TCellArray::const_iterator f2pIter = f2points.begin();
 
-  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++ ){
+//----------------------------------------------------------------------------
+bool 
+VTKViewer_Triangulator
+::Execute(vtkUnstructuredGrid *theInput,
+	  vtkCellData* thInputCD,
+	  vtkIdType theCellId,
+	  int theShowInside,
+	  int theAllVisible,
+	  const char* theCellsVisibility,
+	  vtkPolyData *theOutput,
+	  vtkCellData* theOutputCD,
+	  int theStoreMapping,
+	  std::vector<vtkIdType>& theVTK2ObjIds,
+	  bool theIsCheckConvex)
+{
+  myInput = theInput;
+  myCellId = theCellId;
+  myShowInside = theShowInside;
+  myAllVisible = theAllVisible;
+  myCellsVisibility = theCellsVisibility;
 
-      id2 = face_points[i];
+  vtkPoints *aPoints = InitPoints();
+  vtkIdType aNumPts = GetNbOfPoints();
+  //cout<<"Triangulator - aNumPts = "<<aNumPts<<"\n";
 
-      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(aNumPts == 0)
+    return true;
 
-    if(MYDEBUG){
-      cout << "fId[" << faceId <<"]: ";
-      std::copy(output_faces.begin(), output_faces.end(), std::ostream_iterator<vtkIdType>(cout, " "));
-      cout << endl;
+  // To calculate the bary center of the cell
+  float aCellCenter[3] = {0.0, 0.0, 0.0};
+  {
+    float aPntCoord[3];
+    for (int aPntId = 0; aPntId < aNumPts; aPntId++) {
+      aPoints->GetPoint(GetPointId(aPntId),aPntCoord);
+      //cout<<"\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}\n";
+      aCellCenter[0] += aPntCoord[0];
+      aCellCenter[1] += aPntCoord[1];
+      aCellCenter[2] += aPntCoord[2];
     }
-    
-    face2face_output[faceId] = output_faces;
+    aCellCenter[0] /= aNumPts;
+    aCellCenter[1] /= aNumPts;
+    aCellCenter[2] /= aNumPts;
   }
-}
 
-/*! \fn bool IsConnectedFacesOnOnePlane( TInput* theGrid,vtkIdType theFId1, vtkIdType theFId2,TUIDS FpIds1, TUIDS FpIds2 )
- * \brief Check is connected faces on one plane.
- * \param theGrid - TInput
- * \param theFId1 - id of first face
- * \param theFId2 - id of second face
- * \param FpIds1  - first face points ids.
- * \param FpIds2  - second face points ids.
- * \return TRUE if two faces on one plane, else FALSE.
- */
-bool IsConnectedFacesOnOnePlane( TInput* theGrid,
-                                 vtkIdType theFId1, vtkIdType theFId2,
-				 TUIDS FpIds1, TUIDS FpIds2 )
-{
-  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);
+  float aCellLength = GetCellLength();
+  int aNumFaces = GetNumFaces();
 
-    float b1 = vtkMath::Norm(vec_b1);
+  static float EPS = 1.0E-2;
+  float aDistEps = aCellLength * EPS;
+  //cout<<"\taCellLength = "<<aCellLength<<"; aDistEps = "<<aDistEps<<"\n";
 
-    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] << ")   ";
-      }
-      if(status) 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
-  }
+  // 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));
   
-  return status;
-}
-
-/*! \fn void GetAllFacesOnOnePlane( TPTOIDS theFaces, vtkIdType faceId,TUIDS &new_faces, TCell &new_faces_v2 )
- * \brief Calculate faces which on one plane.
- * \param theFaces - 
- * \param faceId - 
- * \param new_faces - 
- * \param new_faces_v2 - 
- */
-void GetAllFacesOnOnePlane( TPTOIDS theFaces, vtkIdType faceId, 
-                            TUIDS &new_faces, TCell &new_faces_v2 )
-{
-  if (new_faces.find(faceId) != new_faces.end()) return;
+  // To initialize set of points by face that belong to the cell and backward
+  typedef std::set<vtkIdType> TFace2Visibility;
+  TFace2Visibility aFace2Visibility;
   
-  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
-    }
-  }
-  return;
-}
+  typedef std::set<TPointIds> TFace2PointIds;
+  TFace2PointIds aFace2PointIds;
 
-/*! \fn void GetSumm(TCell v1,TCell v2,TCell &output)
- * \brief Gluing two faces (gluing points ids)
- * \param v1 - first face
- * \param v2 - second face
- * \param output - output face.
- */
-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;
-  
-  TCell::iterator v1_iter = v1.begin();
-  
-  for(;v1_iter!=v1.end();v1_iter++){
-    
-    vtkIdType curr_id = *v1_iter;
+  for (int aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
+    vtkCell* aFace = GetFace(aFaceId);
     
-    output.push_back(curr_id);
-    
-    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(theCellId, aFace, myCellIds);
+    if((!myAllVisible && !myCellsVisibility[myCellIds->GetId(0)]) || 
+       myCellIds->GetNumberOfIds() <= 0 || myShowInside)
+    {
+      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;
-}
-
-static void GetAndRemoveIdsOnOneLine(vtkPoints* points,
-				     TUIDS input_points_ids,
-				     TUIDS input_two_points_ids,
-				     TUIDS& out_two_points_ids,
-				     TUIDS& removed_points_ids){
-  if (MYDEBUG_REMOVE) cout << EPS <<endl;
-  float P[3][3];
-  vtkIdType current_points_ids[2];
-  if(MYDEBUG_REMOVE) 
-    if(input_two_points_ids.size()!=2) cout << "Error. Must be two ids in variable input_two_points_ids="<<input_two_points_ids.size()<<endl;
-  TUIDS::const_iterator aInPointsIter = input_two_points_ids.begin();
-  for(int i=0;i<2 && aInPointsIter!=input_two_points_ids.end();aInPointsIter++,i++){
-    current_points_ids[i] = *aInPointsIter;
-    if (MYDEBUG_REMOVE) cout << " " << *aInPointsIter;
-  }
-  if (MYDEBUG_REMOVE) cout << endl;
-  bool iscurrent_points_changed = false;
-  points->GetPoint(current_points_ids[0],P[0]);
-  points->GetPoint(current_points_ids[1],P[1]);
-  TUIDS::iterator aPointsIter = input_points_ids.begin();
-  for(;aPointsIter!=input_points_ids.end();aPointsIter++){
-    if(iscurrent_points_changed){
-      points->GetPoint(current_points_ids[0],P[0]);
-      points->GetPoint(current_points_ids[1],P[1]);
-      iscurrent_points_changed = false;
-      if (MYDEBUG_REMOVE) 
-	cout << " " << current_points_ids[0] << " " << current_points_ids[1] << endl;
-    }
-    // check: is point on line input_two_points_ids
-    points->GetPoint(*aPointsIter,P[2]);
-    if (MYDEBUG_REMOVE) {
-      cout << "\t" << current_points_ids[0] << ":"<<P[0][0]<<","<<P[0][1]<<","<<P[0][2]<<endl;
-      cout << "\t" << current_points_ids[1] << ":"<<P[1][0]<<","<<P[1][1]<<","<<P[1][2]<<endl;
-      cout << "\t" << *aPointsIter << ":"<<P[2][0]<<","<<P[2][1]<<","<<P[2][2]<<endl;
-    }
-  
-    // x-x1=(x2-x1)*t -> coeff[0][0] = (x-x1), coeff[0][1] = x2-x1
-    // y-y1=(y2-y1)*t -> coeff[1][0] = (y-y1), coeff[1][1] = y2-y1
-    // z-z1=(z2-z1)*t -> coeff[2][0] = (z-z1), coeff[2][1] = z2-z1
-    float coeff[3][2];
-    for(int i=0;i<3;i++){
-      coeff[i][0] = P[2][i]-P[0][i];
-      coeff[i][1] = P[1][i]-P[0][i];
-    }
-    bool isok_coord[3];
-    bool isok = true;
-    float t[3];
-    for(int i=0;i<3;i++){
-      isok_coord[i] = false;
-      if( fabs(coeff[i][0]) <= EPS && fabs(coeff[i][1]) <= EPS) {
-	isok_coord[i] = true;
-	continue;
-      }
-      if( fabs(coeff[i][1]) <= EPS && fabs(coeff[i][0]) > EPS) {isok = false;t[i]=1.0/EPS;break;}
-      t[i] = (coeff[i][0])/(coeff[i][1]);
-    }
-    for(int i=0;i<3;i++)
-      if (MYDEBUG_REMOVE) 
-	cout << __LINE__ << "          " 
-	     << coeff[i][0] << ","<<coeff[i][1]
-	     <<"   t="<<t[i]<<" isok_coord="<<isok_coord[i]<<endl;
-    if(!isok) continue;
-
-    if (!isok_coord[0] && !isok_coord[1]){
-      if (fabs(t[0]-t[1]) <= EPS_T) isok = true;
-      else isok = false;
-    }
-    if (MYDEBUG_REMOVE) cout << __LINE__ << "          1000 " << isok << endl;
-    if(!isok) continue;
-    if (!isok_coord[1] && !isok_coord[2]){
-      if (fabs(t[1] - t[2]) <= EPS_T) isok = true;
-      else isok = false;
-    }
-    if (MYDEBUG_REMOVE) cout << __LINE__ << "          2000 " << isok << endl;
-    if(!isok) continue;
-    if (!isok_coord[0] && !isok_coord[2]){
-      if (fabs(t[0] - t[2]) <= EPS_T) isok = true;
-      else isok = false;
-    }
-    if (MYDEBUG_REMOVE) cout << __LINE__ << "          3000 " << isok<<endl;
-    if(!isok) continue;
-    
-    float param[3]; // parametric coord for P[0],P[1],P[2] <--->t[0],t[1],t[2]
-    // anilize bounds of line
-    for(int i=0;i<3;i++){
-      for(int j=0;j<3;j++)
-	if(!isok_coord[j]) param[i] = (P[i][j]-P[0][j])/(P[1][j]-P[0][j]);
-    }
-    if (MYDEBUG_REMOVE) cout << "Params: " << param[0] << "," << param[1] << "," << param[2] << endl;
-    vtkIdType imax,imin;
-    float min,max;
-    for(vtkIdType i=0;i<3;i++)
-      if(!isok_coord[i]){
-	min = param[0];imin=0;
-	max = param[0];imax=0;
-	break;
-      }
-    for(vtkIdType i=0;i<3;i++){
-	if(min > param[i]) {min = param[i]; imin=i;}
-	if(max < param[i]) {max = param[i]; imax=i;}
-    }
-    if (MYDEBUG_REMOVE) 
-      cout << "\t min="<<min<<"  max="<<max<<" - "<<"imin="<<imin<<"  imax="<<imax<<endl;
-    // imin - index of left point
-    // imax - index of right point
-    
-    // add id to removed point
-    vtkIdType rem_id,id1,id2;
-    for(vtkIdType i=0;i<3;i++)
-      if (i!=imin && i!=imax) {rem_id = i; break;}
-    
-    if(rem_id == 0) {
-      rem_id = current_points_ids[0];
-      id1=current_points_ids[1];
-      id2=*aPointsIter;
-    }
-    else if (rem_id == 1) {
-      rem_id = current_points_ids[1];
-      id1=current_points_ids[0];
-      id2=*aPointsIter;
-    }
-    else if (rem_id == 2) {
-      rem_id = *aPointsIter;
-      id1=current_points_ids[0];
-      id2=current_points_ids[1];
-    }
-    if (MYDEBUG_REMOVE) 
-      cout << " " << current_points_ids[0] <<","<<current_points_ids[1]<<"---->"<<id1<<","<<id2<<endl;
-    if((current_points_ids[0] == id1 && current_points_ids[1] == id2) ||
-       (current_points_ids[0] == id2 && current_points_ids[1] == id1))
-      {}
-    else {
-      iscurrent_points_changed = true;
-      current_points_ids[0] = id1;
-      current_points_ids[1] = id2;
-    }
-    
-    removed_points_ids.insert(rem_id);
-  }
-  out_two_points_ids.insert(current_points_ids[0]);
-  out_two_points_ids.insert(current_points_ids[1]);
-}
-
-static vtkSmartPointer<vtkConvexPointSet> RemoveAllUnneededPoints(vtkConvexPointSet* convex){
-  vtkSmartPointer<vtkConvexPointSet> out = vtkConvexPointSet::New();
-  
-  TUIDS two_points,input_points,out_two_points_ids,removed_points_ids,loc_removed_points_ids;
-  vtkIdList* aPointIds = convex->GetPointIds();
-  int numIds = aPointIds->GetNumberOfIds();
-  if (numIds<2) return out;
-  TSet good_point_ids;
-  TSet aLists[numIds-2];
-  for(int i=0;i<numIds-2;i++){
-    for(int j=i+1;j<numIds-1;j++){
-      TPair aPair(i,j);
-      aLists[i].insert(aPair);
-    }
-  }
-  for(vtkIdType i=0;i<numIds-2;i++){
-    TUIDS::iterator aRemIter=removed_points_ids.find(i);
-    if(aRemIter!=removed_points_ids.end()) continue;
-    TSet::iterator aPairIter=aLists[i].begin();
-    loc_removed_points_ids.clear();
-    out_two_points_ids.clear();
-    input_points.clear();
-    two_points.clear();
-    for(;aPairIter!=aLists[i].end();aPairIter++){
-      vtkIdType aFirId = (*aPairIter).first;
-      vtkIdType aSecId = (*aPairIter).second;
-      aRemIter=removed_points_ids.find(aSecId);
-      if(aRemIter!=removed_points_ids.end()) continue;
-      TPair aPair1(aFirId,aSecId);
-      TPair aPair2(aSecId,aFirId);
-      TSet::iterator aGoodIter=good_point_ids.find(aPair1);
-      if(aGoodIter!=good_point_ids.end()) continue;
-      aGoodIter=good_point_ids.find(aPair2);
-      if(aGoodIter!=good_point_ids.end()) continue;
-      two_points.insert(aFirId);
-      two_points.insert(aSecId);
-      if (MYDEBUG_REMOVE) 
-	cout << "\nInput: " << aFirId<<":"<<aPointIds->GetId(aFirId) << "," << aSecId <<":"<<aPointIds->GetId(aSecId)<< "  --- ";
-      for(vtkIdType k=aSecId+1;k<numIds;k++) {
-	input_points.insert(k);
-	if (MYDEBUG_REMOVE) cout << k<<":"<<aPointIds->GetId(k) << ",";
-      }
-      if (MYDEBUG_REMOVE) {
-	cout << endl;
-	cout << "\t";
-	for(TUIDS::iterator aDelIter = loc_removed_points_ids.begin();aDelIter!=loc_removed_points_ids.end();aDelIter++) 
-	  cout << *aDelIter<<",";
-	cout << endl;
-      }
-      GetAndRemoveIdsOnOneLine(convex->Points,
-			       input_points,
-			       two_points,
-			       out_two_points_ids,
-			       loc_removed_points_ids);
-      TUIDS::iterator aOutIter = out_two_points_ids.begin();
-      vtkIdType aFirst=*aOutIter;aOutIter++;vtkIdType aSecond=*aOutIter;
-      TPair aPair(aFirst,aSecond);
-      good_point_ids.insert(aPair);
-      if (MYDEBUG_REMOVE){
-	cout << "Output: ";
-	TUIDS::iterator aIter = out_two_points_ids.begin();
-	for(;aIter!=out_two_points_ids.end();aIter++)
-	  cout << *aIter << ",";
-	cout << " --- ";
-      }
-      TUIDS::iterator aDelIter = loc_removed_points_ids.begin();
-      for(;aDelIter!=loc_removed_points_ids.end();aDelIter++){
-	removed_points_ids.insert(*aDelIter);
-	if (MYDEBUG_REMOVE) cout << *aDelIter << ",";
-      }
-      if (MYDEBUG_REMOVE) cout << endl;
-    }
-  }
-  if (MYDEBUG_REMOVE) {
-    cout << "============ Resultat ================" <<endl;
-    cout << "Removed:";
-    for(TUIDS::iterator aIter=removed_points_ids.begin();aIter!=removed_points_ids.end();aIter++)
-      cout << *aIter << ",";
-    cout << endl;
-  }
-  
-  TUIDS result_ids,all_ids;
-  for(vtkIdType i=0;i<numIds;i++) all_ids.insert(i);
-  std::set_difference(all_ids.begin(),
-		      all_ids.end(),
-		      removed_points_ids.begin(),
-		      removed_points_ids.end(),
-		      std::inserter(result_ids,result_ids.begin()));
-
-  out->Points->SetNumberOfPoints(result_ids.size());
-  out->PointIds->SetNumberOfIds(result_ids.size());
-  int aId=0;
-  if(MYDEBUG_REMOVE) cout << "Result out:";
-  for(TUIDS::iterator aIter=result_ids.begin();aIter!=result_ids.end();aIter++,aId++){
-    float P[3];
-    convex->Points->GetPoint(*aIter,P);
-    out->Points->SetPoint(aId,P);
-    out->PointIds->SetId(aId,aPointIds->GetId(*aIter));
-    if (MYDEBUG_REMOVE) cout << *aIter << ":" << aPointIds->GetId(*aIter) << " , ";
-  }
-  if(MYDEBUG_REMOVE) cout << endl;
-  out->Modified();
-  out->Initialize();
-  
-  return out;
-}
 
-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_in = static_cast<vtkConvexPointSet*>(theGrid->GetCell(cellId))){
-      vtkSmartPointer<vtkConvexPointSet> convex = RemoveAllUnneededPoints(convex_in);
-      TCellArray f2points;
-      float convex_center[3]; // convex center point coorinat
-      int aNbFaces = convex->GetNumberOfFaces();
-      int numPts = convex->GetNumberOfPoints();
-      if(MYDEBUG_REMOVE) cout << "aNbFaces="<<aNbFaces<<endl;
-      if(MYDEBUG_REMOVE) cout << "numPts="<<numPts<<endl;
-      TPTOIDS p2faces; // key=pointId , value facesIds set
-
-      GetCenter(convex->Points,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];
-	}
+  ::TPolygons aPolygons;
 
- 	if (vtkMath::Determinant3x3(v_a,v_b,v_convex2face) < 0){
-	  ReverseIds(aIds);
- 	}
+  for (int aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
+    if(aFace2Visibility.find(aFaceId) == aFace2Visibility.end())
+      continue;
 
-	for(i=0;i<(int)aIds.size();i++){
-	  TUIDS &acell = p2faces[aIds[i]];
-	  acell.insert(faceId);
-	}
-	
-	f2points[faceId] = aIds;
+    vtkCell* aFace = GetFace(aFaceId);
 
-      }
-      
-      TPTOIDS face2face;
-      GetFriends(p2faces,f2points,face2face);
+    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();
+    //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
+      float aCoord[3][3];
+      aPoints->GetPoint(aNewPts[0],aCoord[0]);
+      aPoints->GetPoint(aNewPts[1],aCoord[1]);
+      aPoints->GetPoint(aNewPts[2],aCoord[2]);
       
-      TPTOIDS face2points;
+      // To calculate plane normal
+      float aVector01[3] = { aCoord[1][0] - aCoord[0][0],
+			     aCoord[1][1] - aCoord[0][1],
+			     aCoord[1][2] - aCoord[0][2] };
       
-      // 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
-	
+      float aVector02[3] = { aCoord[2][0] - aCoord[0][0],
+			     aCoord[2][1] - aCoord[0][1],
+			     aCoord[2][2] - aCoord[0][2] };
       
-      TPTOIDS new_face2faces; // which connected and in one plane
+      // To calculate the normal for the triangle
+      float aNormal[3];
+      vtkMath::Cross(aVector02,aVector01,aNormal);
       
-      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
+      vtkMath::Normalize(aNormal);
       
-      if(MYDEBUG)
+      // To calculate what points belong to the plane
+      // To calculate bounds of the point set
+      float aCenter[3] = {0.0, 0.0, 0.0};
       {
-	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;
+	TPointIds::const_iterator anIter = anInitialPointIds.begin();
+	TPointIds::const_iterator anEndIter = anInitialPointIds.end();
+	for(; anIter != anEndIter; anIter++){
+	  float aPntCoord[3];
+	  vtkIdType aPntId = *anIter;
+	  aPoints->GetPoint(aPntId,aPntCoord);
+	  float aDist = vtkPlane::DistanceToPlane(aPntCoord,aNormal,aCoord[0]);
+	  //cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}; aDist = "<<aDist<<"\n";
+	  if(fabs(aDist) < aDistEps){
+	    aPointIds.insert(aPntId);
+	    aCenter[0] += aPntCoord[0];
+	    aCenter[1] += aPntCoord[1];
+	    aCenter[2] += aPntCoord[2];
+	  }
 	}
+	int aNbPoints = aPointIds.size();
+	aCenter[0] /= aNbPoints;
+	aCenter[1] /= aNbPoints;
+	aCenter[2] /= aNbPoints;
+      }
+      
+      //To sinchronize orientation of the cell and its face
+      float aVectorC[3] = { aCenter[0] - aCellCenter[0],
+			    aCenter[1] - aCellCenter[1],
+			    aCenter[2] - aCellCenter[2] };
+      vtkMath::Normalize(aVectorC);
+      
+      float aDot = vtkMath::Dot(aNormal,aVectorC);
+      //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 output_newid2face;
-      TCellArray output_newid2face_v2;
+      // To calculate the primary direction for point set
+      float aVector0[3] = { aCoord[0][0] - aCenter[0],
+			    aCoord[0][1] - aCenter[1],
+			    aCoord[0][2] - aCenter[2] };
+      vtkMath::Normalize(aVector0);
+      
+      //cout<<"\t\taCenter = {"<<aCenter[0]<<", "<<aCenter[1]<<", "<<aCenter[2]<<"}";
+      //cout<<"; aVector0 = {"<<aVector0[0]<<", "<<aVector0[1]<<", "<<aVector0[2]<<"}\n";
+      
+      // 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;
+	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()));
 	  
-	  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);
-	      }
+	  if(anIntersection == anIds){
+	    aRemoveFace2PointIds.insert(anIds);
 	  }
-	  k++;
 	}
       }
       
-      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++){
+	typedef std::map<float,vtkIdType> TSortedPointIds;
+	TSortedPointIds aSortedPointIds;
+	
+	TPointIds::const_iterator anIter = aPointIds.begin();
+	TPointIds::const_iterator anEndIter = aPointIds.end();
+	for(; anIter != anEndIter; anIter++){
+	  float aPntCoord[3];
+	  vtkIdType aPntId = *anIter;
+	  aPoints->GetPoint(aPntId,aPntCoord);
+	  float aVector[3] = { aPntCoord[0] - aCenter[0],
+			       aPntCoord[1] - aCenter[1],
+			       aPntCoord[2] - aCenter[2] };
+	  vtkMath::Normalize(aVector);
 	  
-	  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;
+	  float aCross[3];
+	  vtkMath::Cross(aVector,aVector0,aCross);
+	  bool aGreaterThanPi = vtkMath::Dot(aCross,aNormal) < 0;
+	  float aCosinus = vtkMath::Dot(aVector,aVector0);
+	  if(aCosinus > 1.0)
+	    aCosinus = 1.0;
+	  if(aCosinus < -1.0)
+	    aCosinus = -1.0;
+	  static float a2Pi = 2.0 * vtkMath::Pi();
+	  float anAngle = acos(aCosinus);
+	  //cout<<"\t\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}";
+	  //cout<<"; aGreaterThanPi = "<<aGreaterThanPi<<"; aCosinus = "<<aCosinus<<"; anAngle = "<<anAngle<<"\n";
+	  if(aGreaterThanPi){
+	    anAngle = a2Pi - anAngle;
+	    //cout<<"\t\t\t\tanAngle = "<<anAngle<<"\n";
 	  }
-	  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;
+	  aSortedPointIds[anAngle] = aPntId;
+	}
 
-	  } // end new_faces_iter
-	  
-	} // new_face2face_iter
+	if(!aSortedPointIds.empty()){
+	  int aNumFacePts = aSortedPointIds.size();
+	  ::TConnectivities aConnectivities(aNumFacePts);
+	  TSortedPointIds::const_iterator anIter = aSortedPointIds.begin();
+	  TSortedPointIds::const_iterator anEndIter = aSortedPointIds.end();
+	  for(vtkIdType anId = 0; anIter != anEndIter; anIter++, anId++){
+	    vtkIdType aPntId = anIter->second;
+	    aConnectivities[anId] = GetConnectivity(aPntId);
+	  }
+	  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){
+    int aNbPolygons = aPolygons.size();
+    for (int aPolygonId = 0; aPolygonId < aNbPolygons; aPolygonId++) {
+      ::TPolygon& aPolygon = aPolygons[aPolygonId];
+      float* aNormal = aPolygon.myNormal;
+      float* anOrigin = aPolygon.myOrigin;
+      //cout<<"\taPolygonId = "<<aPolygonId<<"\n";
+      //cout<<"\t\taNormal = {"<<aNormal[0]<<", "<<aNormal[1]<<", "<<aNormal[2]<<"}";
+      //cout<<"; anOrigin = {"<<anOrigin[0]<<", "<<anOrigin[1]<<", "<<anOrigin[2]<<"}\n";
+      for(vtkIdType aPntId = 0; aPntId < aNumPts; aPntId++){
+	float aPntCoord[3];
+	vtkIdType anId = GetPointId(aPntId);
+	aPoints->GetPoint(anId,aPntCoord);
+	float aDist = vtkPlane::Evaluate(aNormal,anOrigin,aPntCoord);
+	//cout<<"\t\taPntId = "<<anId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}; aDist = "<<aDist<<"\n";
+	if(aDist < -aDistEps)
+	  return false;
       }
-      outputCellArray = output_new_face2ids;//f2points;
     }
-  } else {
-    // not implemented
   }
+
+
+  // To pass resulting set of the polygons to the output
+  {
+    int aNbPolygons = aPolygons.size();
+    for (int aPolygonId = 0; aPolygonId < aNbPolygons; aPolygonId++) {
+      ::TPolygon& aPolygon = aPolygons[aPolygonId];
+      TConnectivities& aConnectivities = aPolygon.myConnectivities;
+      int aNbPoints = aConnectivities.size();
+      vtkIdType aNewCellId = theOutput->InsertNextCell(VTK_POLYGON,aNbPoints,&aConnectivities[0]);
+      if(theStoreMapping)
+	theVTK2ObjIds.push_back(theCellId);
+      theOutputCD->CopyData(thInputCD,theCellId,aNewCellId);
+    }
+  }
+
+  //cout<<"\tTriangulator - Ok\n";
+  return true;
+}
+
+//----------------------------------------------------------------------------
+VTKViewer_OrderedTriangulator
+::VTKViewer_OrderedTriangulator():
+  myCell(vtkGenericCell::New())
+{}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_OrderedTriangulator
+::~VTKViewer_OrderedTriangulator()
+{
+  myCell->Delete();
+}
+
+//----------------------------------------------------------------------------
+vtkPoints*
+VTKViewer_OrderedTriangulator
+::InitPoints()
+{
+  myInput->GetCell(myCellId,myCell);
+  return myInput->GetPoints();
+}
+
+vtkIdType 
+VTKViewer_OrderedTriangulator
+::GetNbOfPoints()
+{
+  return myCell->GetNumberOfPoints();
+}
+
+vtkIdType 
+VTKViewer_OrderedTriangulator
+::GetPointId(vtkIdType thePointId)
+{
+  return myCell->GetPointId(thePointId);
 }
+
+float 
+VTKViewer_OrderedTriangulator
+::GetCellLength()
+{
+  return sqrt(myCell->GetLength2());
+}
+
+vtkIdType 
+VTKViewer_OrderedTriangulator
+::GetNumFaces()
+{
+  return myCell->GetNumberOfFaces();
+}
+
+vtkCell*
+VTKViewer_OrderedTriangulator
+::GetFace(vtkIdType theFaceId)
+{
+  return myCell->GetFace(theFaceId);
+}
+
+void 
+VTKViewer_OrderedTriangulator
+::GetCellNeighbors(vtkIdType theCellId,
+		   vtkCell* theFace,
+		   vtkIdList* theCellIds)
+{
+  vtkIdList *anIdList = theFace->PointIds;  
+  myInput->GetCellNeighbors(theCellId, anIdList, theCellIds);
+}
+
+vtkIdType 
+VTKViewer_OrderedTriangulator
+::GetConnectivity(vtkIdType thePntId)
+{
+  return thePntId;
+}
+
+//----------------------------------------------------------------------------
+VTKViewer_DelaunayTriangulator
+::VTKViewer_DelaunayTriangulator():
+  myUnstructuredGrid(vtkUnstructuredGrid::New()),
+  myGeometryFilter(vtkGeometryFilter::New()),
+  myDelaunay3D(vtkDelaunay3D::New()),
+  myFaceIds(vtkIdList::New()),
+  myPoints(vtkPoints::New()),
+  myPolyData(NULL),
+  myPointIds(NULL)
+{
+  myDelaunay3D->SetInput(myUnstructuredGrid);
+  myGeometryFilter->SetInput(myDelaunay3D->GetOutput());
+}
+
+
+//----------------------------------------------------------------------------
+VTKViewer_DelaunayTriangulator
+::~VTKViewer_DelaunayTriangulator()
+{
+  myUnstructuredGrid->Delete();
+  myGeometryFilter->Delete();
+  myDelaunay3D->Delete();
+  myFaceIds->Delete();
+  myPoints->Delete();
+}
+
+
+//----------------------------------------------------------------------------
+vtkPoints* 
+VTKViewer_DelaunayTriangulator
+::InitPoints()
+{
+  myUnstructuredGrid->Initialize();
+  myUnstructuredGrid->Allocate();
+  myUnstructuredGrid->SetPoints(myPoints);
+
+  vtkIdType aNumPts;
+  myInput->GetCellPoints(myCellId,aNumPts,myPointIds); 
+  {
+    float aPntCoord[3];
+    myPoints->SetNumberOfPoints(aNumPts);
+    vtkPoints *anInputPoints = myInput->GetPoints();
+    for (int aPntId = 0; aPntId < aNumPts; aPntId++) {
+      anInputPoints->GetPoint(myPointIds[aPntId],aPntCoord);
+      myPoints->SetPoint(aPntId,aPntCoord);
+    }
+  }
+
+  myGeometryFilter->Update();
+  myPolyData = myGeometryFilter->GetOutput();
+
+  return myPoints;
+}
+
+vtkIdType 
+VTKViewer_DelaunayTriangulator
+::GetNbOfPoints()
+{
+  return myPoints->GetNumberOfPoints();
+}
+
+vtkIdType 
+VTKViewer_DelaunayTriangulator
+::GetPointId(vtkIdType thePointId)
+{
+  return thePointId;
+}
+
+float 
+VTKViewer_DelaunayTriangulator
+::GetCellLength()
+{
+  return myPolyData->GetLength();
+}
+
+vtkIdType 
+VTKViewer_DelaunayTriangulator
+::GetNumFaces()
+{
+  return myPolyData->GetNumberOfCells();
+}
+
+vtkCell*
+VTKViewer_DelaunayTriangulator
+::GetFace(vtkIdType theFaceId)
+{
+  return myPolyData->GetCell(theFaceId);
+}
+
+void 
+VTKViewer_DelaunayTriangulator
+::GetCellNeighbors(vtkIdType theCellId,
+		   vtkCell* theFace,
+		   vtkIdList* theCellIds)
+{
+  myFaceIds->Reset();
+  vtkIdList *anIdList = theFace->PointIds;  
+  myFaceIds->InsertNextId(myPointIds[anIdList->GetId(0)]);
+  myFaceIds->InsertNextId(myPointIds[anIdList->GetId(1)]);
+  myFaceIds->InsertNextId(myPointIds[anIdList->GetId(2)]);
+
+  myInput->GetCellNeighbors(theCellId, myFaceIds, theCellIds);
+}
+
+
+vtkIdType 
+VTKViewer_DelaunayTriangulator
+::GetConnectivity(vtkIdType thePntId)
+{
+  return myPointIds[thePntId];
 }
diff --git a/src/VTKViewer/VTKViewer_ConvexTool.h b/src/VTKViewer/VTKViewer_ConvexTool.h
index ed303af03..f56a502b4 100644
--- a/src/VTKViewer/VTKViewer_ConvexTool.h
+++ b/src/VTKViewer/VTKViewer_ConvexTool.h
@@ -21,28 +21,186 @@
 #ifndef _VTKViewer_ConvexTool_H
 #define _VTKViewer_ConvexTool_H
 
-#include <vtkUnstructuredGrid.h>
+#include "VTKViewer.h"
+
 #include <vector>
-#include <map>
 
-typedef std::vector<vtkIdType> TCell; // ptsIds
-typedef std::map<vtkIdType,TCell> TCellArray; // CellId, TCell
+#include <vtkSystemIncludes.h>
+
+class vtkUnstructuredGrid;
+class vtkGeometryFilter;
+class vtkGenericCell;
+class vtkDelaunay3D;
+class vtkPolyData;
+class vtkCellData;
+class vtkPoints;
+class vtkIdList;
+class vtkCell;
+
+//----------------------------------------------------------------------------
+class VTKVIEWER_EXPORT VTKViewer_Triangulator
+{
+ public:
+  VTKViewer_Triangulator();
+
+  ~VTKViewer_Triangulator();
+
+  bool 
+  Execute(vtkUnstructuredGrid *theInput,
+	  vtkCellData* thInputCD,
+	  vtkIdType theCellId,
+	  int theShowInside,
+	  int theAllVisible,
+	  const char* theCellsVisibility,
+	  vtkPolyData *theOutput,
+	  vtkCellData* theOutputCD,
+	  int theStoreMapping,
+	  std::vector<vtkIdType>& theVTK2ObjIds,
+	  bool theIsCheckConvex);
+
+ protected:
+  vtkIdList* myCellIds;
+
+  vtkUnstructuredGrid *myInput;
+  vtkIdType myCellId;
+  int myShowInside;
+  int myAllVisible;
+  const char* myCellsVisibility;
+
+  virtual
+  vtkPoints* 
+  InitPoints() = 0;
+
+  virtual
+  vtkIdType 
+  GetNbOfPoints() = 0;
+
+  virtual
+  vtkIdType 
+  GetPointId(vtkIdType thePointId) = 0;
 
-/*! This package \namespace CONVEX_TOOL used for: 
- *  calculation of VTK_POLYGON cell array from VTK_TRIANGLE (triangulation)
- *  of VTK_CONVEX_POINT_SET cell type.
- */
-namespace CONVEX_TOOL
+  virtual
+  float 
+  GetCellLength() = 0;
+
+  virtual
+  vtkIdType 
+  GetNumFaces() = 0;
+
+  virtual
+  vtkCell* 
+  GetFace(vtkIdType theFaceId) = 0;
+
+  virtual
+  void 
+  GetCellNeighbors(vtkIdType theCellId,
+		   vtkCell* theFace,
+		   vtkIdList* theCellIds) = 0;
+
+  virtual
+  vtkIdType 
+  GetConnectivity(vtkIdType thePntId) = 0;
+};
+
+
+//----------------------------------------------------------------------------
+class VTKVIEWER_EXPORT VTKViewer_OrderedTriangulator : public VTKViewer_Triangulator
 {
-  /*! \fn void CONVEX_TOOL::GetPolygonalFaces(vtkUnstructuredGrid* theCell,int cellId,TCellArray &outputCellArray)
-   *  \brief Main function.
-   *  \param theCell - vtkUnstructuredGrid cell pointer
-   *  \param cellId  - id of cell type VTK_CONVEX_POINT_SET
-   *  \retval outputCellArray - output array with new cells types VTK_POLYGON
-   */
-  void
-    WriteToFile(vtkUnstructuredGrid* theDataSet, const std::string& theFileName);
-  void GetPolygonalFaces(vtkUnstructuredGrid* theCell,int cellId,TCellArray &outputCellArray);
-}
+ public:
+
+  VTKViewer_OrderedTriangulator();
+
+  ~VTKViewer_OrderedTriangulator();
+
+ protected:
+  vtkGenericCell *myCell;
+
+  virtual
+  vtkPoints* 
+  InitPoints();
+
+  virtual
+  vtkIdType 
+  GetNbOfPoints();
+
+  vtkIdType 
+  GetPointId(vtkIdType thePointId);
+
+  virtual
+  float 
+  GetCellLength();
+
+  virtual
+  vtkIdType 
+  GetNumFaces();
+
+  virtual
+  vtkCell* 
+  GetFace(vtkIdType theFaceId);
+
+  virtual
+  void 
+  GetCellNeighbors(vtkIdType theCellId,
+		   vtkCell* theFace,
+		   vtkIdList* theCellIds);
+
+  virtual
+  vtkIdType 
+  GetConnectivity(vtkIdType thePntId);
+};
+
+
+//----------------------------------------------------------------------------
+class VTKVIEWER_EXPORT VTKViewer_DelaunayTriangulator : public VTKViewer_Triangulator
+{
+ public:
+
+  VTKViewer_DelaunayTriangulator();
+
+  ~VTKViewer_DelaunayTriangulator();
+
+ protected:
+  vtkUnstructuredGrid* myUnstructuredGrid;
+  vtkGeometryFilter* myGeometryFilter;
+  vtkDelaunay3D* myDelaunay3D;
+  vtkPolyData* myPolyData;
+  vtkIdType *myPointIds;
+  vtkIdList* myFaceIds;
+  vtkPoints* myPoints;
+
+  virtual
+  vtkPoints* 
+  InitPoints();
+
+  virtual
+  vtkIdType 
+  GetNbOfPoints();
+
+  vtkIdType 
+  GetPointId(vtkIdType thePointId);
+
+  virtual
+  float 
+  GetCellLength();
+
+  virtual
+  vtkIdType 
+  GetNumFaces();
+
+  virtual
+  vtkCell* 
+  GetFace(vtkIdType theFaceId);
+
+  virtual
+  void 
+  GetCellNeighbors(vtkIdType theCellId,
+		   vtkCell* theFace,
+		   vtkIdList* theCellIds);
+
+  virtual
+  vtkIdType 
+  GetConnectivity(vtkIdType thePntId);
+};
+
 
 #endif // _VTKViewer_ConvexTool_H
diff --git a/src/VTKViewer/VTKViewer_GeometryFilter.cxx b/src/VTKViewer/VTKViewer_GeometryFilter.cxx
index df096cf7a..83d6093ce 100755
--- a/src/VTKViewer/VTKViewer_GeometryFilter.cxx
+++ b/src/VTKViewer/VTKViewer_GeometryFilter.cxx
@@ -27,6 +27,7 @@
 //  $Header$
 
 #include "VTKViewer_GeometryFilter.h"
+#include "VTKViewer_ConvexTool.h"
 
 #include <vtkSmartPointer.h>
 #include <vtkCellArray.h>
@@ -45,11 +46,6 @@
 #include <vtkVoxel.h>
 #include <vtkWedge.h>
 
-#include <vtkMath.h>
-#include <vtkPlane.h>
-#include <vtkDelaunay3D.h>
-#include <vtkGeometryFilter.h>
-
 #include <algorithm>
 #include <iterator>
 #include <vector>
@@ -62,7 +58,7 @@
   #endif
 #endif
 
-#define USE_ROBUST_TRIANGULATION
+//#define USE_ROBUST_TRIANGULATION
 
 //----------------------------------------------------------------------------
 vtkCxxRevisionMacro(VTKViewer_GeometryFilter, "$Revision$");
@@ -128,21 +124,12 @@ VTKViewer_GeometryFilter
   vtkPolyData *output = this->GetOutput();
   vtkPointData *outputPD = output->GetPointData();
   
-#ifdef USE_ROBUST_TRIANGULATION
-  vtkUnstructuredGrid* anUnstructuredGrid = vtkUnstructuredGrid::New();
-  vtkPoints* aDelaunayPoints = vtkPoints::New();
-
-  vtkDelaunay3D* aDelaunay3D = vtkDelaunay3D::New();
-  aDelaunay3D->SetInput(anUnstructuredGrid);
-  
-  vtkGeometryFilter* aGeometryFilter = vtkGeometryFilter::New();
-  aGeometryFilter->SetInput(aDelaunay3D->GetOutput());
-#endif
+  VTKViewer_OrderedTriangulator anOrderedTriangulator;
+  VTKViewer_DelaunayTriangulator aDelaunayTriangulator;
 
   vtkCellData *outputCD = output->GetCellData();
   vtkGenericCell *cell = vtkGenericCell::New();
 
-
   vtkIdList *cellIds = vtkIdList::New();
   vtkIdList *faceIds = vtkIdList::New();
 
@@ -306,299 +293,30 @@ VTKViewer_GeometryFilter
           break;
 	  
 	case VTK_CONVEX_POINT_SET: {
-	  //cout<<"cellId = "<<cellId<<"\n";
-
-	  vtkIdType aNumPts;
-	  vtkPoints *aPoints;
-#ifdef USE_ROBUST_TRIANGULATION
-	  aPoints = aDelaunayPoints;
-	  anUnstructuredGrid->Initialize();
-	  anUnstructuredGrid->Allocate();
-	  anUnstructuredGrid->SetPoints(aDelaunayPoints);
-
-	  vtkIdType *aPts;
-	  input->GetCellPoints(cellId,aNumPts,aPts); 
-	  {
-	    float aPntCoord[3];
-	    aDelaunayPoints->SetNumberOfPoints(aNumPts);
-	    vtkPoints *anInputPoints = input->GetPoints();
-	    for (int aPntId = 0; aPntId < aNumPts; aPntId++) {
-	      anInputPoints->GetPoint(aPts[aPntId],aPntCoord);
-	      aDelaunayPoints->SetPoint(aPntId,aPntCoord);
-	    }
-	  }
-#else
-	  input->GetCell(cellId,cell);
-	  aPoints = input->GetPoints();
-	  aNumPts = cell->GetNumberOfPoints();
-#endif
-	  // To calculate the bary center of the cell
-	  float aCellCenter[3] = {0.0, 0.0, 0.0};
-	  {
-	    float aPntCoord[3];
-	    for (int aPntId = 0; aPntId < aNumPts; aPntId++) {
-#ifdef USE_ROBUST_TRIANGULATION
-	      aPoints->GetPoint(aPntId,aPntCoord);
-#else
-	      aPoints->GetPoint(cell->GetPointId(aPntId),aPntCoord);
-#endif
-	      //cout<<"\t\taPntId = "<<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;
-	  }
-
-#ifdef USE_ROBUST_TRIANGULATION
-	  aGeometryFilter->Update();
-	  vtkPolyData* aPolyData = aGeometryFilter->GetOutput();
-
-	  float aCellLength = aPolyData->GetLength();
-          int aNumFaces = aPolyData->GetNumberOfCells();
-#else
-	  float aCellLength = sqrt(cell->GetLength2());
-          int aNumFaces = cell->GetNumberOfFaces();
-#endif
-	  
-	  static float EPS = 1.0E-5;
-	  float aDistEps = aCellLength * EPS;
-
-	  // To initialize set of points that belong to the cell
-	  typedef std::set<vtkIdType> TPointIds;
-	  TPointIds anInitialPointIds;
-	  for(vtkIdType aPntId = 0; aPntId < aNumPts; aPntId++){
-#ifdef USE_ROBUST_TRIANGULATION
-	    anInitialPointIds.insert(aPntId);
-#else
-	    anInitialPointIds.insert(cell->GetPointId(aPntId));
-#endif
-	  }
-
-	  // To initialize set of points by face that belong to the cell and backward
-	  typedef std::set<vtkIdType> TFace2Visibility;
-	  TFace2Visibility aFace2Visibility;
-
-	  typedef std::set<TPointIds> TFace2PointIds;
-	  TFace2PointIds aFace2PointIds;
-
-          for (int aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
-#ifdef USE_ROBUST_TRIANGULATION
-            vtkCell* aFace = aPolyData->GetCell(aFaceId);
-#else
-            vtkCell* aFace = cell->GetFace(aFaceId);
-#endif
-	    vtkIdList *anIdList = aFace->PointIds;  
-	    aNewPts[0] = anIdList->GetId(0);
-	    aNewPts[1] = anIdList->GetId(1);
-	    aNewPts[2] = anIdList->GetId(2);
-
-#ifdef USE_ROBUST_TRIANGULATION
-            faceIds->Reset();
-            faceIds->InsertNextId(aPts[aNewPts[0]]);
-            faceIds->InsertNextId(aPts[aNewPts[1]]);
-            faceIds->InsertNextId(aPts[aNewPts[2]]);
-            input->GetCellNeighbors(cellId, faceIds, cellIds);
-#else
-            input->GetCellNeighbors(cellId, anIdList, cellIds);
-#endif
-	    if((!allVisible && !cellVis[cellIds->GetId(0)]) || 
-	       cellIds->GetNumberOfIds() <= 0 ||
-	       myShowInside)
-	    {
-	      TPointIds aPointIds;
-	      aPointIds.insert(aNewPts[0]);
-	      aPointIds.insert(aNewPts[1]);
-	      aPointIds.insert(aNewPts[2]);
-
-	      aFace2PointIds.insert(aPointIds);
-	      aFace2Visibility.insert(aFaceId);
-	    }
-	  }
-
-          for (int aFaceId = 0; aFaceId < aNumFaces; aFaceId++) {
-            if(aFace2Visibility.find(aFaceId) == aFace2Visibility.end())
-	      continue;
-
-#ifdef USE_ROBUST_TRIANGULATION
-            vtkCell* aFace = aPolyData->GetCell(aFaceId);
-#else
-            vtkCell* aFace = cell->GetFace(aFaceId);
-#endif
-	    vtkIdList *anIdList = aFace->PointIds;
-	    aNewPts[0] = anIdList->GetId(0);
-	    aNewPts[1] = anIdList->GetId(1);
-	    aNewPts[2] = 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();
-	    //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
-	      float 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
-	      float aVector01[3] = { aCoord[1][0] - aCoord[0][0],
-				     aCoord[1][1] - aCoord[0][1],
-				     aCoord[1][2] - aCoord[0][2] };
-	      
-	      float aVector02[3] = { aCoord[2][0] - aCoord[0][0],
-				     aCoord[2][1] - aCoord[0][1],
-				     aCoord[2][2] - aCoord[0][2] };
+	  bool anIsOk = anOrderedTriangulator.Execute(input,
+						      cd,
+						      cellId,
+						      myShowInside,
+						      allVisible,
+						      cellVis,
+						      output,
+						      outputCD,
+						      myStoreMapping,
+						      myVTK2ObjIds,
+						      true);
+	  if(!anIsOk)
+	    aDelaunayTriangulator.Execute(input,
+					  cd,
+					  cellId,
+					  myShowInside,
+					  allVisible,
+					  cellVis,
+					  output,
+					  outputCD,
+					  myStoreMapping,
+					  myVTK2ObjIds,
+					  false);
 	      
-	      float aCross21[3];
-	      vtkMath::Cross(aVector02,aVector01,aCross21);
-	      
-	      vtkMath::Normalize(aCross21);
-	      
-	      // To calculate what points belong to the plane
-	      // To calculate bounds of the point set
-	      float aCenter[3] = {0.0, 0.0, 0.0};
-	      {
-		TPointIds::const_iterator anIter = anInitialPointIds.begin();
-		TPointIds::const_iterator anEndIter = anInitialPointIds.end();
-		for(; anIter != anEndIter; anIter++){
-		  float aPntCoord[3];
-		  vtkIdType aPntId = *anIter;
-		  aPoints->GetPoint(aPntId,aPntCoord);
-		  float aDist = vtkPlane::DistanceToPlane(aPntCoord,aCross21,aCoord[0]);
-		  //cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}; aDist = "<<aDist<<"\n";
-		  if(fabs(aDist) < aDistEps){
-		    aPointIds.insert(aPntId);
-		    aCenter[0] += aPntCoord[0];
-		    aCenter[1] += aPntCoord[1];
-		    aCenter[2] += aPntCoord[2];
-		  }
-		}
-		int aNbPoints = aPointIds.size();
-		aCenter[0] /= aNbPoints;
-		aCenter[1] /= aNbPoints;
-		aCenter[2] /= aNbPoints;
-	      }
-	      
-	      // To calculate the primary direction for point set
-	      float aVector0[3] = { aCoord[0][0] - aCenter[0],
-				    aCoord[0][1] - aCenter[1],
-				    aCoord[0][2] - aCenter[2] };
-
-	      //To sinchronize orientation of the cell and its face
-	      float aVectorC[3] = { aCenter[0] - aCellCenter[0],
-				    aCenter[1] - aCellCenter[1],
-				    aCenter[2] - aCellCenter[2] };
-	      vtkMath::Normalize(aVectorC);
-
-	      float aDot = vtkMath::Dot(aCross21,aVectorC);
-	      //cout<<"\t\taCross21 = {"<<aCross21[0]<<", "<<aCross21[1]<<", "<<aCross21[2]<<"}";
-	      //cout<<"; aVectorC = {"<<aVectorC[0]<<", "<<aVectorC[1]<<", "<<aVectorC[2]<<"}\n";
-	      //cout<<"\t\taDot = "<<aDot<<"\n";
-	      if(aDot > 0){
-		aCross21[0] = -aCross21[0];
-		aCross21[1] = -aCross21[1];
-		aCross21[2] = -aCross21[2];
-	      }
-		
-	      vtkMath::Normalize(aVector0);
-	      
-	      //cout<<"\t\taCenter = {"<<aCenter[0]<<", "<<aCenter[1]<<", "<<aCenter[2]<<"}";
-	      //cout<<"; aVector0 = {"<<aVector0[0]<<", "<<aVector0[1]<<", "<<aVector0[2]<<"}\n";
-
-	      // To calculate the set of points by face those that belong to the plane
-	      TFace2PointIds aRemoveFace2PointIds;
-	      {
-		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(anIntersection == anIds){
-		    aRemoveFace2PointIds.insert(anIds);
-		  }
-		}
-	      }
-
-	      // 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);
-		}
-	      }
-
-	      // To sort the planar set of the points accrding to the angle
-	      {
-		typedef std::map<float,vtkIdType> TSortedPointIds;
-		TSortedPointIds aSortedPointIds;
-
-		TPointIds::const_iterator anIter = aPointIds.begin();
-		TPointIds::const_iterator anEndIter = aPointIds.end();
-		for(; anIter != anEndIter; anIter++){
-		  float aPntCoord[3];
-		  vtkIdType aPntId = *anIter;
-		  aPoints->GetPoint(aPntId,aPntCoord);
-		  float aVector[3] = { aPntCoord[0] - aCenter[0],
-				       aPntCoord[1] - aCenter[1],
-				       aPntCoord[2] - aCenter[2] };
-		  vtkMath::Normalize(aVector);
-
-		  float aCross[3];
-		  vtkMath::Cross(aVector,aVector0,aCross);
-		  bool aGreaterThanPi = vtkMath::Dot(aCross,aCross21) < 0;
-		  float aCosinus = vtkMath::Dot(aVector,aVector0);
-		  if(aCosinus > 1.0)
-		    aCosinus = 1.0;
-		  if(aCosinus < -1.0)
-		    aCosinus = -1.0;
-		  static float a2Pi = 2.0 * vtkMath::Pi();
-		  float anAngle = acos(aCosinus);
-		  //cout<<"\t\taPntId = "<<aPntId<<" {"<<aPntCoord[0]<<", "<<aPntCoord[1]<<", "<<aPntCoord[2]<<"}";
-		  //cout<<"; aGreaterThanPi = "<<aGreaterThanPi<<"; aCosinus = "<<aCosinus<<"; anAngle = "<<anAngle<<"\n";
-		  if(aGreaterThanPi)
-		    anAngle = a2Pi - anAngle;
-		  aSortedPointIds[anAngle] = aPntId;
-		  //cout<<"\t\t\tanAngle = "<<anAngle<<"\n";
-		}
-		if(!aSortedPointIds.empty()){
-		  aCellType = VTK_POLYGON;
-		  int numFacePts = aSortedPointIds.size();
-		  std::vector<vtkIdType> aConnectivities(numFacePts);
-		  TSortedPointIds::const_iterator anIter = aSortedPointIds.begin();
-		  TSortedPointIds::const_iterator anEndIter = aSortedPointIds.end();
-		  for(vtkIdType anId = 0; anIter != anEndIter; anIter++, anId++){
-		    vtkIdType aPntId = anIter->second;
-#ifdef USE_ROBUST_TRIANGULATION
-		    aConnectivities[anId] = aPts[aPntId];
-#else
-		    aConnectivities[anId] = aPntId;
-#endif
-		  }
-		  newCellId = output->InsertNextCell(aCellType,numFacePts,&aConnectivities[0]);
-		  if(myStoreMapping)
-		    myVTK2ObjIds.push_back(cellId);
-		  outputCD->CopyData(cd,cellId,newCellId);
-		}
-	      }
-	    }
-	  }
-
 	  break;
 	}
         case VTK_TETRA: {
@@ -1029,14 +747,6 @@ VTKViewer_GeometryFilter
   vtkDebugMacro(<<"Extracted " << input->GetNumberOfPoints() << " points,"
   << output->GetNumberOfCells() << " cells.");
 
-#ifdef USE_ROBUST_TRIANGULATION
-  anUnstructuredGrid->Delete();
-  aDelaunayPoints->Delete();
-
-  aDelaunay3D->Delete();
-  aGeometryFilter->Delete();
-#endif
-
   cell->Delete();
 
   cellIds->Delete();