From: ouv Date: Mon, 20 Sep 2010 12:58:06 +0000 (+0000) Subject: Issue 0020226: [CEA 329] Invalid glyphs position in vector fields on polyhedrons... X-Git-Tag: V5_1_5a1~6 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=17b2585cee3ff488988ead537679056d75d2fc5e;p=modules%2Fgui.git Issue 0020226: [CEA 329] Invalid glyphs position in vector fields on polyhedrons cell. --- diff --git a/src/VTKViewer/Makefile.am b/src/VTKViewer/Makefile.am index b92eb34b8..56da23db4 100755 --- a/src/VTKViewer/Makefile.am +++ b/src/VTKViewer/Makefile.am @@ -59,7 +59,8 @@ salomeinclude_HEADERS = \ VTKViewer_MarkerWidget.h \ VTKViewer_MarkerDlg.h \ VTKViewer_PolyDataMapper.h \ - VTKViewer_DataSetMapper.h + VTKViewer_DataSetMapper.h \ + VTKViewer_CellCenters.h dist_libVTKViewer_la_SOURCES = \ VTKViewer_CellLocationsArray.cxx \ @@ -87,7 +88,8 @@ dist_libVTKViewer_la_SOURCES = \ VTKViewer_MarkerWidget.cxx \ VTKViewer_MarkerDlg.cxx \ VTKViewer_PolyDataMapper.cxx \ - VTKViewer_DataSetMapper.cxx + VTKViewer_DataSetMapper.cxx \ + VTKViewer_CellCenters.cxx MOC_FILES = \ VTKViewer_RenderWindow_moc.cxx \ diff --git a/src/VTKViewer/VTKViewer_CellCenters.cxx b/src/VTKViewer/VTKViewer_CellCenters.cxx new file mode 100644 index 000000000..91c998856 --- /dev/null +++ b/src/VTKViewer/VTKViewer_CellCenters.cxx @@ -0,0 +1,184 @@ +// Copyright (C) 2007-2010 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 +// + +#include "VTKViewer_CellCenters.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +vtkCxxRevisionMacro(VTKViewer_CellCenters, "$Revision$"); +vtkStandardNewMacro(VTKViewer_CellCenters); + +/*! + * Class : VTKViewer_CellCenters + * Description : Filter computing geometrical centers of given cells + * (differs from native vtk filter by small fix for VTK_CONVEX_POINT_SET cells) + */ + +/*! + Constructor +*/ +VTKViewer_CellCenters::VTKViewer_CellCenters() +{ +} + +/*! + Redefined main method +*/ +int VTKViewer_CellCenters::RequestData( + vtkInformation *vtkNotUsed(request), + vtkInformationVector **inputVector, + vtkInformationVector *outputVector) +{ + // get the info objects + vtkInformation *inInfo = inputVector[0]->GetInformationObject(0); + vtkInformation *outInfo = outputVector->GetInformationObject(0); + + // get the input and ouptut + vtkDataSet *input = vtkDataSet::SafeDownCast( + inInfo->Get(vtkDataObject::DATA_OBJECT())); + vtkPolyData *output = vtkPolyData::SafeDownCast( + outInfo->Get(vtkDataObject::DATA_OBJECT())); + + vtkIdType cellId, numCells; + int subId; + vtkCellData *inCD; + vtkPointData *outPD; + vtkPoints *newPts; + vtkCell *cell; + double x[3], pcoords[3]; + double *weights; + + inCD=input->GetCellData(); + outPD=output->GetPointData(); + + if ( (numCells = input->GetNumberOfCells()) < 1 ) + { + vtkWarningMacro(<<"No cells to generate center points for"); + return 1; + } + + newPts = vtkPoints::New(); + newPts->SetNumberOfPoints(numCells); + weights = new double [input->GetMaxCellSize()]; + + int abort=0; + vtkIdType progressInterval = numCells/10 + 1; + int hasEmptyCells = 0; + for (cellId=0; cellId < numCells && !abort; cellId++) + { + if ( ! (cellId % progressInterval) ) + { + vtkDebugMacro(<<"Processing #" << cellId); + this->UpdateProgress (0.5*cellId/numCells); + abort = this->GetAbortExecute(); + } + + cell = input->GetCell(cellId); + if (cell->GetCellType() != VTK_EMPTY_CELL) + { + // fix for VTK_CONVEX_POINT_SET cells + if (cell->GetCellType() == VTK_CONVEX_POINT_SET ) + { + x[0] = x[1] = x[2] = 0; + vtkPoints* aPoints = cell->GetPoints(); + int aNbPoints = aPoints->GetNumberOfPoints(); + for( int i = 0; i < aNbPoints; i++ ) + { + double aCoord[3]; + aPoints->GetPoint( i, aCoord ); + x[0] += aCoord[0]; + x[1] += aCoord[1]; + x[2] += aCoord[2]; + } + x[0] /= aNbPoints; + x[1] /= aNbPoints; + x[2] /= aNbPoints; + } + else + { + subId = cell->GetParametricCenter(pcoords); + cell->EvaluateLocation(subId, pcoords, x, weights); + } + newPts->SetPoint(cellId,x); + } + else + { + hasEmptyCells = 1; + } + } + + if ( this->VertexCells ) + { + vtkIdType pts[1]; + vtkCellData *outCD=output->GetCellData(); + vtkCellArray *verts = vtkCellArray::New(); + verts->Allocate(verts->EstimateSize(1,numCells),1); + + for (cellId=0; cellId < numCells && !abort; cellId++) + { + if ( ! (cellId % progressInterval) ) + { + vtkDebugMacro(<<"Processing #" << cellId); + this->UpdateProgress (0.5+0.5*cellId/numCells); + abort = this->GetAbortExecute(); + } + + cell = input->GetCell(cellId); + if (cell->GetCellType() != VTK_EMPTY_CELL) + { + pts[0] = cellId; + verts->InsertNextCell(1,pts); + } + } + + output->SetVerts(verts); + verts->Delete(); + if (!hasEmptyCells) + { + outCD->PassData(inCD); //only if verts are generated + } + } + + // clean up and update output + output->SetPoints(newPts); + newPts->Delete(); + + if (!hasEmptyCells) + { + outPD->PassData(inCD); //because number of points = number of cells + } + if (weights) + { + delete [] weights; + } + + return 1; +} diff --git a/src/VTKViewer/VTKViewer_CellCenters.h b/src/VTKViewer/VTKViewer_CellCenters.h new file mode 100644 index 000000000..7057d7961 --- /dev/null +++ b/src/VTKViewer/VTKViewer_CellCenters.h @@ -0,0 +1,60 @@ +// Copyright (C) 2007-2010 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_CELLCENTERS_H +#define VTKVIEWER_CELLCENTERS_H + +#include "VTKViewer.h" + +#include + +#ifdef WIN32 +#pragma warning ( disable:4251 ) +#endif + +/*! + * Class : VTKViewer_CellCenters + * Description : Filter computing geometrical centers of given cells + * (differs from native vtk filter by small fix for VTK_CONVEX_POINT_SET cells) + */ +class VTKVIEWER_EXPORT VTKViewer_CellCenters : public vtkCellCenters +{ +public: + vtkTypeRevisionMacro(VTKViewer_CellCenters,vtkCellCenters); + + static VTKViewer_CellCenters *New(); + +protected: + VTKViewer_CellCenters(); + + virtual int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *); + +private: + VTKViewer_CellCenters(const VTKViewer_CellCenters&); // Not implemented. + void operator=(const VTKViewer_CellCenters&); // Not implemented. +}; + +#ifdef WIN32 +#pragma warning ( default:4251 ) +#endif + +#endif