Back to index

salome-smesh  6.5.0
SMESH_FaceOrientationFilter.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 #include "SMESH_FaceOrientationFilter.h"
00021 #include "SMESH_ActorUtils.h"
00022 
00023 #include "SUIT_Session.h"
00024 #include "SUIT_ResourceMgr.h"
00025 
00026 #include <VTKViewer_CellCenters.h>
00027 
00028 #include <vtkCellData.h>
00029 #include <vtkDataSet.h>
00030 #include <vtkPolyData.h>
00031 #include <vtkObjectFactory.h>
00032 #include <vtkInformation.h>
00033 #include <vtkInformationVector.h>
00034 
00035 #include <vtkFloatArray.h>
00036 #include <vtkCellArray.h>
00037 #include <vtkMaskPoints.h>
00038 #include <vtkGlyph3D.h>
00039 #include <vtkGlyphSource2D.h>
00040 
00041 #include <QColor>
00042 
00043 #define PI   3.14159265359
00044 
00045 vtkCxxRevisionMacro(SMESH_FaceOrientationFilter, "$Revision: 1.2.2.2.4.3.8.1 $");
00046 vtkStandardNewMacro(SMESH_FaceOrientationFilter);
00047 
00053 SMESH_FaceOrientationFilter::SMESH_FaceOrientationFilter()
00054 {
00055   SUIT_ResourceMgr* mgr = SUIT_Session::session()->resourceMgr();
00056   myOrientationScale = mgr->doubleValue( "SMESH", "orientation_scale", 0.1 );
00057   my3dVectors = mgr->booleanValue( "SMESH", "orientation_3d_vectors", false );
00058 
00059   myArrowPolyData = CreateArrowPolyData();
00060 
00061   myFacePolyData = vtkPolyData::New();
00062 
00063   myFaceCenters = VTKViewer_CellCenters::New();
00064   myFaceCenters->SetInput(myFacePolyData);
00065 
00066   myFaceMaskPoints = vtkMaskPoints::New();
00067   myFaceMaskPoints->SetInput(myFaceCenters->GetOutput());
00068   myFaceMaskPoints->SetOnRatio(1);
00069 
00070   myGlyphSource = vtkGlyphSource2D::New();
00071   myGlyphSource->SetGlyphTypeToThickArrow();
00072   myGlyphSource->SetFilled(0);
00073   myGlyphSource->SetCenter(0.5, 0.0, 0.0);
00074 
00075   myBaseGlyph = vtkGlyph3D::New();
00076   myBaseGlyph->SetInput(myFaceMaskPoints->GetOutput());
00077   myBaseGlyph->SetVectorModeToUseVector();
00078   myBaseGlyph->SetScaleModeToDataScalingOff();
00079   myBaseGlyph->SetColorModeToColorByScalar();
00080   myBaseGlyph->SetSource(my3dVectors ? myArrowPolyData : myGlyphSource->GetOutput());
00081 }
00082 
00083 SMESH_FaceOrientationFilter::~SMESH_FaceOrientationFilter()
00084 {
00085   myArrowPolyData->Delete();
00086   myFacePolyData->Delete();
00087   myFaceCenters->Delete();
00088   myFaceMaskPoints->Delete();
00089   myGlyphSource->Delete();
00090   myBaseGlyph->Delete();
00091 }
00092 
00093 void SMESH_FaceOrientationFilter::SetOrientationScale( vtkFloatingPointType theScale )
00094 {
00095   myOrientationScale = theScale;
00096   Modified();
00097 }
00098 
00099 void SMESH_FaceOrientationFilter::Set3dVectors( bool theState )
00100 {
00101   my3dVectors = theState;
00102   myBaseGlyph->SetSource(my3dVectors ? myArrowPolyData : myGlyphSource->GetOutput());
00103   Modified();
00104 }
00105 
00106 vtkPolyData* SMESH_FaceOrientationFilter::CreateArrowPolyData()
00107 {
00108   vtkPoints* points = vtkPoints::New();
00109   vtkCellArray* polys = vtkCellArray::New();
00110 
00111   float l1 = 0.8;
00112   float l2 = 1.0;
00113   int n = 16;
00114   float r1 = 0.04;
00115   float r2 = 0.08;
00116   float angle = 2. * PI / n;
00117   float p[3];
00118   vtkIdType c3[3];
00119   vtkIdType c4[4];
00120 
00121   float p0[3] = { 0.0, 0.0, 0.0 };
00122   float p1[3] = {  l1, 0.0, 0.0 };
00123   float p2[3] = {  l2, 0.0, 0.0 };
00124 
00125   points->InsertPoint( 0, p0 );
00126   points->InsertPoint( 1, p1 );
00127   points->InsertPoint( 2, p2 );
00128 
00129   // shaft
00130   for( int i = 0; i < n; i++ )
00131   {
00132     p[0] = 0;
00133     p[1] = r1 * sin( i * angle );
00134     p[2] = r1 * cos( i * angle );
00135     points->InsertPoint( i + 3, p );
00136 
00137     p[0] = l1;
00138     points->InsertPoint( i + 3 + n, p );
00139   }
00140 
00141   // insert the last cells outside a loop
00142   {
00143     c3[0] = 0;
00144     c3[1] = 3;
00145     c3[2] = 3 + n - 1;
00146     polys->InsertNextCell( 3, c3 );
00147 
00148     c4[0] = 3;
00149     c4[1] = 3 + n - 1;
00150     c4[2] = 3 + 2 * n - 1;
00151     c4[3] = 3 + n;
00152     polys->InsertNextCell( 4, c4 );
00153   }
00154   for( int i = 0; i < n - 1; i++ )
00155   {
00156     c3[0] = 0;
00157     c3[1] = i + 3;
00158     c3[2] = i + 4;
00159     polys->InsertNextCell( 3, c3 );
00160 
00161     c4[0] = i + 3;
00162     c4[1] = i + 4;
00163     c4[2] = i + 4 + n;
00164     c4[3] = i + 3 + n;
00165     polys->InsertNextCell( 4, c4 );
00166   }
00167 
00168   // cone
00169   for( int i = 0; i < n; i++ )
00170   {
00171     p[0] = l1;
00172     p[1] = r2 * sin( i * angle );
00173     p[2] = r2 * cos( i * angle );
00174     points->InsertPoint( i + 3 + 2 * n, p );
00175   }
00176 
00177   // insert the last cells outside a loop
00178   {
00179     c3[0] = 1;
00180     c3[1] = 3 + 2 * n;
00181     c3[2] = 3 + 2 * n + n - 1;
00182     polys->InsertNextCell( 3, c3 );
00183 
00184     c3[0] = 2;
00185     polys->InsertNextCell( 3, c3 );
00186   }
00187   for( int i = 0; i < n - 1; i++ )
00188   {
00189     c3[0] = 1;
00190     c3[1] = 3 + i + 2 * n;
00191     c3[2] = 3 + i + 2 * n + 1;
00192     polys->InsertNextCell( 3, c3 );
00193 
00194     c3[0] = 2;
00195     polys->InsertNextCell( 3, c3 );
00196   }
00197 
00198   vtkPolyData* aPolyData = vtkPolyData::New();
00199 
00200   aPolyData->SetPoints(points);
00201   points->Delete();
00202 
00203   aPolyData->SetPolys(polys);
00204   polys->Delete();
00205 
00206   return aPolyData;
00207 }
00208 
00209 void GetFaceParams( vtkCell* theFace, double theNormal[3], double& theSize ) 
00210 {
00211   vtkPoints* aPoints = theFace->GetPoints();
00212 
00213   // here we get first 3 points from the face and calculate the normal as a cross-product of vectors
00214   double x0 = aPoints->GetPoint(0)[0], y0 = aPoints->GetPoint(0)[1], z0 = aPoints->GetPoint(0)[2];
00215   double x1 = aPoints->GetPoint(1)[0], y1 = aPoints->GetPoint(1)[1], z1 = aPoints->GetPoint(1)[2];
00216   double x2 = aPoints->GetPoint(2)[0], y2 = aPoints->GetPoint(2)[1], z2 = aPoints->GetPoint(2)[2];
00217 
00218   theNormal[0] = ( y1 - y0 ) * ( z2 - z0 ) - ( z1 - z0 ) * ( y2 - y0 );
00219   theNormal[1] = ( z1 - z0 ) * ( x2 - x0 ) - ( x1 - x0 ) * ( z2 - z0 );
00220   theNormal[2] = ( x1 - x0 ) * ( y2 - y0 ) - ( y1 - y0 ) * ( x2 - x0 );
00221 
00222   double* aBounds = theFace->GetBounds();
00223   theSize = pow( pow( aBounds[1] - aBounds[0], 2 ) +
00224                  pow( aBounds[3] - aBounds[2], 2 ) +
00225                  pow( aBounds[5] - aBounds[4], 2 ), 0.5 );
00226 }
00227 
00231 int SMESH_FaceOrientationFilter::RequestData(
00232   vtkInformation *request,
00233   vtkInformationVector **inputVector,
00234   vtkInformationVector *outputVector)
00235 {
00236   // get the info objects
00237   vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
00238   vtkInformation *outInfo = outputVector->GetInformationObject(0);
00239 
00240   // get the input and ouptut
00241   vtkDataSet *input = vtkDataSet::SafeDownCast(
00242     inInfo->Get(vtkDataObject::DATA_OBJECT()));
00243   vtkPolyData *output = vtkPolyData::SafeDownCast(
00244     outInfo->Get(vtkDataObject::DATA_OBJECT()));
00245 
00246   myFacePolyData->Initialize();
00247   myFacePolyData->ShallowCopy(input);
00248 
00249   vtkCellArray* aFaces = vtkCellArray::New();
00250 
00251   vtkFloatArray* aVectors = vtkFloatArray::New();
00252   aVectors->SetNumberOfComponents(3);
00253 
00254   int anAllFaces = 0;
00255   double anAverageSize = 0;
00256 
00257   vtkIdList* aNeighborIds = vtkIdList::New();
00258 
00259   for(int aCellId = 0, aNbCells = input->GetNumberOfCells(); aCellId < aNbCells; aCellId++)
00260   {
00261     vtkCell* aCell = input->GetCell(aCellId);
00262 
00263     if( aCell->GetNumberOfFaces() == 0 && aCell->GetNumberOfPoints() > 2 ) // cell is a face
00264     {
00265       double aSize, aNormal[3];
00266       GetFaceParams( aCell, aNormal, aSize );
00267 
00268       aFaces->InsertNextCell(aCell);
00269       aVectors->InsertNextTuple(aNormal);
00270 
00271       anAllFaces++;
00272       anAverageSize += aSize;
00273 
00274       continue;
00275     }
00276 
00277     for(int aFaceId = 0, aNbFaces = aCell->GetNumberOfFaces(); aFaceId < aNbFaces; aFaceId++)
00278     {
00279       vtkCell* aFace = aCell->GetFace(aFaceId);
00280 
00281       input->GetCellNeighbors( aCellId, aFace->PointIds, aNeighborIds );
00282       if( aNeighborIds->GetNumberOfIds() > 0 )
00283         continue;
00284 
00285       double aSize, aNormal[3];
00286       GetFaceParams( aFace, aNormal, aSize );
00287 
00288       aFaces->InsertNextCell(aFace->GetPointIds());
00289       aVectors->InsertNextTuple(aNormal);
00290 
00291       anAllFaces++;
00292       anAverageSize += aSize;
00293     }
00294   }
00295   aNeighborIds->Delete();
00296 
00297   myFacePolyData->SetPolys(aFaces);
00298   aFaces->Delete();
00299 
00300   myFacePolyData->GetCellData()->SetScalars(0);
00301   myFacePolyData->GetCellData()->SetVectors(aVectors);
00302   aVectors->Delete();
00303 
00304   if( anAllFaces == 0 )
00305     return 0;
00306 
00307   anAverageSize /= anAllFaces;
00308   anAverageSize *= myOrientationScale;
00309 
00310   myBaseGlyph->SetScaleFactor( anAverageSize );
00311   myBaseGlyph->Update();
00312 
00313   output->ShallowCopy( myBaseGlyph->GetOutput() );
00314 
00315   return 1;
00316 }
00317 
00318 int SMESH_FaceOrientationFilter::FillInputPortInformation(int, vtkInformation *info)
00319 {
00320   info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataSet");
00321   return 1;
00322 }