src/Plugins/MedReader/IO/vtkMedDriver30.cxx

   1 // Copyright (C) 2010-2012  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 "vtkMedDriver30.h"
  21
  22 #include "vtkObjectFactory.h"
  23 #include "vtkStringArray.h"
  24 #include "vtkDataArray.h"
  25 #include "vtkIdTypeArray.h"
  26 #include "vtkDoubleArray.h"
  27 #include "vtkIntArray.h"
  28 #include "vtkCharArray.h"
  29
  30 #include "vtkMedFile.h"
  31 #include "vtkMedCartesianGrid.h"
  32 #include "vtkMedPolarGrid.h"
  33 #include "vtkMedCurvilinearGrid.h"
  34 #include "vtkMedUnstructuredGrid.h"
  35 #include "vtkMedField.h"
  36 #include "vtkMedMesh.h"
  37 #include "vtkMedFamily.h"
  38 #include "vtkMedUtilities.h"
  39 #include "vtkMedEntityArray.h"
  40 #include "vtkMedLocalization.h"
  41 #include "vtkMedProfile.h"
  42 #include "vtkMedFieldOverEntity.h"
  43 #include "vtkMedFieldStep.h"
  44 #include "vtkMedGroup.h"
  45 #include "vtkMedIntArray.h"
  46 #include "vtkMedInterpolation.h"
  47 #include "vtkMedFieldOnProfile.h"
  48 #include "vtkMedFraction.h"
  49 #include "vtkMedLink.h"
  50 #include "vtkMedStructElement.h"
  51 #include "vtkMedConstantAttribute.h"
  52 #include "vtkMedVariableAttribute.h"
  53
  54 #include "vtkMultiProcessController.h"
  55
  56 #include <string>
  57 #include <vector>
  58 #include <algorithm>
  59 using namespace std;
  60
  61 vtkCxxRevisionMacro(vtkMedDriver30, "$Revision$")
  62 vtkStandardNewMacro(vtkMedDriver30)
  63
  64 vtkMedDriver30::vtkMedDriver30() : vtkMedDriver()
  65 {
  66 }
  67
  68 vtkMedDriver30::~vtkMedDriver30()
  69 {
  70 }
  71
  72 // Description:
  73 // load all Information data associated with this cartesian grid.
  74 void vtkMedDriver30::ReadRegularGridInformation(vtkMedRegularGrid* grid)
  75 {
  76   FileOpen open(this);
  77
  78   grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
  79
  80   for(int axis=0; axis < grid->GetDimension(); axis++)
  81     {
  82     med_int size;
  83     med_bool coordinatechangement, geotransformation;
  84
  85     if ((size = MEDmeshnEntity(
  86                   this->FileId,
  87                   grid->GetParentMesh()->GetName(),
  88                   grid->GetComputeStep().TimeIt,
  89                   grid->GetComputeStep().IterationIt,
  90                   MED_NODE,
  91                   MED_NONE,
  92                   (med_data_type)(MED_COORDINATE_AXIS1 + axis),
  93                   MED_NO_CMODE,
  94                   &coordinatechangement,
  95                   &geotransformation)) < 0)
  96       {
  97       vtkErrorMacro("ERROR : number of coordinates on X axis ...");
  98       }
  99
 100     grid->SetAxisSize(axis, size);
 101     }
 102
 103   med_int ncell = 1;
 104   if(grid->GetAxisSize(0) > 1)
 105     ncell *= grid->GetAxisSize(0)-1;
 106   if(grid->GetAxisSize(1) > 1)
 107     ncell *= grid->GetAxisSize(1)-1;
 108   if(grid->GetAxisSize(2) > 1)
 109     ncell *= grid->GetAxisSize(2)-1;
 110
 111   vtkMedEntity entity;
 112   entity.EntityType = MED_CELL;
 113
 114   switch(grid->GetDimension())
 115     {
 116     case 0 :
 117       entity.GeometryType = MED_POINT1;
 118       break;
 119     case 1 :
 120       entity.GeometryType = MED_SEG2;
 121       break;
 122     case 2 :
 123       entity.GeometryType = MED_QUAD4;
 124       break;
 125     case 3 :
 126       entity.GeometryType = MED_HEXA8;
 127       break;
 128     default :
 129         vtkErrorMacro("Unsupported dimension for curvilinear grid : "
 130                       << grid->GetDimension());
 131     return;
 132     }
 133
 134   vtkMedEntityArray* array = vtkMedEntityArray::New();
 135   array->SetParentGrid(grid);
 136   array->SetNumberOfEntity(ncell);
 137   array->SetEntity(entity);
 138   array->SetConnectivity(MED_NODAL);
 139   grid->AppendEntityArray(array);
 140   array->Delete();
 141   // this triggers the creation of undefined families
 142   this->LoadFamilyIds(array);
 143
 144 }
 145
 146 void  vtkMedDriver30::LoadRegularGridCoordinates(vtkMedRegularGrid* grid)
 147 {
 148   FileOpen open(this);
 149
 150   for(int axis=0; axis < grid->GetParentMesh()->GetNumberOfAxis(); axis++)
 151     {
 152
 153     vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
 154     grid->SetAxisCoordinate(axis, coords);
 155     coords->Delete();
 156
 157     coords->SetNumberOfComponents(1);
 158     coords->SetNumberOfTuples(grid->GetAxisSize(axis));
 159
 160     if (MEDmeshGridIndexCoordinateRd(
 161           this->FileId,
 162           grid->GetParentMesh()->GetName(),
 163           grid->GetComputeStep().TimeIt,
 164           grid->GetComputeStep().IterationIt,
 165           axis+1,
 166           (med_float*)coords->GetVoidPointer(0)) < 0)
 167       {
 168       vtkErrorMacro("ERROR : read axis " << axis << " coordinates ...");
 169       grid->SetAxisCoordinate(axis, NULL);
 170       return;
 171       }
 172     }
 173 }
 174
 175 // Description:
 176 // load all Information data associated with this standard grid.
 177 void vtkMedDriver30::ReadCurvilinearGridInformation(vtkMedCurvilinearGrid* grid)
 178 {
 179   FileOpen open(this);
 180
 181   grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
 182
 183   med_int size;
 184   med_bool coordinatechangement, geotransformation;
 185
 186   if ((size = MEDmeshnEntity(
 187                 this->FileId,
 188                 grid->GetParentMesh()->GetName(),
 189                 grid->GetComputeStep().TimeIt,
 190                 grid->GetComputeStep().IterationIt,
 191                 MED_NODE,
 192                 MED_NONE,
 193                 MED_COORDINATE,
 194                 MED_NO_CMODE,
 195                 &coordinatechangement,
 196                 &geotransformation)) < 0)
 197     {
 198     vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshnEntity");
 199     }
 200
 201   grid->SetNumberOfPoints(size);
 202
 203   med_int axissize[3];
 204   if(MEDmeshGridStructRd(
 205       this->FileId,
 206       grid->GetParentMesh()->GetName(),
 207       grid->GetComputeStep().TimeIt,
 208       grid->GetComputeStep().IterationIt,
 209       axissize) < 0)
 210     {
 211     vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshGridStructRd");
 212     }
 213   grid->SetAxisSize(0, (axissize[0] <= 0 ? 1: axissize[0]));
 214   grid->SetAxisSize(1, (axissize[1] <= 0 ? 1: axissize[1]));
 215   grid->SetAxisSize(2, (axissize[2] <= 0 ? 1: axissize[2]));
 216
 217   // A test to verify the number of points : total number of
 218   // points should be equal to the product of each axis size
 219   if(size != grid->GetAxisSize(0)*grid->GetAxisSize(1)*grid->GetAxisSize(2))
 220     {
 221     vtkErrorMacro("The total number of points of a Curvilinear grid should "
 222                   << "be the product of each axis size!");
 223     }
 224
 225   med_int ncell = 1;
 226   if(grid->GetAxisSize(0) > 1)
 227     ncell *= grid->GetAxisSize(0)-1;
 228   if(grid->GetAxisSize(1) > 1)
 229     ncell *= grid->GetAxisSize(1)-1;
 230   if(grid->GetAxisSize(2) > 1)
 231     ncell *= grid->GetAxisSize(2)-1;
 232
 233   vtkMedEntity entity;
 234   entity.EntityType = MED_CELL;
 235
 236   switch(grid->GetDimension())
 237     {
 238     case 0 :
 239       entity.GeometryType = MED_POINT1;
 240       break;
 241     case 1 :
 242       entity.GeometryType = MED_SEG2;
 243       break;
 244     case 2 :
 245       entity.GeometryType = MED_QUAD4;
 246       break;
 247     case 3 :
 248       entity.GeometryType = MED_HEXA8;
 249       break;
 250     default :
 251         vtkErrorMacro("Unsupported dimension for curvilinear grid : "
 252                       << grid->GetDimension());
 253     return;
 254     }
 255
 256   vtkMedEntityArray* array = vtkMedEntityArray::New();
 257   array->SetParentGrid(grid);
 258   array->SetNumberOfEntity(ncell);
 259   array->SetEntity(entity);
 260   array->SetConnectivity(MED_NODAL);
 261   grid->AppendEntityArray(array);
 262   array->Delete();
 263   // this triggers the creation of undefined families
 264   this->LoadFamilyIds(array);
 265
 266 }
 267
 268 // Description : read the number of entity of all geometry type
 269 // for a given entity type and a given connectivity mode
 270 void vtkMedDriver30::ReadNumberOfEntity(
 271     vtkMedUnstructuredGrid* grid,
 272     med_entity_type entityType,
 273     med_connectivity_mode connectivity)
 274 {
 275   FileOpen open(this);
 276
 277   med_bool changement, transformation;
 278
 279   const char* meshName = grid->GetParentMesh()->GetName();
 280
 281   const vtkMedComputeStep& cs = grid->GetComputeStep();
 282
 283   med_int nentity = MEDmeshnEntity(
 284                         this->FileId,
 285                         meshName,
 286                         cs.TimeIt,
 287                         cs.IterationIt,
 288                         entityType,
 289                         MED_GEO_ALL,
 290                         MED_UNDEF_DATATYPE ,
 291                         connectivity,
 292                         &changement,
 293                         &transformation );
 294
 295   for(med_int geotypeit = 1; geotypeit <= nentity; geotypeit++)
 296     {
 297     // read cell informations
 298     vtkMedEntity entity;
 299     entity.EntityType = entityType;
 300
 301     char geometryName[MED_NAME_SIZE+1] = "";
 302
 303     // this gives us the med_geometry_type
 304     if( MEDmeshEntityInfo( FileId, meshName,
 305                            cs.TimeIt,
 306                            cs.IterationIt,
 307                            entityType,
 308                            geotypeit,
 309                            geometryName,
 310                            &entity.GeometryType) < 0)
 311       {
 312       vtkErrorMacro("MEDmeshEntityInfo");
 313       continue;
 314       }
 315
 316     entity.GeometryName = geometryName;
 317     med_int ncell = 0;
 318
 319     if(entity.GeometryType == MED_POLYGON)
 320       {
 321       // read the number of cells of this type
 322       ncell = MEDmeshnEntity(this->FileId,
 323                              meshName,
 324                              cs.TimeIt,
 325                              cs.IterationIt,
 326                              entity.EntityType,
 327                              entity.GeometryType,
 328                              MED_INDEX_NODE,
 329                              connectivity,
 330                              &changement,
 331                              &transformation ) - 1;
 332       }
 333     else if(entity.GeometryType == MED_POLYHEDRON)
 334       {
 335       // read the number of cells of this type
 336       ncell = MEDmeshnEntity(this->FileId,
 337                              meshName,
 338                              cs.TimeIt,
 339                              cs.IterationIt,
 340                              entity.EntityType,
 341                              entity.GeometryType,
 342                              MED_INDEX_FACE,
 343                              connectivity,
 344                              &changement,
 345                              &transformation  ) - 1;
 346       }
 347     else
 348       {
 349       ncell = MEDmeshnEntity(this->FileId,
 350                              meshName,
 351                              cs.TimeIt,
 352                              cs.IterationIt,
 353                              entity.EntityType,
 354                              entity.GeometryType,
 355                              MED_CONNECTIVITY,
 356                              connectivity,
 357                              &changement,
 358                              &transformation  );
 359       }
 360
 361     if(ncell > 0)
 362       {
 363       vtkMedEntityArray* array = vtkMedEntityArray::New();
 364       array->SetParentGrid(grid);
 365       array->SetNumberOfEntity(ncell);
 366       array->SetEntity(entity);
 367       array->SetConnectivity(connectivity);
 368       grid->AppendEntityArray(array);
 369       array->Delete();
 370       // this triggers the creation of undefined families
 371       this->LoadFamilyIds(array);
 372       }
 373     }
 374 }
 375
 376 // Description:
 377 // load all Information data associated with this unstructured grid.
 378 void vtkMedDriver30::ReadUnstructuredGridInformation(
 379     vtkMedUnstructuredGrid* grid)
 380 {
 381   FileOpen open(this);
 382
 383   vtkMedMesh *mesh = grid->GetParentMesh();
 384
 385   const char *meshName = mesh->GetName();
 386   med_connectivity_mode connectivity;
 387
 388   med_bool changement;
 389   med_bool transformation;
 390   med_int profilesize;
 391
 392   char profilename[MED_NAME_SIZE+1];
 393   memset(profilename, '\0', MED_NAME_SIZE+1);
 394
 395   vtkMedComputeStep cs = grid->GetComputeStep();
 396
 397   // first check if we have points
 398   vtkIdType npoints = MEDmeshnEntityWithProfile(
 399                         this->FileId,
 400                         meshName,
 401                         cs.TimeIt,
 402                         cs.IterationIt,
 403                         MED_NODE,
 404                         MED_NONE,
 405                         MED_COORDINATE,
 406                         MED_NODAL,
 407                         MED_COMPACT_PFLMODE,
 408                         profilename,
 409                         &profilesize,
 410                         &changement,
 411                         &transformation);
 412
 413   if(npoints > 0)
 414     {
 415     grid->SetNumberOfPoints(npoints);
 416     }
 417   else
 418     {
 419     if(grid->GetPreviousGrid() == NULL)
 420       {
 421       vtkErrorMacro("No point and no previous grid");
 422       }
 423     grid->SetUsePreviousCoordinates(true);
 424     grid->SetNumberOfPoints(grid->GetPreviousGrid()->GetNumberOfPoints());
 425     return;
 426     }
 427
 428   this->ReadNumberOfEntity(grid, MED_CELL, MED_NODAL);
 429   this->ReadNumberOfEntity(grid, MED_CELL, MED_DESCENDING);
 430   this->ReadNumberOfEntity(grid, MED_DESCENDING_FACE, MED_NODAL);
 431   this->ReadNumberOfEntity(grid, MED_DESCENDING_FACE, MED_DESCENDING);
 432   this->ReadNumberOfEntity(grid, MED_DESCENDING_EDGE, MED_NODAL);
 433   this->ReadNumberOfEntity(grid, MED_DESCENDING_EDGE, MED_DESCENDING);
 434   this->ReadNumberOfEntity(grid, MED_STRUCT_ELEMENT, MED_NODAL);
 435   this->ReadNumberOfEntity(grid, MED_STRUCT_ELEMENT, MED_DESCENDING);
 436
 437   // create the point vtkMedEntityArray support
 438   vtkMedEntity entity;
 439   entity.EntityType = MED_NODE;
 440   entity.GeometryType = MED_POINT1;
 441   vtkMedEntityArray* pea = vtkMedEntityArray::New();
 442   pea->SetEntity(entity);
 443   pea->SetParentGrid(grid);
 444   pea->SetNumberOfEntity(grid->GetNumberOfPoints());
 445   grid->AppendEntityArray(pea);
 446   pea->Delete();
 447
 448   this->LoadFamilyIds(pea);
 449 }
 450
 451 void vtkMedDriver30::ReadFamilyInformation(vtkMedMesh* mesh, vtkMedFamily* family)
 452 {
 453   FileOpen open(this);
 454
 455   med_int familyid;
 456   med_int ngroup;
 457   char* groupNames = NULL;
 458   const  char* meshName = mesh->GetName();
 459
 460   ngroup = MEDnFamilyGroup(FileId, meshName, family->GetMedIterator());
 461
 462   bool has_no_group = false;
 463   if(ngroup <= 0)
 464     {
 465     if(ngroup < 0)
 466       {
 467       vtkErrorMacro("Error while reading the number of groups");
 468       }
 469     ngroup = 1;
 470     has_no_group = true;
 471     }
 472
 473   groupNames = new char[ngroup * MED_LNAME_SIZE + 1];
 474   memset(groupNames, '\0', ngroup * MED_LNAME_SIZE + 1);
 475
 476   // special case for files written by med < 3,
 477   // I have to use the 23 interface
 478   if(mesh->GetParentFile()->GetVersionMajor() < 3)
 479     {
 480     med_int *attributenumber = NULL;
 481     med_int *attributevalue = NULL;
 482     char *attributedes = NULL;
 483
 484     med_int nattr = MEDnFamily23Attribute(
 485                       this->FileId,
 486                       meshName,
 487                       family->GetMedIterator());
 488
 489     if(nattr < 0)
 490       {
 491       vtkErrorMacro("MEDnFamily23Attribute");
 492       }
 493
 494     if(nattr > 0)
 495       {
 496       attributenumber = new med_int[nattr];
 497       attributevalue = new med_int[nattr];
 498       attributedes = new char[nattr*MED_COMMENT_SIZE+1];
 499       memset(attributedes, '\0', nattr*MED_COMMENT_SIZE+1);
 500       }
 501
 502     char familyName[MED_NAME_SIZE+1] = "";
 503
 504     if(MEDfamily23Info (this->FileId,
 505                         meshName,
 506                         family->GetMedIterator(),
 507                         familyName,
 508                         attributenumber,
 509                         attributevalue,
 510                         attributedes,
 511                         &familyid,
 512                         groupNames ) < 0)
 513       {
 514       vtkDebugMacro("MEDfamily23Info");
 515       }
 516
 517     family->SetName(familyName);
 518
 519     if(attributenumber != NULL)
 520       delete[] attributenumber;
 521     if(attributevalue != NULL)
 522       delete[] attributevalue;
 523     if(attributedes != NULL)
 524       delete[] attributedes;
 525     }
 526   else
 527     {
 528     char familyName[MED_NAME_SIZE+1] = "";
 529     if(MEDfamilyInfo( this->FileId,
 530                       meshName,
 531                       family->GetMedIterator(),
 532                       familyName,
 533                       &familyid,
 534                       groupNames ) < 0)
 535       {
 536       vtkErrorMacro(
 537           "vtkMedDriver30::ReadInformation(vtkMedFamily* family)"
 538           << " cannot read family informations.");
 539       return;
 540       }
 541     family->SetName(familyName);
 542     }
 543
 544   family->SetId(familyid);
 545
 546   if( familyid == 0 ) {
 547     family->SetPointOrCell( vtkMedUtilities::OnCell );
 548     mesh->AppendCellFamily( family );
 549   } else {
 550     //rnv: improve algorithm to determine entity of the family:
 551     //     1) Read Nb nodes
 552     //     2) Read families of the nodes
 553     //     3) If list of the families of the nodes contains familyid => vtkMedUtilities::OnPoint
 554     //     otherwise => vtkMedUtilities::OnCell
 555     med_bool  v1, v2;
 556     med_int size = MEDmeshnEntity(
 557                                   this->FileId,
 558                                   meshName,
 559                                   MED_NO_DT,
 560                                   MED_NO_IT,
 561                                   MED_NODE,
 562                                   MED_NO_GEOTYPE,
 563                                   MED_COORDINATE,
 564                                   MED_NO_CMODE,
 565                                   &v1,
 566                                   &v2 );
 567     if( size < 0 ) {
 568       vtkErrorMacro( "vtkMedDriver30::ReadInformation(vtkMedFamily* family)"
 569                      <<" cannot read nb Nodes" );
 570       return;
 571     }
 572
 573     vector<med_int> n_fams;
 574     n_fams.resize( size );
 575
 576     med_int ret_val = MEDmeshEntityFamilyNumberRd( this->FileId,
 577                                                    meshName,
 578                                                    MED_NO_DT,
 579                                                    MED_NO_IT,
 580                                                    MED_NODE,
 581                                                    MED_NO_GEOTYPE ,
 582                                                    &n_fams[0] );
 583     // Remove ZERO FAMILY
 584     remove( n_fams.begin(),n_fams.end(), 0 );
 585
 586     bool isOnPoints = ( ret_val >= 0) && (find( n_fams.begin(), n_fams.end(), familyid ) != n_fams.end() );
 587     if( isOnPoints ) {
 588       family->SetPointOrCell(vtkMedUtilities::OnPoint);
 589       mesh->AppendPointFamily(family);
 590     } else {
 591       family->SetPointOrCell(vtkMedUtilities::OnCell);
 592       mesh->AppendCellFamily(family);
 593     }
 594   }
 595
 596   family->AllocateNumberOfGroup(ngroup);
 597   // if there where no group, set the name to the default value
 598   if(has_no_group)
 599     {
 600     memcpy(groupNames, vtkMedUtilities::NoGroupName,
 601            strlen(vtkMedUtilities::NoGroupName));
 602     }
 603
 604   for(int index = 0; index < ngroup; index++)
 605     {
 606     char realGroupName[MED_LNAME_SIZE + 1];
 607     memset(realGroupName, '\0', MED_LNAME_SIZE + 1);
 608     memcpy(realGroupName, groupNames + index * MED_LNAME_SIZE,
 609         MED_LNAME_SIZE * sizeof(char));
 610     vtkMedGroup* group = mesh->GetOrCreateGroup(family->GetPointOrCell(),
 611         realGroupName);
 612
 613     family->SetGroup(index, group);
 614     }
 615
 616   delete[] groupNames;
 617
 618   if(familyid == 0)
 619     {
 620     vtkMedFamily* famzero = vtkMedFamily::New();
 621     mesh->AppendPointFamily(famzero);
 622     famzero->Delete();
 623
 624     famzero->SetName(family->GetName());
 625     famzero->SetMedIterator(family->GetMedIterator());
 626     famzero->SetId(family->GetId());
 627     famzero->SetPointOrCell(vtkMedUtilities::OnPoint);
 628     famzero->AllocateNumberOfGroup(family->GetNumberOfGroup());
 629     for(int gid=0; gid<family->GetNumberOfGroup(); gid++)
 630       {
 631       vtkMedGroup* group = mesh->GetOrCreateGroup(
 632           vtkMedUtilities::OnPoint,
 633           family->GetGroup(gid)->GetName());
 634       famzero->SetGroup(gid, group);
 635       mesh->AppendPointGroup(group);
 636       }
 637     }
 638 }
 639
 640 void  vtkMedDriver30::ReadLinkInformation(vtkMedLink* link)
 641 {
 642   med_int size;
 643   char linkMeshName[MED_NAME_SIZE+1] = "";
 644   if(MEDlinkInfo(this->FileId,
 645                  link->GetMedIterator(),
 646                  linkMeshName,
 647                  &size) < 0)
 648     {
 649     vtkErrorMacro("MEDlinkInfo");
 650     return;
 651     }
 652   link->SetMeshName(linkMeshName);
 653   if(size <= 0)
 654     return;
 655
 656   char* path = new char[size + 1];
 657   memset(path, '\0', size+1);
 658   if(MEDlinkRd(this->FileId, link->GetMeshName(), path) < 0)
 659     {
 660     vtkErrorMacro("MEDlinkRd");
 661     memset(path, '\0', size+1);
 662     }
 663
 664   link->SetLink(path);
 665
 666   delete[] path;
 667 }
 668
 669 void vtkMedDriver30::ReadFileInformation(vtkMedFile* file)
 670 {
 671   FileOpen open(this);
 672
 673   char comment[MED_COMMENT_SIZE+1] = "";
 674
 675   MEDfileCommentRd(this->FileId,
 676                   comment);
 677
 678   file->SetComment(comment);
 679
 680   med_int major, minor, release;
 681   MEDfileNumVersionRd(this->FileId, &major, &minor, &release);
 682   file->SetVersionMajor(major);
 683   file->SetVersionMinor(minor);
 684   file->SetVersionRelease(release);
 685
 686   int nlink = MEDnLink(this->FileId);
 687   file->AllocateNumberOfLink(nlink);
 688   for(int linkid=0; linkid<nlink; linkid++)
 689     {
 690     vtkMedLink* link = file->GetLink(linkid);
 691     link->SetMedIterator(linkid+1);
 692     this->ReadLinkInformation(link);
 693     }
 694
 695   int nprof = MEDnProfile(FileId);
 696   // Reading id s not possible in parallel if the file contains Profiles
 697   vtkMultiProcessController* controller = vtkMultiProcessController::GetGlobalController();
 698   if (controller != NULL)
 699     if ((nprof != 0) && (controller->GetNumberOfProcesses() > 1))
 700     {
 701       vtkWarningMacro("ATTENTION: The MED Reader cannot read profiles when used in parallel");
 702     return;
 703     }
 704   file->AllocateNumberOfProfile(nprof);
 705   for(int i = 0; i < nprof; i++)
 706     {
 707     vtkMedProfile* profile = file->GetProfile(i);
 708     profile->SetMedIterator(i + 1);
 709     profile->SetParentFile(file);
 710     this->ReadProfileInformation(profile);
 711     }
 712
 713   int nloc = MEDnLocalization(this->FileId);
 714   file->AllocateNumberOfLocalization(nloc);
 715   for(int i = 0; i < nloc; i++)
 716     {
 717     vtkMedLocalization* loc = file->GetLocalization(i);
 718     loc->SetMedIterator(i + 1);
 719     loc->SetParentFile(file);
 720     this->ReadLocalizationInformation(loc);
 721     }
 722
 723   int nsupportmesh = MEDnSupportMesh(this->FileId);
 724   file->AllocateNumberOfSupportMesh(nsupportmesh);
 725   for(int i = 0; i < nsupportmesh; i++)
 726     {
 727     vtkMedMesh* supportmesh = file->GetSupportMesh(i);
 728     supportmesh->SetMedIterator(i + 1);
 729     supportmesh->SetParentFile(file);
 730     this->ReadSupportMeshInformation(supportmesh);
 731     }
 732
 733   int nmesh = MEDnMesh(this->FileId);
 734   file->AllocateNumberOfMesh(nmesh);
 735   for(int i = 0; i < nmesh; i++)
 736     {
 737     vtkMedMesh* mesh = file->GetMesh(i);
 738     mesh->SetMedIterator(i + 1);
 739     mesh->SetParentFile(file);
 740     this->ReadMeshInformation(mesh);
 741     }
 742
 743   int nfields = MEDnField(this->FileId);
 744   file->AllocateNumberOfField(nfields);
 745   for(int i = 0; i < nfields; i++)
 746     {
 747     vtkMedField* field = file->GetField(i);
 748     field->SetMedIterator(i + 1);
 749     field->SetParentFile(file);
 750     this->ReadFieldInformation(field);
 751     field->ComputeFieldType();
 752     while(field->HasManyFieldTypes())
 753       {
 754       vtkMedField* newfield = vtkMedField::New();
 755       int type = field->GetFirstType();
 756       newfield->ExtractFieldType(field, type);
 757       file->AppendField(newfield);
 758       newfield->Delete();
 759       }
 760     }
 761
 762   int nstruct = MEDnStructElement(this->FileId);
 763
 764   file->AllocateNumberOfStructElement(nstruct);
 765   for(int i = 0; i < nstruct; i++)
 766     {
 767     vtkMedStructElement* structelem = file->GetStructElement(i);
 768     structelem->SetMedIterator(i+1);
 769     structelem->SetParentFile(file);
 770     this->ReadStructElementInformation(structelem);
 771     }
 772
 773 }
 774
 775 void  vtkMedDriver30::ReadStructElementInformation(
 776     vtkMedStructElement* structelem)
 777 {
 778
 779   FileOpen open(this);
 780
 781   char modelname[MED_NAME_SIZE+1] = "";
 782   med_geometry_type mgeotype;
 783   med_int modeldim;
 784   char supportmeshname[MED_NAME_SIZE+1] = "";
 785   med_entity_type sentitytype;
 786   med_int snbofnode;
 787   med_int snbofcell;
 788   med_geometry_type sgeotype;
 789   med_int nbofconstantattribute;
 790   med_bool anyprofile;
 791   med_int nbofvariableattribute;
 792
 793   if(MEDstructElementInfo (this->FileId,
 794                            structelem->GetMedIterator(),
 795                            modelname,
 796                            &mgeotype,
 797                            &modeldim,
 798                            supportmeshname,
 799                            &sentitytype,
 800                            &snbofnode,
 801                            &snbofcell,
 802                            &sgeotype,
 803                            &nbofconstantattribute,
 804                            &anyprofile,
 805                            &nbofvariableattribute) < 0)
 806     {
 807     vtkErrorMacro("Error in MEDstructElementInfo");
 808     return;
 809     }
 810   structelem->SetName(modelname);
 811   structelem->SetGeometryType(mgeotype);
 812   structelem->SetModelDimension(modeldim);
 813   structelem->SetSupportMeshName(supportmeshname);
 814   structelem->SetSupportEntityType(sentitytype);
 815   structelem->SetSupportNumberOfNode(snbofnode);
 816   structelem->SetSupportNumberOfCell(snbofcell);
 817   structelem->SetSupportGeometryType(sgeotype);
 818   structelem->AllocateNumberOfConstantAttribute(nbofconstantattribute);
 819   structelem->AllocateNumberOfVariableAttribute(nbofvariableattribute);
 820   structelem->SetAnyProfile(anyprofile);
 821
 822   for(int attit = 0; attit < nbofconstantattribute; attit ++)
 823     {
 824     vtkMedConstantAttribute* constatt = structelem->GetConstantAttribute(attit);
 825     constatt->SetMedIterator(attit+1);
 826     constatt->SetParentStructElement(structelem);
 827     this->ReadConstantAttributeInformation(constatt);
 828     }
 829
 830   for(int attit = 0; attit < nbofvariableattribute; attit ++)
 831     {
 832     vtkMedVariableAttribute* varatt = structelem->GetVariableAttribute(attit);
 833     varatt->SetMedIterator(attit+1);
 834     varatt->SetParentStructElement(structelem);
 835     this->ReadVariableAttributeInformation(varatt);
 836     }
 837 }
 838
 839 void vtkMedDriver30::ReadConstantAttributeInformation(
 840     vtkMedConstantAttribute* constAttr)
 841 {
 842
 843   FileOpen open(this);
 844
 845   char constattname[MED_NAME_SIZE+1] = "";
 846   med_attribute_type constatttype;
 847   med_int nbofcomponent;
 848   med_entity_type sentitytype;
 849   char profilename[MED_NAME_SIZE+1] = "";
 850   med_int profilesize;
 851
 852   if(MEDstructElementConstAttInfo(
 853       this->FileId,
 854       constAttr->GetParentStructElement()->GetName(),
 855       constAttr->GetMedIterator(),
 856       constattname,
 857       &constatttype,
 858       &nbofcomponent,
 859       &sentitytype,
 860       profilename,
 861       &profilesize)   < 0)
 862     {
 863     vtkErrorMacro("MEDstructElementConstAttInfo error");
 864     return;
 865     }
 866
 867   constAttr->SetName(constattname);
 868   constAttr->SetAttributeType(constatttype);
 869   constAttr->SetNumberOfComponent(nbofcomponent);
 870   constAttr->SetSupportEntityType(sentitytype);
 871   constAttr->SetProfileName(profilename);
 872   constAttr->SetProfileSize(profilesize);
 873
 874   vtkAbstractArray* values = vtkMedUtilities::NewArray(constatttype);
 875   if(values == NULL)
 876     return;
 877   constAttr->SetValues(values);
 878   values->Delete();
 879
 880   values->SetNumberOfComponents(nbofcomponent);
 881   vtkIdType ntuple = 0;
 882   if((strcmp(profilename, MED_NO_PROFILE) != 0) &&
 883      (strcmp(profilename, "\0") != 0))
 884     {
 885     ntuple = profilesize;
 886     }
 887   else if(constAttr->GetSupportEntityType() == MED_CELL)
 888     {
 889     ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfCell();
 890     }
 891   else
 892     {
 893     ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfNode();
 894     }
 895   values->SetNumberOfTuples(ntuple);
 896
 897   void* ptr = NULL;
 898   vtkSmartPointer<vtkCharArray> buffer = vtkSmartPointer<vtkCharArray>::New();
 899   if(constatttype != MED_ATT_NAME)
 900     {
 901     ptr = values->GetVoidPointer(0);
 902     }
 903   else
 904     {
 905     buffer->SetNumberOfValues(MED_NAME_SIZE*nbofcomponent*ntuple);
 906     ptr = buffer->GetVoidPointer(0);
 907     }
 908
 909   if(MEDstructElementConstAttRd (this->FileId,
 910         constAttr->GetParentStructElement()->GetName(),
 911         constAttr->GetName(), ptr) < 0)
 912     {
 913     vtkErrorMacro("MEDstructElementConstAttRd");
 914     return;
 915     }
 916
 917   if(constatttype == MED_ATT_NAME)
 918     {
 919     char name[MED_NAME_SIZE+1] = "";
 920     char* nameptr = (char*) ptr;
 921     vtkStringArray* names = vtkStringArray::SafeDownCast(values);
 922     for(vtkIdType id = 0; id < nbofcomponent*ntuple; id++)
 923       {
 924       memset(name, '\0', MED_NAME_SIZE+1);
 925       strncpy(name, nameptr + id * MED_NAME_SIZE, MED_NAME_SIZE);
 926       names->SetValue(id, name);
 927       }
 928     }
 929
 930   return;
 931 }
 932
 933 void vtkMedDriver30::ReadVariableAttributeInformation(
 934     vtkMedVariableAttribute* varAttr)
 935 {
 936
 937   FileOpen open(this);
 938
 939   char varattname[MED_NAME_SIZE+1] = "";
 940   med_attribute_type varatttype;
 941   med_int nbofcomponent;
 942
 943   if(MEDstructElementVarAttInfo (
 944       this->FileId,
 945       varAttr->GetParentStructElement()->GetName(),
 946       varAttr->GetMedIterator(),
 947       varattname,
 948       &varatttype,
 949       &nbofcomponent) < 0)
 950     {
 951     vtkErrorMacro("MEDstructElementVarAttInfo");
 952     return;
 953     }
 954
 955   varAttr->SetName(varattname);
 956   varAttr->SetAttributeType(varatttype);
 957   varAttr->SetNumberOfComponent(nbofcomponent);
 958
 959   return;
 960 }
 961
 962 void vtkMedDriver30::ReadProfileInformation(vtkMedProfile* profile)
 963 {
 964   FileOpen open(this);
 965
 966   med_int nelem;
 967   char profileName[MED_NAME_SIZE+1] = "";
 968
 969   if(MEDprofileInfo(this->FileId,
 970                     profile->GetMedIterator(),
 971                     profileName,
 972                     &nelem) < 0)
 973     {
 974     vtkErrorMacro("cannot read information on profile"
 975         << profile->GetMedIterator());
 976     }
 977   profile->SetName(profileName);
 978   profile->SetNumberOfElement(nelem);
 979 }
 980
 981 void vtkMedDriver30::ReadFieldInformation(vtkMedField* field)
 982 {
 983   FileOpen open(this);
 984
 985   if (field->GetMedIterator() == 0)
 986     return;
 987
 988   int ncomp = MEDfieldnComponent(FileId, field->GetMedIterator());
 989
 990   if(ncomp < 0)
 991     {
 992     field->SetNumberOfComponent(-1);
 993     vtkErrorMacro("cannot read the number of component of field "
 994         << field->GetMedIterator())
 995     return;
 996     }
 997
 998   field->SetNumberOfComponent(ncomp);
 999
1000   char* units = new char[MED_SNAME_SIZE * ncomp + 1];
1001   char* componentNames = new char[MED_SNAME_SIZE * ncomp + 1];
1002   memset(units, '\0', MED_SNAME_SIZE * ncomp + 1);
1003   memset(componentNames, '\0', MED_SNAME_SIZE * ncomp + 1);
1004
1005   //med_type_champ dataType;
1006   med_field_type dataType;
1007   med_int nstep;
1008   med_bool localmesh;
1009
1010   char name[MED_NAME_SIZE+1] = "";
1011   char meshName[MED_NAME_SIZE+1] = "";
1012   char unit[MED_SNAME_SIZE+1] = "";
1013
1014   if( MEDfieldInfo( FileId,
1015                     field->GetMedIterator(),
1016                     name,
1017                     meshName,
1018                     &localmesh,
1019                     &dataType,
1020                     componentNames,
1021                     units,
1022                     unit,
1023                     &nstep) < 0)
1024     {
1025     vtkErrorMacro("cannot read the informations on field "
1026         << field->GetMedIterator())
1027     return;
1028     }
1029
1030   field->SetName(name);
1031   field->SetMeshName(meshName);
1032   field->SetTimeUnit(unit);
1033   field->SetDataType(dataType);
1034   field->SetLocal(localmesh);
1035
1036   for(int comp = 0; comp < ncomp; comp++)
1037     {
1038     char unit[MED_NAME_SIZE + 1] = "";
1039     memcpy(unit, units + MED_SNAME_SIZE * comp, MED_SNAME_SIZE * sizeof(char));
1040     field->GetUnit()->SetValue(comp, unit);
1041
1042     char compName[MED_SNAME_SIZE + 1] = "";
1043     memcpy(compName, componentNames + MED_SNAME_SIZE * comp, MED_SNAME_SIZE
1044         * sizeof(char));
1045     field->GetComponentName()->SetValue(comp, compName);
1046     }
1047
1048   delete[] units;
1049   delete[] componentNames;
1050
1051   med_int ninterp = MEDfieldnInterp(FileId, field->GetName());
1052   if(ninterp < 0)
1053     {
1054     vtkErrorMacro("Error in MEDfieldnInterp");
1055     return;
1056     }
1057
1058   field->AllocateNumberOfInterpolation(ninterp);
1059
1060   for(med_int interpit=0; interpit<ninterp; interpit++)
1061     {
1062     vtkMedInterpolation* interp = field->GetInterpolation(interpit);
1063     interp->SetMedIterator(interpit + 1);
1064     this->ReadInterpolationInformation(interp);
1065     }
1066
1067   vtkMedFieldStep* previousStep = NULL;
1068
1069   for(med_int csit = 0; csit < nstep; csit++)
1070     {
1071     vtkMedFieldStep* step = vtkMedFieldStep::New();
1072     step->SetMedIterator(csit + 1);
1073     step->SetParentField(field);
1074     this->ReadFieldStepInformation(step, true);
1075     field->AddFieldStep(step);
1076     step->SetPreviousStep(previousStep);
1077     previousStep = step;
1078     step->Delete();
1079     }
1080 }
1081
1082 void vtkMedDriver30::ReadFieldStepInformation(vtkMedFieldStep* step, bool readAllEntityInfo)
1083 {
1084   vtkMedComputeStep cs;
1085   vtkMedComputeStep meshcs;
1086   vtkMedField* field = step->GetParentField();
1087
1088   FileOpen open(this);
1089
1090   if( MEDfieldComputingStepMeshInfo(
1091         FileId,
1092         field->GetName(),
1093         step->GetMedIterator(),
1094         &cs.TimeIt,
1095         &cs.IterationIt,
1096         &cs.TimeOrFrequency,
1097         &meshcs.TimeIt,
1098         &meshcs.IterationIt) < 0)
1099     {
1100     vtkErrorMacro("Error in MEDfieldComputingStepMeshInfo");
1101     return;
1102     }
1103
1104   step->SetComputeStep(cs);
1105   step->SetMeshComputeStep(meshcs);
1106
1107   if(!readAllEntityInfo || step->GetEntityInfoLoaded())
1108     return;
1109
1110   step->SetEntityInfoLoaded(1);
1111
1112   vtkMedFile* file = field->GetParentFile();
1113   vtkMedMesh* mesh = file->GetMesh(field->GetMeshName());
1114
1115   if(mesh == NULL)
1116     return;
1117
1118   std::set<vtkMedEntity> tmp_entities;
1119   std::set<vtkMedEntity> entities;
1120   mesh->GatherMedEntities(tmp_entities);
1121
1122   std::set<vtkMedEntity>::iterator tmp_entity_it = tmp_entities.begin();
1123   while(tmp_entity_it != tmp_entities.end())
1124     {
1125     vtkMedEntity entity = *tmp_entity_it;
1126     tmp_entity_it++;
1127     entities.insert(entity);
1128     if(entity.EntityType == MED_CELL)
1129       {
1130       vtkMedEntity newEntity;
1131       newEntity.EntityType = MED_NODE_ELEMENT;
1132       newEntity.GeometryType = entity.GeometryType;
1133       newEntity.GeometryName = entity.GeometryName;
1134       entities.insert(newEntity);
1135       }
1136     }
1137
1138   vtkMedEntity pointEntity;
1139   pointEntity.EntityType = MED_NODE;
1140   pointEntity.GeometryType = MED_NONE;
1141   pointEntity.GeometryName = MED_NAME_BLANK;
1142   entities.insert(pointEntity);
1143
1144   std::set<vtkMedEntity>::iterator entity_it = entities.begin();
1145   while(entity_it != entities.end())
1146     {
1147     vtkMedEntity entity = *entity_it;
1148     entity_it++;
1149
1150     med_int nvalues = 0;
1151     med_int nprofile;
1152     char profilename[MED_NAME_SIZE+1] = "";
1153     char localizationname[MED_NAME_SIZE+1] = "";
1154
1155     nprofile = MEDfieldnProfile(
1156         this->FileId,
1157         field->GetName(),
1158         step->GetComputeStep().TimeIt,
1159         step->GetComputeStep().IterationIt,
1160         entity.EntityType,
1161         entity.GeometryType,
1162         profilename,
1163         localizationname);
1164     if(nprofile < 0)
1165       {
1166       vtkErrorMacro("MEDfieldnProfile");
1167       continue;
1168       }
1169
1170     bool hasprofile = (nprofile > 0);
1171     if(!hasprofile)
1172       {
1173       nprofile = 1;
1174       }
1175
1176     med_int profilesize;
1177     med_int nintegrationpoint;
1178
1179     for(int pid=0; pid<nprofile; pid++)
1180       {
1181       med_int medid = (hasprofile ? pid+1 : -1);
1182       nvalues = MEDfieldnValueWithProfile(
1183                   this->FileId,
1184                   field->GetName(),
1185                   step->GetComputeStep().TimeIt,
1186                   step->GetComputeStep().IterationIt,
1187                   entity.EntityType,
1188                   entity.GeometryType,
1189                   medid,
1190                   MED_COMPACT_PFLMODE,
1191                   profilename,
1192                   &profilesize,
1193                   localizationname,
1194                   &nintegrationpoint);
1195
1196       if(nvalues < 0)
1197         {
1198         vtkErrorMacro("MEDfieldnValueWithProfile");
1199         continue;
1200         }
1201       else if(nvalues > 0)
1202         {
1203         // I have found a profile with values, stop the loop here
1204         break;
1205         }
1206       }
1207
1208     if(nvalues > 0)
1209       {
1210       vtkMedFieldOverEntity* fieldOverEntity = vtkMedFieldOverEntity::New();
1211       step->AppendFieldOverEntity(fieldOverEntity);
1212       fieldOverEntity->Delete();
1213
1214       fieldOverEntity->SetParentStep(step);
1215       fieldOverEntity->SetEntity(entity);
1216
1217       this->ReadFieldOverEntityInformation(fieldOverEntity);
1218       }
1219     }
1220 }
1221
1222 void vtkMedDriver30::ReadFieldOverEntityInformation(vtkMedFieldOverEntity* fieldOverEntity)
1223 {
1224   FileOpen open(this);
1225
1226   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
1227   vtkMedField* field = step->GetParentField();
1228   vtkMedEntity entity = fieldOverEntity->GetEntity();
1229
1230   const char* fieldName = field->GetName();
1231   const vtkMedComputeStep& cs = step->GetComputeStep();
1232
1233   char profilename[MED_NAME_SIZE+1] = "";
1234   char localizationname[MED_NAME_SIZE+1] = "";
1235
1236   med_int nProfiles = MEDfieldnProfile(
1237                         this->FileId,
1238                         fieldName,
1239                         cs.TimeIt,
1240                         cs.IterationIt,
1241                         entity.EntityType,
1242                         entity.GeometryType,
1243                         profilename,
1244                         localizationname);
1245
1246   if(nProfiles < 0)
1247     {
1248     vtkErrorMacro("MEDfieldnProfile");
1249     }
1250   else if(nProfiles == 0)
1251     {
1252     fieldOverEntity->SetHasProfile(0);
1253     nProfiles = 1;
1254     }
1255   else
1256     {
1257     fieldOverEntity->SetHasProfile(1);
1258     }
1259   fieldOverEntity->AllocateNumberOfFieldOnProfile(nProfiles);
1260   for(int profit = 0; profit < nProfiles; profit++)
1261     {
1262     vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(profit);
1263     med_int medid = (fieldOverEntity->GetHasProfile()? profit+1: -1);
1264     fop->SetMedIterator(medid);
1265     fop->SetParentFieldOverEntity(fieldOverEntity);
1266     this->ReadFieldOnProfileInformation(fop);
1267     }
1268 }
1269
1270 void vtkMedDriver30::ReadFieldOnProfileInformation(vtkMedFieldOnProfile* fop)
1271 {
1272   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
1273   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
1274   vtkMedField* field = step->GetParentField();
1275
1276   const vtkMedComputeStep& cs = step->GetComputeStep();
1277   med_int profilesize;
1278   med_int nbofintegrationpoint;
1279
1280   char profileName[MED_NAME_SIZE+1] = "";
1281   char localizationName[MED_NAME_SIZE+1] = "";
1282
1283   med_int nvalue = MEDfieldnValueWithProfile(FileId,
1284                     field->GetName(),
1285                     cs.TimeIt,
1286                     cs.IterationIt,
1287                     fieldOverEntity->GetEntity().EntityType,
1288                     fieldOverEntity->GetEntity().GeometryType,
1289                     fop->GetMedIterator(),
1290                     MED_COMPACT_STMODE,
1291                     profileName,
1292                     &profilesize,
1293                     localizationName,
1294                     &nbofintegrationpoint);
1295
1296   if(nvalue < 0)
1297     {
1298     vtkErrorMacro("Error while reading MEDfieldnValueWithProfile");
1299     }
1300
1301   fop->SetProfileName(profileName);
1302   fop->SetLocalizationName(localizationName);
1303   fop->SetNumberOfValues(nvalue);
1304   fop->SetNumberOfIntegrationPoint(nbofintegrationpoint);
1305   fop->SetProfileSize(profilesize);
1306 }
1307
1308 void vtkMedDriver30::ReadMeshInformation(vtkMedMesh* mesh)
1309 {
1310   FileOpen open(this);
1311
1312   med_int mdim = 0;
1313   med_int sdim = 0;
1314   med_mesh_type meshtype;
1315
1316   med_sorting_type sortingtype;
1317   med_int nstep = 0;
1318   med_axis_type axistype;
1319   med_int naxis;
1320
1321   if ( (naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator())) <0 )
1322     {
1323     vtkDebugMacro("Error reading mesh axis number");
1324     }
1325
1326   if(naxis == 0)
1327     naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator());
1328
1329   char axisname[3*MED_SNAME_SIZE+1] = "";
1330   char axisunit[3*MED_SNAME_SIZE+1] = "";
1331   char name[MED_NAME_SIZE+1] = "";
1332   char description[MED_COMMENT_SIZE+1] = "";
1333   char timeUnit[MED_SNAME_SIZE+1] = "";
1334
1335   if( MEDmeshInfo( this->FileId,
1336         mesh->GetMedIterator(),
1337         name,
1338         &sdim,
1339         &mdim,
1340         &meshtype,
1341         description,
1342         timeUnit,
1343         &sortingtype,
1344         &nstep,
1345         &axistype,
1346         axisname,
1347         axisunit ) )
1348     {
1349     vtkErrorMacro("Error reading mesh");
1350     }
1351
1352   mesh->SetName(name);
1353   mesh->SetDescription(description);
1354   mesh->SetTimeUnit(timeUnit);
1355   mesh->SetSpaceDimension(sdim);
1356   mesh->SetMeshDimension(mdim);
1357   mesh->SetMeshType(meshtype);
1358   mesh->SetSortingType(sortingtype);
1359   mesh->SetAxisType(axistype);
1360   mesh->SetNumberOfAxis(naxis);
1361
1362   for(int axis = 0; axis < naxis; axis++)
1363     {
1364     char theaxisname[MED_SNAME_SIZE+1] = "";
1365     char theaxisunit[MED_SNAME_SIZE+1] = "";
1366     strncpy(theaxisname, axisname + axis*MED_SNAME_SIZE, MED_SNAME_SIZE);
1367     strncpy(theaxisunit, axisunit + axis*MED_SNAME_SIZE, MED_SNAME_SIZE);
1368     mesh->GetAxisName()->SetValue(axis, theaxisname);
1369     mesh->GetAxisUnit()->SetValue(axis, theaxisunit);
1370     }
1371
1372   char universalName[MED_LNAME_SIZE+1] = "";
1373
1374   if(MEDmeshUniversalNameRd(this->FileId, name,
1375       universalName) < 0)
1376     {
1377     vtkDebugMacro("MEDmeshUniversalNameRd < 0");
1378     }
1379   mesh->SetUniversalName(universalName);
1380
1381   // read the Information data of all families.
1382   // writing the family 0 is optional,
1383   // but I need it, so add it if it is not here.
1384
1385   med_int nfam = MEDnFamily(this->FileId, name);
1386
1387   for(int index = 0; index < nfam; index++)
1388     {
1389     vtkMedFamily* family = vtkMedFamily::New();
1390     family->SetMedIterator(index + 1);
1391     this->ReadFamilyInformation(mesh, family);
1392     family->Delete();
1393     }
1394
1395   // this creates a family 0 if none has been read
1396   vtkMedFamily* familyZeroOnCell = mesh->GetOrCreateCellFamilyById(0);
1397   vtkMedFamily* familyZeroOnPoint = mesh->GetOrCreatePointFamilyById(0);
1398
1399   // Load Information regarding the grid type
1400   if(meshtype == MED_STRUCTURED_MESH)
1401     {
1402     // Do it for structured data
1403     med_grid_type mtg;
1404     if(MEDmeshGridTypeRd(FileId, name, &mtg) < 0)
1405       {
1406       vtkErrorMacro("Error during structured grid Information loading.");
1407       return;
1408       }
1409     mesh->SetStructuredGridType(mtg);
1410     }
1411
1412   vtkMedGrid* previousGrid = NULL;
1413   for(int gid=1; gid <= nstep; gid++)
1414     {
1415     vtkMedComputeStep cs;
1416     if(MEDmeshComputationStepInfo(FileId,
1417                                   name,
1418                                   gid,
1419                                   &cs.TimeIt,
1420                                   &cs.IterationIt,
1421                                   &cs.TimeOrFrequency) < 0)
1422       {
1423       vtkErrorMacro("MEDmeshComputationStepInfo error");
1424       }
1425     // Load Information regarding the grid type
1426     vtkMedGrid* grid = NULL;
1427     if(meshtype == MED_STRUCTURED_MESH)
1428       {
1429       switch(mesh->GetStructuredGridType())
1430         {
1431         case MED_CARTESIAN_GRID:
1432           grid = vtkMedCartesianGrid::New();
1433           break;
1434         case MED_POLAR_GRID:
1435           grid = vtkMedPolarGrid::New();
1436           break;
1437         case MED_CURVILINEAR_GRID:
1438           grid = vtkMedCurvilinearGrid::New();
1439           break;
1440         default:
1441           vtkErrorMacro("Unknown structured grid type " << mesh->GetStructuredGridType());
1442           return;
1443         }
1444       }
1445     else //(mesh->GetType() == MED_STRUCTURED_MESH)
1446       {
1447       grid = vtkMedUnstructuredGrid::New();
1448       }
1449     grid->SetParentMesh(mesh);
1450     grid->SetComputeStep(cs);
1451     this->ReadGridInformation(grid);
1452     mesh->AddGridStep(grid);
1453     grid->Delete();
1454     grid->SetPreviousGrid(previousGrid);
1455     previousGrid = grid;
1456     }
1457 }
1458
1459 void vtkMedDriver30::ReadLocalizationInformation(vtkMedLocalization* loc)
1460 {
1461   FileOpen open(this);
1462
1463   med_int ngp;
1464   med_int spaceDimension;
1465   med_geometry_type type_geo;
1466   med_geometry_type sectiongeotype;
1467   med_int nsectionmeshcell;
1468
1469   char name[MED_NAME_SIZE+1] = "";
1470   char interpolationName[MED_NAME_SIZE+1] = "";
1471   char sectionName[MED_NAME_SIZE+1] = "";
1472
1473   if(MEDlocalizationInfo(
1474       this->FileId,
1475       loc->GetMedIterator(),
1476       name,
1477       &type_geo,
1478       &spaceDimension,
1479       &ngp,
1480       interpolationName,
1481       sectionName,
1482       &nsectionmeshcell,
1483       &sectiongeotype ) < 0)
1484     {
1485     vtkErrorMacro("Reading information on quadrature points definition : "
1486         << loc->GetMedIterator());
1487     }
1488
1489   loc->SetName(name);
1490   loc->SetInterpolationName(interpolationName);
1491   loc->SetSectionName(sectionName);
1492   loc->SetNumberOfQuadraturePoint(ngp);
1493   loc->SetGeometryType(type_geo);
1494   loc->SetSpaceDimension(spaceDimension);
1495   loc->SetNumberOfCellInSection(nsectionmeshcell);
1496   loc->SetSectionGeometryType(sectiongeotype);
1497
1498   med_float *localCoordinates = new med_float[loc->GetSizeOfPointLocalCoordinates()];
1499   med_float *pqLocalCoordinates = new med_float[loc->GetSizeOfQuadraturePointLocalCoordinates()];
1500   med_float *weights = new med_float[loc->GetSizeOfWeights()];
1501
1502   if(MEDlocalizationRd(FileId,
1503       loc->GetName(),
1504       MED_FULL_INTERLACE,
1505       localCoordinates,
1506       pqLocalCoordinates,
1507       weights) < 0)
1508     {
1509     vtkErrorMacro("MEDlocalizationRd : " << loc->GetName());
1510     }
1511
1512   vtkDoubleArray* lc = loc->GetPointLocalCoordinates();
1513   vtkDoubleArray *pqlc = loc->GetQuadraturePointLocalCoordinates();
1514   vtkDoubleArray *w = loc->GetWeights();
1515
1516   lc->SetNumberOfValues(loc->GetSizeOfPointLocalCoordinates());
1517   for(int i=0; i<loc->GetSizeOfPointLocalCoordinates(); i++)
1518     {
1519     lc->SetValue(i, localCoordinates[i]);
1520     }
1521
1522   pqlc->SetNumberOfValues(loc->GetSizeOfQuadraturePointLocalCoordinates());
1523   for(int i=0; i<loc->GetSizeOfQuadraturePointLocalCoordinates(); i++)
1524     {
1525     pqlc->SetValue(i, pqLocalCoordinates[i]);
1526     }
1527
1528   w->SetNumberOfValues(loc->GetSizeOfWeights());
1529   for(int i=0; i<loc->GetSizeOfWeights(); i++)
1530     {
1531     w->SetValue(i, weights[i]);
1532     }
1533 }
1534
1535 void vtkMedDriver30::ReadInterpolationInformation(vtkMedInterpolation* interp)
1536 {
1537
1538   med_geometry_type geotype;
1539   med_bool cellnode;
1540   med_int nbofbasisfunc;
1541   med_int nbofvariable;
1542   med_int maxdegree;
1543   med_int nmaxcoef;
1544
1545   char name[MED_NAME_SIZE+1] = "";
1546
1547   if(MEDinterpInfo (this->FileId,
1548                     interp->GetMedIterator(),
1549                     name,
1550                     &geotype, &cellnode, &nbofbasisfunc,
1551                     &nbofvariable, &maxdegree, &nmaxcoef) < 0)
1552     {
1553     vtkErrorMacro("MEDinterpInfo");
1554     return;
1555     }
1556
1557   interp->SetName(name);
1558   interp->SetGeometryType(geotype);
1559   interp->SetIsCellNode(cellnode);
1560   interp->SetNumberOfVariable(nbofvariable);
1561   interp->SetMaximumDegree(maxdegree);
1562   interp->SetMaximumNumberOfCoefficient(nmaxcoef);
1563   interp->AllocateNumberOfBasisFunction(nbofbasisfunc);
1564
1565   for(int basisid=0; basisid < nbofbasisfunc; basisid++)
1566     {
1567     vtkMedFraction* func = interp->GetBasisFunction(basisid);
1568     func->SetNumberOfVariable(nbofvariable);
1569
1570     med_int ncoef = MEDinterpBaseFunctionCoefSize (
1571         this->FileId,
1572         interp->GetName(),
1573         basisid+1);
1574     func->SetNumberOfCoefficients(ncoef);
1575
1576     if(ncoef <= 0 || nbofvariable <= 0)
1577       continue;
1578
1579     med_int *power = new med_int[nbofvariable * ncoef];
1580     med_float *coefficient = new med_float[ncoef];
1581
1582     if(MEDinterpBaseFunctionRd  (
1583         this->FileId,
1584         interp->GetName(),
1585         basisid+1,
1586         &ncoef,
1587         power,
1588         coefficient) < 0)
1589       {
1590       vtkErrorMacro("MEDinterpBaseFunctionRd");
1591       continue;
1592       }
1593     vtkDoubleArray* coeffs = func->GetCoefficients();
1594     for(int cid=0; cid < ncoef; cid++)
1595       {
1596       coeffs->SetValue(cid, coefficient[cid]);
1597       }
1598     vtkIntArray* powers = func->GetPowers();
1599     for(int pid=0; pid < ncoef*nbofvariable; pid++)
1600       {
1601       powers->SetValue(pid, power[pid]);
1602       }
1603
1604     delete[] power;
1605     delete[] coefficient;
1606     }
1607 }
1608
1609 void vtkMedDriver30::ReadSupportMeshInformation(
1610     vtkMedMesh* supportMesh)
1611 {
1612   FileOpen open(this);
1613
1614   char supportmeshname[MED_NAME_SIZE+1] = "";
1615   char description[MED_COMMENT_SIZE+1] = "";
1616   med_int spacedim;
1617   med_int meshdim;
1618   med_axis_type axistype;
1619   char axisname[3*MED_SNAME_SIZE+1] = "";
1620   char axisunit[3*MED_SNAME_SIZE+1] = "";
1621
1622   if(MEDsupportMeshInfo (this->FileId,
1623                          supportMesh->GetMedIterator(),
1624                          supportmeshname,
1625                          &spacedim,
1626                          &meshdim,
1627                          description,
1628                          &axistype,
1629                          axisname,
1630                          axisunit) < 0)
1631     {
1632     vtkErrorMacro("MEDsupportMeshInfo");
1633     }
1634
1635   supportMesh->SetName(supportmeshname);
1636   supportMesh->SetDescription(description);
1637   supportMesh->SetSpaceDimension(spacedim);
1638   supportMesh->SetMeshDimension(meshdim);
1639   supportMesh->SetAxisType(axistype);
1640   for(int dim = 0; dim < 3; dim++)
1641     {
1642     char axisname_dim[MED_SNAME_SIZE+1] = "";
1643     char axisunit_dim[MED_SNAME_SIZE+1] = "";
1644
1645     strncpy(axisname_dim, axisname+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
1646     strncpy(axisunit_dim, axisunit+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
1647
1648     supportMesh->GetAxisName()->SetValue(dim, axisname_dim);
1649     supportMesh->GetAxisUnit()->SetValue(dim, axisunit_dim);
1650     }
1651
1652   return;
1653 }
1654
1655 void vtkMedDriver30::LoadFamilyIds(vtkMedEntityArray* array)
1656 {
1657   if(array->IsFamilyIdsLoaded())
1658     return;
1659
1660   FileOpen open(this);
1661
1662   vtkMedGrid* grid = array->GetParentGrid();
1663
1664   vtkMedComputeStep cs = grid->GetComputeStep();
1665
1666   // first, find if the family ids are implicit or explicit
1667   med_bool changement, transformation;
1668
1669   med_int nfamid = MEDmeshnEntity(this->FileId,
1670                       grid->GetParentMesh()->GetName(),
1671                       cs.TimeIt,
1672                       cs.IterationIt,
1673                       array->GetEntity().EntityType,
1674                       array->GetEntity().GeometryType,
1675                       MED_FAMILY_NUMBER,
1676                       MED_NO_CMODE,
1677                       &changement,
1678                       &transformation);
1679
1680   if(nfamid == array->GetNumberOfEntity())
1681     {
1682
1683     vtkMedIntArray* famIds = vtkMedIntArray::New();
1684     array->SetFamilyIds(famIds);
1685     famIds->Delete();
1686
1687     famIds->SetNumberOfTuples(nfamid);
1688
1689     if ( MEDmeshEntityFamilyNumberRd(
1690             this->FileId,
1691             grid->GetParentMesh()->GetName(),
1692             cs.TimeIt,
1693             cs.IterationIt,
1694             array->GetEntity().EntityType,
1695             array->GetEntity().GeometryType,
1696             famIds->GetPointer(0) ) < 0)
1697       {
1698       vtkWarningMacro("Error loading the family ids of entity "
1699         << array->GetEntity().EntityType
1700         << " " << array->GetEntity().GeometryType
1701         << " on mesh " << grid->GetParentMesh()->GetName());
1702       array->SetFamilyIds(NULL);
1703       }
1704     }
1705   else
1706     {
1707     vtkDebugMacro("NumberOfEntity != Number of family ids");
1708     array->SetFamilyIds(NULL);
1709     }
1710
1711   array->ComputeFamilies();
1712 }
1713
1714 void vtkMedDriver30::LoadPointGlobalIds(vtkMedGrid* grid)
1715 {
1716   if(grid->IsPointGlobalIdsLoaded())
1717     return;
1718
1719   FileOpen open(this);
1720
1721   vtkMedIntArray* globalIds = vtkMedIntArray::New();
1722   grid->SetPointGlobalIds(globalIds);
1723   globalIds->Delete();
1724
1725   globalIds->SetNumberOfTuples(grid->GetNumberOfPoints());
1726
1727   if( MEDmeshEntityNumberRd (
1728         this->FileId,
1729         grid->GetParentMesh()->GetName(),
1730         grid->GetComputeStep().TimeIt,
1731         grid->GetComputeStep().IterationIt,
1732         MED_NODE,
1733         MED_NONE,
1734         globalIds->GetPointer(0) ) < 0)
1735     {
1736     grid->SetPointGlobalIds(NULL);
1737     }
1738 }
1739
1740 void vtkMedDriver30::LoadCoordinates(vtkMedGrid* grid)
1741 {
1742   if(grid->IsCoordinatesLoaded())
1743     return;
1744
1745   vtkMedRegularGrid* rgrid = vtkMedRegularGrid::SafeDownCast(grid);
1746   if(rgrid != NULL)
1747     {
1748     this->LoadRegularGridCoordinates(rgrid);
1749     return;
1750     }
1751
1752   vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
1753   vtkMedCurvilinearGrid* cgrid = vtkMedCurvilinearGrid::SafeDownCast(grid);
1754   if(ugrid == NULL && cgrid == NULL)
1755     {
1756     //TODO : deal with structured grids
1757     vtkWarningMacro("this kind of grid is not yet supported");
1758     return;
1759     }
1760
1761   if(grid->GetUsePreviousCoordinates())
1762     {
1763     vtkMedGrid* previousgrid = grid->GetPreviousGrid();
1764     if(previousgrid == NULL)
1765       {
1766       vtkErrorMacro("coordiantes have not changed, "
1767                     << "but there is no previous grid!");
1768       return;
1769       }
1770
1771     this->LoadCoordinates(previousgrid);
1772     if(ugrid != NULL)
1773       ugrid->SetCoordinates(vtkMedUnstructuredGrid::SafeDownCast(previousgrid)
1774                             ->GetCoordinates());
1775     if(cgrid != NULL)
1776       cgrid->SetCoordinates(vtkMedCurvilinearGrid::SafeDownCast(previousgrid)
1777                             ->GetCoordinates());
1778     }
1779   else
1780     {
1781
1782     FileOpen open(this);
1783
1784     vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
1785     if(ugrid != NULL)
1786       ugrid->SetCoordinates(coords);
1787     if(cgrid != NULL)
1788       cgrid->SetCoordinates(coords);
1789     coords->Delete();
1790
1791     vtkMedComputeStep cs = grid->GetComputeStep();
1792
1793     coords->SetNumberOfComponents(grid->GetParentMesh()->GetSpaceDimension());
1794     coords->SetNumberOfTuples(grid->GetNumberOfPoints());
1795
1796     if ( MEDmeshNodeCoordinateRd( this->FileId,
1797                                   grid->GetParentMesh()->GetName(),
1798                                   cs.TimeIt,
1799                                   cs.IterationIt,
1800                                   MED_FULL_INTERLACE,
1801                                   (med_float*) coords->GetVoidPointer(0) ) < 0)
1802       {
1803       vtkErrorMacro("Load Coordinates for mesh "
1804           << grid->GetParentMesh()->GetName());
1805       }
1806     }
1807 }
1808
1809 void vtkMedDriver30::LoadProfile(vtkMedProfile* profile)
1810 {
1811   if(!profile || profile->IsLoaded())
1812     return;
1813
1814   FileOpen open(this);
1815
1816   vtkMedIntArray* indices = vtkMedIntArray::New();
1817   profile->SetIds(indices);
1818   indices->Delete();
1819
1820   indices->SetNumberOfTuples(profile->GetNumberOfElement());
1821
1822   char name[MED_NAME_SIZE+1] = "";
1823
1824   if( MEDprofileRd(this->FileId,
1825                    profile->GetName(),
1826                    indices->GetPointer(0) ) < 0)
1827     {
1828     vtkErrorMacro("Reading profile indices ");
1829     }
1830 }
1831
1832 void vtkMedDriver30::LoadConnectivity(vtkMedEntityArray* array)
1833 {
1834   if(array->IsConnectivityLoaded())
1835     return;
1836
1837   FileOpen open(this);
1838
1839   vtkMedGrid* grid = array->GetParentGrid();
1840
1841   grid = array->GetParentGrid();
1842
1843   const char* meshName = grid->GetParentMesh()->GetName();
1844
1845   vtkMedIntArray* conn = vtkMedIntArray::New();
1846   array->SetConnectivityArray(conn);
1847   conn->Delete();
1848
1849   vtkMedComputeStep cs = grid->GetComputeStep();
1850
1851   med_bool change;
1852   med_bool transformation;
1853
1854   if(array->GetEntity().GeometryType == MED_POLYGON)
1855     {
1856     // first check if we have points
1857     med_int connSize = MEDmeshnEntity(
1858                             this->FileId,
1859                             meshName,
1860                             cs.TimeIt,
1861                             cs.IterationIt,
1862                             array->GetEntity().EntityType,
1863                             MED_POLYGON,
1864                             MED_CONNECTIVITY,
1865                             array->GetConnectivity(),
1866                             &change,
1867                             &transformation);
1868
1869     if (connSize < 0)
1870       {
1871       vtkErrorMacro(<< "Error while reading polygons connectivity size."
1872                     << endl );
1873             return;
1874       }
1875
1876     conn->SetNumberOfTuples(connSize);
1877
1878     // How many polygon in the mesh in nodal connectivity mode
1879     // For the polygons, we get the size of array index
1880     med_int indexsize;
1881     if ((indexsize = MEDmeshnEntity(this->FileId,
1882                                     meshName,
1883                                     cs.TimeIt,
1884                                     cs.IterationIt,
1885                                     array->GetEntity().EntityType,
1886                                     MED_POLYGON,
1887                                     MED_INDEX_NODE,
1888                                     array->GetConnectivity(),
1889                                     &change,
1890                                     &transformation )) < 0)
1891       {
1892       vtkErrorMacro(<< "Error while reading polygons array index." << endl );
1893             return;
1894       }
1895
1896     vtkMedIntArray* index = vtkMedIntArray::New();
1897     array->SetFaceIndex(index);
1898     index->Delete();
1899
1900     index->SetNumberOfTuples( indexsize );
1901
1902     if ( MEDmeshPolygonRd(this->FileId,
1903                           meshName,
1904                           cs.TimeIt,
1905                           cs.IterationIt,
1906                           array->GetEntity().EntityType,
1907                           array->GetConnectivity(),
1908                           index->GetPointer(0),
1909                           conn->GetPointer(0) ) < 0)
1910       {
1911       vtkErrorMacro(<< "MEDmeshPolygonRd");
1912       return;
1913       }
1914     }
1915   else if(array->GetEntity().GeometryType == MED_POLYHEDRON)
1916     {
1917
1918     vtkIdType connSize = MEDmeshnEntity(this->FileId,
1919                                         meshName,
1920                                         grid->GetComputeStep().TimeIt,
1921                                         grid->GetComputeStep().IterationIt,
1922                                         array->GetEntity().EntityType,
1923                                         MED_POLYHEDRON,
1924                                         MED_CONNECTIVITY,
1925                                         array->GetConnectivity(),
1926                                         &change,
1927                                         &transformation);
1928     if (connSize < 0)
1929       {
1930       vtkErrorMacro(<< "Error while reading polyhedrons connectivity size."
1931                     << endl );
1932             return;
1933       }
1934
1935     conn->SetNumberOfTuples(connSize);
1936
1937     vtkMedIntArray* faceIndex = vtkMedIntArray::New();
1938     array->SetFaceIndex(faceIndex);
1939     faceIndex->Delete();
1940
1941     vtkMedIntArray* nodeIndex = vtkMedIntArray::New();
1942     array->SetNodeIndex(nodeIndex);
1943     nodeIndex->Delete();
1944
1945     vtkIdType np = array->GetNumberOfEntity() + 1;
1946     faceIndex->SetNumberOfTuples(np);
1947
1948     med_int nodeIndexSize;
1949
1950     if ((nodeIndexSize = MEDmeshnEntity(this->FileId,
1951                                         meshName,
1952                                         grid->GetComputeStep().TimeIt,
1953                                         grid->GetComputeStep().IterationIt,
1954                                         array->GetEntity().EntityType,
1955                                         MED_POLYHEDRON,
1956                                         MED_INDEX_NODE,
1957                                         array->GetConnectivity(),
1958                                         &change,
1959                                         &transformation )) < 0)
1960       {
1961       vtkErrorMacro(<< "Error while reading polygons array index." << endl );
1962             return;
1963       }
1964
1965     nodeIndex->SetNumberOfTuples(nodeIndexSize);
1966
1967     if (MEDmeshPolyhedronRd(this->FileId,
1968                             meshName,
1969                             cs.TimeIt,
1970                             cs.IterationIt,
1971                             array->GetEntity().EntityType,
1972                             array->GetConnectivity(),
1973                             faceIndex->GetPointer(0),
1974                             nodeIndex->GetPointer(0),
1975                             conn->GetPointer(0) ) < 0)
1976       {
1977       vtkErrorMacro("Error while reading connectivity of polyhedrons");
1978       return;
1979       }
1980
1981     }
1982   else
1983     {
1984     bool doReadConnectivity = true;
1985     if(array->GetConnectivity() == MED_NODAL)
1986       {
1987       if(array->GetEntity().EntityType == MED_STRUCT_ELEMENT)
1988         {
1989         if(array->GetStructElement() == NULL)
1990           {
1991           vtkErrorMacro("Entity type = MED_STRUCT_ELEMENT, but StructElement is not set!");
1992           return;
1993           }
1994         vtkIdType ntuple = array->GetNumberOfEntity()
1995                            * array->GetStructElement()->GetConnectivitySize();
1996
1997         conn->SetNumberOfTuples(ntuple);
1998         // particles are special : connectivity is not stored in the med file
1999         if(strcmp(array->GetStructElement()->GetName(), MED_PARTICLE_NAME) == 0 )
2000           {
2001           for(vtkIdType cellId = 0; cellId < ntuple; cellId++)
2002             {
2003             conn->SetValue(cellId, cellId+1);
2004             }
2005           doReadConnectivity = false;
2006           }
2007         }
2008       else
2009         {
2010         conn->SetNumberOfTuples(array->GetNumberOfEntity()
2011             * vtkMedUtilities::GetNumberOfPoint(
2012                 array->GetEntity().GeometryType));
2013         }
2014       }
2015     else
2016       {
2017       conn->SetNumberOfTuples(array->GetNumberOfEntity()
2018           * vtkMedUtilities::GetNumberOfSubEntity(
2019               array->GetEntity().GeometryType));
2020       }
2021
2022     if  (this->ParallelFileId == -1) // also (array->GetFilter() == NULL)
2023       {
2024       if ( (MEDmeshElementConnectivityRd(
2025             this->FileId,
2026             meshName,
2027             cs.TimeIt,
2028             cs.IterationIt,
2029             array->GetEntity().EntityType,
2030             array->GetEntity().GeometryType,
2031             array->GetConnectivity(),
2032             MED_FULL_INTERLACE,
2033             conn->GetPointer(0)) ) < 0)
2034         {
2035         vtkErrorMacro("Error while load connectivity of cells "
2036             << array->GetEntity().GeometryType);
2037         }
2038       }
2039     else
2040       {
2041       med_filter filter = MED_FILTER_INIT;
2042
2043       int    start;
2044       int    stride;
2045       int    count;
2046       int    blocksize;
2047       int    lastblocksize;
2048       array->GetFilter()->GetFilterSizes(start, stride, count,
2049                                    blocksize, lastblocksize );
2050
2051       med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
2052                                         array->GetEntity().GeometryType);
2053
2054       if ( MEDfilterBlockOfEntityCr( this->ParallelFileId,
2055               array->GetNumberOfEntity(),
2056             1, // one is for mesh elements, more than 1 is for fields
2057               nbofconstituentpervalue,
2058             MED_ALL_CONSTITUENT,
2059             MED_FULL_INTERLACE,
2060             MED_COMPACT_STMODE,
2061             MED_NO_PROFILE,
2062             (med_size)start,
2063             (med_size)stride,
2064             (med_size)count,
2065             (med_size)blocksize,
2066             (med_size)lastblocksize,
2067             &filter ) < 0 )
2068         {
2069         vtkErrorMacro("Filter creation ");
2070         }
2071
2072         if ( (MEDmeshElementConnectivityAdvancedRd(
2073               this->ParallelFileId,
2074               meshName,
2075               cs.TimeIt,
2076               cs.IterationIt,
2077               array->GetEntity().EntityType,
2078               array->GetEntity().GeometryType,
2079               array->GetConnectivity(),
2080               &filter,
2081               conn->GetPointer(0)) ) < 0)
2082           {
2083           vtkErrorMacro("Error while load connectivity of cells "
2084               << array->GetEntity().GeometryType);
2085           }
2086
2087       if ( MEDfilterClose( &filter ) < 0)
2088           {
2089         vtkErrorMacro("ERROR : filter closing ...");
2090           }
2091
2092       }
2093     }
2094 }
2095
2096 void vtkMedDriver30::LoadCellGlobalIds(vtkMedEntityArray* array)
2097 {
2098   if(array->IsGlobalIdsLoaded())
2099     return;
2100
2101   FileOpen open(this);
2102
2103   vtkMedIntArray* globalIds = vtkMedIntArray::New();
2104   array->SetGlobalIds(globalIds);
2105   globalIds->Delete();
2106
2107   globalIds->SetNumberOfTuples(array->GetNumberOfEntity());
2108
2109   vtkMedGrid* grid = array->GetParentGrid();
2110   vtkMedComputeStep cs = grid->GetComputeStep();
2111
2112   if( MEDmeshEntityNumberRd (
2113         this->FileId,
2114         grid->GetParentMesh()->GetName(),
2115         cs.TimeIt,
2116         cs.IterationIt,
2117         array->GetEntity().EntityType,
2118         array->GetEntity().GeometryType,
2119         globalIds->GetPointer(0) ) < 0)
2120     {
2121     array->SetGlobalIds(NULL);
2122     }
2123 }
2124
2125 void vtkMedDriver30::LoadField(vtkMedFieldOnProfile* fop, med_storage_mode mode)
2126 {
2127   FileOpen open(this);
2128
2129   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2130   vtkMedFieldStep *step = fieldOverEntity->GetParentStep();
2131   vtkMedField* field = step->GetParentField();
2132   const vtkMedComputeStep& cs = step->GetComputeStep();
2133
2134   vtkDataArray* data = vtkMedUtilities::NewArray(field->GetDataType());
2135   fop->SetData(data);
2136   data->Delete();
2137
2138   med_int size;
2139   if(mode == MED_COMPACT_STMODE)
2140     {
2141     size = fop->GetNumberOfValues();
2142     }
2143   else
2144     {
2145     med_int profilesize;
2146     med_int nbofintegrationpoint;
2147     char profileName[MED_NAME_SIZE+1] = "";
2148     char localizationName[MED_NAME_SIZE+1] = "";
2149     size = MEDfieldnValueWithProfile(this->FileId,
2150                 field->GetName(),
2151                 cs.TimeIt,
2152                 cs.IterationIt,
2153                 fieldOverEntity->GetEntity().EntityType,
2154                 fieldOverEntity->GetEntity().GeometryType,
2155                 fop->GetMedIterator(),
2156                 MED_GLOBAL_STMODE,
2157                 profileName,
2158                 &profilesize,
2159                 localizationName,
2160                 &nbofintegrationpoint);
2161     }
2162
2163   if(fop->GetNumberOfIntegrationPoint() > 1)
2164     {
2165     size *= fop->GetNumberOfIntegrationPoint();
2166     }
2167
2168   data->SetNumberOfComponents(field->GetNumberOfComponent());
2169   data->SetNumberOfTuples(size);
2170   if  (this->ParallelFileId == -1)
2171     {
2172     if ( MEDfieldValueWithProfileRd(
2173           this->FileId,
2174           field->GetName(),
2175           cs.TimeIt,
2176           cs.IterationIt,
2177           fieldOverEntity->GetEntity().EntityType,
2178           fieldOverEntity->GetEntity().GeometryType,
2179           mode,
2180           fop->GetProfileName(),
2181           MED_FULL_INTERLACE,
2182           MED_ALL_CONSTITUENT,
2183           (unsigned char*) data->GetVoidPointer(0) ) < 0)
2184       {
2185       vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
2186       }
2187     }
2188   else
2189     {
2190   if  (field->GetFieldType() == vtkMedField::CellField)
2191     {
2192     med_filter filter = MED_FILTER_INIT;
2193
2194     int    start;
2195     int    stride;
2196     int    count;
2197     int    blocksize;
2198     int    lastblocksize;
2199     fop->GetFilter()->GetFilterSizes(start, stride, count,
2200                                  blocksize, lastblocksize );
2201
2202     if ( MEDfilterBlockOfEntityCr( this->ParallelFileId,
2203         fop->GetNumberOfValues(),
2204           1, // one is for mesh elements, more than 1 is for fields
2205           field->GetNumberOfComponent(),
2206           MED_ALL_CONSTITUENT,
2207           MED_FULL_INTERLACE,
2208           MED_COMPACT_STMODE,
2209           MED_NO_PROFILE,
2210           (med_size)start,
2211           (med_size)stride,
2212           (med_size)count,
2213           (med_size)blocksize,
2214           (med_size)lastblocksize,
2215           &filter ) < 0 )
2216       {
2217       vtkErrorMacro("Filter creation ");
2218       }
2219
2220     if ( MEDfieldValueAdvancedRd(
2221             this->ParallelFileId,
2222             field->GetName(),
2223             cs.TimeIt,
2224             cs.IterationIt,
2225             fieldOverEntity->GetEntity().EntityType,
2226             fieldOverEntity->GetEntity().GeometryType,
2227             &filter,
2228             (unsigned char*) data->GetVoidPointer(0) ) < 0)
2229         {
2230         vtkErrorMacro("Error on MEDfieldValueAdvancedRd");
2231         }
2232
2233     if ( MEDfilterClose( &filter ) < 0)
2234         {
2235       vtkErrorMacro("ERROR : filter closing ...");
2236         }
2237       }
2238   else
2239     {//TODO : option utilisateur pour desactiver ou non les champs avec profile en //
2240     if ( MEDfieldValueWithProfileRd(
2241               this->FileId,
2242               field->GetName(),
2243               cs.TimeIt,
2244               cs.IterationIt,
2245               fieldOverEntity->GetEntity().EntityType,
2246               fieldOverEntity->GetEntity().GeometryType,
2247               mode,
2248               fop->GetProfileName(),
2249               MED_FULL_INTERLACE,
2250               MED_ALL_CONSTITUENT,
2251               (unsigned char*) data->GetVoidPointer(0) ) < 0)
2252           {
2253           vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
2254           }
2255     }
2256     }
2257 }
2258
2259 void vtkMedDriver30::LoadVariableAttribute(vtkMedVariableAttribute* varatt,
2260                                            vtkMedEntityArray* array)
2261 {
2262   FileOpen open(this);
2263
2264   void  *value = NULL;
2265
2266   vtkAbstractArray* valuearray = array->GetVariableAttributeValue(varatt);
2267   // first test if this is already loaded
2268   if(valuearray != NULL && valuearray->GetNumberOfTuples() > 0)
2269     return;
2270
2271   if(valuearray == NULL)
2272     {
2273     valuearray = vtkMedUtilities::NewArray(varatt->GetAttributeType());
2274     array->SetVariableAttributeValues(varatt, valuearray);
2275     valuearray->Delete();
2276     }
2277
2278   valuearray->SetNumberOfComponents(varatt->GetNumberOfComponent());
2279   valuearray->SetNumberOfTuples(array->GetNumberOfEntity());
2280   valuearray->SetName(varatt->GetName());
2281
2282   vtkSmartPointer<vtkCharArray> chararray = vtkSmartPointer<vtkCharArray>::New();
2283
2284   if(varatt->GetAttributeType() != MED_ATT_NAME)
2285     {
2286     value = valuearray->GetVoidPointer(0);
2287     }
2288   else
2289     {
2290     chararray->SetNumberOfValues(varatt->GetNumberOfComponent() *
2291                                   array->GetNumberOfEntity() *
2292                                   MED_NAME_SIZE);
2293
2294     value = chararray->GetVoidPointer(0);
2295     }
2296
2297   vtkMedComputeStep cs = array->GetParentGrid()->GetComputeStep();
2298
2299   if(MEDmeshStructElementVarAttRd(
2300       this->FileId,
2301       array->GetParentGrid()->GetParentMesh()->GetName(),
2302       cs.TimeIt,
2303       cs.IterationIt,
2304       varatt->GetParentStructElement()->GetGeometryType(),
2305       varatt->GetName(),
2306       value) < 0)
2307     {
2308
2309     if(cs.IterationIt == MED_NO_IT && cs.TimeIt == MED_NO_DT && cs.TimeOrFrequency == MED_UNDEF_DT)
2310       {
2311       vtkErrorMacro("MEDmeshStructElementVarAttRd");
2312       return;
2313       }
2314     // try to see if I can reuse
2315     // the variable attributes of the NO_DT, NO_IT compute step
2316     vtkMedComputeStep nocs;
2317     nocs.IterationIt = MED_NO_IT;
2318     nocs.TimeIt = MED_NO_DT;
2319     nocs.TimeOrFrequency = MED_UNDEF_DT;
2320     vtkMedEntityArray* nocs_array =
2321         array->GetParentGrid()->GetParentMesh()->GetGridStep(nocs)->GetEntityArray(array->GetEntity());
2322     if(nocs_array == NULL)
2323       {
2324       nocs_array = array->GetParentGrid()->GetParentMesh()->GetGridStep(0)->GetEntityArray(array->GetEntity());
2325       }
2326
2327     if(nocs_array == NULL || nocs_array == array)
2328       {
2329       // try to force load the default compute step.
2330       if(MEDmeshStructElementVarAttRd(
2331           this->FileId,
2332           array->GetParentGrid()->GetParentMesh()->GetName(),
2333           nocs.TimeIt,
2334           nocs.IterationIt,
2335           varatt->GetParentStructElement()->GetGeometryType(),
2336           varatt->GetName(),
2337           value) < 0)
2338         {
2339         vtkErrorMacro("MEDmeshStructElementVarAttRd");
2340         return;
2341         }
2342       }
2343     else
2344       {
2345       this->LoadVariableAttribute(varatt, nocs_array);
2346       array->SetVariableAttributeValues(varatt, nocs_array->GetVariableAttributeValue(varatt));
2347       return;
2348       }
2349     }
2350
2351   // If I am here, it means that I read the values
2352   if(varatt->GetAttributeType() == MED_ATT_NAME)
2353     {
2354     char current_name[MED_NAME_SIZE+1] = "";
2355     vtkStringArray* sarray = vtkStringArray::SafeDownCast(valuearray);
2356     for(vtkIdType id = 0; id < varatt->GetNumberOfComponent() *
2357                        array->GetNumberOfEntity(); id++)
2358       {
2359       memset(current_name, '\0', MED_NAME_SIZE+1);
2360       strncpy(current_name, ((char*)value) + id*MED_NAME_SIZE, MED_NAME_SIZE);
2361       sarray->SetValue(id, current_name);
2362       }
2363     }
2364
2365   return;
2366 }
2367
2368 void vtkMedDriver30::PrintSelf(ostream& os, vtkIndent indent)
2369 {
2370   this->Superclass::PrintSelf(os, indent);
2371 }