Back to index

salome-paravis  6.5.0
vtkTableTo3D.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 "vtkTableTo3D.h"
00021 
00022 #include "vtkSmartPointer.h"
00023 #include "vtkDoubleArray.h"
00024 #include "vtkVariantArray.h"
00025 #include "vtkObjectFactory.h"
00026 #include "vtkPointData.h"
00027 #include "vtkPoints.h"
00028 #include "vtkPolyData.h"
00029 #include "vtkTable.h"
00030 #include "vtkInformation.h"
00031 #include "vtkStructuredGrid.h"
00032 #include "vtkStructuredGridGeometryFilter.h"
00033 #include "vtkWarpScalar.h"
00034 #include "vtkContourFilter.h"
00035 
00036 vtkStandardNewMacro(vtkTableTo3D);
00037 vtkCxxRevisionMacro(vtkTableTo3D, "$Revision: 1.1.2.2.12.1 $");
00038 
00039 
00040 vtkTableTo3D::vtkTableTo3D()
00041 {
00042   this->ScaleFactor = 1.0;
00043   this->UseOptimusScale = true;
00044   this->PresentationType = TABLETO3D_SURFACE;
00045   this->NumberOfContours = 32;
00046 }
00047 
00048 vtkTableTo3D::~vtkTableTo3D()
00049 {
00050 }
00051 
00052 int vtkTableTo3D::FillInputPortInformation(
00053   int vtkNotUsed(port), vtkInformation* info)
00054 {
00055   info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkTable");
00056   return 1;
00057 }
00058 
00059 int vtkTableTo3D::RequestData(vtkInformation* vtkNotUsed(request),
00060   vtkInformationVector** inputVector, vtkInformationVector* outputVector)
00061 {
00062   vtkTable* input = vtkTable::GetData(inputVector[0], 0);
00063   vtkPolyData* output = vtkPolyData::GetData(outputVector, 0);
00064 
00065   if (input->GetNumberOfRows() == 0 ||input->GetNumberOfColumns() < 2)
00066     {
00067       return 1;
00068     }
00069   
00070   vtkIdType xSize = input->GetNumberOfRows();
00071   vtkIdType ySize = input->GetNumberOfColumns() - 1; 
00072   vtkIdType nbPoints = xSize * ySize;
00073 
00074   vtkDataArray* xAxis = vtkDataArray::SafeDownCast(input->GetColumn(0));
00075   if (!xAxis)
00076     {
00077       vtkErrorMacro("The first column is not numeric.");
00078       return 1;
00079     }
00080     
00081   double xRange = xAxis->GetTuple1(xSize - 1) - xAxis->GetTuple1(0);
00082   double yDelta = xRange / ySize;
00083   
00084   vtkSmartPointer<vtkDoubleArray> yAxis = 
00085     vtkSmartPointer<vtkDoubleArray>::New();
00086   yAxis->SetNumberOfValues(ySize);
00087   for (vtkIdType i = 0; i < ySize; i++ )
00088     {
00089       yAxis->SetValue(i, i*yDelta);
00090     }
00091 
00092   vtkSmartPointer<vtkPoints> points = 
00093     vtkSmartPointer<vtkPoints>::New();
00094   points->SetNumberOfPoints(nbPoints);
00095 
00096   vtkSmartPointer<vtkIntArray> pointsIdMapper = 
00097     vtkSmartPointer<vtkIntArray>::New();
00098   pointsIdMapper->SetName("POINTS_ID_MAPPER");
00099   pointsIdMapper->SetNumberOfComponents(2);
00100   pointsIdMapper->SetNumberOfTuples(nbPoints);
00101   int *pointsIdMapperPtr = pointsIdMapper->GetPointer(0);
00102 
00103   for (vtkIdType i = 0, pntId = 0; i < ySize; i++) 
00104     {
00105       for (vtkIdType j = 0; j < xSize; j++, pntId++) 
00106        {
00107          points->SetPoint(pntId, xAxis->GetTuple1(j), 
00108                                yAxis->GetValue(i), 
00109                                   0.0);
00110 
00111          *pointsIdMapperPtr++ = pntId;
00112          *pointsIdMapperPtr++ = 0;
00113        }
00114     }
00115 
00116   vtkSmartPointer<vtkDoubleArray> scalars = 
00117     vtkSmartPointer<vtkDoubleArray>::New();
00118   scalars->SetNumberOfComponents(1);
00119   scalars->SetNumberOfTuples(nbPoints);
00120   double *scalarsPtr = scalars->GetPointer(0);
00121   for (vtkIdType i = 0; i < ySize; i++) 
00122     {
00123       vtkDataArray* col = 
00124        vtkDataArray::SafeDownCast(input->GetColumn(i + 1));
00125       
00126       if (!col)
00127        {
00128          vtkErrorMacro("Column "<< i <<"is not numeric.");
00129          return 1;
00130        }
00131       
00132       for ( vtkIdType j = 0; j < xSize; j++ ) 
00133        {
00134          double value = col->GetTuple1(j);
00135          *scalarsPtr++ = value;
00136        }
00137     }
00138 
00139   vtkSmartPointer<vtkStructuredGrid> structuredGrid = 
00140     vtkSmartPointer<vtkStructuredGrid>::New();
00141   structuredGrid->SetPoints(points);
00142 
00143   structuredGrid->SetDimensions(xSize, ySize, 1);
00144 
00145   // structuredGrid->GetPointData()->AddArray(pointsIdMapper);
00146   if (input->GetInformation()->Has(vtkDataObject::FIELD_NAME()))
00147     {
00148       scalars->SetName(input->GetInformation()->Get(vtkDataObject::FIELD_NAME()));
00149     }
00150   else
00151     {
00152       scalars->SetName("Table");
00153     }
00154   structuredGrid->GetPointData()->SetScalars(scalars);
00155 
00156   vtkSmartPointer<vtkStructuredGridGeometryFilter> geomFilter = 
00157     vtkSmartPointer<vtkStructuredGridGeometryFilter>::New();
00158   geomFilter->SetInput(structuredGrid);
00159   geomFilter->Update();
00160   
00161   vtkSmartPointer<vtkWarpScalar> warpScalar = 
00162     vtkSmartPointer<vtkWarpScalar>::New();
00163   
00164   double scaleFactor = this->ScaleFactor;
00165   if (this->UseOptimusScale)
00166     {
00167       double range[2];
00168       geomFilter->GetOutput()->GetScalarRange(range);
00169       double length = geomFilter->GetOutput()->GetLength();
00170       if (range[1] > 0)
00171        {
00172          scaleFactor = length / range[1] * 0.3;
00173        }
00174       else
00175        {
00176          scaleFactor = 0;
00177        }
00178     }
00179 
00180   if (this->PresentationType == TABLETO3D_SURFACE)
00181     {
00182       warpScalar->SetInput(geomFilter->GetOutput());
00183       warpScalar->SetScaleFactor(scaleFactor);
00184     }
00185   else
00186     {
00187       vtkSmartPointer<vtkContourFilter> contourFilter = 
00188        vtkSmartPointer<vtkContourFilter>::New();
00189       contourFilter->SetInput(geomFilter->GetOutput());
00190       contourFilter->GenerateValues(this->NumberOfContours, 
00191                                 geomFilter->GetOutput()->GetScalarRange());
00192       warpScalar->SetInput(contourFilter->GetOutput());
00193       warpScalar->SetScaleFactor(scaleFactor);
00194     }
00195 
00196   warpScalar->Update();
00197   output->ShallowCopy(warpScalar->GetPolyDataOutput());
00198   
00199   return 1;
00200 }
00201 
00202 //----------------------------------------------------------------------------
00203 void vtkTableTo3D::PrintSelf(ostream& os, vtkIndent indent)
00204 {
00205   this->Superclass::PrintSelf(os, indent);
00206 
00207   os << indent << "ScaleFactor: " << this->ScaleFactor << endl;
00208   os << indent << "UseOptimusScale: "
00209      << (this->UseOptimusScale? "true" : "false") << endl;
00210   os << indent << "PresentationType: " 
00211      << ((this->PresentationType == TABLETO3D_SURFACE)? "Surface" : "Contour")
00212      << endl;
00213   os << indent << "NumberOfContours: " << this->NumberOfContours << endl;
00214 }