Back to index

salome-paravis  6.5.0
vtkELNOFilter.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2010-2012  CEA/DEN, EDF R&D
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 "vtkELNOFilter.h"
00021 #include "vtkInformation.h"
00022 #include "vtkInformationVector.h"
00023 #include "vtkObjectFactory.h"
00024 #include "vtkPolyDataAlgorithm.h"
00025 #include "vtkPolyData.h"
00026 #include "vtkIdTypeArray.h"
00027 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
00028 #include "vtkQuadratureSchemeDefinition.h"
00029 #include "vtkUnstructuredGrid.h"
00030 
00031 vtkCxxRevisionMacro(vtkELNOFilter, "$Revision: 1.1.4.3.12.1 $");
00032 vtkStandardNewMacro(vtkELNOFilter);
00033 
00034 vtkELNOFilter::vtkELNOFilter()
00035 {
00036   this->ShrinkFactor = 0.5;
00037 }
00038 
00039 vtkELNOFilter::~vtkELNOFilter()
00040 {
00041 }
00042 
00043 int vtkELNOFilter::RequestData(vtkInformation *request,
00044     vtkInformationVector **input, vtkInformationVector *output)
00045 {
00046   vtkUnstructuredGrid *usgIn = vtkUnstructuredGrid::SafeDownCast(
00047       input[0]->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00048 
00049   vtkPolyData *pdOut = vtkPolyData::SafeDownCast(
00050       output->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00051 
00052   vtkDataArray* array = this->GetInputArrayToProcess(0, input);
00053   vtkIdTypeArray* offsets = vtkIdTypeArray::SafeDownCast(
00054       this->GetInputArrayToProcess(0, input));
00055 
00056   if(usgIn == NULL || offsets == NULL || pdOut == NULL)
00057     {
00058     vtkDebugMacro("vtkELNOFilter no correctly configured : offsets = " << offsets);
00059     return 1;
00060     }
00061 
00062   vtkInformation *info = offsets->GetInformation();
00063   vtkInformationQuadratureSchemeDefinitionVectorKey *key =
00064       vtkQuadratureSchemeDefinition::DICTIONARY();
00065   if(!key->Has(info))
00066     {
00067     vtkDebugMacro("Dictionary is not present in array " << offsets->GetName() << " " << offsets << " Aborting." );
00068     return 1;
00069     }
00070 
00071   int res = this->Superclass::RequestData(request, input, output);
00072   if(res == 0)
00073     {
00074     return 0;
00075     }
00076 
00077   int dictSize = key->Size(info);
00078   vtkQuadratureSchemeDefinition **dict =
00079       new vtkQuadratureSchemeDefinition *[dictSize];
00080   key->GetRange(info, dict, 0, 0, dictSize);
00081 
00082   vtkIdType ncell = usgIn->GetNumberOfCells();
00083   vtkPoints *points = pdOut->GetPoints();
00084   vtkIdType start = 0;
00085   for(vtkIdType cellId = 0; cellId < ncell; cellId++)
00086     {
00087     vtkIdType offset = offsets->GetValue(cellId);
00088     int cellType = usgIn->GetCellType(cellId);
00089     // a simple check to see if a scheme really exists for this cell type.
00090     // should not happen if the cell type has not been modified.
00091     if(dict[cellType] == NULL)
00092       continue;
00093     int np = dict[cellType]->GetNumberOfQuadraturePoints();
00094     double center[3] = {0, 0, 0};
00095     for(int id = start; id < start + np; id++)
00096       {
00097       double *position = points->GetPoint(id);
00098       center[0] += position[0];
00099       center[1] += position[1];
00100       center[2] += position[2];
00101       }
00102     center[0] /= np;
00103     center[1] /= np;
00104     center[2] /= np;
00105     for(int id = start; id < start + np; id++)
00106       {
00107       double *position = points->GetPoint(id);
00108       double newpos[3];
00109       newpos[0] = position[0] * this->ShrinkFactor + center[0] * (1
00110           - this->ShrinkFactor);
00111       newpos[1] = position[1] * this->ShrinkFactor + center[1] * (1
00112           - this->ShrinkFactor);
00113       newpos[2] = position[2] * this->ShrinkFactor + center[2] * (1
00114           - this->ShrinkFactor);
00115       points->SetPoint(id, newpos);
00116       }
00117     start += np;
00118     }
00119 
00120   return 1;
00121 }
00122 
00123 void vtkELNOFilter::PrintSelf(ostream& os, vtkIndent indent)
00124 {
00125   this->Superclass::PrintSelf(os, indent);
00126 
00127   os << indent << "ShrinkFactor : " << this->ShrinkFactor << endl;
00128 }