From 2f38498d1a2b63b6895ed9fa0600e8b1311ac9d2 Mon Sep 17 00:00:00 2001 From: rnv Date: Wed, 27 May 2009 11:18:35 +0000 Subject: [PATCH] Implementation of the issue 20115: [CEA 308] Quadratic elements visualization. --- src/SVTK/SVTK_DeviceActor.cxx | 27 + src/SVTK/SVTK_DeviceActor.h | 8 + src/VTKViewer/Makefile.am | 6 +- src/VTKViewer/VTKViewer_Actor.cxx | 30 + src/VTKViewer/VTKViewer_Actor.h | 9 + src/VTKViewer/VTKViewer_ArcBuilder.cxx | 601 +++++++++++++++++++++ src/VTKViewer/VTKViewer_ArcBuilder.h | 179 ++++++ src/VTKViewer/VTKViewer_GeometryFilter.cxx | 441 ++++++++++++--- src/VTKViewer/VTKViewer_GeometryFilter.h | 15 + 9 files changed, 1239 insertions(+), 77 deletions(-) create mode 100644 src/VTKViewer/VTKViewer_ArcBuilder.cxx create mode 100644 src/VTKViewer/VTKViewer_ArcBuilder.h diff --git a/src/SVTK/SVTK_DeviceActor.cxx b/src/SVTK/SVTK_DeviceActor.cxx index 786add600..30e5e8678 100644 --- a/src/SVTK/SVTK_DeviceActor.cxx +++ b/src/SVTK/SVTK_DeviceActor.cxx @@ -656,3 +656,30 @@ vtkDataSetMapper* SVTK_DeviceActor::GetDataSetMapper() { return myMapper; } + +/*! + * On/Off representation 2D quadratic element as arked polygon + */ +void SVTK_DeviceActor::SetQuadraticArcMode(bool theFlag){ + myGeomFilter->SetQuadraticArcMode(theFlag); +} + +/*! + * Return true if 2D quadratic element displayed as arked polygon + */ +bool SVTK_DeviceActor::GetQuadraticArcMode(){ + return myGeomFilter->GetQuadraticArcMode(); +} +/*! + * Set Max angle for representation 2D quadratic element as arked polygon + */ +void SVTK_DeviceActor::SetQuadraticArcAngle(vtkFloatingPointType theMaxAngle){ + myGeomFilter->SetQuadraticArcAngle(theMaxAngle); +} + +/*! + * Return Max angle of the representation 2D quadratic element as arked polygon + */ +vtkFloatingPointType SVTK_DeviceActor::GetQuadraticArcAngle(){ + return myGeomFilter->GetQuadraticArcAngle(); +} diff --git a/src/SVTK/SVTK_DeviceActor.h b/src/SVTK/SVTK_DeviceActor.h index 95949133b..18c8cee97 100644 --- a/src/SVTK/SVTK_DeviceActor.h +++ b/src/SVTK/SVTK_DeviceActor.h @@ -225,6 +225,14 @@ class SVTK_EXPORT SVTK_DeviceActor: public vtkLODActor vtkDataSetMapper* GetDataSetMapper(); + //---------------------------------------------------------------------------- + //! Setting for displaying quadratic elements + virtual void SetQuadraticArcMode(bool theFlag); + virtual bool GetQuadraticArcMode(); + + virtual void SetQuadraticArcAngle(vtkFloatingPointType theMaxAngle); + virtual vtkFloatingPointType GetQuadraticArcAngle(); + protected: SVTK::Representation::Type myRepresentation; vtkProperty *myProperty; diff --git a/src/VTKViewer/Makefile.am b/src/VTKViewer/Makefile.am index cb3e1078d..e89b79601 100755 --- a/src/VTKViewer/Makefile.am +++ b/src/VTKViewer/Makefile.am @@ -49,7 +49,8 @@ salomeinclude_HEADERS = \ VTKViewer_ViewManager.h \ VTKViewer_ViewModel.h \ VTKViewer_ViewWindow.h \ - VTKViewer_Functor.h + VTKViewer_Functor.h \ + VTKViewer_ArcBuilder.h dist_libVTKViewer_la_SOURCES = \ VTKViewer_CellLocationsArray.cxx \ @@ -70,7 +71,8 @@ dist_libVTKViewer_la_SOURCES = \ VTKViewer_ViewManager.cxx \ VTKViewer_ViewModel.cxx \ VTKViewer_ConvexTool.cxx \ - VTKViewer_ViewWindow.cxx + VTKViewer_ViewWindow.cxx \ + VTKViewer_ArcBuilder.cxx MOC_FILES = \ VTKViewer_RenderWindow_moc.cxx \ diff --git a/src/VTKViewer/VTKViewer_Actor.cxx b/src/VTKViewer/VTKViewer_Actor.cxx index d08612771..920e3b08c 100755 --- a/src/VTKViewer/VTKViewer_Actor.cxx +++ b/src/VTKViewer/VTKViewer_Actor.cxx @@ -673,4 +673,34 @@ VTKViewer_Actor myIsHighlighted = theIsHighlight; } +/*! + * On/Off representation 2D quadratic element as arked polygon + */ +void VTKViewer_Actor::SetQuadraticArcMode(bool theFlag){ + myGeomFilter->SetQuadraticArcMode(theFlag); +} + +/*! + * Return true if 2D quadratic element displayed as arked polygon + */ +bool VTKViewer_Actor::GetQuadraticArcMode() const{ + return myGeomFilter->GetQuadraticArcMode(); +} +/*! + * Set Max angle for representation 2D quadratic element as arked polygon + */ +void VTKViewer_Actor::SetQuadraticArcAngle(vtkFloatingPointType theMaxAngle){ + myGeomFilter->SetQuadraticArcAngle(theMaxAngle); +} + +/*! + * Return Max angle of the representation 2D quadratic element as arked polygon + */ +vtkFloatingPointType VTKViewer_Actor::GetQuadraticArcAngle() const{ + return myGeomFilter->GetQuadraticArcAngle(); +} + + + + vtkCxxSetObjectMacro(VTKViewer_Actor,PreviewProperty,vtkProperty); diff --git a/src/VTKViewer/VTKViewer_Actor.h b/src/VTKViewer/VTKViewer_Actor.h index 77e54f245..439e45dbf 100755 --- a/src/VTKViewer/VTKViewer_Actor.h +++ b/src/VTKViewer/VTKViewer_Actor.h @@ -301,6 +301,15 @@ class VTKVIEWER_EXPORT VTKViewer_Actor : public vtkLODActor void SetPreviewProperty(vtkProperty* theProperty); + //---------------------------------------------------------------------------- + //! Setting for displaying quadratic elements + virtual void SetQuadraticArcMode(bool theFlag); + virtual bool GetQuadraticArcMode() const; + + virtual void SetQuadraticArcAngle(vtkFloatingPointType theMaxAngle); + virtual vtkFloatingPointType GetQuadraticArcAngle() const; + + protected: //---------------------------------------------------------------------------- bool myIsResolveCoincidentTopology; diff --git a/src/VTKViewer/VTKViewer_ArcBuilder.cxx b/src/VTKViewer/VTKViewer_ArcBuilder.cxx new file mode 100644 index 000000000..0c077bbdc --- /dev/null +++ b/src/VTKViewer/VTKViewer_ArcBuilder.cxx @@ -0,0 +1,601 @@ +// Copyright (C) 2007-2008 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. +// +// 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 +// +// File : VTKViewer_ArcBuilder.cxx +// Author : Roman NIKOLAEV +// Module : GUI +// +#include "VTKViewer_ArcBuilder.h" + +#include + +//VTK includes +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#define PRECISION 10e-4 +#define ANGLE_PRECISION 0.5 +//#define _MY_DEBUG_ + +#ifdef _MY_DEBUG_ +#include +#endif + + +bool CheckAngle(const double compare, const double angle){ + if((angle <= compare - ANGLE_PRECISION) || (angle >= compare + ANGLE_PRECISION)) + return true; + else + return false; +} + +//------------------------------------------------------------------------ +/*! + *Class XYZ + *Constructor + */ +XYZ::XYZ(double X, double Y, double Z): + x(X),y(Y),z(Z) {} + +/*! + * Empty Constructor + */ +XYZ::XYZ(){} + +/*! + * Destructor + */ +XYZ::~XYZ() +{} + +/* + * Calculate modulus + */ +double XYZ::Modulus () const { + return sqrt (x * x + y * y + z * z); +} + +//------------------------------------------------------------------------ +/*! + * Class Pnt + * Constructor + */ +Pnt::Pnt(double X, + double Y, + double Z): + coord(X,Y,Z) {} + + +Pnt::Pnt(){} +/*! + * Destructor + */ +Pnt::~Pnt() +{} + +//------------------------------------------------------------------------ +/*! + * Class Vec + * Constructor + */ +Vec::Vec (const double Xv, + const double Yv, + const double Zv) +{ + double D = sqrt (Xv * Xv + Yv * Yv + Zv * Zv); + if(D != 0) { + coord.SetX(Xv / D); + coord.SetY(Yv / D); + coord.SetZ(Zv / D); + } +} + +/*! + * Destructor + */ +Vec::~Vec(){} + + +/*! + * Calculate angle between vectors in radians + */ +double Vec::AngleBetween(const Vec & Other) +{ + double res; + double numerator = GetXYZ().X()*Other.GetXYZ().X()+GetXYZ().Y()*Other.GetXYZ().Y()+GetXYZ().Z()*Other.GetXYZ().Z(); + double denumerator = GetXYZ().Modulus()*Other.GetXYZ().Modulus(); + res = acos(numerator/denumerator); + return res; +} + +/*! + * Calculate angle between vectors in degrees + */ +double Vec::AngleBetweenInGrad(const Vec & Other){ + return AngleBetween(Other)*vtkMath::DoubleRadiansToDegrees(); +} + +/* + * Calculate vector multiplication +*/ +Vec Vec::VectMultiplication(const Vec & Other) const{ + double x = GetXYZ().Y()*Other.GetXYZ().Z() - GetXYZ().Z()*Other.GetXYZ().Y(); + double y = GetXYZ().Z()*Other.GetXYZ().X() - GetXYZ().X()*Other.GetXYZ().Z(); + double z = GetXYZ().X()*Other.GetXYZ().Y() - GetXYZ().Y()*Other.GetXYZ().X(); + Vec *aRes = new Vec(x,y,z); + return *aRes; +} + +/*---------------------Class Plane --------------------------------*/ +/*! + * Constructor + */ +Plane::Plane(const Pnt& thePnt1, const Pnt& thePnt2, const Pnt& thePnt3) +{ + CalculatePlane(thePnt1,thePnt2,thePnt3); +} + +/*! + * Destructor + */ +Plane::~Plane() +{} + +/* + * Return plane normale + */ +Vec Plane::GetNormale() const{ + return Vec(myA,myB,myC); +} + +/* + * Calculate A,B,C coeefs of plane + */ +void Plane::CalculatePlane(const Pnt& thePnt1, const Pnt& thePnt2, const Pnt& thePnt3){ + + double x1 = thePnt1.GetXYZ().X(); + double x2 = thePnt2.GetXYZ().X(); + double x3 = thePnt3.GetXYZ().X(); + + double y1 = thePnt1.GetXYZ().Y(); + double y2 = thePnt2.GetXYZ().Y(); + double y3 = thePnt3.GetXYZ().Y(); + + double z1 = thePnt1.GetXYZ().Z(); + double z2 = thePnt2.GetXYZ().Z(); + double z3 = thePnt3.GetXYZ().Z(); + + myA = y1*(z2 - z3) + y2*(z3 - z1) + y3*(z1 - z2); + myB = z1*(x2 - x3) + z2*(x3 - x1) + z3*(x1 - x2); + myC = x1*(y2 - y3) + x2*(y3 - y1) + x3*(y1 - y2); + +#ifdef _MY_DEBUG_ + cout<<"Plane A: "<Update(); +#ifdef _MY_DEBUG_ + cout<<"Need Rotation!!!"<GetPoint(0,coords); + myPnt1 = Pnt(coords[0],coords[1],coords[2]); + aTransformedGrid->GetPoint(1,coords); + myPnt2 = Pnt(coords[0],coords[1],coords[2]); + aTransformedGrid->GetPoint(2,coords); + myPnt3 = Pnt(coords[0],coords[1],coords[2]); + vtkUnstructuredGrid* anArc = BuildArc(); + vtkUnstructuredGrid* anTransArc; + if(needRotation) { + anTransArc = TransformGrid(anArc,aAxis,-anAngle); + anTransArc->Update(); + } + else + anTransArc = anArc; + + myPoints = anTransArc->GetPoints(); + myStatus = Arc_Done; + } + } + else{ +#ifdef _MY_DEBUG_ + std::cout<<"Points lay on the one line !"<Update(); + myPoints = aGrid->GetPoints(); + } +} + +/*! + * Destructor + */ +VTKViewer_ArcBuilder::~VTKViewer_ArcBuilder() +{} + + +/* + * Add to the vtkUnstructuredGrid points from input list + */ +vtkUnstructuredGrid* VTKViewer_ArcBuilder::BuildGrid(const PntList& theList) const +{ + int aListsize = theList.size(); + vtkUnstructuredGrid* aGrid = NULL; + + if(aListsize != 0) { + aGrid = vtkUnstructuredGrid::New(); + vtkPoints* aPoints = vtkPoints::New(); + aPoints->SetNumberOfPoints(aListsize); + + aGrid->Allocate(aListsize); + + PntList::const_iterator it = theList.begin(); + int aCounter = 0; + for(;it != theList.end();it++) { + Pnt aPnt = *it; + aPoints->InsertPoint(aCounter, + aPnt.GetXYZ().X(), + aPnt.GetXYZ().Y(), + aPnt.GetXYZ().Z()); + vtkVertex* aVertex = vtkVertex::New(); + aVertex->GetPointIds()->SetId(0, aCounter); + aGrid->InsertNextCell(aVertex->GetCellType(), aVertex->GetPointIds()); + aCounter++; + aVertex->Delete(); + } + aGrid->SetPoints(aPoints); + aPoints->Delete(); + } + return aGrid; +} + + +vtkUnstructuredGrid* +VTKViewer_ArcBuilder::TransformGrid(vtkUnstructuredGrid* theGrid, + const Vec& theAxis, const double angle) const +{ + vtkTransform *aTransform = vtkTransform::New(); + aTransform->RotateWXYZ(angle, theAxis.GetXYZ().X(), theAxis.GetXYZ().Y(), theAxis.GetXYZ().Z()); + vtkTransformFilter* aTransformFilter = vtkTransformFilter::New(); + aTransformFilter->SetTransform(aTransform); + aTransformFilter->SetInput(theGrid); + aTransform->Delete(); + return aTransformFilter->GetUnstructuredGridOutput(); +} + + +void VTKViewer_ArcBuilder::GetAngle(const double theAngle){ + myAngle = theAngle; +} + +vtkUnstructuredGrid* VTKViewer_ArcBuilder::BuildArc(){ + double x1 = myPnt1.GetXYZ().X(); double x2 = myPnt2.GetXYZ().X(); double x3 = myPnt3.GetXYZ().X(); + double y1 = myPnt1.GetXYZ().Y(); double y2 = myPnt2.GetXYZ().Y(); double y3 = myPnt3.GetXYZ().Y(); + double z = myPnt1.GetXYZ().Z(); //Points on plane || XOY + + + double K1; + double K2; + K1 = (y2 - y1)/(x2 - x1); + K2 = (y3 - y2)/(x3 - x2); +#ifdef _MY_DEBUG_ + std::cout<<"K1 : "<< K1 < infinity + xCenter = (K1*(y1-y3) + (x1+x2))/2.0; + + else if(x2 == x1) //K1 --> infinity + xCenter = (K2*(y1-y3) - (x2+x3))/(-2.0); + + else + xCenter = (K1*K2*(y1-y3) + K2*(x1+x2) - K1*(x2+x3))/ (2.0*(K2-K1)); + + double yCenter; + + if ( K1 == 0 ) + yCenter = (-1/K2)*(xCenter - (x2+x3)/2.0) + (y2 + y3)/2.0; + else + yCenter = (-1/K1)*(xCenter - (x1+x2)/2.0) + (y1 + y2)/2.0; + + double zCenter = z; +#ifdef _MY_DEBUG_ + std::cout<<"xCenter : "< vtkMath::DoublePi()) + aTotalAngle = 2*vtkMath::DoublePi()-aTotalAngle; + */ + + double aTotalAngle = 0; + IncOrder aOrder = GetArcAngle(angle1,angle2,angle3,&aTotalAngle); + + vtkUnstructuredGrid *aC = NULL; + + if(aTotalAngle > aMaxAngle) { + int nbSteps = int(aTotalAngle/aMaxAngle)+1; + double anIncrementAngle = aTotalAngle/nbSteps; + double aCurrentAngle = angle1; + if(aOrder == VTKViewer_ArcBuilder::MINUS) + aCurrentAngle-=anIncrementAngle; + else + aCurrentAngle+=anIncrementAngle; +#ifdef _MY_DEBUG_ + cout<<"Total angle :"<GetCell(cellId); + //Get All points from input cell + aCell->GetPoints()->GetPoint(0,coord); + Pnt P0(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(1,coord); + Pnt P1(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(2,coord); + Pnt P2(coord[0],coord[1],coord[2]); + + VTKViewer_ArcBuilder aBuilder(P0,P2,P1,myMaxArcAngle); + if (aBuilder.GetStatus() != VTKViewer_ArcBuilder::Arc_Done) { + return aResult; + } + else{ + vtkPoints* aPoints = aBuilder.GetPoints(); + vtkIdType aNbPts = aPoints->GetNumberOfPoints(); + aNewPoints = new vtkIdType[aNbPts]; + vtkIdType curID; + vtkIdType aCellType = VTK_POLY_LINE; + + aNewPoints[0] = pts[0]; + for(vtkIdType idx = 1; idx < aNbPts-1;idx++) { + curID = output->GetPoints()->InsertNextPoint(aPoints->GetPoint(idx)); + aNewPoints[idx] = curID; + } + aNewPoints[aNbPts-1] = pts[1]; + + aResult = output->InsertNextCell(aCellType,aNbPts,aNewPoints); + return aResult; + } +} + +/*! + * Add all points from the input vector theCollection into thePoints. + * Array theIds - it is array with ids of added points. + */ +vtkIdType MergevtkPoints(const std::vector& theCollection, vtkPoints* thePoints, vtkIdType* &theIds){ + vtkIdType aNbPoints = 0; + vtkIdType anIdCounter = 0; + vtkIdType aNewPntId = 0; + + //Compute number of points + std::vector::const_iterator it = theCollection.begin(); + for(;it != theCollection.end();it++){ + vtkPoints* aPoints = *it; + if(aPoints) { + aNbPoints += aPoints->GetNumberOfPoints()-1; //Because we add all points except last + } + } + it = theCollection.begin(); + theIds = new vtkIdType[aNbPoints]; + // ..and add all points + for(;it != theCollection.end();it++){ + vtkPoints* aPoints = *it; + + if(aPoints){ + for(vtkIdType idx = 0;idx < aPoints->GetNumberOfPoints()-1;idx++){ + aNewPntId = thePoints->InsertNextPoint(aPoints->GetPoint(idx)); + theIds[anIdCounter] = aNewPntId; + anIdCounter++; + } + } + } + return aNbPoints; +} diff --git a/src/VTKViewer/VTKViewer_ArcBuilder.h b/src/VTKViewer/VTKViewer_ArcBuilder.h new file mode 100644 index 000000000..3ee5acaa1 --- /dev/null +++ b/src/VTKViewer/VTKViewer_ArcBuilder.h @@ -0,0 +1,179 @@ +// Copyright (C) 2007-2008 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. +// +// 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 +// +#ifndef VTKVIEWER_ARCBUILDER_H +#define VTKVIEWER_ARCBUILDER_H + +#include "VTKViewer.h" +#include +#include + +class vtkUnstructuredGrid; +class vtkPoints; +class vtkPolyData; + +class Pnt; + +typedef std::list PntList; + +vtkIdType MergevtkPoints(const std::vector& theCollection, + vtkPoints* thePoints, + vtkIdType* &theIds); + +vtkIdType Build1DArc(vtkIdType cellId, + vtkUnstructuredGrid* input, + vtkPolyData *output, + vtkIdType *pts, + vtkFloatingPointType myMaxArcAngle); + + +/*! + * Class for represenation coordinates X,Y,Z + */ +class XYZ{ + public: + + XYZ(); + XYZ(double , double , double); + ~XYZ(); + + double X()const {return x;} + double Y()const {return y;} + double Z()const {return z;} + + void SetX(const double X) { x=X; } + void SetY(const double Y) { y=Y; } + void SetZ(const double Z) { z=Z; } + + void Coord (double& X, double& Y, double& Z) const { X = x; Y = y; Z = z; } + + double Modulus () const; + + private: + double x; + double y; + double z; +}; + +/*! + Class for the representation point in the 3D space. +*/ +class Pnt{ + public: + Pnt(); + Pnt(double , double , double ); + ~Pnt(); + + void Coord (double& X, double& Y, double& Z) const {coord.Coord(X,Y,Z);} + XYZ GetXYZ() const {return coord;} + + private: + XYZ coord; +}; + +/*! + Class for the representation Vector in the 3D space. +*/ +class Vec{ + public: + + Vec(const double Xv, const double Yv, const double Zv); + ~Vec(); + + XYZ GetXYZ() const {return coord;} + + double AngleBetween(const Vec & Other); + double AngleBetweenInGrad(const Vec & Other); + + Vec VectMultiplication(const Vec & Other) const; + + private: + XYZ coord; +}; + +/*! + Class for the representation plane in the 3D. +*/ +class Plane{ + + public: + Plane(const Pnt& thePnt1, const Pnt& thePnt2, const Pnt& thePnt3); + ~Plane(); + + double A() const {return myA;} + double B() const {return myB;} + double C() const {return myC;} + + Vec GetNormale() const; + + private: + void CalculatePlane(const Pnt& thePnt1, const Pnt& thePnt2, const Pnt& thePnt3); + + private: + double myA; + double myB; + double myC; +}; + + +class VTKViewer_ArcBuilder{ + public: + enum ArcStatus {Arc_Done=0, Arc_Error}; + VTKViewer_ArcBuilder(const Pnt& thePnt1, + const Pnt& thePnt2, + const Pnt& thePnt3, + double theAngle); + + ~VTKViewer_ArcBuilder(); + + Vec GetNormale(); + + ArcStatus GetStatus(){return myStatus;} + + void GetAngle(const double theAngle); + + static double GetPointAngleOnCircle(const double theXCenter, const double theYCenter, + const double theXPoint, const double theYPoint); + + vtkPoints* GetPoints(); + + private: + + enum IncOrder{MINUS=0,PLUS}; + + vtkUnstructuredGrid* BuildGrid(const PntList& theList) const; + vtkUnstructuredGrid* TransformGrid(vtkUnstructuredGrid* theGrid, const Vec& theAxis, const double angle) const; + vtkUnstructuredGrid* BuildArc(); + IncOrder GetArcAngle( const double& P1, const double& P2, const double& P3, double* Ang); + + + + private: + Pnt myPnt1; + Pnt myPnt2; + Pnt myPnt3; + + double myAngle; + ArcStatus myStatus; + vtkPoints* myPoints; +}; + +#endif //VTKVIEWER_ARCBUILDER_H diff --git a/src/VTKViewer/VTKViewer_GeometryFilter.cxx b/src/VTKViewer/VTKViewer_GeometryFilter.cxx index aa3055f2b..7e8d12e8d 100755 --- a/src/VTKViewer/VTKViewer_GeometryFilter.cxx +++ b/src/VTKViewer/VTKViewer_GeometryFilter.cxx @@ -26,6 +26,7 @@ // #include "VTKViewer_GeometryFilter.h" #include "VTKViewer_ConvexTool.h" +#include "VTKViewer_ArcBuilder.h" #include #include @@ -45,6 +46,7 @@ #include #include #include +#include #include #include @@ -58,6 +60,7 @@ #endif #endif +//#define __MYDEBUG__ //#define USE_ROBUST_TRIANGULATION vtkCxxRevisionMacro(VTKViewer_GeometryFilter, "$Revision$"); @@ -67,7 +70,9 @@ VTKViewer_GeometryFilter ::VTKViewer_GeometryFilter(): myShowInside(0), myStoreMapping(0), - myIsWireframeMode(0) + myIsWireframeMode(0), + myIsBuildArc(false), + myMaxArcAngle(2) {} @@ -479,36 +484,54 @@ VTKViewer_GeometryFilter case VTK_QUADRATIC_PYRAMID: if(!myIsWireframeMode){ input->GetCell(cellId,cell); - vtkIdList *pts = vtkIdList::New(); + vtkIdList *lpts = vtkIdList::New(); vtkPoints *coords = vtkPoints::New(); vtkIdList *cellIds = vtkIdList::New(); vtkIdType newCellId; if ( cell->GetCellDimension() == 1 ) { - aCellType = VTK_LINE; - numFacePts = 2; - cell->Triangulate(0,pts,coords); - for (i=0; i < pts->GetNumberOfIds(); i+=2) { - aNewPts[0] = pts->GetId(i); - aNewPts[1] = pts->GetId(i+1); - newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); - if(myStoreMapping) - myVTK2ObjIds.push_back(cellId); - outputCD->CopyData(cd,cellId,newCellId); + vtkIdType arcResult = -1; + if(myIsBuildArc) { + arcResult = Build1DArc(cellId, input, output, pts, myMaxArcAngle); + newCellId = arcResult; } + + if(!myIsBuildArc || arcResult == -1 ) { + aCellType = VTK_LINE; + numFacePts = 2; + cell->Triangulate(0,lpts,coords); + for (i=0; i < lpts->GetNumberOfIds(); i+=2) { + aNewPts[0] = lpts->GetId(i); + aNewPts[1] = lpts->GetId(i+1); + newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + outputCD->CopyData(cd,cellId,newCellId); + } + } + else { + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + outputCD->CopyData(cd,cellId,newCellId); + } } else if ( cell->GetCellDimension() == 2 ) { - aCellType = VTK_TRIANGLE; - numFacePts = 3; - cell->Triangulate(0,pts,coords); - for (i=0; i < pts->GetNumberOfIds(); i+=3) { - aNewPts[0] = pts->GetId(i); - aNewPts[1] = pts->GetId(i+1); - aNewPts[2] = pts->GetId(i+2); - newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); - if(myStoreMapping) - myVTK2ObjIds.push_back(cellId); - outputCD->CopyData(cd,cellId,newCellId); + if(!myIsBuildArc) { + aCellType = VTK_TRIANGLE; + numFacePts = 3; + cell->Triangulate(0,lpts,coords); + for (i=0; i < lpts->GetNumberOfIds(); i+=3) { + aNewPts[0] = lpts->GetId(i); + aNewPts[1] = lpts->GetId(i+1); + aNewPts[2] = lpts->GetId(i+2); + newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + outputCD->CopyData(cd,cellId,newCellId); + } + } + else{ + BuildArcedPolygon(cellId,input,output,true); } } else //3D nonlinear cell @@ -519,11 +542,11 @@ VTKViewer_GeometryFilter vtkCell *face = cell->GetFace(j); input->GetCellNeighbors(cellId, face->PointIds, cellIds); if ( cellIds->GetNumberOfIds() <= 0 || myShowInside ) { - face->Triangulate(0,pts,coords); - for (i=0; i < pts->GetNumberOfIds(); i+=3) { - aNewPts[0] = pts->GetId(i); - aNewPts[1] = pts->GetId(i+1); - aNewPts[2] = pts->GetId(i+2); + face->Triangulate(0,lpts,coords); + for (i=0; i < lpts->GetNumberOfIds(); i+=3) { + aNewPts[0] = lpts->GetId(i); + aNewPts[1] = lpts->GetId(i+1); + aNewPts[2] = lpts->GetId(i+2); newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); if(myStoreMapping) myVTK2ObjIds.push_back(cellId); @@ -531,65 +554,72 @@ VTKViewer_GeometryFilter } } } - } //3d cell + } //3d nonlinear cell cellIds->Delete(); coords->Delete(); - pts->Delete(); + lpts->Delete(); break; }else{ switch(aCellType){ case VTK_QUADRATIC_EDGE: { - aCellType = VTK_POLY_LINE; - numFacePts = 3; + vtkIdType arcResult =-1; + if(myIsBuildArc) { + arcResult = Build1DArc(cellId, input, output, pts,myMaxArcAngle); + newCellId = arcResult; + } + if(!myIsBuildArc || arcResult == -1) { + aCellType = VTK_POLY_LINE; + numFacePts = 3; - aNewPts[0] = pts[0]; - aNewPts[2] = pts[1]; - aNewPts[1] = pts[2]; + aNewPts[0] = pts[0]; + aNewPts[2] = pts[1]; + aNewPts[1] = pts[2]; - newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); - if(myStoreMapping) - myVTK2ObjIds.push_back(cellId); + newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + } + + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); outputCD->CopyData(cd,cellId,newCellId); break; } - case VTK_QUADRATIC_TRIANGLE: { - aCellType = VTK_POLYGON; - numFacePts = 6; - - aNewPts[0] = pts[0]; - aNewPts[1] = pts[3]; - aNewPts[2] = pts[1]; - aNewPts[3] = pts[4]; - aNewPts[4] = pts[2]; - aNewPts[5] = pts[5]; - - newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); - if(myStoreMapping) - myVTK2ObjIds.push_back(cellId); - - outputCD->CopyData(cd,cellId,newCellId); - break; + case VTK_QUADRATIC_TRIANGLE: { + if(!myIsBuildArc) { + aCellType = VTK_POLYGON; + numFacePts = 6; + + aNewPts[0] = pts[0]; + aNewPts[1] = pts[3]; + aNewPts[2] = pts[1]; + aNewPts[3] = pts[4]; + aNewPts[4] = pts[2]; + aNewPts[5] = pts[5]; + + newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + } + else + BuildArcedPolygon(cellId,input,output); + break; } - case VTK_QUADRATIC_QUAD: { - aCellType = VTK_POLYGON; - numFacePts = 8; - - aNewPts[0] = pts[0]; - aNewPts[1] = pts[4]; - aNewPts[2] = pts[1]; - aNewPts[3] = pts[5]; - aNewPts[4] = pts[2]; - aNewPts[5] = pts[6]; - aNewPts[6] = pts[3]; - aNewPts[7] = pts[7]; - - newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); - if(myStoreMapping) - myVTK2ObjIds.push_back(cellId); - - outputCD->CopyData(cd,cellId,newCellId); - break; + case VTK_QUADRATIC_QUAD: { + if(!myIsBuildArc) { + aCellType = VTK_POLYGON; + numFacePts = 8; + + aNewPts[0] = pts[0]; + aNewPts[1] = pts[4]; + aNewPts[2] = pts[1]; + aNewPts[3] = pts[5]; + aNewPts[4] = pts[2]; + aNewPts[5] = pts[6]; + aNewPts[6] = pts[3]; + aNewPts[7] = pts[7]; + newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + } + else + BuildArcedPolygon(cellId,input,output); + break; } case VTK_QUADRATIC_TETRA: { aCellType = VTK_POLYGON; @@ -914,7 +944,7 @@ VTKViewer_GeometryFilter } //switch } //if visible } //for all cells - + output->Squeeze(); vtkDebugMacro(<<"Extracted " << input->GetNumberOfPoints() << " points," @@ -997,3 +1027,264 @@ vtkIdType VTKViewer_GeometryFilter::GetElemObjId( int theVtkID ) return -1; return myVTK2ObjIds[theVtkID]; } + + +void VTKViewer_GeometryFilter::BuildArcedPolygon(vtkIdType cellId, vtkUnstructuredGrid* input, vtkPolyData *output, bool triangulate){ + vtkFloatingPointType coord[3]; + vtkIdType aCellType = VTK_POLYGON; + vtkIdType *aNewPoints; + vtkIdType aNbPoints = 0; + vtkIdType newCellId; + + //Input and output cell data + vtkCellData *cd = input->GetCellData(); + vtkCellData *outputCD = output->GetCellData(); + + vtkCell* aCell = input->GetCell(cellId); + switch(aCell->GetCellType()) { + case VTK_QUADRATIC_TRIANGLE: + { + //Get All points from input cell + aCell->GetPoints()->GetPoint(0,coord); + Pnt P0(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(1,coord); + Pnt P1(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(2,coord); + Pnt P2(coord[0],coord[1],coord[2]); + + aCell->GetPoints()->GetPoint(3,coord); + Pnt P3(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(4,coord); + Pnt P4(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(5,coord); + Pnt P5(coord[0],coord[1],coord[2]); + + //Build arc usinf 0, 3 and 1 points + VTKViewer_ArcBuilder aBuilder1(P0,P3,P1,myMaxArcAngle); +#ifdef __MYDEBUG__ + if(aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done) { + cout<<"Triangle arc 1 done !!!"< aCollection; + aCollection.push_back(aBuilder1.GetPoints()); + aCollection.push_back(aBuilder2.GetPoints()); + aCollection.push_back(aBuilder3.GetPoints()); + + + //----------------------------------------------------------------------------------------- + if(triangulate){ + vtkIdType numFacePts = 3; + vtkIdList *pts = vtkIdList::New(); + vtkPoints *coords = vtkPoints::New(); + aCellType = VTK_TRIANGLE; + vtkIdType aNewPts[numFacePts]; + vtkIdType aTriangleId; + + vtkPolygon *aPlg = vtkPolygon::New(); + aNbPoints = MergevtkPoints(aCollection, aPlg->GetPoints(), aNewPoints); + aPlg->GetPointIds()->SetNumberOfIds(aNbPoints); + + for(vtkIdType i = 0; i < aNbPoints;i++) { + aPlg->GetPointIds()->SetId(i, aNewPoints[i]); + } + + aPlg->Triangulate(0,pts,coords); + + for (vtkIdType i=0; i < pts->GetNumberOfIds(); i+=3) { + aNewPts[0] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i)); + aNewPts[1] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+1)); + aNewPts[2] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+2)); + + aTriangleId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + outputCD->CopyData(cd,cellId,aTriangleId); + } + pts->Delete(); + coords->Delete(); + aPlg->Delete(); + } + //--------------------------------------------------------------------------------------------------- + else { + aNbPoints = MergevtkPoints(aCollection, output->GetPoints(), aNewPoints); + newCellId = output->InsertNextCell(aCellType,aNbPoints,aNewPoints); + outputCD->CopyData(cd,cellId,newCellId); + + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + } + break; + } //VTK_QUADRATIC_TRIANGLE + + case VTK_QUADRATIC_QUAD: + { + //Get All points from input cell + aCell->GetPoints()->GetPoint(0,coord); + Pnt P0(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(1,coord); + Pnt P1(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(2,coord); + Pnt P2(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(3,coord); + Pnt P3(coord[0],coord[1],coord[2]); + + aCell->GetPoints()->GetPoint(4,coord); + Pnt P4(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(5,coord); + Pnt P5(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(6,coord); + Pnt P6(coord[0],coord[1],coord[2]); + aCell->GetPoints()->GetPoint(7,coord); + Pnt P7(coord[0],coord[1],coord[2]); + + //Build arc usinf 0, 4 and 1 points + VTKViewer_ArcBuilder aBuilder1(P0,P4,P1,myMaxArcAngle); +#ifdef __MYDEBUG__ + if(aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done) { + cout<<"Quadrangle arc 1 done !!!"< aCollection; + aCollection.push_back(aBuilder1.GetPoints()); + aCollection.push_back(aBuilder2.GetPoints()); + aCollection.push_back(aBuilder3.GetPoints()); + aCollection.push_back(aBuilder4.GetPoints()); + + //----------------------------------------------------------------------------------------- + if(triangulate){ + vtkIdType numFacePts = 3; + vtkIdList *pts = vtkIdList::New(); + vtkPoints *coords = vtkPoints::New(); + aCellType = VTK_TRIANGLE; + vtkIdType aNewPts[numFacePts]; + vtkIdType aTriangleId; + + vtkPolygon *aPlg = vtkPolygon::New(); + aNbPoints = MergevtkPoints(aCollection, aPlg->GetPoints(), aNewPoints); + aPlg->GetPointIds()->SetNumberOfIds(aNbPoints); + + for(vtkIdType i = 0; i < aNbPoints;i++) { + aPlg->GetPointIds()->SetId(i, aNewPoints[i]); + } + + aPlg->Triangulate(0,pts,coords); + + for (vtkIdType i=0; i < pts->GetNumberOfIds(); i+=3) { + aNewPts[0] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i)); + aNewPts[1] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+1)); + aNewPts[2] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+2)); + + aTriangleId = output->InsertNextCell(aCellType,numFacePts,aNewPts); + + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + outputCD->CopyData(cd,cellId,aTriangleId); + } + pts->Delete(); + coords->Delete(); + aPlg->Delete(); + } + else { + aNbPoints = MergevtkPoints(aCollection, output->GetPoints(), aNewPoints); + newCellId = output->InsertNextCell(aCellType,aNbPoints,aNewPoints); + outputCD->CopyData(cd,cellId,newCellId); + + if(myStoreMapping) + myVTK2ObjIds.push_back(cellId); + } + break; + } //VTK_QUADRATIC_QUAD + + //Unsupported cell type + default: + break; + } //switch + + if (aNewPoints) + delete [] aNewPoints; +} + + +void VTKViewer_GeometryFilter::SetQuadraticArcMode(bool theFlag) +{ + if(myIsBuildArc != theFlag) { + myIsBuildArc = theFlag; + this->Modified(); + } +} +bool VTKViewer_GeometryFilter::GetQuadraticArcMode() const +{ + return myIsBuildArc; +} + +void VTKViewer_GeometryFilter::SetQuadraticArcAngle(vtkFloatingPointType theMaxAngle) +{ + if(myMaxArcAngle != theMaxAngle) { + myMaxArcAngle = theMaxAngle; + this->Modified(); + } +} + +vtkFloatingPointType VTKViewer_GeometryFilter:: GetQuadraticArcAngle() const +{ + return myMaxArcAngle; +} diff --git a/src/VTKViewer/VTKViewer_GeometryFilter.h b/src/VTKViewer/VTKViewer_GeometryFilter.h index b482eb5af..ae50f4af0 100755 --- a/src/VTKViewer/VTKViewer_GeometryFilter.h +++ b/src/VTKViewer/VTKViewer_GeometryFilter.h @@ -32,6 +32,8 @@ #pragma warning ( disable:4251 ) #endif +class vtkUnstructuredGrid; + /*! \brief This class used same as vtkGeometryFilter. See documentation on VTK for more information. */ class VTKVIEWER_EXPORT VTKViewer_GeometryFilter : public vtkGeometryFilter @@ -86,6 +88,13 @@ public: */ virtual vtkIdType GetElemObjId(int theVtkID); + virtual void SetQuadraticArcMode(bool theFlag); + virtual bool GetQuadraticArcMode() const; + + virtual void SetQuadraticArcAngle(vtkFloatingPointType theMaxAngle); + virtual vtkFloatingPointType GetQuadraticArcAngle() const; + + protected: /*! \fn VTKViewer_GeometryFilter(); * \brief Constructor which sets \a myShowInside = 0 and \a myStoreMapping = 0 @@ -104,6 +113,9 @@ protected: * \brief Filter culculation method for data object type is VTK_UNSTRUCTURED_GRID. */ int UnstructuredGridExecute (vtkDataSet *, vtkPolyData *, vtkInformation *); + + + void BuildArcedPolygon(vtkIdType cellId, vtkUnstructuredGrid* input, vtkPolyData *output, bool triangulate = false); private: typedef std::vector TVectorId; @@ -113,6 +125,9 @@ private: int myShowInside; int myStoreMapping; int myIsWireframeMode; + + vtkFloatingPointType myMaxArcAngle; // define max angle for mesh 2D quadratic element in the degrees + bool myIsBuildArc; // flag for representation 2D quadratic element as arked polygon }; #ifdef WIN32 -- 2.39.2