Back to index

salome-gui  6.5.0
VTKViewer_ArcBuilder.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 //  File   : VTKViewer_ArcBuilder.cxx
00021 //  Author : Roman NIKOLAEV
00022 //  Module : GUI
00023 //
00024 #include "VTKViewer_ArcBuilder.h"
00025 
00026 #include <math.h>
00027 #include <float.h>
00028 
00029 //VTK includes
00030 #include <vtkMath.h>
00031 #include <vtkUnstructuredGrid.h>
00032 #include <vtkTransformFilter.h>
00033 #include <vtkTransform.h>
00034 #include <vtkPoints.h>
00035 #include <vtkVertex.h>
00036 #include <vtkCellArray.h>
00037 #include <vtkTriangle.h>
00038 #include <vtkPolyData.h>
00039 #include <vtkPointData.h>
00040 
00041 #define PRECISION 10e-4
00042 #define ANGLE_PRECISION 0.5
00043 //#define _MY_DEBUG_
00044 
00045 #ifdef _MY_DEBUG_
00046 #include <iostream>
00047 #endif
00048 
00049 
00050 bool CheckAngle(const double compare, const double angle){
00051   if((angle <= compare - ANGLE_PRECISION) || (angle >= compare + ANGLE_PRECISION))
00052     return true;
00053   else
00054     return false;
00055 }
00056 
00057 //------------------------------------------------------------------------
00062 XYZ::XYZ(double X, double Y, double Z):
00063  x(X),y(Y),z(Z) {}
00064 
00068 XYZ::XYZ(){}
00069 
00073 XYZ::~XYZ()
00074 {}
00075 
00076 /*
00077  * Calculate modulus
00078  */
00079 double XYZ::Modulus () const {
00080   return sqrt (x * x + y * y + z * z);
00081 }
00082 
00083 //------------------------------------------------------------------------
00088 Pnt::Pnt(double X, 
00089          double Y, 
00090          double Z,
00091          double ScalarValue):
00092   coord(X,Y,Z),
00093   scalarValue(ScalarValue)
00094 {
00095 }
00096 
00097 Pnt::Pnt()
00098 {
00099 }
00100 
00104 Pnt::~Pnt()
00105 {
00106 }
00107 
00108 //------------------------------------------------------------------------
00113 Vec::Vec (const double Xv,
00114           const double Yv,
00115           const double Zv)
00116 {
00117   double D = sqrt (Xv * Xv + Yv * Yv + Zv * Zv);
00118   if(D != 0) {
00119     coord.SetX(Xv / D);
00120     coord.SetY(Yv / D);
00121     coord.SetZ(Zv / D);
00122   }
00123 }
00124 
00128 Vec::~Vec(){}
00129 
00130 
00134 double Vec::AngleBetween(const Vec & Other)
00135 {
00136   double res;
00137   double numerator = GetXYZ().X()*Other.GetXYZ().X()+GetXYZ().Y()*Other.GetXYZ().Y()+GetXYZ().Z()*Other.GetXYZ().Z();
00138   double denumerator = GetXYZ().Modulus()*Other.GetXYZ().Modulus();
00139   double d = numerator/denumerator;
00140   if( d < -1 && d > -1 - PRECISION )
00141     d = -1;
00142   else if( d > 1 && d < 1 + PRECISION )
00143     d = 1;
00144   res = acos( d );
00145   return res;
00146 }
00147 
00151 double Vec::AngleBetweenInGrad(const Vec & Other){
00152   return AngleBetween(Other)*vtkMath::DegreesFromRadians(1.);
00153 }
00154 
00155 /*
00156  * Calculate vector multiplication
00157 */
00158 Vec Vec::VectMultiplication(const Vec & Other) const{
00159   double x = GetXYZ().Y()*Other.GetXYZ().Z() - GetXYZ().Z()*Other.GetXYZ().Y();
00160   double y = GetXYZ().Z()*Other.GetXYZ().X() - GetXYZ().X()*Other.GetXYZ().Z();
00161   double z = GetXYZ().X()*Other.GetXYZ().Y() - GetXYZ().Y()*Other.GetXYZ().X();
00162   Vec *aRes  = new Vec(x,y,z);
00163   return *aRes;
00164 }
00165 
00166 /*---------------------Class Plane --------------------------------*/
00170 Plane::Plane(const Pnt& thePnt1, const Pnt& thePnt2, const Pnt& thePnt3)
00171 {
00172   CalculatePlane(thePnt1,thePnt2,thePnt3);
00173 }
00174 
00178 Plane::~Plane()
00179 {}
00180 
00181 /*
00182  * Return plane normale
00183  */
00184 Vec Plane::GetNormale() const{
00185   return Vec(myA,myB,myC);
00186 }
00187 
00188 /*
00189  * Calculate A,B,C coeefs of plane
00190  */
00191 void Plane::CalculatePlane(const Pnt& thePnt1, const Pnt& thePnt2, const Pnt& thePnt3){
00192   
00193   double x1 = thePnt1.GetXYZ().X();
00194   double x2 = thePnt2.GetXYZ().X();
00195   double x3 = thePnt3.GetXYZ().X();
00196   
00197   double y1 = thePnt1.GetXYZ().Y();
00198   double y2 = thePnt2.GetXYZ().Y();
00199   double y3 = thePnt3.GetXYZ().Y();
00200   
00201   double z1 = thePnt1.GetXYZ().Z();
00202   double z2 = thePnt2.GetXYZ().Z();
00203   double z3 = thePnt3.GetXYZ().Z();
00204   
00205   myA = y1*(z2 - z3) + y2*(z3 - z1) + y3*(z1 - z2);
00206   myB = z1*(x2 - x3) + z2*(x3 - x1) + z3*(x1 - x2);
00207   myC = x1*(y2 - y3) + x2*(y3 - y1) + x3*(y1 - y2);
00208   
00209 #ifdef _MY_DEBUG_
00210   cout<<"Plane A: "<<myA<<endl;
00211   cout<<"Plane B: "<<myB<<endl;
00212   cout<<"Plane C: "<<myC<<endl;
00213 #endif  
00214 
00215 } 
00216 
00217 
00218 /*---------------------Class VTKViewer_ArcBuilder --------------------------------*/
00222 VTKViewer_ArcBuilder::VTKViewer_ArcBuilder(const Pnt& thePnt1,
00223                                            const Pnt& thePnt2,
00224                                            const Pnt& thePnt3,
00225                                            double theAngle):
00226   myStatus(Arc_Error),
00227   myAngle(theAngle)
00228 {
00229   Vec V1(thePnt2.GetXYZ().X()-thePnt1.GetXYZ().X(),
00230          thePnt2.GetXYZ().Y()-thePnt1.GetXYZ().Y(),
00231          thePnt2.GetXYZ().Z()-thePnt1.GetXYZ().Z());
00232   
00233   Vec V2(thePnt2.GetXYZ().X()-thePnt3.GetXYZ().X(),
00234          thePnt2.GetXYZ().Y()-thePnt3.GetXYZ().Y(),
00235          thePnt2.GetXYZ().Z()-thePnt3.GetXYZ().Z());
00236 
00237   double angle = V1.AngleBetweenInGrad(V2);
00238   
00239   //Check that points are not belong one line
00240 #ifdef _MY_DEBUG_
00241   cout<<"Angle for check: "<<angle<<endl;
00242 #endif
00243   if(CheckAngle(180,angle)) {
00244     
00245     // Build plane by three points
00246     Plane aPlane(thePnt1,thePnt2,thePnt3);
00247     
00248     //Plane normales
00249     Vec aPlaneNormale = aPlane.GetNormale();
00250 #ifdef _MY_DEBUG_
00251     std::cout<<"X Normale: "<<aPlaneNormale.GetXYZ().X()<<std::endl;
00252     std::cout<<"Y Normale: "<<aPlaneNormale.GetXYZ().Y()<<std::endl;
00253     std::cout<<"Z Normale: "<<aPlaneNormale.GetXYZ().Z()<<std::endl;
00254 #endif
00255     //OZ vector
00256     Vec OZ(0,0,1);
00257     
00258     //Projection plane normale on XOY
00259     Vec aNormaleProjection (aPlaneNormale.GetXYZ().X(),
00260                             aPlaneNormale.GetXYZ().Y(),
00261                             0);
00262 #ifdef _MY_DEBUG_
00263     std::cout<<"X Normale Projection: "<<aNormaleProjection.GetXYZ().X()<<std::endl;
00264     std::cout<<"Y Normale Projection: "<<aNormaleProjection.GetXYZ().Y()<<std::endl;
00265     std::cout<<"Z Normale Projection: "<<aNormaleProjection.GetXYZ().Z()<<std::endl;
00266 #endif
00267     
00268     //Rotation axis
00269     Vec aAxis = aNormaleProjection.VectMultiplication(OZ);
00270 #ifdef _MY_DEBUG_
00271     std::cout<<"X Axis: "<<aAxis.GetXYZ().X()<<std::endl;
00272     std::cout<<"Y Axis: "<<aAxis.GetXYZ().Y()<<std::endl;
00273     std::cout<<"Z Axis: "<<aAxis.GetXYZ().Z()<<std::endl;
00274 #endif   
00275     //Rotation angle
00276     double anAngle = OZ.AngleBetweenInGrad(aPlaneNormale);
00277 #ifdef _MY_DEBUG_
00278     std::cout<<"Rotation Angle :"<<anAngle<<endl;
00279 #endif
00280     PntList aInputPnts;
00281     aInputPnts.push_back(thePnt1);
00282     aInputPnts.push_back(thePnt2);
00283     aInputPnts.push_back(thePnt3);
00284     
00285     vtkUnstructuredGrid* aGrid = BuildGrid(aInputPnts);
00286     
00287     bool needRotation = true;
00288     if(anAngle  == 0 || anAngle == 180)
00289       needRotation = false;
00290     
00291     if(aGrid) {
00292       vtkUnstructuredGrid* aTransformedGrid;
00293       if(needRotation) {
00294         aTransformedGrid = TransformGrid(aGrid,aAxis,anAngle);    
00295         aTransformedGrid->Update();
00296 #ifdef _MY_DEBUG_
00297         cout<<"Need Rotation!!!"<<endl;
00298 #endif
00299       }
00300       else {
00301         aTransformedGrid = aGrid;
00302 #ifdef _MY_DEBUG_
00303         cout<<"Rotation does not need!!!"<<endl;
00304 #endif
00305       }
00306       
00307       double coords[3];
00308       aTransformedGrid->GetPoint(0,coords);
00309       myPnt1 = Pnt(coords[0],coords[1],coords[2], thePnt1.GetScalarValue());
00310       aTransformedGrid->GetPoint(1,coords);
00311       myPnt2 = Pnt(coords[0],coords[1],coords[2], thePnt2.GetScalarValue());
00312       aTransformedGrid->GetPoint(2,coords);
00313       myPnt3 = Pnt(coords[0],coords[1],coords[2], thePnt3.GetScalarValue());
00314       std::vector<double> aScalarValues;
00315       vtkUnstructuredGrid* anArc = BuildArc(aScalarValues);
00316       vtkUnstructuredGrid* anTransArc;
00317       if(needRotation) {
00318         anTransArc = TransformGrid(anArc,aAxis,-anAngle);
00319         anTransArc->Update();
00320       }
00321       else
00322         anTransArc = anArc;
00323       
00324       myPoints = anTransArc->GetPoints();
00325       myScalarValues = aScalarValues;
00326       myStatus = Arc_Done;
00327     }
00328   }
00329   else{
00330 #ifdef _MY_DEBUG_
00331     std::cout<<"Points lay on the one line !"<<endl;
00332 #endif           
00333     PntList aList;
00334     aList.push_back(thePnt1);
00335     aList.push_back(thePnt2);
00336     aList.push_back(thePnt3);
00337     vtkUnstructuredGrid* aGrid = BuildGrid(aList);
00338     aGrid->Update();
00339     myPoints = aGrid->GetPoints();
00340 
00341     myScalarValues.clear();
00342     myScalarValues.push_back(thePnt1.GetScalarValue());
00343     myScalarValues.push_back(thePnt2.GetScalarValue());
00344     myScalarValues.push_back(thePnt3.GetScalarValue());
00345     myStatus = Arc_Done;
00346   }
00347 }
00348 
00352 VTKViewer_ArcBuilder::~VTKViewer_ArcBuilder()
00353 {}
00354 
00355 
00356 /*
00357  * Add to the vtkUnstructuredGrid points from input list
00358  */
00359 vtkUnstructuredGrid* VTKViewer_ArcBuilder::BuildGrid(const PntList& theList) const
00360 {
00361   int aListsize = theList.size();  
00362   vtkUnstructuredGrid* aGrid = NULL;
00363   
00364   if(aListsize != 0) {
00365     aGrid = vtkUnstructuredGrid::New();
00366     vtkPoints* aPoints = vtkPoints::New();
00367     aPoints->SetNumberOfPoints(aListsize);
00368     
00369     aGrid->Allocate(aListsize);
00370     
00371     PntList::const_iterator it = theList.begin();
00372     int aCounter = 0;
00373     for(;it != theList.end();it++) {
00374       Pnt aPnt =  *it;
00375       aPoints->InsertPoint(aCounter, 
00376                            aPnt.GetXYZ().X(), 
00377                            aPnt.GetXYZ().Y(),
00378                            aPnt.GetXYZ().Z());
00379       vtkVertex* aVertex = vtkVertex::New();
00380       aVertex->GetPointIds()->SetId(0, aCounter);
00381       aGrid->InsertNextCell(aVertex->GetCellType(), aVertex->GetPointIds());
00382       aCounter++;
00383       aVertex->Delete();
00384     }
00385     aGrid->SetPoints(aPoints);
00386     aPoints->Delete();
00387   }
00388   return aGrid;
00389 }
00390 
00391 
00392 vtkUnstructuredGrid* 
00393 VTKViewer_ArcBuilder::TransformGrid(vtkUnstructuredGrid* theGrid, 
00394                                     const Vec& theAxis, const double angle) const
00395 {
00396   vtkTransform *aTransform = vtkTransform::New();
00397   aTransform->RotateWXYZ(angle, theAxis.GetXYZ().X(), theAxis.GetXYZ().Y(), theAxis.GetXYZ().Z());
00398   vtkTransformFilter* aTransformFilter  = vtkTransformFilter::New();
00399   aTransformFilter->SetTransform(aTransform);
00400   aTransformFilter->SetInput(theGrid);
00401   aTransform->Delete();
00402   return aTransformFilter->GetUnstructuredGridOutput();
00403 }
00404 
00405 
00406 void VTKViewer_ArcBuilder::GetAngle(const double theAngle){
00407   myAngle = theAngle;
00408 }
00409 
00410 double InterpolateScalarValue(int index, int count, double firstValue, double middleValue, double lastValue)
00411 {
00412   bool isFirstHalf = index <= count / 2;
00413   double first = isFirstHalf ? firstValue : lastValue;
00414   double last = isFirstHalf ? middleValue : middleValue;
00415   double ratio = (double)index / (double)count;
00416   double position = isFirstHalf ? ratio * 2 : ( 1 - ratio ) * 2;
00417   double value = first + (last - first) * position;
00418   return value;
00419 }
00420 
00421 vtkUnstructuredGrid* VTKViewer_ArcBuilder::BuildArc(std::vector<double>& theScalarValues){
00422   double x1 = myPnt1.GetXYZ().X(); double x2 = myPnt2.GetXYZ().X(); double x3 = myPnt3.GetXYZ().X();
00423   double y1 = myPnt1.GetXYZ().Y(); double y2 = myPnt2.GetXYZ().Y(); double y3 = myPnt3.GetXYZ().Y();
00424   double z =  myPnt1.GetXYZ().Z();  //Points on plane || XOY
00425   
00426 
00427   theScalarValues.clear();
00428 
00429   double K1 = 0;
00430   double K2 = 0;
00431   bool okK1 = false, okK2 = false;
00432   if ( fabs(x2 - x1) > DBL_MIN )
00433     K1 = (y2 - y1)/(x2 - x1), okK1 = true;
00434   if ( fabs(x3 - x2) > DBL_MIN )
00435     K2 = (y3 - y2)/(x3 - x2), okK2 = true;
00436   
00437 #ifdef _MY_DEBUG_   
00438   std::cout<<"K1 : "<< K1 <<endl; 
00439   std::cout<<"K2 : "<< K2 <<endl; 
00440 #endif
00441   double xCenter;
00442   if( !okK2 ) //K2 --> infinity
00443     xCenter = (K1*(y1-y3) + (x1+x2))/2.0;
00444   
00445   else if( !okK1 ) //K1 --> infinity
00446     xCenter = (K2*(y1-y3) - (x2+x3))/(-2.0);
00447   
00448   else 
00449     xCenter = (K1*K2*(y1-y3) + K2*(x1+x2) - K1*(x2+x3))/ (2.0*(K2-K1));
00450   
00451   double yCenter;
00452   
00453   if ( K1 == 0 )
00454     yCenter =  (-1/K2)*(xCenter - (x2+x3)/2.0) + (y2 + y3)/2.0;
00455   else 
00456     yCenter =  (-1/K1)*(xCenter - (x1+x2)/2.0) + (y1 + y2)/2.0;
00457   
00458 #ifdef _MY_DEBUG_   
00459   double zCenter = z;
00460   std::cout<<"xCenter : "<<xCenter<<endl;
00461   std::cout<<"yCenter : "<<yCenter<<endl;
00462   std::cout<<"zCenter : "<<zCenter<<endl;
00463 #endif
00464   double aRadius = sqrt((x1 - xCenter)*(x1 - xCenter) + (y1 - yCenter)*(y1 - yCenter));
00465   
00466   double angle1 = GetPointAngleOnCircle(xCenter,yCenter,x1,y1);
00467   double angle2 = GetPointAngleOnCircle(xCenter,yCenter,x2,y2);
00468   double angle3 = GetPointAngleOnCircle(xCenter,yCenter,x3,y3);
00469   
00470   
00471   double aMaxAngle = vtkMath::RadiansFromDegrees(1.)*myAngle*2;
00472   
00473   /*  double aTotalAngle =  fabs(angle3 - angle1);
00474   
00475   if (aTotalAngle > vtkMath::DoublePi())
00476     aTotalAngle = 2*vtkMath::DoublePi()-aTotalAngle;
00477   */
00478   
00479   double aTotalAngle = 0;
00480   IncOrder aOrder = GetArcAngle(angle1,angle2,angle3,&aTotalAngle);
00481   
00482   vtkUnstructuredGrid *aC = NULL;
00483 
00484   if(aTotalAngle > aMaxAngle) {
00485     int nbSteps = int(aTotalAngle/aMaxAngle)+1;
00486     double anIncrementAngle = aTotalAngle/nbSteps;
00487     double aCurrentAngle = angle1;
00488     if(aOrder == VTKViewer_ArcBuilder::MINUS)
00489       aCurrentAngle-=anIncrementAngle;
00490     else
00491       aCurrentAngle+=anIncrementAngle;
00492 #ifdef _MY_DEBUG_
00493     cout<<"Total angle :"<<aTotalAngle<<endl;
00494     cout<<"Max Increment Angle :"<<aMaxAngle<<endl;
00495     cout<<"Real Increment angle :"<<anIncrementAngle<<endl;
00496     cout<<"Nb Steps : "<<nbSteps<<endl;
00497 #endif
00498   
00499     PntList aList;
00500     aList.push_back(myPnt1);
00501     theScalarValues.push_back(myPnt1.GetScalarValue());
00502     for(int i=1;i<=nbSteps-1;i++){
00503       double x = xCenter + aRadius*cos(aCurrentAngle);
00504       double y = yCenter + aRadius*sin(aCurrentAngle);
00505       double value = InterpolateScalarValue(i, nbSteps, myPnt1.GetScalarValue(), myPnt2.GetScalarValue(), myPnt3.GetScalarValue());
00506       Pnt aPnt(x,y,z,value);
00507 
00508       aList.push_back(aPnt);
00509       theScalarValues.push_back(aPnt.GetScalarValue());
00510       if(aOrder == VTKViewer_ArcBuilder::MINUS)
00511         aCurrentAngle-=anIncrementAngle;
00512       else
00513         aCurrentAngle+=anIncrementAngle;
00514     }
00515     aList.push_back(myPnt3);
00516     theScalarValues.push_back(myPnt3.GetScalarValue());
00517     
00518     aC = BuildGrid(aList);
00519 #ifdef _MY_DEBUG_
00520   cout<<"angle1 : "<<angle1*vtkMath::DoubleRadiansToDegrees()<<endl;
00521   cout<<"angle2 : "<<angle2*vtkMath::DoubleRadiansToDegrees()<<endl;
00522   cout<<"angle3 : "<<angle3*vtkMath::DoubleRadiansToDegrees()<<endl;
00523 #endif
00524   }
00525   else{
00526     PntList aList;
00527     aList.push_back(myPnt1);
00528     aList.push_back(myPnt2);
00529     aList.push_back(myPnt3);
00530     aC = BuildGrid(aList);
00531 
00532     theScalarValues.push_back(myPnt1.GetScalarValue());
00533     theScalarValues.push_back(myPnt2.GetScalarValue());
00534     theScalarValues.push_back(myPnt3.GetScalarValue());
00535   }
00536   return aC;
00537 }
00538 
00539 double VTKViewer_ArcBuilder::
00540 GetPointAngleOnCircle(const double theXCenter, const double theYCenter,
00541                       const double theXPoint, const double theYPoint){
00542   double result = atan2(theYCenter - theYPoint, theXPoint - theXCenter);
00543   if(result < 0 )
00544     result = result+vtkMath::DoublePi()*2;
00545   return vtkMath::DoublePi()*2-result;
00546   return result;
00547 }
00548 
00549 vtkPoints* VTKViewer_ArcBuilder::GetPoints(){
00550   return myPoints;
00551 }
00552 
00553 const std::vector<double>& VTKViewer_ArcBuilder::GetScalarValues()
00554 {
00555   return myScalarValues;
00556 }
00557 
00558 VTKViewer_ArcBuilder::IncOrder VTKViewer_ArcBuilder::GetArcAngle( const double& P1, const double& P2, const double& P3,double* Ang){
00559   IncOrder aResult;
00560   if(P1 < P2 && P2 < P3){
00561     *Ang = P3 - P1;
00562     aResult = VTKViewer_ArcBuilder::PLUS;
00563   }
00564   else if((P1 < P3 && P3 < P2) || (P2 < P1 && P1 < P3)){
00565     *Ang = 2*vtkMath::DoublePi() - P3 + P1;
00566     aResult = VTKViewer_ArcBuilder::MINUS;
00567   }
00568   else if((P2 < P3 && P3 < P1) || (P3 < P1 && P1 < P2)){
00569     *Ang = 2*vtkMath::DoublePi() - P1 + P3;
00570     aResult = VTKViewer_ArcBuilder::PLUS;
00571   }
00572   else if(P3 < P2 && P2 < P1){
00573     *Ang = P1 - P3;
00574     aResult = VTKViewer_ArcBuilder::MINUS;
00575   }
00576   return aResult;
00577 }
00578 
00579 //------------------------------------------------------------------------
00580 Pnt CreatePnt(vtkCell* cell, vtkDataArray* scalars, vtkIdType index)
00581 {
00582   vtkFloatingPointType coord[3];
00583   cell->GetPoints()->GetPoint(index, coord);
00584   vtkIdType pointId = cell->GetPointId(index);
00585   double scalarValue = scalars ? scalars->GetTuple1(pointId) : 0;
00586   Pnt point(coord[0], coord[1], coord[2], scalarValue);
00587   return point;
00588 }
00589 
00590 //------------------------------------------------------------------------
00591 vtkIdType Build1DArc(vtkIdType cellId, vtkUnstructuredGrid* input, 
00592                      vtkPolyData *output,vtkIdType *pts, 
00593                      vtkFloatingPointType myMaxArcAngle){
00594   
00595   vtkIdType aResult = -1;
00596   vtkIdType *aNewPoints;
00597 
00598   vtkDataArray* inputScalars = input->GetPointData()->GetScalars();
00599   vtkDataArray* outputScalars = output->GetPointData()->GetScalars();
00600 
00601   vtkCell* aCell = input->GetCell(cellId);
00602   //Get All points from input cell
00603   Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
00604   Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
00605   Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
00606 
00607   VTKViewer_ArcBuilder aBuilder(P0,P2,P1,myMaxArcAngle);
00608   if (aBuilder.GetStatus() != VTKViewer_ArcBuilder::Arc_Done) {
00609     return aResult;
00610   }
00611   else{
00612     vtkPoints* aPoints = aBuilder.GetPoints();
00613     std::vector<double> aScalarValues = aBuilder.GetScalarValues();
00614     vtkIdType aNbPts = aPoints->GetNumberOfPoints();
00615     aNewPoints = new vtkIdType[aNbPts];
00616     vtkIdType curID;
00617     vtkIdType aCellType = VTK_POLY_LINE;
00618     
00619     aNewPoints[0] = pts[0];
00620     for(vtkIdType idx = 1; idx < aNbPts-1;idx++) {
00621       curID = output->GetPoints()->InsertNextPoint(aPoints->GetPoint(idx));
00622       if( outputScalars )
00623         outputScalars->InsertNextTuple1(aScalarValues[idx]);
00624       aNewPoints[idx] = curID;
00625     }
00626     aNewPoints[aNbPts-1] = pts[1];
00627     
00628     aResult = output->InsertNextCell(aCellType,aNbPts,aNewPoints);
00629     return aResult;
00630   } 
00631 }
00632 
00637 vtkIdType MergevtkPoints(const std::vector<vtkPoints*>& theCollection,
00638                          const std::vector< std::vector<double> >& theScalarCollection,
00639                          vtkPoints* thePoints,
00640                          std::map<int, double>& thePntId2ScalarValue,
00641                          vtkIdType* &theIds){
00642   vtkIdType aNbPoints = 0;
00643   vtkIdType anIdCounter = 0;
00644   vtkIdType aNewPntId = 0;
00645   
00646   //Compute number of points
00647   std::vector<vtkPoints*>::const_iterator it = theCollection.begin();
00648   for(;it != theCollection.end();it++){
00649     vtkPoints* aPoints = *it;
00650     if(aPoints) { 
00651       aNbPoints += aPoints->GetNumberOfPoints()-1; //Because we add all points except last
00652     }
00653   }
00654   it = theCollection.begin();
00655   std::vector< std::vector<double> >::const_iterator itScalar = theScalarCollection.begin();
00656   theIds = new vtkIdType[aNbPoints];
00657   // ..and add all points
00658   for(;it != theCollection.end() && itScalar != theScalarCollection.end(); it++, itScalar++){
00659     vtkPoints* aPoints = *it;
00660     std::vector<double> aScalarValues = *itScalar;
00661     
00662     if(aPoints){
00663       for(vtkIdType idx = 0;idx < aPoints->GetNumberOfPoints()-1;idx++){
00664         aNewPntId = thePoints->InsertNextPoint(aPoints->GetPoint(idx));
00665         theIds[anIdCounter] = aNewPntId;
00666         thePntId2ScalarValue[ aNewPntId ] = aScalarValues[idx];
00667         anIdCounter++;
00668       }
00669     }
00670   }
00671   return aNbPoints;
00672 }