Back to index

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