src/Plugins/MedReader/IO/vtkMedRegularGrid.cxx

   1 // Copyright (C) 2010-2011  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 #include "vtkMedRegularGrid.h"
  21
  22 #include "vtkObjectFactory.h"
  23 #include "vtkDataArray.h"
  24 #include "vtkStructuredGrid.h"
  25
  26 #include "vtkMedUtilities.h"
  27 #include "vtkMedFile.h"
  28 #include "vtkMedDriver.h"
  29 #include "vtkMedMesh.h"
  30 #include "vtkMedFamilyOnEntityOnProfile.h"
  31 #include "vtkMedFamilyOnEntity.h"
  32 #include "vtkMedEntityArray.h"
  33 #include "vtkMedProfile.h"
  34 #include "vtkMedFamily.h"
  35
  36 vtkCxxGetObjectVectorMacro(vtkMedRegularGrid, AxisCoordinate, vtkDataArray);
  37 vtkCxxSetAbstractObjectVectorMacro(vtkMedRegularGrid, AxisCoordinate, vtkDataArray);
  38
  39 // vtkCxxRevisionMacro(vtkMedRegularGrid, "$Revision$")
  40 vtkStandardNewMacro(vtkMedRegularGrid)
  41
  42 vtkMedRegularGrid::vtkMedRegularGrid()
  43 {
  44   this->AxisCoordinate = new vtkObjectVector<vtkDataArray>();
  45 }
  46
  47 vtkMedRegularGrid::~vtkMedRegularGrid()
  48 {
  49   delete this->AxisCoordinate;
  50 }
  51
  52 void vtkMedRegularGrid::SetDimension(med_int dim)
  53 {
  54   if(dim < 0)
  55     dim = 0;
  56   this->AxisSize.resize(dim);
  57   this->SetNumberOfAxisCoordinate(dim);
  58 }
  59
  60 int vtkMedRegularGrid::GetDimension()
  61 {
  62   return this->AxisSize.size();
  63 }
  64
  65 void  vtkMedRegularGrid::SetAxisSize(int axis, med_int size)
  66 {
  67   if(axis < 0)
  68     return;
  69
  70   if(axis >= this->GetDimension())
  71     {
  72     this->SetDimension(axis+1);
  73     }
  74
  75   this->AxisSize[axis] = size;
  76 }
  77
  78 med_int vtkMedRegularGrid::GetAxisSize(int dim)
  79 {
  80   if(dim < 0 || dim >= this->AxisSize.size())
  81     return 0;
  82   return this->AxisSize[dim];
  83 }
  84
  85 med_int vtkMedRegularGrid::GetNumberOfPoints()
  86 {
  87   med_int npts = 1;
  88   for(int dim = 0; dim < this->AxisSize.size(); dim++)
  89     {
  90     npts *= this->AxisSize[dim];
  91     }
  92   return npts;
  93 }
  94
  95 void  vtkMedRegularGrid::LoadCoordinates()
  96 {
  97   this->GetParentMesh()->GetParentFile()->GetMedDriver()->LoadCoordinates(this);
  98 }
  99
 100 int vtkMedRegularGrid::IsCoordinatesLoaded()
 101 {
 102   bool res =  this->GetDimension() == this->AxisSize.size() &&
 103       this->GetDimension() == this->AxisCoordinate->size();
 104
 105   if(!res)
 106     return 0;
 107
 108   med_int nloadedcoords = 1;
 109   for(int axis=0; axis < this->GetDimension(); axis++)
 110     {
 111     vtkDataArray* axiscoords = this->GetAxisCoordinate(axis);
 112     if(axiscoords != NULL)
 113       nloadedcoords *= axiscoords->GetNumberOfTuples();
 114     else
 115       nloadedcoords = 0;
 116     }
 117   return nloadedcoords == this->GetNumberOfPoints();
 118 }
 119
 120 double* vtkMedRegularGrid::GetCoordTuple(med_int index)
 121 {
 122   this->CoordTuple[0] = 0;
 123   this->CoordTuple[1] = 0;
 124   this->CoordTuple[2] = 0;
 125   med_int prevmod = 1;
 126   for(int axis=0; axis < this->GetDimension(); axis++)
 127     {
 128     med_int modulo = prevmod * this->AxisSize[axis];
 129     med_int axisindex = (index % modulo) / prevmod;
 130     prevmod = modulo;
 131     vtkDataArray* coords = this->GetAxisCoordinate(axis);
 132     this->CoordTuple[axis] = coords->
 133                              GetTuple1(axisindex);
 134     }
 135
 136   return this->CoordTuple;
 137 }
 138
 139 vtkDataSet* vtkMedRegularGrid::CreateVTKDataSet(vtkMedFamilyOnEntityOnProfile* foep)
 140 {
 141   vtkStructuredGrid* vtkgrid = vtkStructuredGrid::New();
 142
 143   vtkPoints* points = vtkPoints::New();
 144   vtkgrid->SetPoints(points);
 145   points->Delete();
 146
 147   vtkIdType dims[3] = {this->GetAxisSize(0),
 148                     this->GetAxisSize(1),
 149                     this->GetAxisSize(2)};
 150
 151   for(int dim=0; dim<3; dim++)
 152     dims[dim] = (dims[dim] >= 1 ? dims[dim] : 1);
 153
 154   vtkgrid->SetDimensions(dims[0], dims[1], dims[2]);
 155
 156   this->LoadCoordinates();
 157
 158   if(this->GetAxisCoordinate(0) == NULL)
 159     {
 160     vtkgrid->Delete();
 161     return NULL;
 162     }
 163
 164   vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
 165   coords->SetNumberOfComponents(3);
 166   coords->SetNumberOfTuples(this->GetNumberOfPoints());
 167   vtkgrid->GetPoints()->SetData(coords);
 168   coords->Delete();
 169
 170   med_int npts;
 171   double coord[3] = {0, 0, 0};
 172
 173   npts = this->GetNumberOfPoints();
 174   for(med_int id=0; id<npts; id++)
 175     {
 176     double * tuple = this->GetCoordTuple(id);
 177     for(int dim=0; dim<this->GetDimension(); dim++)
 178       {
 179       coord[dim] = tuple[dim];
 180       }
 181     coords->SetTuple(id, coord);
 182     }
 183
 184   if(foep->GetProfile() != NULL)
 185     {
 186     foep->GetProfile()->Load();
 187     vtkMedIntArray* pids = foep->GetProfile()->GetIds();
 188     med_int previd = -1;
 189     for(med_int pid=0; pid<pids->GetNumberOfTuples(); pid++)
 190       {
 191       med_int id = pids->GetValue(pid) - 1;
 192       for(med_int theid=previd+1; theid<id; theid++)
 193         {
 194         vtkgrid->BlankCell(theid);
 195         }
 196
 197       previd = id;
 198       }
 199     }
 200
 201   if(foep->GetFamilyOnEntity()->GetEntityArray()->GetNumberOfFamilyOnEntity() > 1)
 202     {
 203     med_int famid = foep->GetFamilyOnEntity()->GetFamily()->GetId();
 204     vtkMedEntityArray* ea = foep->GetFamilyOnEntity()->GetEntityArray();
 205     for(med_int id=0; id<vtkgrid->GetNumberOfCells(); id++)
 206       {
 207       if(ea->GetFamilyId(id) != famid)
 208         vtkgrid->BlankCell(id);
 209       }
 210     }
 211   return vtkgrid;
 212 }
 213
 214 void vtkMedRegularGrid::PrintSelf(ostream& os, vtkIndent indent)
 215 {
 216   this->Superclass::PrintSelf(os, indent);
 217 }