Back to index

salome-gui  6.5.0
VTKViewer_ExtractUnstructuredGrid.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
00004 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
00005 //
00006 // This library is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 2.1 of the License.
00010 //
00011 // This library is distributed in the hope that it will be useful,
00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014 // Lesser General Public License for more details.
00015 //
00016 // You should have received a copy of the GNU Lesser General Public
00017 // License along with this library; if not, write to the Free Software
00018 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00019 //
00020 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00021 //
00022 
00023 // File:    VTKViewer_ExtractUnstructuredGrid.cxx
00024 // Author:  Alexey PETROV
00025 
00026 #include "VTKViewer_ExtractUnstructuredGrid.h"
00027 #include "VTKViewer_CellLocationsArray.h"
00028 
00029 #include <vtkUnsignedCharArray.h>
00030 #include <vtkUnstructuredGrid.h>
00031 #include <vtkObjectFactory.h>
00032 #include <vtkCellArray.h>
00033 #include <vtkIdList.h>
00034 #include <vtkCell.h>
00035 #include <vtkInformation.h>
00036 #include <vtkInformationVector.h>
00037 #include <vtkVersion.h>
00038 
00039 #include "utilities.h"
00040 
00041 #if defined __GNUC__
00042   #if __GNUC__ == 2
00043     #define __GNUC_2__
00044   #endif
00045 #endif
00046 
00047 #define VTK_XVERSION (VTK_MAJOR_VERSION*10000+VTK_MINOR_VERSION*100+VTK_BUILD_VERSION)
00048 
00049 vtkStandardNewMacro(VTKViewer_ExtractUnstructuredGrid);
00050 
00051 
00052 VTKViewer_ExtractUnstructuredGrid::VTKViewer_ExtractUnstructuredGrid():
00053   myExtractionMode(eCells), myChangeMode(ePassAll)
00054 {}
00055 
00056 
00057 VTKViewer_ExtractUnstructuredGrid::~VTKViewer_ExtractUnstructuredGrid(){}
00058 
00059 
00060 void VTKViewer_ExtractUnstructuredGrid::RegisterCell(vtkIdType theCellId){
00061 //  if(0 && MYDEBUG) MESSAGE("RegisterCell - theCellId = "<<theCellId);
00062   myCellIds.insert(theCellId);
00063   Modified();
00064 }
00065 
00066 
00067 void VTKViewer_ExtractUnstructuredGrid::RegisterCellsWithType(vtkIdType theCellType){
00068 //  if(0 && MYDEBUG) MESSAGE("RegisterCellsWithType - theCellType = "<<theCellType);
00069   myCellTypes.insert(theCellType);
00070   //MESSAGE("myCellTypes.insert " << theCellType);
00071   Modified();
00072 }
00073 
00074 
00075 void VTKViewer_ExtractUnstructuredGrid::SetStoreMapping(int theStoreMapping){
00076   myStoreMapping = theStoreMapping != 0;
00077   this->Modified();
00078 }
00079 
00080 vtkIdType VTKViewer_ExtractUnstructuredGrid::GetInputId(int theOutId) const
00081 {
00082   if ( myCellIds.empty() && myCellTypes.empty() )
00083     return theOutId;
00084 
00085   if ( theOutId<0 || theOutId >= (int)myOut2InId.size() )
00086     return -1;
00087   return myOut2InId[theOutId];
00088 }
00089 
00090 vtkIdType VTKViewer_ExtractUnstructuredGrid::GetOutputId(int theInId) const{
00091   if(myCellIds.empty() && myCellTypes.empty()) return theInId;
00092   TMapId::const_iterator anIter = myIn2OutId.find(theInId);
00093   if(anIter == myIn2OutId.end()) return -1;
00094   return anIter->second;
00095 }
00096 
00097 
00098 inline void InsertCell(vtkUnstructuredGrid *theInput,
00099                        vtkCellArray *theConnectivity, 
00100                        vtkUnsignedCharArray* theCellTypesArray,
00101                        vtkIdTypeArray*& theFaces,
00102                        vtkIdTypeArray*& theFaceLocations,
00103                        vtkIdType theCellId, 
00104                        vtkIdList *theIdList,
00105                        bool theStoreMapping,
00106                        vtkIdType theOutId, 
00107                        VTKViewer_ExtractUnstructuredGrid::TVectorId& theOut2InId,
00108                        VTKViewer_ExtractUnstructuredGrid::TMapId& theIn2OutId)
00109 {
00110   vtkCell *aCell = theInput->GetCell(theCellId);
00111   vtkIdList *aPntIds = aCell->GetPointIds();
00112   vtkIdType aNbIds = aPntIds->GetNumberOfIds();
00113   theIdList->SetNumberOfIds(aNbIds);
00114   for(vtkIdType i = 0; i < aNbIds; i++){
00115     theIdList->SetId(i,aPntIds->GetId(i));
00116   }
00117   vtkIdType aCellType = aCell->GetCellType();
00118 #if VTK_XVERSION > 50700
00119   if (aCellType != VTK_POLYHEDRON)
00120     {
00121 #endif
00122       theConnectivity->InsertNextCell(theIdList);
00123       if (theFaceLocations)
00124         theFaceLocations->InsertNextValue(-1);
00125 #if VTK_XVERSION > 50700
00126     }
00127   else
00128     {
00129       //MESSAGE("InsertCell type VTK_POLYHEDRON " << theStoreMapping);
00130       if (!theFaces)
00131         {
00132           theFaces = vtkIdTypeArray::New();
00133           theFaces->Allocate(theCellTypesArray->GetSize());
00134           theFaceLocations = vtkIdTypeArray::New();
00135           theFaceLocations->Allocate(theCellTypesArray->GetSize());
00136           // FaceLocations must be padded until the current position
00137           for (vtkIdType i = 0; i <= theCellTypesArray->GetMaxId(); i++)
00138             {
00139               theFaceLocations->InsertNextValue(-1);
00140             }
00141         }
00142       // insert face location
00143       theFaceLocations->InsertNextValue(theFaces->GetMaxId() + 1);
00144 
00145       // insert cell connectivity and faces stream
00146       vtkIdType nfaces;
00147       vtkIdType* face;
00148       vtkIdType realnpts;
00149       theInput->GetFaceStream(theCellId, nfaces, face);
00150       vtkUnstructuredGrid::DecomposeAPolyhedronCell(
00151           nfaces, face, realnpts, theConnectivity, theFaces);
00152     }
00153 #endif
00154 
00155   theCellTypesArray->InsertNextValue(aCellType);
00156   if(theStoreMapping){
00157     theOut2InId.push_back(theCellId);
00158     theIn2OutId[theCellId] = theOutId;
00159   }
00160 }
00161 
00162 inline void InsertPointCell(vtkCellArray *theConnectivity, 
00163                             vtkUnsignedCharArray* theCellTypesArray,
00164                             vtkIdType theCellId,
00165                             vtkIdList *theIdList,
00166                             bool theStoreMapping,
00167                             vtkIdType theOutId, 
00168                             VTKViewer_ExtractUnstructuredGrid::TVectorId& theOut2InId,
00169                             VTKViewer_ExtractUnstructuredGrid::TMapId& theIn2OutId)
00170 {
00171   theIdList->SetId(0,theCellId);
00172   theConnectivity->InsertNextCell(theIdList);
00173   theCellTypesArray->InsertNextValue(VTK_VERTEX);
00174   if(theStoreMapping){
00175     theOut2InId.push_back(theCellId);
00176     theIn2OutId[theCellId] = theOutId;
00177   }
00178 }
00179 
00180 
00181 // int VTKViewer_ExtractUnstructuredGrid::RequestData(
00182 //   vtkInformation *vtkNotUsed(request),
00183 //   vtkInformationVector **inputVector,
00184 //   vtkInformationVector *outputVector)
00185 void VTKViewer_ExtractUnstructuredGrid::Execute()
00186 {
00187   /*
00188   not ported yet to the new executive-based pipeline architecture.
00189 
00190   // get the info objects
00191   vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
00192   vtkInformation *outInfo = outputVector->GetInformationObject(0);
00193 
00194   // get the input and ouptut
00195   vtkUnstructuredGrid *anInput = vtkUnstructuredGrid::SafeDownCast(
00196     inInfo->Get(vtkDataObject::DATA_OBJECT()));
00197   vtkUnstructuredGrid *anOutput = vtkUnstructuredGrid::SafeDownCast(
00198     outInfo->Get(vtkDataObject::DATA_OBJECT()));
00199   */
00200   vtkUnstructuredGrid *anInput = this->GetInput();
00201   vtkUnstructuredGrid *anOutput = this->GetOutput();
00202   
00203   myOut2InId.clear();  myIn2OutId.clear();
00204 
00205 /*  if(MYDEBUG){
00206     MESSAGE("Execute - anInput->GetNumberOfCells() = "<<anInput->GetNumberOfCells());
00207     MESSAGE("Execute - myCellTypes.size() = "<<myCellTypes.size());
00208     MESSAGE("Execute - myCellIds.size() = "<<myCellIds.size());
00209     MESSAGE("Execute - myExtractionMode = "<<myExtractionMode);
00210     MESSAGE("Execute - myChangeMode = "<<myChangeMode);
00211   }*/
00212   if(myExtractionMode == eCells){
00213     if(myChangeMode == ePassAll || (myCellIds.empty() && myCellTypes.empty() && myChangeMode == eRemoving)){
00214       if(vtkIdType aNbElems = anInput->GetNumberOfCells()){
00215         if(myStoreMapping) myOut2InId.reserve(aNbElems);
00216         anOutput->ShallowCopy(anInput);
00217         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00218           if(myStoreMapping){
00219             myOut2InId.push_back(aCellId);
00220             myIn2OutId[aCellId] = anOutId;
00221           }
00222         }
00223       }
00224     }else{
00225       vtkIdList *anIdList = vtkIdList::New();
00226       vtkCellArray *aConnectivity = vtkCellArray::New();
00227       vtkIdType aNbElems = anInput->GetNumberOfCells();
00228       aConnectivity->Allocate(2*aNbElems,0);
00229       vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
00230       aCellTypesArray->SetNumberOfComponents(1);
00231       aCellTypesArray->Allocate(aNbElems*aCellTypesArray->GetNumberOfComponents());
00232 
00233       vtkIdTypeArray *newFaces = 0;
00234       vtkIdTypeArray *newFaceLocations = 0;
00235 
00236       if(!myCellIds.empty() && myCellTypes.empty()){
00237         if(myStoreMapping) myOut2InId.reserve(myCellIds.size());
00238         if(myChangeMode == eAdding){
00239           for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00240             if(myCellIds.find(aCellId) != myCellIds.end()){
00241               InsertCell(anInput,aConnectivity,aCellTypesArray,newFaces,newFaceLocations,aCellId,anIdList,
00242                          myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00243             }
00244           }
00245         }else{
00246           for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00247             if(myCellIds.find(aCellId) == myCellIds.end()){
00248               InsertCell(anInput,aConnectivity,aCellTypesArray,newFaces,newFaceLocations,aCellId,anIdList,
00249                          myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00250             }
00251           }
00252         }
00253       }else if(myCellIds.empty() && !myCellTypes.empty()){
00254         if(myChangeMode == eAdding){
00255           for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00256             vtkIdType aType = anInput->GetCellType(aCellId);
00257             if(myCellTypes.find(aType) != myCellTypes.end()){
00258               InsertCell(anInput,aConnectivity,aCellTypesArray,newFaces,newFaceLocations,aCellId,anIdList,
00259                          myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00260             }
00261           }
00262         }else{
00263           for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00264             vtkIdType aType = anInput->GetCellType(aCellId);
00265             if(myCellTypes.find(aType) == myCellTypes.end()){
00266               InsertCell(anInput,aConnectivity,aCellTypesArray,newFaces,newFaceLocations,aCellId,anIdList,
00267                          myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00268             }
00269           }
00270         }
00271       }else if(!myCellIds.empty() && !myCellTypes.empty()){
00272         if(myChangeMode == eAdding){
00273           for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00274             vtkIdType aType = anInput->GetCellType(aCellId);
00275             if(myCellTypes.find(aType) != myCellTypes.end()){
00276               if(myCellIds.find(aCellId) != myCellIds.end()){
00277                 InsertCell(anInput,aConnectivity,aCellTypesArray,newFaces,newFaceLocations,aCellId,anIdList,
00278                            myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00279               }
00280             }
00281           }
00282         }else{
00283           for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00284             vtkIdType aType = anInput->GetCellType(aCellId);
00285             if(myCellTypes.find(aType) == myCellTypes.end()){
00286               if(myCellIds.find(aCellId) == myCellIds.end()){
00287                 InsertCell(anInput,aConnectivity,aCellTypesArray,newFaces,newFaceLocations,aCellId,anIdList,
00288                            myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00289               }
00290             }
00291           }
00292         }
00293       }
00294       if((aNbElems = aConnectivity->GetNumberOfCells())){
00295         VTKViewer_CellLocationsArray* aCellLocationsArray = VTKViewer_CellLocationsArray::New();
00296         aCellLocationsArray->SetNumberOfComponents(1);
00297         aCellLocationsArray->SetNumberOfTuples(aNbElems);
00298         aConnectivity->InitTraversal();
00299         for(vtkIdType i = 0, *pts, npts; aConnectivity->GetNextCell(npts,pts); i++){
00300           aCellLocationsArray->SetValue(i,aConnectivity->GetTraversalLocation(npts));
00301         }
00302 #if VTK_XVERSION > 50700
00303         anOutput->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity,newFaceLocations,newFaces);
00304 #else
00305         anOutput->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity);
00306 #endif
00307         anOutput->SetPoints(anInput->GetPoints());
00308         aCellLocationsArray->Delete();
00309       }
00310       aCellTypesArray->Delete();
00311       aConnectivity->Delete();
00312       anIdList->Delete();
00313     }
00314   }else{
00315     vtkIdList *anIdList = vtkIdList::New();
00316     anIdList->SetNumberOfIds(1);
00317     vtkCellArray *aConnectivity = vtkCellArray::New();
00318     vtkIdType aNbElems = anInput->GetNumberOfPoints();
00319     aConnectivity->Allocate(2*aNbElems,0);
00320     vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
00321     aCellTypesArray->SetNumberOfComponents(1);
00322     aCellTypesArray->Allocate(aNbElems*aCellTypesArray->GetNumberOfComponents());
00323     // additional condition has been added to treat a case described in IPAL21372
00324     // note that it is significant only when myExtractionMode == ePoints
00325     if(myChangeMode == ePassAll || (myCellIds.empty() && myCellTypes.empty() && myChangeMode == eRemoving) ||
00326        !anInput->GetCellTypesArray()){
00327       if(myStoreMapping) myOut2InId.reserve(aNbElems);
00328       for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00329         InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00330                         myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00331       }
00332     }else if(!myCellIds.empty() && myCellTypes.empty()){
00333       if(myStoreMapping) myOut2InId.reserve(myCellIds.size());
00334       if(myChangeMode == eAdding){
00335         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00336           if(myCellIds.find(aCellId) != myCellIds.end()){
00337             InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00338                             myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00339           }
00340         }
00341       }else{
00342         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00343           if(myCellIds.find(aCellId) == myCellIds.end()){
00344             InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00345                             myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00346           }
00347         }
00348       }
00349     }else if(myCellIds.empty() && !myCellTypes.empty()){
00350       if(myChangeMode == eAdding){
00351         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00352           vtkIdType aType = anInput->GetCellType(aCellId);
00353           if(myCellTypes.find(aType) != myCellTypes.end()){
00354             InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00355                             myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00356           }
00357         }
00358       }else{
00359         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00360           vtkIdType aType = anInput->GetCellType(aCellId);
00361           if(myCellTypes.find(aType) == myCellTypes.end()){
00362             InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00363                             myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00364           }
00365         }
00366       }
00367     }else if(!myCellIds.empty() && !myCellTypes.empty()){
00368       if(myChangeMode == eAdding){
00369         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00370           vtkIdType aType = anInput->GetCellType(aCellId);
00371           if(myCellTypes.find(aType) != myCellTypes.end()){
00372             if(myCellIds.find(aCellId) != myCellIds.end()){
00373               InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00374                               myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00375             }
00376           }
00377         }
00378       }else{
00379         for(vtkIdType aCellId = 0, anOutId = 0; aCellId < aNbElems; aCellId++,anOutId++){
00380           vtkIdType aType = anInput->GetCellType(aCellId);
00381           if(myCellTypes.find(aType) == myCellTypes.end()){
00382             if(myCellIds.find(aCellId) == myCellIds.end()){
00383               InsertPointCell(aConnectivity,aCellTypesArray,aCellId,anIdList,
00384                               myStoreMapping,anOutId,myOut2InId,myIn2OutId);
00385             }
00386           }
00387         }
00388       }
00389     }
00390     if((aNbElems = aConnectivity->GetNumberOfCells())){
00391       VTKViewer_CellLocationsArray* aCellLocationsArray = VTKViewer_CellLocationsArray::New();
00392       aCellLocationsArray->SetNumberOfComponents(1);
00393       aCellLocationsArray->SetNumberOfTuples(aNbElems);
00394       aConnectivity->InitTraversal();
00395       for(vtkIdType i = 0, *pts, npts; aConnectivity->GetNextCell(npts,pts); i++){
00396         aCellLocationsArray->SetValue(i,aConnectivity->GetTraversalLocation(npts));
00397       }
00398 #if VTK_XVERSION > 50700
00399       anOutput->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity,0, 0);
00400 #else
00401       anOutput->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity);
00402 #endif
00403       anOutput->SetPoints(anInput->GetPoints());
00404       aCellLocationsArray->Delete();
00405     }
00406     aCellTypesArray->Delete();
00407     aConnectivity->Delete();
00408     anIdList->Delete();
00409   }
00410 /*  if(MYDEBUG){
00411     MESSAGE("Execute - anOutput->GetNumberOfCells() = "<<anOutput->GetNumberOfCells());
00412     if(myStoreMapping){
00413       MESSAGE("Execute - myOut2InId.size() = "<<myOut2InId.size());
00414       MESSAGE("Execute - myIn2OutId.size() = "<<myIn2OutId.size());
00415     }
00416   }*/
00417 //  return 1;
00418 }