Back to index

salome-paravis  6.5.0
vtkElevationSurfaceFilter.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 "vtkElevationSurfaceFilter.h"
00021 #include "vtkInformation.h"
00022 #include "vtkInformationVector.h"
00023 #include "vtkObjectFactory.h"
00024 #include "vtkElevationSurfaceFilter.h"
00025 #include "vtkPolyData.h"
00026 #include "vtkIdTypeArray.h"
00027 #include "vtkPolyData.h"
00028 #include "vtkUnstructuredGrid.h"
00029 #include "vtkDemandDrivenPipeline.h"
00030 #include "vtkStreamingDemandDrivenPipeline.h"
00031 #include "vtkGenericCell.h"
00032 #include "vtkSmartPointer.h"
00033 #include "vtkPoints.h"
00034 #include "vtkCellArray.h"
00035 #include "vtkPointData.h"
00036 #include "vtkCellData.h"
00037 //#include "vtkDataSetSurfaceFilter.h"
00038 
00039 #include <math.h>
00040 
00041 vtkCxxRevisionMacro(vtkElevationSurfaceFilter, "$Revision: 1.1.4.3.2.1 $");
00042 vtkStandardNewMacro(vtkElevationSurfaceFilter);
00043 
00044 vtkElevationSurfaceFilter::vtkElevationSurfaceFilter()
00045 {
00046   this->SetNumberOfInputPorts(1);
00047   this->SetNumberOfOutputPorts(1);
00048 
00049   this->ScaleFactor = 0.5;
00050   this->Direction[0] = 0.0;
00051   this->Direction[1] = 0.0;
00052   this->Direction[2] = 1.0;
00053   this->AutoDetectDirection = 1;
00054 }
00055 
00056 vtkElevationSurfaceFilter::~vtkElevationSurfaceFilter()
00057 {
00058 }
00059 
00060 //----------------------------------------------------------------------------
00061 int vtkElevationSurfaceFilter::ProcessRequest(vtkInformation* request,
00062                                          vtkInformationVector** inputVector,
00063                                          vtkInformationVector* outputVector)
00064 {
00065   // generate the data
00066   if(request->Has(vtkDemandDrivenPipeline::REQUEST_DATA()))
00067     {
00068     return this->RequestData(request, inputVector, outputVector);
00069     }
00070 
00071   if(request->Has(vtkStreamingDemandDrivenPipeline::REQUEST_UPDATE_EXTENT()))
00072     {
00073     return this->RequestUpdateExtent(request, inputVector, outputVector);
00074     }
00075 
00076   // execute information
00077   if(request->Has(vtkDemandDrivenPipeline::REQUEST_INFORMATION()))
00078     {
00079     return this->RequestInformation(request, inputVector, outputVector);
00080     }
00081 
00082   return this->Superclass::ProcessRequest(request, inputVector, outputVector);
00083 }
00084 
00085 //----------------------------------------------------------------------------
00086 int vtkElevationSurfaceFilter::FillOutputPortInformation(
00087   int vtkNotUsed(port), vtkInformation* info)
00088 {
00089   // now add our info
00090   info->Set(vtkDataObject::DATA_TYPE_NAME(), "vtkUnstructuredGrid");
00091   return 1;
00092 }
00093 
00094 //----------------------------------------------------------------------------
00095 int vtkElevationSurfaceFilter::FillInputPortInformation(
00096   int vtkNotUsed(port), vtkInformation* info)
00097 {
00098   info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkPointSet");
00099   return 1;
00100 }
00101 
00102 //----------------------------------------------------------------------------
00103 int vtkElevationSurfaceFilter::RequestUpdateExtent(
00104   vtkInformation* vtkNotUsed(request),
00105   vtkInformationVector** inputVector,
00106   vtkInformationVector* vtkNotUsed(outputVector))
00107 {
00108   int numInputPorts = this->GetNumberOfInputPorts();
00109   for (int i=0; i<numInputPorts; i++)
00110     {
00111     int numInputConnections = this->GetNumberOfInputConnections(i);
00112     for (int j=0; j<numInputConnections; j++)
00113       {
00114       vtkInformation* inputInfo = inputVector[i]->GetInformationObject(j);
00115       inputInfo->Set(vtkStreamingDemandDrivenPipeline::EXACT_EXTENT(), 1);
00116       }
00117     }
00118   return 1;
00119 }
00120 
00121 int vtkElevationSurfaceFilter::RequestInformation(vtkInformation *request,
00122     vtkInformationVector **input, vtkInformationVector *output)
00123 {
00124   return 1;
00125 }
00126 
00127 int vtkElevationSurfaceFilter::RequestData(vtkInformation *request,
00128     vtkInformationVector **input, vtkInformationVector *output)
00129 {
00130        vtkPointSet *psIn = vtkPointSet::SafeDownCast(
00131       input[0]->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00132 
00133   //vtkPolyData *psIn = vtkPolyData::SafeDownCast(
00134   //    input[0]->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00135 
00136   /*vtkSmartPointer<vtkDataSetSurfaceFilter> surfaceFilter = vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
00137   surfaceFilter->SetInput(psIn);
00138   vtkPolyData* psIn = vtkPolyData::SafeDownCast(surfaceFilter->GetOutput());*/
00139 
00140   vtkUnstructuredGrid *usgOut = vtkUnstructuredGrid::SafeDownCast(
00141       output->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
00142 
00143   vtkDataArray* array = this->GetInputArrayToProcess(0, input);
00144 
00145   if(psIn == NULL || array == NULL || usgOut == NULL
00146      || array->GetNumberOfComponents() != 1)
00147     {
00148     vtkDebugMacro("vtkElevationSurfaceFilter no correctly configured");
00149     return 1;
00150     }
00151 
00152   double dir[3];
00153   if(this->AutoDetectDirection)
00154     {
00155     this->ComputeDirection(psIn, dir);
00156     }
00157   else
00158     {
00159     dir[0] = this->Direction[0];
00160     dir[1] = this->Direction[1];
00161     dir[2] = this->Direction[2];
00162     }
00163 
00164   double len = dir[0]*dir[0] + dir[1]*dir[1] + dir[2]*dir[2];
00165 
00166   if(len == 0)
00167     {
00168     dir[2] = 1;
00169     len = 1.0;
00170     }
00171 
00172   len = sqrt(len);
00173 
00174   dir[0] /= len;
00175   dir[1] /= len;
00176   dir[2] /= len;
00177 
00178   dir[0] *= this->GetScaleFactor();
00179   dir[1] *= this->GetScaleFactor();
00180   dir[2] *= this->GetScaleFactor();
00181 
00182   usgOut->Allocate(psIn->GetNumberOfCells());
00183 
00184   vtkSmartPointer<vtkPoints> newPts = vtkSmartPointer<vtkPoints>::New();
00185   usgOut->SetPoints(newPts);
00186 
00187   usgOut->GetPointData()->CopyAllocate(psIn->GetPointData(),
00188           2*psIn->GetNumberOfPoints());
00189   usgOut->GetCellData()->CopyAllocate(psIn->GetCellData(),
00190                 psIn->GetNumberOfCells());
00191 
00192   vtkIdType ncell = psIn->GetNumberOfCells();
00193   vtkSmartPointer<vtkIdList> newIds = vtkSmartPointer<vtkIdList>::New();
00194   vtkSmartPointer<vtkIdList> polyhedronIds = vtkSmartPointer<vtkIdList>::New();
00195   vtkSmartPointer<vtkIdList> neighbors = vtkSmartPointer<vtkIdList>::New();
00196   for(vtkIdType cellId=0; cellId < ncell; cellId++)
00197     {
00198     vtkCell* cell = psIn->GetCell(cellId);
00199     if(cell->GetCellDimension() != 2)
00200       continue;
00201 
00202     unsigned char newCellType = VTK_EMPTY_CELL;
00203     unsigned char oldCellType = cell->GetCellType();
00204     switch(oldCellType)
00205       {
00206       case VTK_TRIANGLE :
00207         newCellType = VTK_WEDGE;
00208         break;
00209       case VTK_QUAD :
00210         newCellType = VTK_HEXAHEDRON;
00211         break;
00212       case VTK_POLYGON :
00213         newCellType = VTK_POLYHEDRON;
00214       // default : add new cell types to extrude here
00215       }
00216     if(newCellType == VTK_EMPTY_CELL)
00217       continue;
00218 
00219     double cellScalar = array->GetTuple1(cellId);
00220 
00221     vtkIdList* oldIds = cell->GetPointIds();
00222     int oldPtsNumber = oldIds->GetNumberOfIds();
00223     int newPtsNumber = oldPtsNumber * 2;
00224     newIds->SetNumberOfIds(newPtsNumber);
00225     double coords[VTK_CELL_SIZE*3];
00226     for(int ptid = 0; ptid < oldPtsNumber; ptid++)
00227       {
00228        psIn->GetPoint(oldIds->GetId(ptid), coords + 3*ptid);
00229       }
00230     for(int ptid = 0; ptid < oldPtsNumber; ptid++)
00231       {
00232       coords[(ptid+oldPtsNumber)*3+0] = coords[ptid*3+0] + cellScalar*dir[0];
00233       coords[(ptid+oldPtsNumber)*3+1] = coords[ptid*3+1] + cellScalar*dir[1];
00234       coords[(ptid+oldPtsNumber)*3+2] = coords[ptid*3+2] + cellScalar*dir[2];
00235       }
00236     double minScalar;
00237     bool minInitialized = false;
00238     for(int ptid = 0; ptid < oldPtsNumber; ptid++)
00239       {
00240       neighbors->Initialize();
00241       psIn->GetPointCells(oldIds->GetId(ptid), neighbors);
00242       for(int neiCellIt = 0; neiCellIt < neighbors->GetNumberOfIds(); neiCellIt++)
00243         {
00244         vtkIdType  neigCellId = neighbors->GetId(neiCellIt);
00245         if(neigCellId == cellId)
00246           continue;
00247         double neighborScalar = array->GetTuple1(neigCellId);
00248         if(neighborScalar != 0.0)
00249                 {
00250           if(!minInitialized)
00251               minScalar = neighborScalar;
00252           else
00253                   minScalar = (neighborScalar < minScalar ? neighborScalar : minScalar);
00254                 minInitialized = true;
00255                 }
00256         }
00257       if(!minInitialized)
00258        minScalar = 0.0;
00259       }
00260     for(int ptid = 0; ptid < oldPtsNumber; ptid++)
00261       {
00262          if(cellScalar != 0)
00263            {
00264         coords[(ptid)*3+0] = coords[ptid*3+0] + minScalar*dir[0];
00265         coords[(ptid)*3+1] = coords[ptid*3+1] + minScalar*dir[1];
00266         coords[(ptid)*3+2] = coords[ptid*3+2] + minScalar*dir[2];
00267        }
00268       else
00269        {
00270         coords[(ptid+oldPtsNumber)*3+0] = coords[ptid*3+0] + minScalar*dir[0];
00271         coords[(ptid+oldPtsNumber)*3+1] = coords[ptid*3+1] + minScalar*dir[1];
00272         coords[(ptid+oldPtsNumber)*3+2] = coords[ptid*3+2] + minScalar*dir[2];
00273        }
00274       }
00275     for(int ptid=0; ptid<newPtsNumber; ptid++)
00276       {
00277       vtkIdType newId = newPts->InsertNextPoint(coords + 3*ptid);
00278       newIds->SetId(ptid, newId);
00279       usgOut->GetPointData()->CopyData(psIn->GetPointData(),
00280                                        oldIds->GetId(ptid % oldPtsNumber),
00281                                        newIds->GetId(ptid));
00282       }
00283     vtkIdType newCellId;
00284     if(newCellType == VTK_POLYHEDRON)
00285       {
00286       polyhedronIds->Initialize();
00287       // in the polyhedron case, I will generate a quad for each edge
00288       // of the input, and two capping faces
00289       polyhedronIds->InsertNextId(2+oldPtsNumber);
00290       // insert the bottom face
00291       polyhedronIds->InsertNextId(oldPtsNumber);
00292       for(int ptid = 0; ptid < oldPtsNumber; ptid++)
00293         {
00294         polyhedronIds->InsertNextId(newIds->GetId(ptid));
00295         }
00296       // insert the top face
00297       polyhedronIds->InsertNextId(oldPtsNumber);
00298       for(int ptid = oldPtsNumber; ptid < 2*oldPtsNumber; ptid++)
00299         {
00300         polyhedronIds->InsertNextId(newIds->GetId(ptid));
00301         }
00302       // insert the bording quads
00303       for(int qid = 0; qid < oldPtsNumber; qid++)
00304         {
00305         polyhedronIds->InsertNextId(4);
00306         polyhedronIds->InsertNextId(newIds->GetId(qid));
00307         polyhedronIds->InsertNextId(newIds->GetId(qid+oldPtsNumber));
00308         polyhedronIds->InsertNextId(newIds->GetId(qid+((oldPtsNumber+1)%oldPtsNumber)));
00309         polyhedronIds->InsertNextId(newIds->GetId((qid+1)%oldPtsNumber));
00310         }
00311       newIds->Initialize();
00312       for(int jj=0; jj<polyhedronIds->GetNumberOfIds(); jj++)
00313         {
00314         newIds->InsertNextId(polyhedronIds->GetId(jj));
00315         }
00316       }
00317     newCellId = usgOut->InsertNextCell(newCellType, newIds);
00318     usgOut->GetCellData()->CopyData(psIn->GetCellData(),
00319                                    cellId,
00320                                    newCellId);
00321     }
00322 
00323   usgOut->GetFieldData()->ShallowCopy(psIn->GetFieldData());
00324 
00325   usgOut->Squeeze();
00326 
00327   return 1;
00328 }
00329 
00330 void  vtkElevationSurfaceFilter::ComputeDirection(vtkPointSet* psIn, double *outDir)
00331 {
00332   double tmp[2][3] = {{0, 0, 0}, {0, 0, 0}};
00333   outDir[0] = outDir[1] = outDir[2] = 0;
00334 
00335   vtkPoints* pts = psIn->GetPoints();
00336   vtkSmartPointer<vtkGenericCell> cell = vtkSmartPointer<vtkGenericCell>::New();
00337 
00338   for(vtkIdType cellId = 0; cellId < psIn->GetNumberOfCells(); cellId++)
00339     {
00340     psIn->GetCell(cellId, cell);
00341     if(cell->GetCellDimension() != 2)
00342       continue;
00343 
00344     vtkIdList* ptIds = cell->GetPointIds();
00345     for(int i=0; i<ptIds->GetNumberOfIds(); i++)
00346       {
00347       vtkIdType firstId = ptIds->GetId(i);
00348       vtkIdType secondId = ptIds->GetId((i+1)%ptIds->GetNumberOfIds());
00349       pts->GetPoint(firstId, tmp[0]);
00350       pts->GetPoint(secondId, tmp[1]);
00351       outDir[0] += tmp[0][1]*tmp[1][2] - tmp[0][2]*tmp[1][1];
00352       outDir[1] += tmp[0][2]*tmp[1][0] - tmp[0][0]*tmp[1][2];
00353       outDir[2] += tmp[0][0]*tmp[1][1] - tmp[0][1]*tmp[1][0];
00354       }
00355     }
00356 }
00357 
00358 void vtkElevationSurfaceFilter::PrintSelf(ostream& os, vtkIndent indent)
00359 {
00360   this->Superclass::PrintSelf(os, indent);
00361 
00362   os << indent << "ScaleFactor : " << this->ScaleFactor << endl;
00363 }