Back to index

salome-paravis  6.5.0
vtkMedDriver30.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 "vtkMedDriver30.h"
00021 
00022 #include "vtkObjectFactory.h"
00023 #include "vtkStringArray.h"
00024 #include "vtkDataArray.h"
00025 #include "vtkIdTypeArray.h"
00026 #include "vtkDoubleArray.h"
00027 #include "vtkIntArray.h"
00028 #include "vtkCharArray.h"
00029 
00030 #include "vtkMedFile.h"
00031 #include "vtkMedCartesianGrid.h"
00032 #include "vtkMedPolarGrid.h"
00033 #include "vtkMedCurvilinearGrid.h"
00034 #include "vtkMedUnstructuredGrid.h"
00035 #include "vtkMedField.h"
00036 #include "vtkMedMesh.h"
00037 #include "vtkMedFamily.h"
00038 #include "vtkMedUtilities.h"
00039 #include "vtkMedEntityArray.h"
00040 #include "vtkMedLocalization.h"
00041 #include "vtkMedProfile.h"
00042 #include "vtkMedFieldOverEntity.h"
00043 #include "vtkMedFieldStep.h"
00044 #include "vtkMedGroup.h"
00045 #include "vtkMedIntArray.h"
00046 #include "vtkMedInterpolation.h"
00047 #include "vtkMedFieldOnProfile.h"
00048 #include "vtkMedFraction.h"
00049 #include "vtkMedLink.h"
00050 #include "vtkMedStructElement.h"
00051 #include "vtkMedConstantAttribute.h"
00052 #include "vtkMedVariableAttribute.h"
00053 
00054 #include "vtkMultiProcessController.h"
00055 
00056 #include <string>
00057 using namespace std;
00058 
00059 vtkCxxRevisionMacro(vtkMedDriver30, "$Revision: 1.1.4.9.2.3 $")
00060 vtkStandardNewMacro(vtkMedDriver30)
00061 
00062 vtkMedDriver30::vtkMedDriver30() : vtkMedDriver()
00063 {
00064 }
00065 
00066 vtkMedDriver30::~vtkMedDriver30()
00067 {
00068 }
00069 
00070 // Description:
00071 // load all Information data associated with this cartesian grid.
00072 void vtkMedDriver30::ReadRegularGridInformation(vtkMedRegularGrid* grid)
00073 {
00074   FileOpen open(this);
00075 
00076   grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
00077 
00078   for(int axis=0; axis < grid->GetDimension(); axis++)
00079     {
00080     med_int size;
00081     med_bool coordinatechangement, geotransformation;
00082 
00083     if ((size = MEDmeshnEntity(
00084                   this->FileId,
00085                   grid->GetParentMesh()->GetName(),
00086                   grid->GetComputeStep().TimeIt,
00087                   grid->GetComputeStep().IterationIt,
00088                   MED_NODE,
00089                   MED_NONE,
00090                   (med_data_type)(MED_COORDINATE_AXIS1 + axis),
00091                   MED_NO_CMODE,
00092                   &coordinatechangement,
00093                   &geotransformation)) < 0)
00094       {
00095       vtkErrorMacro("ERROR : number of coordinates on X axis ...");
00096       }
00097 
00098     grid->SetAxisSize(axis, size);
00099     }
00100 
00101   med_int ncell = 1;
00102   if(grid->GetAxisSize(0) > 1)
00103     ncell *= grid->GetAxisSize(0)-1;
00104   if(grid->GetAxisSize(1) > 1)
00105     ncell *= grid->GetAxisSize(1)-1;
00106   if(grid->GetAxisSize(2) > 1)
00107     ncell *= grid->GetAxisSize(2)-1;
00108 
00109   vtkMedEntity entity;
00110   entity.EntityType = MED_CELL;
00111 
00112   switch(grid->GetDimension())
00113     {
00114     case 0 :
00115       entity.GeometryType = MED_POINT1;
00116       break;
00117     case 1 :
00118       entity.GeometryType = MED_SEG2;
00119       break;
00120     case 2 :
00121       entity.GeometryType = MED_QUAD4;
00122       break;
00123     case 3 :
00124       entity.GeometryType = MED_HEXA8;
00125       break;
00126     default :
00127         vtkErrorMacro("Unsupported dimension for curvilinear grid : "
00128                       << grid->GetDimension());
00129     return;
00130     }
00131 
00132   vtkMedEntityArray* array = vtkMedEntityArray::New();
00133   array->SetParentGrid(grid);
00134   array->SetNumberOfEntity(ncell);
00135   array->SetEntity(entity);
00136   array->SetConnectivity(MED_NODAL);
00137   grid->AppendEntityArray(array);
00138   array->Delete();
00139   // this triggers the creation of undefined families
00140   this->LoadFamilyIds(array);
00141 
00142 }
00143 
00144 void  vtkMedDriver30::LoadRegularGridCoordinates(vtkMedRegularGrid* grid)
00145 {
00146   FileOpen open(this);
00147 
00148   for(int axis=0; axis < grid->GetParentMesh()->GetNumberOfAxis(); axis++)
00149     {
00150 
00151     vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
00152     grid->SetAxisCoordinate(axis, coords);
00153     coords->Delete();
00154 
00155     coords->SetNumberOfComponents(1);
00156     coords->SetNumberOfTuples(grid->GetAxisSize(axis));
00157 
00158     if (MEDmeshGridIndexCoordinateRd(
00159           this->FileId,
00160           grid->GetParentMesh()->GetName(),
00161           grid->GetComputeStep().TimeIt,
00162           grid->GetComputeStep().IterationIt,
00163           axis+1,
00164           (med_float*)coords->GetVoidPointer(0)) < 0)
00165       {
00166       vtkErrorMacro("ERROR : read axis " << axis << " coordinates ...");
00167       grid->SetAxisCoordinate(axis, NULL);
00168       return;
00169       }
00170     }
00171 }
00172 
00173 // Description:
00174 // load all Information data associated with this standard grid.
00175 void vtkMedDriver30::ReadCurvilinearGridInformation(vtkMedCurvilinearGrid* grid)
00176 {
00177   FileOpen open(this);
00178 
00179   grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
00180 
00181   med_int size;
00182   med_bool coordinatechangement, geotransformation;
00183 
00184   if ((size = MEDmeshnEntity(
00185                 this->FileId,
00186                 grid->GetParentMesh()->GetName(),
00187                 grid->GetComputeStep().TimeIt,
00188                 grid->GetComputeStep().IterationIt,
00189                 MED_NODE,
00190                 MED_NONE,
00191                 MED_COORDINATE,
00192                 MED_NO_CMODE,
00193                 &coordinatechangement,
00194                 &geotransformation)) < 0)
00195     {
00196     vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshnEntity");
00197     }
00198 
00199   grid->SetNumberOfPoints(size);
00200 
00201   med_int axissize[3];
00202   if(MEDmeshGridStructRd(
00203       this->FileId,
00204       grid->GetParentMesh()->GetName(),
00205       grid->GetComputeStep().TimeIt,
00206       grid->GetComputeStep().IterationIt,
00207       axissize) < 0)
00208     {
00209     vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshGridStructRd");
00210     }
00211   grid->SetAxisSize(0, (axissize[0] <= 0 ? 1: axissize[0]));
00212   grid->SetAxisSize(1, (axissize[1] <= 0 ? 1: axissize[1]));
00213   grid->SetAxisSize(2, (axissize[2] <= 0 ? 1: axissize[2]));
00214 
00215   // A test to verify the number of points : total number of
00216   // points should be equal to the product of each axis size
00217   if(size != grid->GetAxisSize(0)*grid->GetAxisSize(1)*grid->GetAxisSize(2))
00218     {
00219     vtkErrorMacro("The total number of points of a Curvilinear grid should "
00220                   << "be the product of each axis size!");
00221     }
00222 
00223   med_int ncell = 1;
00224   if(grid->GetAxisSize(0) > 1)
00225     ncell *= grid->GetAxisSize(0)-1;
00226   if(grid->GetAxisSize(1) > 1)
00227     ncell *= grid->GetAxisSize(1)-1;
00228   if(grid->GetAxisSize(2) > 1)
00229     ncell *= grid->GetAxisSize(2)-1;
00230 
00231   vtkMedEntity entity;
00232   entity.EntityType = MED_CELL;
00233 
00234   switch(grid->GetDimension())
00235     {
00236     case 0 :
00237       entity.GeometryType = MED_POINT1;
00238       break;
00239     case 1 :
00240       entity.GeometryType = MED_SEG2;
00241       break;
00242     case 2 :
00243       entity.GeometryType = MED_QUAD4;
00244       break;
00245     case 3 :
00246       entity.GeometryType = MED_HEXA8;
00247       break;
00248     default :
00249         vtkErrorMacro("Unsupported dimension for curvilinear grid : "
00250                       << grid->GetDimension());
00251     return;
00252     }
00253 
00254   vtkMedEntityArray* array = vtkMedEntityArray::New();
00255   array->SetParentGrid(grid);
00256   array->SetNumberOfEntity(ncell);
00257   array->SetEntity(entity);
00258   array->SetConnectivity(MED_NODAL);
00259   grid->AppendEntityArray(array);
00260   array->Delete();
00261   // this triggers the creation of undefined families
00262   this->LoadFamilyIds(array);
00263 
00264 }
00265 
00266 // Description : read the number of entity of all geometry type
00267 // for a given entity type and a given connectivity mode
00268 void vtkMedDriver30::ReadNumberOfEntity(
00269     vtkMedUnstructuredGrid* grid,
00270     med_entity_type entityType,
00271     med_connectivity_mode connectivity)
00272 {
00273   FileOpen open(this);
00274 
00275   med_bool changement, transformation;
00276 
00277   const char* meshName = grid->GetParentMesh()->GetName();
00278 
00279   const vtkMedComputeStep& cs = grid->GetComputeStep();
00280 
00281   med_int nentity = MEDmeshnEntity(
00282                         this->FileId,
00283                         meshName,
00284                         cs.TimeIt,
00285                         cs.IterationIt,
00286                         entityType,
00287                         MED_GEO_ALL,
00288                         MED_UNDEF_DATATYPE ,
00289                         connectivity,
00290                         &changement,
00291                         &transformation );
00292 
00293   for(med_int geotypeit = 1; geotypeit <= nentity; geotypeit++)
00294     {
00295     // read cell informations
00296     vtkMedEntity entity;
00297     entity.EntityType = entityType;
00298 
00299     char geometryName[MED_NAME_SIZE+1] = "";
00300 
00301     // this gives us the med_geometry_type
00302     if( MEDmeshEntityInfo( FileId, meshName,
00303                            cs.TimeIt,
00304                            cs.IterationIt,
00305                            entityType,
00306                            geotypeit,
00307                            geometryName,
00308                            &entity.GeometryType) < 0)
00309       {
00310       vtkErrorMacro("MEDmeshEntityInfo");
00311       continue;
00312       }
00313 
00314     entity.GeometryName = geometryName;
00315     med_int ncell = 0;
00316 
00317     if(entity.GeometryType == MED_POLYGON)
00318       {
00319       // read the number of cells of this type
00320       ncell = MEDmeshnEntity(this->FileId,
00321                              meshName,
00322                              cs.TimeIt,
00323                              cs.IterationIt,
00324                              entity.EntityType,
00325                              entity.GeometryType,
00326                              MED_INDEX_NODE,
00327                              connectivity,
00328                              &changement,
00329                              &transformation ) - 1;
00330       }
00331     else if(entity.GeometryType == MED_POLYHEDRON)
00332       {
00333       // read the number of cells of this type
00334       ncell = MEDmeshnEntity(this->FileId,
00335                              meshName,
00336                              cs.TimeIt,
00337                              cs.IterationIt,
00338                              entity.EntityType,
00339                              entity.GeometryType,
00340                              MED_INDEX_FACE,
00341                              connectivity,
00342                              &changement,
00343                              &transformation  ) - 1;
00344       }
00345     else
00346       {
00347       ncell = MEDmeshnEntity(this->FileId,
00348                              meshName,
00349                              cs.TimeIt,
00350                              cs.IterationIt,
00351                              entity.EntityType,
00352                              entity.GeometryType,
00353                              MED_CONNECTIVITY,
00354                              connectivity,
00355                              &changement,
00356                              &transformation  );
00357       }
00358 
00359     if(ncell > 0)
00360       {
00361       vtkMedEntityArray* array = vtkMedEntityArray::New();
00362       array->SetParentGrid(grid);
00363       array->SetNumberOfEntity(ncell);
00364       array->SetEntity(entity);
00365       array->SetConnectivity(connectivity);
00366       grid->AppendEntityArray(array);
00367       array->Delete();
00368       // this triggers the creation of undefined families
00369       this->LoadFamilyIds(array);
00370       }
00371     }
00372 }
00373 
00374 // Description:
00375 // load all Information data associated with this unstructured grid.
00376 void vtkMedDriver30::ReadUnstructuredGridInformation(
00377     vtkMedUnstructuredGrid* grid)
00378 {
00379   FileOpen open(this);
00380 
00381   vtkMedMesh *mesh = grid->GetParentMesh();
00382 
00383   const char *meshName = mesh->GetName();
00384   med_connectivity_mode connectivity;
00385 
00386   med_bool changement;
00387   med_bool transformation;
00388   med_int profilesize;
00389 
00390   char profilename[MED_NAME_SIZE+1];
00391   memset(profilename, '\0', MED_NAME_SIZE+1);
00392 
00393   vtkMedComputeStep cs = grid->GetComputeStep();
00394 
00395   // first check if we have points
00396   vtkIdType npoints = MEDmeshnEntityWithProfile(
00397                         this->FileId,
00398                         meshName,
00399                         cs.TimeIt,
00400                         cs.IterationIt,
00401                         MED_NODE,
00402                         MED_NONE,
00403                         MED_COORDINATE,
00404                         MED_NODAL,
00405                         MED_COMPACT_PFLMODE,
00406                         profilename,
00407                         &profilesize,
00408                         &changement,
00409                         &transformation);
00410 
00411   if(npoints > 0)
00412     {
00413     grid->SetNumberOfPoints(npoints);
00414     }
00415   else
00416     {
00417     if(grid->GetPreviousGrid() == NULL)
00418       {
00419       vtkErrorMacro("No point and no previous grid");
00420       }
00421     grid->SetUsePreviousCoordinates(true);
00422     grid->SetNumberOfPoints(grid->GetPreviousGrid()->GetNumberOfPoints());
00423     return;
00424     }
00425 
00426   this->ReadNumberOfEntity(grid, MED_CELL, MED_NODAL);
00427   this->ReadNumberOfEntity(grid, MED_CELL, MED_DESCENDING);
00428   this->ReadNumberOfEntity(grid, MED_DESCENDING_FACE, MED_NODAL);
00429   this->ReadNumberOfEntity(grid, MED_DESCENDING_FACE, MED_DESCENDING);
00430   this->ReadNumberOfEntity(grid, MED_DESCENDING_EDGE, MED_NODAL);
00431   this->ReadNumberOfEntity(grid, MED_DESCENDING_EDGE, MED_DESCENDING);
00432   this->ReadNumberOfEntity(grid, MED_STRUCT_ELEMENT, MED_NODAL);
00433   this->ReadNumberOfEntity(grid, MED_STRUCT_ELEMENT, MED_DESCENDING);
00434 
00435   // create the point vtkMedEntityArray support
00436   vtkMedEntity entity;
00437   entity.EntityType = MED_NODE;
00438   entity.GeometryType = MED_POINT1;
00439   vtkMedEntityArray* pea = vtkMedEntityArray::New();
00440   pea->SetEntity(entity);
00441   pea->SetParentGrid(grid);
00442   pea->SetNumberOfEntity(grid->GetNumberOfPoints());
00443   grid->AppendEntityArray(pea);
00444   pea->Delete();
00445 
00446   this->LoadFamilyIds(pea);
00447 }
00448 
00449 void vtkMedDriver30::ReadFamilyInformation(vtkMedMesh* mesh, vtkMedFamily* family)
00450 {
00451   FileOpen open(this);
00452 
00453   med_int familyid;
00454   med_int ngroup;
00455   char* groupNames = NULL;
00456   const  char* meshName = mesh->GetName();
00457 
00458   ngroup = MEDnFamilyGroup(FileId, meshName, family->GetMedIterator());
00459 
00460   bool has_no_group = false;
00461   if(ngroup <= 0)
00462     {
00463     if(ngroup < 0)
00464       {
00465       vtkErrorMacro("Error while reading the number of groups");
00466       }
00467     ngroup = 1;
00468     has_no_group = true;
00469     }
00470 
00471   groupNames = new char[ngroup * MED_LNAME_SIZE + 1];
00472   memset(groupNames, '\0', ngroup * MED_LNAME_SIZE + 1);
00473 
00474   // special case for files written by med < 3,
00475   // I have to use the 23 interface
00476   if(mesh->GetParentFile()->GetVersionMajor() < 3)
00477     {
00478     med_int *attributenumber = NULL;
00479     med_int *attributevalue = NULL;
00480     char *attributedes = NULL;
00481 
00482     med_int nattr = MEDnFamily23Attribute(
00483                       this->FileId,
00484                       meshName,
00485                       family->GetMedIterator());
00486 
00487     if(nattr < 0)
00488       {
00489       vtkErrorMacro("MEDnFamily23Attribute");
00490       }
00491 
00492     if(nattr > 0)
00493       {
00494       attributenumber = new med_int[nattr];
00495       attributevalue = new med_int[nattr];
00496       attributedes = new char[nattr*MED_COMMENT_SIZE+1];
00497       memset(attributedes, '\0', nattr*MED_COMMENT_SIZE+1);
00498       }
00499 
00500     char familyName[MED_NAME_SIZE+1] = "";
00501 
00502     if(MEDfamily23Info (this->FileId,
00503                         meshName,
00504                         family->GetMedIterator(),
00505                         familyName,
00506                         attributenumber,
00507                         attributevalue,
00508                         attributedes,
00509                         &familyid,
00510                         groupNames ) < 0)
00511       {
00512       vtkDebugMacro("MEDfamily23Info");
00513       }
00514 
00515     family->SetName(familyName);
00516 
00517     if(attributenumber != NULL)
00518       delete[] attributenumber;
00519     if(attributevalue != NULL)
00520       delete[] attributevalue;
00521     if(attributedes != NULL)
00522       delete[] attributedes;
00523     }
00524   else
00525     {
00526     char familyName[MED_NAME_SIZE+1] = "";
00527     if(MEDfamilyInfo( this->FileId,
00528                       meshName,
00529                       family->GetMedIterator(),
00530                       familyName,
00531                       &familyid,
00532                       groupNames ) < 0)
00533       {
00534       vtkErrorMacro(
00535           "vtkMedDriver30::ReadInformation(vtkMedFamily* family)"
00536           << " cannot read family informations.");
00537       return;
00538       }
00539     family->SetName(familyName);
00540     }
00541 
00542   family->SetId(familyid);
00543 
00544   if(familyid <= 0)
00545     {
00546     family->SetPointOrCell(vtkMedUtilities::OnCell);
00547     mesh->AppendCellFamily(family);
00548     }
00549   else
00550     {
00551     family->SetPointOrCell(vtkMedUtilities::OnPoint);
00552     mesh->AppendPointFamily(family);
00553     }
00554 
00555   family->AllocateNumberOfGroup(ngroup);
00556   // if there where no group, set the name to the default value
00557   if(has_no_group)
00558     {
00559     memcpy(groupNames, vtkMedUtilities::NoGroupName,
00560            strlen(vtkMedUtilities::NoGroupName));
00561     }
00562 
00563   for(int index = 0; index < ngroup; index++)
00564     {
00565     char realGroupName[MED_LNAME_SIZE + 1];
00566     memset(realGroupName, '\0', MED_LNAME_SIZE + 1);
00567     memcpy(realGroupName, groupNames + index * MED_LNAME_SIZE,
00568         MED_LNAME_SIZE * sizeof(char));
00569     vtkMedGroup* group = mesh->GetOrCreateGroup(family->GetPointOrCell(),
00570         realGroupName);
00571 
00572     family->SetGroup(index, group);
00573     }
00574 
00575   delete[] groupNames;
00576 
00577   if(familyid == 0)
00578     {
00579     vtkMedFamily* famzero = vtkMedFamily::New();
00580     mesh->AppendPointFamily(famzero);
00581     famzero->Delete();
00582 
00583     famzero->SetName(family->GetName());
00584     famzero->SetMedIterator(family->GetMedIterator());
00585     famzero->SetId(family->GetId());
00586     famzero->SetPointOrCell(vtkMedUtilities::OnPoint);
00587     famzero->AllocateNumberOfGroup(family->GetNumberOfGroup());
00588     for(int gid=0; gid<family->GetNumberOfGroup(); gid++)
00589       {
00590       vtkMedGroup* group = mesh->GetOrCreateGroup(
00591           vtkMedUtilities::OnPoint,
00592           family->GetGroup(gid)->GetName());
00593       famzero->SetGroup(gid, group);
00594       mesh->AppendPointGroup(group);
00595       }
00596     }
00597 }
00598 
00599 void  vtkMedDriver30::ReadLinkInformation(vtkMedLink* link)
00600 {
00601   med_int size;
00602   char linkMeshName[MED_NAME_SIZE+1] = "";
00603   if(MEDlinkInfo(this->FileId,
00604                  link->GetMedIterator(),
00605                  linkMeshName,
00606                  &size) < 0)
00607     {
00608     vtkErrorMacro("MEDlinkInfo");
00609     return;
00610     }
00611   link->SetMeshName(linkMeshName);
00612   if(size <= 0)
00613     return;
00614 
00615   char* path = new char[size + 1];
00616   memset(path, '\0', size+1);
00617   if(MEDlinkRd(this->FileId, link->GetMeshName(), path) < 0)
00618     {
00619     vtkErrorMacro("MEDlinkRd");
00620     memset(path, '\0', size+1);
00621     }
00622 
00623   link->SetLink(path);
00624 
00625   delete[] path;
00626 }
00627 
00628 void vtkMedDriver30::ReadFileInformation(vtkMedFile* file)
00629 {
00630   FileOpen open(this);
00631 
00632   char comment[MED_COMMENT_SIZE+1] = "";
00633 
00634   MEDfileCommentRd(this->FileId,
00635                   comment);
00636 
00637   file->SetComment(comment);
00638 
00639   med_int major, minor, release;
00640   MEDfileNumVersionRd(this->FileId, &major, &minor, &release);
00641   file->SetVersionMajor(major);
00642   file->SetVersionMinor(minor);
00643   file->SetVersionRelease(release);
00644 
00645   int nlink = MEDnLink(this->FileId);
00646   file->AllocateNumberOfLink(nlink);
00647   for(int linkid=0; linkid<nlink; linkid++)
00648     {
00649     vtkMedLink* link = file->GetLink(linkid);
00650     link->SetMedIterator(linkid+1);
00651     this->ReadLinkInformation(link);
00652     }
00653 
00654   int nprof = MEDnProfile(FileId);
00655   // Reading id s not possible in parallel if the file contains Profiles
00656   vtkMultiProcessController* controller = vtkMultiProcessController::GetGlobalController();
00657   if (controller != NULL)
00658     if ((nprof != 0) && (controller->GetNumberOfProcesses() > 1))
00659     {
00660       vtkWarningMacro("ATTENTION: The MED Reader cannot read profiles when used in parallel");
00661     return;
00662     }
00663   file->AllocateNumberOfProfile(nprof);
00664   for(int i = 0; i < nprof; i++)
00665     {
00666     vtkMedProfile* profile = file->GetProfile(i);
00667     profile->SetMedIterator(i + 1);
00668     profile->SetParentFile(file);
00669     this->ReadProfileInformation(profile);
00670     }
00671 
00672   int nloc = MEDnLocalization(this->FileId);
00673   file->AllocateNumberOfLocalization(nloc);
00674   for(int i = 0; i < nloc; i++)
00675     {
00676     vtkMedLocalization* loc = file->GetLocalization(i);
00677     loc->SetMedIterator(i + 1);
00678     loc->SetParentFile(file);
00679     this->ReadLocalizationInformation(loc);
00680     }
00681 
00682   int nsupportmesh = MEDnSupportMesh(this->FileId);
00683   file->AllocateNumberOfSupportMesh(nsupportmesh);
00684   for(int i = 0; i < nsupportmesh; i++)
00685     {
00686     vtkMedMesh* supportmesh = file->GetSupportMesh(i);
00687     supportmesh->SetMedIterator(i + 1);
00688     supportmesh->SetParentFile(file);
00689     this->ReadSupportMeshInformation(supportmesh);
00690     }
00691 
00692   int nmesh = MEDnMesh(this->FileId);
00693   file->AllocateNumberOfMesh(nmesh);
00694   for(int i = 0; i < nmesh; i++)
00695     {
00696     vtkMedMesh* mesh = file->GetMesh(i);
00697     mesh->SetMedIterator(i + 1);
00698     mesh->SetParentFile(file);
00699     this->ReadMeshInformation(mesh);
00700     }
00701 
00702   int nfields = MEDnField(this->FileId);
00703   file->AllocateNumberOfField(nfields);
00704   for(int i = 0; i < nfields; i++)
00705     {
00706     vtkMedField* field = file->GetField(i);
00707     field->SetMedIterator(i + 1);
00708     field->SetParentFile(file);
00709     this->ReadFieldInformation(field);
00710     field->ComputeFieldType();
00711     while(field->HasManyFieldTypes())
00712       {
00713       vtkMedField* newfield = vtkMedField::New();
00714       int type = field->GetFirstType();
00715       newfield->ExtractFieldType(field, type);
00716       file->AppendField(newfield);
00717       newfield->Delete();
00718       }
00719     }
00720 
00721   int nstruct = MEDnStructElement(this->FileId);
00722 
00723   file->AllocateNumberOfStructElement(nstruct);
00724   for(int i = 0; i < nstruct; i++)
00725     {
00726     vtkMedStructElement* structelem = file->GetStructElement(i);
00727     structelem->SetMedIterator(i+1);
00728     structelem->SetParentFile(file);
00729     this->ReadStructElementInformation(structelem);
00730     }
00731 
00732 }
00733 
00734 void  vtkMedDriver30::ReadStructElementInformation(
00735     vtkMedStructElement* structelem)
00736 {
00737 
00738   FileOpen open(this);
00739 
00740   char modelname[MED_NAME_SIZE+1] = "";
00741   med_geometry_type mgeotype;
00742   med_int modeldim;
00743   char supportmeshname[MED_NAME_SIZE+1] = "";
00744   med_entity_type sentitytype;
00745   med_int snbofnode;
00746   med_int snbofcell;
00747   med_geometry_type sgeotype;
00748   med_int nbofconstantattribute;
00749   med_bool anyprofile;
00750   med_int nbofvariableattribute;
00751 
00752   if(MEDstructElementInfo (this->FileId,
00753                            structelem->GetMedIterator(),
00754                            modelname,
00755                            &mgeotype,
00756                            &modeldim,
00757                            supportmeshname,
00758                            &sentitytype,
00759                            &snbofnode,
00760                            &snbofcell,
00761                            &sgeotype,
00762                            &nbofconstantattribute,
00763                            &anyprofile,
00764                            &nbofvariableattribute) < 0)
00765     {
00766     vtkErrorMacro("Error in MEDstructElementInfo");
00767     return;
00768     }
00769   structelem->SetName(modelname);
00770   structelem->SetGeometryType(mgeotype);
00771   structelem->SetModelDimension(modeldim);
00772   structelem->SetSupportMeshName(supportmeshname);
00773   structelem->SetSupportEntityType(sentitytype);
00774   structelem->SetSupportNumberOfNode(snbofnode);
00775   structelem->SetSupportNumberOfCell(snbofcell);
00776   structelem->SetSupportGeometryType(sgeotype);
00777   structelem->AllocateNumberOfConstantAttribute(nbofconstantattribute);
00778   structelem->AllocateNumberOfVariableAttribute(nbofvariableattribute);
00779   structelem->SetAnyProfile(anyprofile);
00780 
00781   for(int attit = 0; attit < nbofconstantattribute; attit ++)
00782     {
00783     vtkMedConstantAttribute* constatt = structelem->GetConstantAttribute(attit);
00784     constatt->SetMedIterator(attit+1);
00785     constatt->SetParentStructElement(structelem);
00786     this->ReadConstantAttributeInformation(constatt);
00787     }
00788 
00789   for(int attit = 0; attit < nbofvariableattribute; attit ++)
00790     {
00791     vtkMedVariableAttribute* varatt = structelem->GetVariableAttribute(attit);
00792     varatt->SetMedIterator(attit+1);
00793     varatt->SetParentStructElement(structelem);
00794     this->ReadVariableAttributeInformation(varatt);
00795     }
00796 }
00797 
00798 void vtkMedDriver30::ReadConstantAttributeInformation(
00799     vtkMedConstantAttribute* constAttr)
00800 {
00801 
00802   FileOpen open(this);
00803 
00804   char constattname[MED_NAME_SIZE+1] = "";
00805   med_attribute_type constatttype;
00806   med_int nbofcomponent;
00807   med_entity_type sentitytype;
00808   char profilename[MED_NAME_SIZE+1] = "";
00809   med_int profilesize;
00810 
00811   if(MEDstructElementConstAttInfo(
00812       this->FileId,
00813       constAttr->GetParentStructElement()->GetName(),
00814       constAttr->GetMedIterator(),
00815       constattname,
00816       &constatttype,
00817       &nbofcomponent,
00818       &sentitytype,
00819       profilename,
00820       &profilesize)   < 0)
00821     {
00822     vtkErrorMacro("MEDstructElementConstAttInfo error");
00823     return;
00824     }
00825 
00826   constAttr->SetName(constattname);
00827   constAttr->SetAttributeType(constatttype);
00828   constAttr->SetNumberOfComponent(nbofcomponent);
00829   constAttr->SetSupportEntityType(sentitytype);
00830   constAttr->SetProfileName(profilename);
00831   constAttr->SetProfileSize(profilesize);
00832 
00833   vtkAbstractArray* values = vtkMedUtilities::NewArray(constatttype);
00834   if(values == NULL)
00835     return;
00836   constAttr->SetValues(values);
00837   values->Delete();
00838 
00839   values->SetNumberOfComponents(nbofcomponent);
00840   vtkIdType ntuple = 0;
00841   if((strcmp(profilename, MED_NO_PROFILE) != 0) &&
00842      (strcmp(profilename, "\0") != 0))
00843     {
00844     ntuple = profilesize;
00845     }
00846   else if(constAttr->GetSupportEntityType() == MED_CELL)
00847     {
00848     ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfCell();
00849     }
00850   else
00851     {
00852     ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfNode();
00853     }
00854   values->SetNumberOfTuples(ntuple);
00855 
00856   void* ptr = NULL;
00857   vtkSmartPointer<vtkCharArray> buffer = vtkSmartPointer<vtkCharArray>::New();
00858   if(constatttype != MED_ATT_NAME)
00859     {
00860     ptr = values->GetVoidPointer(0);
00861     }
00862   else
00863     {
00864     buffer->SetNumberOfValues(MED_NAME_SIZE*nbofcomponent*ntuple);
00865     ptr = buffer->GetVoidPointer(0);
00866     }
00867 
00868   if(MEDstructElementConstAttRd (this->FileId,
00869         constAttr->GetParentStructElement()->GetName(),
00870         constAttr->GetName(), ptr) < 0)
00871     {
00872     vtkErrorMacro("MEDstructElementConstAttRd");
00873     return;
00874     }
00875 
00876   if(constatttype == MED_ATT_NAME)
00877     {
00878     char name[MED_NAME_SIZE+1] = "";
00879     char* nameptr = (char*) ptr;
00880     vtkStringArray* names = vtkStringArray::SafeDownCast(values);
00881     for(vtkIdType id = 0; id < nbofcomponent*ntuple; id++)
00882       {
00883       memset(name, '\0', MED_NAME_SIZE+1);
00884       strncpy(name, nameptr + id * MED_NAME_SIZE, MED_NAME_SIZE);
00885       names->SetValue(id, name);
00886       }
00887     }
00888 
00889   return;
00890 }
00891 
00892 void vtkMedDriver30::ReadVariableAttributeInformation(
00893     vtkMedVariableAttribute* varAttr)
00894 {
00895 
00896   FileOpen open(this);
00897 
00898   char varattname[MED_NAME_SIZE+1] = "";
00899   med_attribute_type varatttype;
00900   med_int nbofcomponent;
00901 
00902   if(MEDstructElementVarAttInfo (
00903       this->FileId,
00904       varAttr->GetParentStructElement()->GetName(),
00905       varAttr->GetMedIterator(),
00906       varattname,
00907       &varatttype,
00908       &nbofcomponent) < 0)
00909     {
00910     vtkErrorMacro("MEDstructElementVarAttInfo");
00911     return;
00912     }
00913 
00914   varAttr->SetName(varattname);
00915   varAttr->SetAttributeType(varatttype);
00916   varAttr->SetNumberOfComponent(nbofcomponent);
00917 
00918   return;
00919 }
00920 
00921 void vtkMedDriver30::ReadProfileInformation(vtkMedProfile* profile)
00922 {
00923   FileOpen open(this);
00924 
00925   med_int nelem;
00926   char profileName[MED_NAME_SIZE+1] = "";
00927 
00928   if(MEDprofileInfo(this->FileId,
00929                     profile->GetMedIterator(),
00930                     profileName,
00931                     &nelem) < 0)
00932     {
00933     vtkErrorMacro("cannot read information on profile"
00934         << profile->GetMedIterator());
00935     }
00936   profile->SetName(profileName);
00937   profile->SetNumberOfElement(nelem);
00938 }
00939 
00940 void vtkMedDriver30::ReadFieldInformation(vtkMedField* field)
00941 {
00942   FileOpen open(this);
00943 
00944   if (field->GetMedIterator() == 0)
00945     return;
00946 
00947   int ncomp = MEDfieldnComponent(FileId, field->GetMedIterator());
00948 
00949   if(ncomp < 0)
00950     {
00951     field->SetNumberOfComponent(-1);
00952     vtkErrorMacro("cannot read the number of component of field "
00953         << field->GetMedIterator())
00954     return;
00955     }
00956 
00957   field->SetNumberOfComponent(ncomp);
00958 
00959   char* units = new char[MED_SNAME_SIZE * ncomp + 1];
00960   char* componentNames = new char[MED_SNAME_SIZE * ncomp + 1];
00961   memset(units, '\0', MED_SNAME_SIZE * ncomp + 1);
00962   memset(componentNames, '\0', MED_SNAME_SIZE * ncomp + 1);
00963 
00964   //med_type_champ dataType;
00965   med_field_type dataType;
00966   med_int nstep;
00967   med_bool localmesh;
00968 
00969   char name[MED_NAME_SIZE+1] = "";
00970   char meshName[MED_NAME_SIZE+1] = "";
00971   char unit[MED_SNAME_SIZE+1] = "";
00972 
00973   if( MEDfieldInfo( FileId,
00974                     field->GetMedIterator(),
00975                     name,
00976                     meshName,
00977                     &localmesh,
00978                     &dataType,
00979                     componentNames,
00980                     units,
00981                     unit,
00982                     &nstep) < 0)
00983     {
00984     vtkErrorMacro("cannot read the informations on field "
00985         << field->GetMedIterator())
00986     return;
00987     }
00988 
00989   field->SetName(name);
00990   field->SetMeshName(meshName);
00991   field->SetTimeUnit(unit);
00992   field->SetDataType(dataType);
00993   field->SetLocal(localmesh);
00994 
00995   for(int comp = 0; comp < ncomp; comp++)
00996     {
00997     char unit[MED_NAME_SIZE + 1] = "";
00998     memcpy(unit, units + MED_SNAME_SIZE * comp, MED_SNAME_SIZE * sizeof(char));
00999     field->GetUnit()->SetValue(comp, unit);
01000 
01001     char compName[MED_SNAME_SIZE + 1] = "";
01002     memcpy(compName, componentNames + MED_SNAME_SIZE * comp, MED_SNAME_SIZE
01003         * sizeof(char));
01004     field->GetComponentName()->SetValue(comp, compName);
01005     }
01006 
01007   delete[] units;
01008   delete[] componentNames;
01009 
01010   med_int ninterp = MEDfieldnInterp(FileId, field->GetName());
01011   if(ninterp < 0)
01012     {
01013     vtkErrorMacro("Error in MEDfieldnInterp");
01014     return;
01015     }
01016 
01017   field->AllocateNumberOfInterpolation(ninterp);
01018 
01019   for(med_int interpit=0; interpit<ninterp; interpit++)
01020     {
01021     vtkMedInterpolation* interp = field->GetInterpolation(interpit);
01022     interp->SetMedIterator(interpit + 1);
01023     this->ReadInterpolationInformation(interp);
01024     }
01025 
01026   vtkMedFieldStep* previousStep = NULL;
01027 
01028   for(med_int csit = 0; csit < nstep; csit++)
01029     {
01030     vtkMedFieldStep* step = vtkMedFieldStep::New();
01031     step->SetMedIterator(csit + 1);
01032     step->SetParentField(field);
01033     this->ReadFieldStepInformation(step, true);
01034     field->AddFieldStep(step);
01035     step->SetPreviousStep(previousStep);
01036     previousStep = step;
01037     step->Delete();
01038     }
01039 }
01040 
01041 void vtkMedDriver30::ReadFieldStepInformation(vtkMedFieldStep* step, bool readAllEntityInfo)
01042 {
01043   vtkMedComputeStep cs;
01044   vtkMedComputeStep meshcs;
01045   vtkMedField* field = step->GetParentField();
01046 
01047   FileOpen open(this);
01048 
01049   if( MEDfieldComputingStepMeshInfo(
01050         FileId,
01051         field->GetName(),
01052         step->GetMedIterator(),
01053         &cs.TimeIt,
01054         &cs.IterationIt,
01055         &cs.TimeOrFrequency,
01056         &meshcs.TimeIt,
01057         &meshcs.IterationIt) < 0)
01058     {
01059     vtkErrorMacro("Error in MEDfieldComputingStepMeshInfo");
01060     return;
01061     }
01062 
01063   step->SetComputeStep(cs);
01064   step->SetMeshComputeStep(meshcs);
01065 
01066   if(!readAllEntityInfo || step->GetEntityInfoLoaded())
01067     return;
01068 
01069   step->SetEntityInfoLoaded(1);
01070   
01071   vtkMedFile* file = field->GetParentFile();
01072   vtkMedMesh* mesh = file->GetMesh(field->GetMeshName());
01073   
01074   if(mesh == NULL)
01075     return;
01076   
01077   std::set<vtkMedEntity> tmp_entities;
01078   std::set<vtkMedEntity> entities;
01079   mesh->GatherMedEntities(tmp_entities);
01080   
01081   std::set<vtkMedEntity>::iterator tmp_entity_it = tmp_entities.begin();
01082   while(tmp_entity_it != tmp_entities.end())
01083     {
01084     vtkMedEntity entity = *tmp_entity_it;
01085     tmp_entity_it++;
01086     entities.insert(entity);
01087     if(entity.EntityType == MED_CELL)
01088       {
01089       vtkMedEntity newEntity;
01090       newEntity.EntityType = MED_NODE_ELEMENT;
01091       newEntity.GeometryType = entity.GeometryType;
01092       newEntity.GeometryName = entity.GeometryName;
01093       entities.insert(newEntity);
01094       }
01095     }
01096   
01097   vtkMedEntity pointEntity;
01098   pointEntity.EntityType = MED_NODE;
01099   pointEntity.GeometryType = MED_NONE;
01100   pointEntity.GeometryName = MED_NAME_BLANK;
01101   entities.insert(pointEntity);
01102 
01103   std::set<vtkMedEntity>::iterator entity_it = entities.begin();
01104   while(entity_it != entities.end())
01105     {
01106     vtkMedEntity entity = *entity_it;
01107     entity_it++;
01108 
01109     med_int nvalues = 0;
01110     med_int nprofile;
01111     char profilename[MED_NAME_SIZE+1] = "";
01112     char localizationname[MED_NAME_SIZE+1] = "";
01113 
01114     nprofile = MEDfieldnProfile(
01115         this->FileId, 
01116         field->GetName(),
01117         step->GetComputeStep().TimeIt,
01118         step->GetComputeStep().IterationIt,
01119         entity.EntityType,
01120         entity.GeometryType,
01121         profilename,
01122         localizationname);
01123     if(nprofile < 0)
01124       {
01125       vtkErrorMacro("MEDfieldnProfile");
01126       continue;
01127       }
01128 
01129     bool hasprofile = (nprofile > 0);
01130     if(!hasprofile)
01131       {
01132       nprofile = 1;
01133       }
01134 
01135     med_int profilesize;
01136     med_int nintegrationpoint;
01137     
01138     for(int pid=0; pid<nprofile; pid++)
01139       {
01140       med_int medid = (hasprofile ? pid+1 : -1);
01141       nvalues = MEDfieldnValueWithProfile(
01142                   this->FileId, 
01143                   field->GetName(),
01144                   step->GetComputeStep().TimeIt,
01145                   step->GetComputeStep().IterationIt,
01146                   entity.EntityType,
01147                   entity.GeometryType,
01148                   medid,
01149                   MED_COMPACT_PFLMODE,
01150                   profilename,
01151                   &profilesize,
01152                   localizationname,
01153                   &nintegrationpoint);
01154             
01155       if(nvalues < 0)
01156         {
01157         vtkErrorMacro("MEDfieldnValueWithProfile");
01158         continue;
01159         }
01160       else if(nvalues > 0)
01161         {
01162         // I have found a profile with values, stop the loop here
01163         break;
01164         }
01165       }
01166 
01167     if(nvalues > 0)
01168       {
01169       vtkMedFieldOverEntity* fieldOverEntity = vtkMedFieldOverEntity::New();
01170       step->AppendFieldOverEntity(fieldOverEntity);
01171       fieldOverEntity->Delete();
01172 
01173       fieldOverEntity->SetParentStep(step);
01174       fieldOverEntity->SetEntity(entity);
01175 
01176       this->ReadFieldOverEntityInformation(fieldOverEntity);
01177       }
01178     }
01179 }
01180 
01181 void vtkMedDriver30::ReadFieldOverEntityInformation(vtkMedFieldOverEntity* fieldOverEntity)
01182 {
01183   FileOpen open(this);
01184 
01185   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
01186   vtkMedField* field = step->GetParentField();
01187   vtkMedEntity entity = fieldOverEntity->GetEntity();
01188 
01189   const char* fieldName = field->GetName();
01190   const vtkMedComputeStep& cs = step->GetComputeStep();
01191 
01192   char profilename[MED_NAME_SIZE+1] = "";
01193   char localizationname[MED_NAME_SIZE+1] = "";
01194 
01195   med_int nProfiles = MEDfieldnProfile(
01196                         this->FileId,
01197                         fieldName,
01198                         cs.TimeIt,
01199                         cs.IterationIt,
01200                         entity.EntityType,
01201                         entity.GeometryType,
01202                         profilename,
01203                         localizationname);
01204 
01205   if(nProfiles < 0)
01206     {
01207     vtkErrorMacro("MEDfieldnProfile");
01208     }
01209   else if(nProfiles == 0)
01210     {
01211     fieldOverEntity->SetHasProfile(0);
01212     nProfiles = 1;
01213     }
01214   else
01215     {
01216     fieldOverEntity->SetHasProfile(1);
01217     }
01218   fieldOverEntity->AllocateNumberOfFieldOnProfile(nProfiles);
01219   for(int profit = 0; profit < nProfiles; profit++)
01220     {
01221     vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(profit);
01222     med_int medid = (fieldOverEntity->GetHasProfile()? profit+1: -1);
01223     fop->SetMedIterator(medid);
01224     fop->SetParentFieldOverEntity(fieldOverEntity);
01225     this->ReadFieldOnProfileInformation(fop);
01226     }
01227 }
01228 
01229 void vtkMedDriver30::ReadFieldOnProfileInformation(vtkMedFieldOnProfile* fop)
01230 {
01231   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
01232   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
01233   vtkMedField* field = step->GetParentField();
01234 
01235   const vtkMedComputeStep& cs = step->GetComputeStep();
01236   med_int profilesize;
01237   med_int nbofintegrationpoint;
01238 
01239   char profileName[MED_NAME_SIZE+1] = "";
01240   char localizationName[MED_NAME_SIZE+1] = "";
01241 
01242   med_int nvalue = MEDfieldnValueWithProfile(FileId,
01243                     field->GetName(),
01244                     cs.TimeIt,
01245                     cs.IterationIt,
01246                     fieldOverEntity->GetEntity().EntityType,
01247                     fieldOverEntity->GetEntity().GeometryType,
01248                     fop->GetMedIterator(),
01249                     MED_COMPACT_STMODE,
01250                     profileName,
01251                     &profilesize,
01252                     localizationName,
01253                     &nbofintegrationpoint);
01254 
01255   if(nvalue < 0)
01256     {
01257     vtkErrorMacro("Error while reading MEDfieldnValueWithProfile");
01258     }
01259 
01260   fop->SetProfileName(profileName);
01261   fop->SetLocalizationName(localizationName);
01262   fop->SetNumberOfValues(nvalue);
01263   fop->SetNumberOfIntegrationPoint(nbofintegrationpoint);
01264   fop->SetProfileSize(profilesize);
01265 }
01266 
01267 void vtkMedDriver30::ReadMeshInformation(vtkMedMesh* mesh)
01268 {
01269   FileOpen open(this);
01270 
01271   med_int mdim = 0;
01272   med_int sdim = 0;
01273   med_mesh_type meshtype;
01274 
01275   med_sorting_type sortingtype;
01276   med_int nstep = 0;
01277   med_axis_type axistype;
01278   med_int naxis;
01279 
01280   if ( (naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator())) <0 )
01281     {
01282     vtkDebugMacro("Error reading mesh axis number");
01283     }
01284 
01285   if(naxis == 0)
01286     naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator());
01287 
01288   char axisname[3*MED_SNAME_SIZE+1] = "";
01289   char axisunit[3*MED_SNAME_SIZE+1] = "";
01290   char name[MED_NAME_SIZE+1] = "";
01291   char description[MED_COMMENT_SIZE+1] = "";
01292   char timeUnit[MED_SNAME_SIZE+1] = "";
01293 
01294   if( MEDmeshInfo( this->FileId,
01295         mesh->GetMedIterator(),
01296         name,
01297         &sdim,
01298         &mdim,
01299         &meshtype,
01300         description,
01301         timeUnit,
01302         &sortingtype,
01303         &nstep,
01304         &axistype,
01305         axisname,
01306         axisunit ) )
01307     {
01308     vtkErrorMacro("Error reading mesh");
01309     }
01310 
01311   mesh->SetName(name);
01312   mesh->SetDescription(description);
01313   mesh->SetTimeUnit(timeUnit);
01314   mesh->SetSpaceDimension(sdim);
01315   mesh->SetMeshDimension(mdim);
01316   mesh->SetMeshType(meshtype);
01317   mesh->SetSortingType(sortingtype);
01318   mesh->SetAxisType(axistype);
01319   mesh->SetNumberOfAxis(naxis);
01320 
01321   for(int axis = 0; axis < naxis; axis++)
01322     {
01323     char theaxisname[MED_SNAME_SIZE+1] = "";
01324     char theaxisunit[MED_SNAME_SIZE+1] = "";
01325     strncpy(theaxisname, axisname + axis*MED_SNAME_SIZE, MED_SNAME_SIZE);
01326     strncpy(theaxisunit, axisunit + axis*MED_SNAME_SIZE, MED_SNAME_SIZE);
01327     mesh->GetAxisName()->SetValue(axis, theaxisname);
01328     mesh->GetAxisUnit()->SetValue(axis, theaxisunit);
01329     }
01330 
01331   char universalName[MED_LNAME_SIZE+1] = "";
01332 
01333   if(MEDmeshUniversalNameRd(this->FileId, name,
01334       universalName) < 0)
01335     {
01336     vtkDebugMacro("MEDmeshUniversalNameRd < 0");
01337     }
01338   mesh->SetUniversalName(universalName);
01339 
01340   // read the Information data of all families.
01341   // writing the family 0 is optional,
01342   // but I need it, so add it if it is not here.
01343 
01344   med_int nfam = MEDnFamily(this->FileId, name);
01345 
01346   for(int index = 0; index < nfam; index++)
01347     {
01348     vtkMedFamily* family = vtkMedFamily::New();
01349     family->SetMedIterator(index + 1);
01350     this->ReadFamilyInformation(mesh, family);
01351     family->Delete();
01352     }
01353 
01354   // this creates a family 0 if none has been read
01355   vtkMedFamily* familyZeroOnCell = mesh->GetOrCreateCellFamilyById(0);
01356   vtkMedFamily* familyZeroOnPoint = mesh->GetOrCreatePointFamilyById(0);
01357 
01358   // Load Information regarding the grid type
01359   if(meshtype == MED_STRUCTURED_MESH)
01360     {
01361     // Do it for structured data
01362     med_grid_type mtg;
01363     if(MEDmeshGridTypeRd(FileId, name, &mtg) < 0)
01364       {
01365       vtkErrorMacro("Error during structured grid Information loading.");
01366       return;
01367       }
01368     mesh->SetStructuredGridType(mtg);
01369     }
01370 
01371   vtkMedGrid* previousGrid = NULL;
01372   for(int gid=1; gid <= nstep; gid++)
01373     {
01374     vtkMedComputeStep cs;
01375     if(MEDmeshComputationStepInfo(FileId,
01376                                   name,
01377                                   gid,
01378                                   &cs.TimeIt,
01379                                   &cs.IterationIt,
01380                                   &cs.TimeOrFrequency) < 0)
01381       {
01382       vtkErrorMacro("MEDmeshComputationStepInfo error");
01383       }
01384     // Load Information regarding the grid type
01385     vtkMedGrid* grid = NULL;
01386     if(meshtype == MED_STRUCTURED_MESH)
01387       {
01388       switch(mesh->GetStructuredGridType())
01389         {
01390         case MED_CARTESIAN_GRID:
01391           grid = vtkMedCartesianGrid::New();
01392           break;
01393         case MED_POLAR_GRID:
01394           grid = vtkMedPolarGrid::New();
01395           break;
01396         case MED_CURVILINEAR_GRID:
01397           grid = vtkMedCurvilinearGrid::New();
01398           break;
01399         default:
01400           vtkErrorMacro("Unknown structured grid type " << mesh->GetStructuredGridType());
01401           return;
01402         }
01403       }
01404     else //(mesh->GetType() == MED_STRUCTURED_MESH)
01405       {
01406       grid = vtkMedUnstructuredGrid::New();
01407       }
01408     grid->SetParentMesh(mesh);
01409     grid->SetComputeStep(cs);
01410     this->ReadGridInformation(grid);
01411     mesh->AddGridStep(grid);
01412     grid->Delete();
01413     grid->SetPreviousGrid(previousGrid);
01414     previousGrid = grid;
01415     }
01416 }
01417 
01418 void vtkMedDriver30::ReadLocalizationInformation(vtkMedLocalization* loc)
01419 {
01420   FileOpen open(this);
01421 
01422   med_int ngp;
01423   med_int spaceDimension;
01424   med_geometry_type type_geo;
01425   med_geometry_type sectiongeotype;
01426   med_int nsectionmeshcell;
01427 
01428   char name[MED_NAME_SIZE+1] = "";
01429   char interpolationName[MED_NAME_SIZE+1] = "";
01430   char sectionName[MED_NAME_SIZE+1] = "";
01431 
01432   if(MEDlocalizationInfo(
01433       this->FileId,
01434       loc->GetMedIterator(),
01435       name,
01436       &type_geo,
01437       &spaceDimension,
01438       &ngp,
01439       interpolationName,
01440       sectionName,
01441       &nsectionmeshcell,
01442       &sectiongeotype ) < 0)
01443     {
01444     vtkErrorMacro("Reading information on quadrature points definition : "
01445         << loc->GetMedIterator());
01446     }
01447 
01448   loc->SetName(name);
01449   loc->SetInterpolationName(interpolationName);
01450   loc->SetSectionName(sectionName);
01451   loc->SetNumberOfQuadraturePoint(ngp);
01452   loc->SetGeometryType(type_geo);
01453   loc->SetSpaceDimension(spaceDimension);
01454   loc->SetNumberOfCellInSection(nsectionmeshcell);
01455   loc->SetSectionGeometryType(sectiongeotype);
01456 
01457   med_float *localCoordinates = new med_float[loc->GetSizeOfPointLocalCoordinates()];
01458   med_float *pqLocalCoordinates = new med_float[loc->GetSizeOfQuadraturePointLocalCoordinates()];
01459   med_float *weights = new med_float[loc->GetSizeOfWeights()];
01460 
01461   if(MEDlocalizationRd(FileId,
01462       loc->GetName(),
01463       MED_FULL_INTERLACE,
01464       localCoordinates,
01465       pqLocalCoordinates,
01466       weights) < 0)
01467     {
01468     vtkErrorMacro("MEDlocalizationRd : " << loc->GetName());
01469     }
01470 
01471   vtkDoubleArray* lc = loc->GetPointLocalCoordinates();
01472   vtkDoubleArray *pqlc = loc->GetQuadraturePointLocalCoordinates();
01473   vtkDoubleArray *w = loc->GetWeights();
01474 
01475   lc->SetNumberOfValues(loc->GetSizeOfPointLocalCoordinates());
01476   for(int i=0; i<loc->GetSizeOfPointLocalCoordinates(); i++)
01477     {
01478     lc->SetValue(i, localCoordinates[i]);
01479     }
01480 
01481   pqlc->SetNumberOfValues(loc->GetSizeOfQuadraturePointLocalCoordinates());
01482   for(int i=0; i<loc->GetSizeOfQuadraturePointLocalCoordinates(); i++)
01483     {
01484     pqlc->SetValue(i, pqLocalCoordinates[i]);
01485     }
01486 
01487   w->SetNumberOfValues(loc->GetSizeOfWeights());
01488   for(int i=0; i<loc->GetSizeOfWeights(); i++)
01489     {
01490     w->SetValue(i, weights[i]);
01491     }
01492 }
01493 
01494 void vtkMedDriver30::ReadInterpolationInformation(vtkMedInterpolation* interp)
01495 {
01496 
01497   med_geometry_type geotype;
01498   med_bool cellnode;
01499   med_int nbofbasisfunc;
01500   med_int nbofvariable;
01501   med_int maxdegree;
01502   med_int nmaxcoef;
01503 
01504   char name[MED_NAME_SIZE+1] = "";
01505 
01506   if(MEDinterpInfo (this->FileId,
01507                     interp->GetMedIterator(),
01508                     name,
01509                     &geotype, &cellnode, &nbofbasisfunc,
01510                     &nbofvariable, &maxdegree, &nmaxcoef) < 0)
01511     {
01512     vtkErrorMacro("MEDinterpInfo");
01513     return;
01514     }
01515 
01516   interp->SetName(name);
01517   interp->SetGeometryType(geotype);
01518   interp->SetIsCellNode(cellnode);
01519   interp->SetNumberOfVariable(nbofvariable);
01520   interp->SetMaximumDegree(maxdegree);
01521   interp->SetMaximumNumberOfCoefficient(nmaxcoef);
01522   interp->AllocateNumberOfBasisFunction(nbofbasisfunc);
01523 
01524   for(int basisid=0; basisid < nbofbasisfunc; basisid++)
01525     {
01526     vtkMedFraction* func = interp->GetBasisFunction(basisid);
01527     func->SetNumberOfVariable(nbofvariable);
01528 
01529     med_int ncoef = MEDinterpBaseFunctionCoefSize (
01530         this->FileId,
01531         interp->GetName(),
01532         basisid+1);
01533     func->SetNumberOfCoefficients(ncoef);
01534 
01535     if(ncoef <= 0 || nbofvariable <= 0)
01536       continue;
01537 
01538     med_int *power = new med_int[nbofvariable * ncoef];
01539     med_float *coefficient = new med_float[ncoef];
01540 
01541     if(MEDinterpBaseFunctionRd  (
01542         this->FileId,
01543         interp->GetName(),
01544         basisid+1,
01545         &ncoef,
01546         power,
01547         coefficient) < 0)
01548       {
01549       vtkErrorMacro("MEDinterpBaseFunctionRd");
01550       continue;
01551       }
01552     vtkDoubleArray* coeffs = func->GetCoefficients();
01553     for(int cid=0; cid < ncoef; cid++)
01554       {
01555       coeffs->SetValue(cid, coefficient[cid]);
01556       }
01557     vtkIntArray* powers = func->GetPowers();
01558     for(int pid=0; pid < ncoef*nbofvariable; pid++)
01559       {
01560       powers->SetValue(pid, power[pid]);
01561       }
01562 
01563     delete[] power;
01564     delete[] coefficient;
01565     }
01566 }
01567 
01568 void vtkMedDriver30::ReadSupportMeshInformation(
01569     vtkMedMesh* supportMesh)
01570 {
01571   FileOpen open(this);
01572 
01573   char supportmeshname[MED_NAME_SIZE+1] = "";
01574   char description[MED_COMMENT_SIZE+1] = "";
01575   med_int spacedim;
01576   med_int meshdim;
01577   med_axis_type axistype;
01578   char axisname[3*MED_SNAME_SIZE+1] = "";
01579   char axisunit[3*MED_SNAME_SIZE+1] = "";
01580 
01581   if(MEDsupportMeshInfo (this->FileId,
01582                          supportMesh->GetMedIterator(),
01583                          supportmeshname,
01584                          &spacedim,
01585                          &meshdim,
01586                          description,
01587                          &axistype,
01588                          axisname,
01589                          axisunit) < 0)
01590     {
01591     vtkErrorMacro("MEDsupportMeshInfo");
01592     }
01593 
01594   supportMesh->SetName(supportmeshname);
01595   supportMesh->SetDescription(description);
01596   supportMesh->SetSpaceDimension(spacedim);
01597   supportMesh->SetMeshDimension(meshdim);
01598   supportMesh->SetAxisType(axistype);
01599   for(int dim = 0; dim < 3; dim++)
01600     {
01601     char axisname_dim[MED_SNAME_SIZE+1] = "";
01602     char axisunit_dim[MED_SNAME_SIZE+1] = "";
01603 
01604     strncpy(axisname_dim, axisname+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
01605     strncpy(axisunit_dim, axisunit+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
01606 
01607     supportMesh->GetAxisName()->SetValue(dim, axisname_dim);
01608     supportMesh->GetAxisUnit()->SetValue(dim, axisunit_dim);
01609     }
01610 
01611   return;
01612 }
01613 
01614 void vtkMedDriver30::LoadFamilyIds(vtkMedEntityArray* array)
01615 {
01616   if(array->IsFamilyIdsLoaded())
01617     return;
01618 
01619   FileOpen open(this);
01620 
01621   vtkMedGrid* grid = array->GetParentGrid();
01622 
01623   vtkMedComputeStep cs = grid->GetComputeStep();
01624 
01625   // first, find if the family ids are implicit or explicit
01626   med_bool changement, transformation;
01627 
01628   med_int nfamid = MEDmeshnEntity(this->FileId,
01629                       grid->GetParentMesh()->GetName(),
01630                       cs.TimeIt,
01631                       cs.IterationIt,
01632                       array->GetEntity().EntityType,
01633                       array->GetEntity().GeometryType,
01634                       MED_FAMILY_NUMBER,
01635                       MED_NO_CMODE,
01636                       &changement,
01637                       &transformation);
01638 
01639   if(nfamid == array->GetNumberOfEntity())
01640     {
01641 
01642     vtkMedIntArray* famIds = vtkMedIntArray::New();
01643     array->SetFamilyIds(famIds);
01644     famIds->Delete();
01645 
01646     famIds->SetNumberOfTuples(nfamid);
01647 
01648     if ( MEDmeshEntityFamilyNumberRd(
01649             this->FileId,
01650             grid->GetParentMesh()->GetName(),
01651             cs.TimeIt,
01652             cs.IterationIt,
01653             array->GetEntity().EntityType,
01654             array->GetEntity().GeometryType,
01655             famIds->GetPointer(0) ) < 0)
01656       {
01657       vtkWarningMacro("Error loading the family ids of entity "
01658         << array->GetEntity().EntityType
01659         << " " << array->GetEntity().GeometryType
01660         << " on mesh " << grid->GetParentMesh()->GetName());
01661       array->SetFamilyIds(NULL);
01662       }
01663     }
01664   else
01665     {
01666     vtkDebugMacro("NumberOfEntity != Number of family ids");
01667     array->SetFamilyIds(NULL);
01668     }
01669 
01670   array->ComputeFamilies();
01671 }
01672 
01673 void vtkMedDriver30::LoadPointGlobalIds(vtkMedGrid* grid)
01674 {
01675   if(grid->IsPointGlobalIdsLoaded())
01676     return;
01677 
01678   FileOpen open(this);
01679 
01680   vtkMedIntArray* globalIds = vtkMedIntArray::New();
01681   grid->SetPointGlobalIds(globalIds);
01682   globalIds->Delete();
01683 
01684   globalIds->SetNumberOfTuples(grid->GetNumberOfPoints());
01685 
01686   if( MEDmeshEntityNumberRd (
01687         this->FileId,
01688         grid->GetParentMesh()->GetName(),
01689         grid->GetComputeStep().TimeIt,
01690         grid->GetComputeStep().IterationIt,
01691         MED_NODE,
01692         MED_NONE,
01693         globalIds->GetPointer(0) ) < 0)
01694     {
01695     grid->SetPointGlobalIds(NULL);
01696     }
01697 }
01698 
01699 void vtkMedDriver30::LoadCoordinates(vtkMedGrid* grid)
01700 {
01701   if(grid->IsCoordinatesLoaded())
01702     return;
01703 
01704   vtkMedRegularGrid* rgrid = vtkMedRegularGrid::SafeDownCast(grid);
01705   if(rgrid != NULL)
01706     {
01707     this->LoadRegularGridCoordinates(rgrid);
01708     return;
01709     }
01710 
01711   vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
01712   vtkMedCurvilinearGrid* cgrid = vtkMedCurvilinearGrid::SafeDownCast(grid);
01713   if(ugrid == NULL && cgrid == NULL)
01714     {
01715     //TODO : deal with structured grids
01716     vtkWarningMacro("this kind of grid is not yet supported");
01717     return;
01718     }
01719 
01720   if(grid->GetUsePreviousCoordinates())
01721     {
01722     vtkMedGrid* previousgrid = grid->GetPreviousGrid();
01723     if(previousgrid == NULL)
01724       {
01725       vtkErrorMacro("coordiantes have not changed, "
01726                     << "but there is no previous grid!");
01727       return;
01728       }
01729 
01730     this->LoadCoordinates(previousgrid);
01731     if(ugrid != NULL)
01732       ugrid->SetCoordinates(vtkMedUnstructuredGrid::SafeDownCast(previousgrid)
01733                             ->GetCoordinates());
01734     if(cgrid != NULL)
01735       cgrid->SetCoordinates(vtkMedCurvilinearGrid::SafeDownCast(previousgrid)
01736                             ->GetCoordinates());
01737     }
01738   else
01739     {
01740 
01741     FileOpen open(this);
01742 
01743     vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
01744     if(ugrid != NULL)
01745       ugrid->SetCoordinates(coords);
01746     if(cgrid != NULL)
01747       cgrid->SetCoordinates(coords);
01748     coords->Delete();
01749 
01750     vtkMedComputeStep cs = grid->GetComputeStep();
01751 
01752     coords->SetNumberOfComponents(grid->GetParentMesh()->GetSpaceDimension());
01753     coords->SetNumberOfTuples(grid->GetNumberOfPoints());
01754 
01755     if ( MEDmeshNodeCoordinateRd( this->FileId,
01756                                   grid->GetParentMesh()->GetName(),
01757                                   cs.TimeIt,
01758                                   cs.IterationIt,
01759                                   MED_FULL_INTERLACE,
01760                                   (med_float*) coords->GetVoidPointer(0) ) < 0)
01761       {
01762       vtkErrorMacro("Load Coordinates for mesh "
01763           << grid->GetParentMesh()->GetName());
01764       }
01765     }
01766 }
01767 
01768 void vtkMedDriver30::LoadProfile(vtkMedProfile* profile)
01769 {
01770   if(!profile || profile->IsLoaded())
01771     return;
01772 
01773   FileOpen open(this);
01774 
01775   vtkMedIntArray* indices = vtkMedIntArray::New();
01776   profile->SetIds(indices);
01777   indices->Delete();
01778 
01779   indices->SetNumberOfTuples(profile->GetNumberOfElement());
01780 
01781   char name[MED_NAME_SIZE+1] = "";
01782 
01783   if( MEDprofileRd(this->FileId,
01784                    profile->GetName(),
01785                    indices->GetPointer(0) ) < 0)
01786     {
01787     vtkErrorMacro("Reading profile indices ");
01788     }
01789 }
01790 
01791 void vtkMedDriver30::LoadConnectivity(vtkMedEntityArray* array)
01792 {
01793   if(array->IsConnectivityLoaded())
01794     return;
01795 
01796   FileOpen open(this);
01797 
01798   vtkMedGrid* grid = array->GetParentGrid();
01799 
01800   grid = array->GetParentGrid();
01801 
01802   const char* meshName = grid->GetParentMesh()->GetName();
01803 
01804   vtkMedIntArray* conn = vtkMedIntArray::New();
01805   array->SetConnectivityArray(conn);
01806   conn->Delete();
01807 
01808   vtkMedComputeStep cs = grid->GetComputeStep();
01809 
01810   med_bool change;
01811   med_bool transformation;
01812 
01813   if(array->GetEntity().GeometryType == MED_POLYGON)
01814     {
01815     // first check if we have points
01816     med_int connSize = MEDmeshnEntity(
01817                             this->FileId,
01818                             meshName,
01819                             cs.TimeIt,
01820                             cs.IterationIt,
01821                             array->GetEntity().EntityType,
01822                             MED_POLYGON,
01823                             MED_CONNECTIVITY,
01824                             array->GetConnectivity(),
01825                             &change,
01826                             &transformation);
01827 
01828     if (connSize < 0)
01829       {
01830       vtkErrorMacro(<< "Error while reading polygons connectivity size."
01831                     << endl );
01832             return;
01833       }
01834 
01835     conn->SetNumberOfTuples(connSize);
01836 
01837     // How many polygon in the mesh in nodal connectivity mode
01838     // For the polygons, we get the size of array index
01839     med_int indexsize;
01840     if ((indexsize = MEDmeshnEntity(this->FileId,
01841                                     meshName,
01842                                     cs.TimeIt,
01843                                     cs.IterationIt,
01844                                     array->GetEntity().EntityType,
01845                                     MED_POLYGON,
01846                                     MED_INDEX_NODE,
01847                                     array->GetConnectivity(),
01848                                     &change,
01849                                     &transformation )) < 0)
01850       {
01851       vtkErrorMacro(<< "Error while reading polygons array index." << endl );
01852             return;
01853       }
01854 
01855     vtkMedIntArray* index = vtkMedIntArray::New();
01856     array->SetFaceIndex(index);
01857     index->Delete();
01858 
01859     index->SetNumberOfTuples( indexsize );
01860 
01861     if ( MEDmeshPolygonRd(this->FileId,
01862                           meshName,
01863                           cs.TimeIt,
01864                           cs.IterationIt,
01865                           array->GetEntity().EntityType,
01866                           array->GetConnectivity(),
01867                           index->GetPointer(0),
01868                           conn->GetPointer(0) ) < 0)
01869       {
01870       vtkErrorMacro(<< "MEDmeshPolygonRd");
01871       return;
01872       }
01873     }
01874   else if(array->GetEntity().GeometryType == MED_POLYHEDRON)
01875     {
01876 
01877     vtkIdType connSize = MEDmeshnEntity(this->FileId,
01878                                         meshName,
01879                                         grid->GetComputeStep().TimeIt,
01880                                         grid->GetComputeStep().IterationIt,
01881                                         array->GetEntity().EntityType,
01882                                         MED_POLYHEDRON,
01883                                         MED_CONNECTIVITY,
01884                                         array->GetConnectivity(),
01885                                         &change,
01886                                         &transformation);
01887     if (connSize < 0)
01888       {
01889       vtkErrorMacro(<< "Error while reading polyhedrons connectivity size."
01890                     << endl );
01891             return;
01892       }
01893 
01894     conn->SetNumberOfTuples(connSize);
01895 
01896     vtkMedIntArray* faceIndex = vtkMedIntArray::New();
01897     array->SetFaceIndex(faceIndex);
01898     faceIndex->Delete();
01899 
01900     vtkMedIntArray* nodeIndex = vtkMedIntArray::New();
01901     array->SetNodeIndex(nodeIndex);
01902     nodeIndex->Delete();
01903 
01904     vtkIdType np = array->GetNumberOfEntity() + 1;
01905     faceIndex->SetNumberOfTuples(np);
01906 
01907     med_int nodeIndexSize;
01908 
01909     if ((nodeIndexSize = MEDmeshnEntity(this->FileId,
01910                                         meshName,
01911                                         grid->GetComputeStep().TimeIt,
01912                                         grid->GetComputeStep().IterationIt,
01913                                         array->GetEntity().EntityType,
01914                                         MED_POLYHEDRON,
01915                                         MED_INDEX_NODE,
01916                                         array->GetConnectivity(),
01917                                         &change,
01918                                         &transformation )) < 0)
01919       {
01920       vtkErrorMacro(<< "Error while reading polygons array index." << endl );
01921             return;
01922       }
01923 
01924     nodeIndex->SetNumberOfTuples(nodeIndexSize);
01925 
01926     if (MEDmeshPolyhedronRd(this->FileId,
01927                             meshName,
01928                             cs.TimeIt,
01929                             cs.IterationIt,
01930                             array->GetEntity().EntityType,
01931                             array->GetConnectivity(),
01932                             faceIndex->GetPointer(0),
01933                             nodeIndex->GetPointer(0),
01934                             conn->GetPointer(0) ) < 0)
01935       {
01936       vtkErrorMacro("Error while reading connectivity of polyhedrons");
01937       return;
01938       }
01939 
01940     }
01941   else
01942     {
01943     bool doReadConnectivity = true;
01944     if(array->GetConnectivity() == MED_NODAL)
01945       {
01946       if(array->GetEntity().EntityType == MED_STRUCT_ELEMENT)
01947         {
01948         if(array->GetStructElement() == NULL)
01949           {
01950           vtkErrorMacro("Entity type = MED_STRUCT_ELEMENT, but StructElement is not set!");
01951           return;
01952           }
01953         vtkIdType ntuple = array->GetNumberOfEntity()
01954                            * array->GetStructElement()->GetConnectivitySize();
01955 
01956         conn->SetNumberOfTuples(ntuple);
01957         // particles are special : connectivity is not stored in the med file
01958         if(strcmp(array->GetStructElement()->GetName(), MED_PARTICLE_NAME) == 0 )
01959           {
01960           for(vtkIdType cellId = 0; cellId < ntuple; cellId++)
01961             {
01962             conn->SetValue(cellId, cellId+1);
01963             }
01964           doReadConnectivity = false;
01965           }
01966         }
01967       else
01968         {
01969         conn->SetNumberOfTuples(array->GetNumberOfEntity()
01970             * vtkMedUtilities::GetNumberOfPoint(
01971                 array->GetEntity().GeometryType));
01972         }
01973       }
01974     else
01975       {
01976       conn->SetNumberOfTuples(array->GetNumberOfEntity()
01977           * vtkMedUtilities::GetNumberOfSubEntity(
01978               array->GetEntity().GeometryType));
01979       }
01980 
01981     if  (this->ParallelFileId == -1) // also (array->GetFilter() == NULL)
01982       {
01983       if ( (MEDmeshElementConnectivityRd(
01984             this->FileId,
01985             meshName,
01986             cs.TimeIt,
01987             cs.IterationIt,
01988             array->GetEntity().EntityType,
01989             array->GetEntity().GeometryType,
01990             array->GetConnectivity(),
01991             MED_FULL_INTERLACE,
01992             conn->GetPointer(0)) ) < 0)
01993         {
01994         vtkErrorMacro("Error while load connectivity of cells "
01995             << array->GetEntity().GeometryType);
01996         }
01997       }
01998     else
01999       {
02000       med_filter filter = MED_FILTER_INIT;
02001 
02002       int    start;
02003       int    stride;
02004       int    count;
02005       int    blocksize;
02006       int    lastblocksize;
02007       array->GetFilter()->GetFilterSizes(start, stride, count,
02008                                    blocksize, lastblocksize );
02009 
02010       med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
02011                                         array->GetEntity().GeometryType);
02012 
02013       if ( MEDfilterBlockOfEntityCr( this->ParallelFileId,
02014               array->GetNumberOfEntity(),
02015             1, // one is for mesh elements, more than 1 is for fields
02016               nbofconstituentpervalue,
02017             MED_ALL_CONSTITUENT,
02018             MED_FULL_INTERLACE,
02019             MED_COMPACT_STMODE,
02020             MED_NO_PROFILE,
02021             (med_size)start,
02022             (med_size)stride,
02023             (med_size)count,
02024             (med_size)blocksize,
02025             (med_size)lastblocksize,
02026             &filter ) < 0 )
02027         {
02028         vtkErrorMacro("Filter creation ");
02029         }
02030 
02031         if ( (MEDmeshElementConnectivityAdvancedRd(
02032               this->ParallelFileId,
02033               meshName,
02034               cs.TimeIt,
02035               cs.IterationIt,
02036               array->GetEntity().EntityType,
02037               array->GetEntity().GeometryType,
02038               array->GetConnectivity(),
02039               &filter,
02040               conn->GetPointer(0)) ) < 0)
02041           {
02042           vtkErrorMacro("Error while load connectivity of cells "
02043               << array->GetEntity().GeometryType);
02044           }
02045 
02046       if ( MEDfilterClose( &filter ) < 0)
02047           {
02048         vtkErrorMacro("ERROR : filter closing ...");
02049           }
02050 
02051       }
02052     }
02053 }
02054 
02055 void vtkMedDriver30::LoadCellGlobalIds(vtkMedEntityArray* array)
02056 {
02057   if(array->IsGlobalIdsLoaded())
02058     return;
02059 
02060   FileOpen open(this);
02061 
02062   vtkMedIntArray* globalIds = vtkMedIntArray::New();
02063   array->SetGlobalIds(globalIds);
02064   globalIds->Delete();
02065 
02066   globalIds->SetNumberOfTuples(array->GetNumberOfEntity());
02067 
02068   vtkMedGrid* grid = array->GetParentGrid();
02069   vtkMedComputeStep cs = grid->GetComputeStep();
02070 
02071   if( MEDmeshEntityNumberRd (
02072         this->FileId,
02073         grid->GetParentMesh()->GetName(),
02074         cs.TimeIt,
02075         cs.IterationIt,
02076         array->GetEntity().EntityType,
02077         array->GetEntity().GeometryType,
02078         globalIds->GetPointer(0) ) < 0)
02079     {
02080     array->SetGlobalIds(NULL);
02081     }
02082 }
02083 
02084 void vtkMedDriver30::LoadField(vtkMedFieldOnProfile* fop, med_storage_mode mode)
02085 {
02086   FileOpen open(this);
02087 
02088   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
02089   vtkMedFieldStep *step = fieldOverEntity->GetParentStep();
02090   vtkMedField* field = step->GetParentField();
02091   const vtkMedComputeStep& cs = step->GetComputeStep();
02092 
02093   vtkDataArray* data = vtkMedUtilities::NewArray(field->GetDataType());
02094   fop->SetData(data);
02095   data->Delete();
02096 
02097   med_int size;
02098   if(mode == MED_COMPACT_STMODE)
02099     {
02100     size = fop->GetNumberOfValues();
02101     }
02102   else
02103     {
02104     med_int profilesize;
02105     med_int nbofintegrationpoint;
02106     char profileName[MED_NAME_SIZE+1] = "";
02107     char localizationName[MED_NAME_SIZE+1] = "";
02108     size = MEDfieldnValueWithProfile(this->FileId,
02109                 field->GetName(),
02110                 cs.TimeIt,
02111                 cs.IterationIt,
02112                 fieldOverEntity->GetEntity().EntityType,
02113                 fieldOverEntity->GetEntity().GeometryType,
02114                 fop->GetMedIterator(),
02115                 MED_GLOBAL_STMODE,
02116                 profileName,
02117                 &profilesize,
02118                 localizationName,
02119                 &nbofintegrationpoint);
02120     }
02121 
02122   if(fop->GetNumberOfIntegrationPoint() > 1)
02123     {
02124     size *= fop->GetNumberOfIntegrationPoint();
02125     }
02126 
02127   data->SetNumberOfComponents(field->GetNumberOfComponent());
02128   data->SetNumberOfTuples(size);
02129   if  (this->ParallelFileId == -1)
02130     {
02131     if ( MEDfieldValueWithProfileRd(
02132           this->FileId,
02133           field->GetName(),
02134           cs.TimeIt,
02135           cs.IterationIt,
02136           fieldOverEntity->GetEntity().EntityType,
02137           fieldOverEntity->GetEntity().GeometryType,
02138           mode,
02139           fop->GetProfileName(),
02140           MED_FULL_INTERLACE,
02141           MED_ALL_CONSTITUENT,
02142           (unsigned char*) data->GetVoidPointer(0) ) < 0)
02143       {
02144       vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
02145       }
02146     }
02147   else
02148     {
02149   if  (field->GetFieldType() == vtkMedField::CellField)
02150     {
02151     med_filter filter = MED_FILTER_INIT;
02152 
02153     int    start;
02154     int    stride;
02155     int    count;
02156     int    blocksize;
02157     int    lastblocksize;
02158     fop->GetFilter()->GetFilterSizes(start, stride, count,
02159                                  blocksize, lastblocksize );
02160 
02161     if ( MEDfilterBlockOfEntityCr( this->ParallelFileId,
02162         fop->GetNumberOfValues(),
02163           1, // one is for mesh elements, more than 1 is for fields
02164           field->GetNumberOfComponent(),
02165           MED_ALL_CONSTITUENT,
02166           MED_FULL_INTERLACE,
02167           MED_COMPACT_STMODE,
02168           MED_NO_PROFILE,
02169           (med_size)start,
02170           (med_size)stride,
02171           (med_size)count,
02172           (med_size)blocksize,
02173           (med_size)lastblocksize,
02174           &filter ) < 0 )
02175       {
02176       vtkErrorMacro("Filter creation ");
02177       }
02178 
02179     if ( MEDfieldValueAdvancedRd(
02180             this->ParallelFileId,
02181             field->GetName(),
02182             cs.TimeIt,
02183             cs.IterationIt,
02184             fieldOverEntity->GetEntity().EntityType,
02185             fieldOverEntity->GetEntity().GeometryType,
02186             &filter,
02187             (unsigned char*) data->GetVoidPointer(0) ) < 0)
02188         {
02189         vtkErrorMacro("Error on MEDfieldValueAdvancedRd");
02190         }
02191 
02192     if ( MEDfilterClose( &filter ) < 0)
02193         {
02194       vtkErrorMacro("ERROR : filter closing ...");
02195         }
02196       }
02197   else
02198     {//TODO : option utilisateur pour desactiver ou non les champs avec profile en //
02199     if ( MEDfieldValueWithProfileRd(
02200               this->FileId,
02201               field->GetName(),
02202               cs.TimeIt,
02203               cs.IterationIt,
02204               fieldOverEntity->GetEntity().EntityType,
02205               fieldOverEntity->GetEntity().GeometryType,
02206               mode,
02207               fop->GetProfileName(),
02208               MED_FULL_INTERLACE,
02209               MED_ALL_CONSTITUENT,
02210               (unsigned char*) data->GetVoidPointer(0) ) < 0)
02211           {
02212           vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
02213           }
02214     }
02215     }
02216 }
02217 
02218 void vtkMedDriver30::LoadVariableAttribute(vtkMedVariableAttribute* varatt,
02219                                            vtkMedEntityArray* array)
02220 {
02221   FileOpen open(this);
02222 
02223   void  *value = NULL;
02224 
02225   vtkAbstractArray* valuearray = array->GetVariableAttributeValue(varatt);
02226   // first test if this is already loaded
02227   if(valuearray != NULL && valuearray->GetNumberOfTuples() > 0)
02228     return;
02229 
02230   if(valuearray == NULL)
02231     {
02232     valuearray = vtkMedUtilities::NewArray(varatt->GetAttributeType());
02233     array->SetVariableAttributeValues(varatt, valuearray);
02234     valuearray->Delete();
02235     }
02236 
02237   valuearray->SetNumberOfComponents(varatt->GetNumberOfComponent());
02238   valuearray->SetNumberOfTuples(array->GetNumberOfEntity());
02239   valuearray->SetName(varatt->GetName());
02240 
02241   vtkSmartPointer<vtkCharArray> chararray = vtkSmartPointer<vtkCharArray>::New();
02242 
02243   if(varatt->GetAttributeType() != MED_ATT_NAME)
02244     {
02245     value = valuearray->GetVoidPointer(0);
02246     }
02247   else
02248     {
02249     chararray->SetNumberOfValues(varatt->GetNumberOfComponent() *
02250                                   array->GetNumberOfEntity() *
02251                                   MED_NAME_SIZE);
02252 
02253     value = chararray->GetVoidPointer(0);
02254     }
02255 
02256   vtkMedComputeStep cs = array->GetParentGrid()->GetComputeStep();
02257 
02258   if(MEDmeshStructElementVarAttRd(
02259       this->FileId,
02260       array->GetParentGrid()->GetParentMesh()->GetName(),
02261       cs.TimeIt,
02262       cs.IterationIt,
02263       varatt->GetParentStructElement()->GetGeometryType(),
02264       varatt->GetName(),
02265       value) < 0)
02266     {
02267 
02268     if(cs.IterationIt == MED_NO_IT && cs.TimeIt == MED_NO_DT && cs.TimeOrFrequency == MED_UNDEF_DT)
02269       {
02270       vtkErrorMacro("MEDmeshStructElementVarAttRd");
02271       return;
02272       }
02273     // try to see if I can reuse
02274     // the variable attributes of the NO_DT, NO_IT compute step
02275     vtkMedComputeStep nocs;
02276     nocs.IterationIt = MED_NO_IT;
02277     nocs.TimeIt = MED_NO_DT;
02278     nocs.TimeOrFrequency = MED_UNDEF_DT;
02279     vtkMedEntityArray* nocs_array =
02280         array->GetParentGrid()->GetParentMesh()->GetGridStep(nocs)->GetEntityArray(array->GetEntity());
02281     if(nocs_array == NULL)
02282       {
02283       nocs_array = array->GetParentGrid()->GetParentMesh()->GetGridStep(0)->GetEntityArray(array->GetEntity());
02284       }
02285 
02286     if(nocs_array == NULL || nocs_array == array)
02287       {
02288       // try to force load the default compute step.
02289       if(MEDmeshStructElementVarAttRd(
02290           this->FileId,
02291           array->GetParentGrid()->GetParentMesh()->GetName(),
02292           nocs.TimeIt,
02293           nocs.IterationIt,
02294           varatt->GetParentStructElement()->GetGeometryType(),
02295           varatt->GetName(),
02296           value) < 0)
02297         {
02298         vtkErrorMacro("MEDmeshStructElementVarAttRd");
02299         return;
02300         }
02301       }
02302     else
02303       {
02304       this->LoadVariableAttribute(varatt, nocs_array);
02305       array->SetVariableAttributeValues(varatt, nocs_array->GetVariableAttributeValue(varatt));
02306       return;
02307       }
02308     }
02309 
02310   // If I am here, it means that I read the values
02311   if(varatt->GetAttributeType() == MED_ATT_NAME)
02312     {
02313     char current_name[MED_NAME_SIZE+1] = "";
02314     vtkStringArray* sarray = vtkStringArray::SafeDownCast(valuearray);
02315     for(vtkIdType id = 0; id < varatt->GetNumberOfComponent() *
02316                        array->GetNumberOfEntity(); id++)
02317       {
02318       memset(current_name, '\0', MED_NAME_SIZE+1);
02319       strncpy(current_name, ((char*)value) + id*MED_NAME_SIZE, MED_NAME_SIZE);
02320       sarray->SetValue(id, current_name);
02321       }
02322     }
02323 
02324   return;
02325 }
02326 
02327 void vtkMedDriver30::PrintSelf(ostream& os, vtkIndent indent)
02328 {
02329   this->Superclass::PrintSelf(os, indent);
02330 }