Back to index

salome-paravis  6.5.0
vtkELNOMeshFilter.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 "vtkELNOMeshFilter.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 #include "vtkShrinkFilter.h"
00031 #include "vtkSmartPointer.h"
00032 #include "vtkPointData.h"
00033 #include "vtkCellData.h"
00034 #include "vtkIdList.h"
00035 
00036 vtkCxxRevisionMacro(vtkELNOMeshFilter, "$Revision: 1.1.4.3.12.1 $")
00037 ;
00038 vtkStandardNewMacro(vtkELNOMeshFilter)
00039 ;
00040 
00041 vtkELNOMeshFilter::vtkELNOMeshFilter()
00042 {
00043 }
00044 
00045 vtkELNOMeshFilter::~vtkELNOMeshFilter()
00046 {
00047 }
00048 
00049 int vtkELNOMeshFilter::RequestData(vtkInformation *request,
00050     vtkInformationVector **input, vtkInformationVector *output)
00051 {
00052   vtkUnstructuredGrid *usgIn = vtkUnstructuredGrid::SafeDownCast(
00053       input[0]->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00054 
00055   vtkUnstructuredGrid *usgOut = vtkUnstructuredGrid::SafeDownCast(
00056       output->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00057 
00058   vtkIdTypeArray* usg_offsets = vtkIdTypeArray::SafeDownCast(
00059       this->GetInputArrayToProcess(0, input));
00060 
00061   if(usgIn == NULL || usg_offsets == NULL || usgOut == NULL)
00062     {
00063     vtkDebugMacro("vtkELNOMeshFilter no correctly configured : offsets = " << usg_offsets);
00064     return 1;
00065     }
00066 
00067   vtkIdTypeArray* a = vtkIdTypeArray::SafeDownCast(
00068       usgIn->GetCellData()->GetArray(usg_offsets->GetName()));
00069 
00070   vtkInformationVector *inArrayVec =
00071       this->Information->Get(INPUT_ARRAYS_TO_PROCESS());
00072 
00073   // first shrink the input
00074   vtkUnstructuredGrid* usgInClone = usgIn->NewInstance();
00075 
00076   usgInClone->ShallowCopy(usgIn);
00077 
00078   vtkSmartPointer<vtkShrinkFilter> shrink =
00079       vtkSmartPointer<vtkShrinkFilter>::New();
00080   shrink->SetInput(usgInClone);
00081   shrink->SetShrinkFactor(0.9999);
00082   shrink->Update();
00083   vtkUnstructuredGrid* shrinked = shrink->GetOutput();
00084 
00085   usgInClone->Delete();
00086 
00087   usgOut->ShallowCopy(shrinked);
00088 
00089   vtkIdTypeArray* offsets = vtkIdTypeArray::SafeDownCast(
00090       shrinked->GetCellData()->GetArray(usg_offsets->GetName()));
00091 
00092   // now copy ELNO data. Start by verifying if it is possible to
00093   // shallow copy the array.
00094   vtkInformation *info = offsets->GetInformation();
00095   vtkInformationQuadratureSchemeDefinitionVectorKey *key =
00096       vtkQuadratureSchemeDefinition::DICTIONARY();
00097   if(!key->Has(info))
00098     {
00099     vtkDebugMacro("Dictionary is not present in array " << offsets->GetName() << " " << offsets << " Aborting." );
00100     return 0;
00101     }
00102   int dictSize = key->Size(info);
00103   vtkQuadratureSchemeDefinition **dict =
00104       new vtkQuadratureSchemeDefinition *[dictSize];
00105   key->GetRange(info, dict, 0, 0, dictSize);
00106 
00107   vtkIdType nVerts = shrinked->GetNumberOfPoints();
00108   vtkIdType ncell = usgIn->GetNumberOfCells();
00109   // first loop through all cells to check if a shallow copy is possible
00110   bool shallowok = true;
00111   vtkIdType previous = 0;
00112 
00113   for(vtkIdType cellId = 0; cellId < ncell; cellId++)
00114     {
00115     vtkIdType offset = offsets->GetValue(cellId);
00116     if(offset != previous)
00117       {
00118       shallowok = false;
00119       break;
00120       }
00121     int cellType = usgIn->GetCellType(cellId);
00122 
00123     if(dict[cellType] == NULL)
00124       {
00125       previous = offset;
00126       }
00127     else
00128       {
00129       previous = offset + dict[cellType]->GetNumberOfQuadraturePoints();
00130       }
00131     }
00132   if(previous != nVerts)
00133     {
00134     shallowok = false;
00135     }
00136 
00137   vtkFieldData* fielddata = usgIn->GetFieldData();
00138   for(int index = 0; index < fielddata->GetNumberOfArrays(); index++)
00139     {
00140     vtkDataArray* data = fielddata->GetArray(index);
00141     if(data == NULL)
00142       continue;
00143 
00144     vtkInformation* info = data->GetInformation();
00145     const char* arrayOffsetName = info->Get(
00146         vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME());
00147 
00148     if(arrayOffsetName == NULL ||
00149        strcmp(arrayOffsetName, offsets->GetName()) != 0)
00150       {
00151       shrinked->GetFieldData()->AddArray(data);
00152       continue;
00153       }
00154 
00155     if(shallowok)
00156       {
00157       usgOut->GetPointData()->AddArray(data);
00158       }
00159     else
00160       {
00161       vtkDataArray* newArray = data->NewInstance();
00162       newArray->SetName(data->GetName());
00163       usgOut->GetPointData()->AddArray(newArray);
00164       newArray->SetNumberOfComponents(data->GetNumberOfComponents());
00165       newArray->SetNumberOfTuples(usgOut->GetNumberOfPoints());
00166       newArray->CopyComponentNames(data);
00167       newArray->Delete();
00168       vtkIdList *ids = vtkIdList::New();
00169 
00170       for(vtkIdType cellId = 0; cellId < ncell; cellId++)
00171         {
00172         vtkIdType offset = offsets->GetValue(cellId);
00173         int cellType = shrinked->GetCellType(cellId);
00174         shrinked->GetCellPoints(cellId, ids);
00175         for(int id = 0; id < dict[cellType]->GetNumberOfQuadraturePoints(); id++)
00176           {
00177           const double * w = dict[cellType]->GetShapeFunctionWeights(id);
00178           int j;
00179           for(j = 0; j < dict[cellType]->GetNumberOfNodes(); j++)
00180             {
00181             if(w[j] == 1.0)
00182               break;
00183             }
00184           if(j == dict[cellType]->GetNumberOfNodes())
00185             {
00186             j = id;
00187             }
00188           newArray->SetTuple(ids->GetId(id), offset + j, data);
00189           }
00190         }
00191       ids->FastDelete();
00192       }
00193     }
00194 
00195   delete[] dict;
00196 
00197   return 1;
00198 }
00199 
00200 void vtkELNOMeshFilter::PrintSelf(ostream& os, vtkIndent indent)
00201 {
00202   this->Superclass::PrintSelf(os, indent);
00203 }