1 // Copyright (C) 2010-2012 CEA/DEN, EDF R&D
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.
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.
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
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "vtkMedDriver30.h"
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"
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"
54 #include "vtkMultiProcessController.h"
61 vtkCxxRevisionMacro(vtkMedDriver30, "$Revision$")
62 vtkStandardNewMacro(vtkMedDriver30)
64 vtkMedDriver30::vtkMedDriver30() : vtkMedDriver()
68 vtkMedDriver30::~vtkMedDriver30()
73 // load all Information data associated with this cartesian grid.
74 void vtkMedDriver30::ReadRegularGridInformation(vtkMedRegularGrid* grid)
78 grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
80 for(int axis=0; axis < grid->GetDimension(); axis++)
83 med_bool coordinatechangement, geotransformation;
85 if ((size = MEDmeshnEntity(
87 grid->GetParentMesh()->GetName(),
88 grid->GetComputeStep().TimeIt,
89 grid->GetComputeStep().IterationIt,
92 (med_data_type)(MED_COORDINATE_AXIS1 + axis),
94 &coordinatechangement,
95 &geotransformation)) < 0)
97 vtkErrorMacro("ERROR : number of coordinates on X axis ...");
100 grid->SetAxisSize(axis, size);
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;
112 entity.EntityType = MED_CELL;
114 switch(grid->GetDimension())
117 entity.GeometryType = MED_POINT1;
120 entity.GeometryType = MED_SEG2;
123 entity.GeometryType = MED_QUAD4;
126 entity.GeometryType = MED_HEXA8;
129 vtkErrorMacro("Unsupported dimension for curvilinear grid : "
130 << grid->GetDimension());
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);
141 // this triggers the creation of undefined families
142 this->LoadFamilyIds(array);
146 void vtkMedDriver30::LoadRegularGridCoordinates(vtkMedRegularGrid* grid)
150 for(int axis=0; axis < grid->GetParentMesh()->GetNumberOfAxis(); axis++)
153 vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
154 grid->SetAxisCoordinate(axis, coords);
157 coords->SetNumberOfComponents(1);
158 coords->SetNumberOfTuples(grid->GetAxisSize(axis));
160 if (MEDmeshGridIndexCoordinateRd(
162 grid->GetParentMesh()->GetName(),
163 grid->GetComputeStep().TimeIt,
164 grid->GetComputeStep().IterationIt,
166 (med_float*)coords->GetVoidPointer(0)) < 0)
168 vtkErrorMacro("ERROR : read axis " << axis << " coordinates ...");
169 grid->SetAxisCoordinate(axis, NULL);
176 // load all Information data associated with this standard grid.
177 void vtkMedDriver30::ReadCurvilinearGridInformation(vtkMedCurvilinearGrid* grid)
181 grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
184 med_bool coordinatechangement, geotransformation;
186 if ((size = MEDmeshnEntity(
188 grid->GetParentMesh()->GetName(),
189 grid->GetComputeStep().TimeIt,
190 grid->GetComputeStep().IterationIt,
195 &coordinatechangement,
196 &geotransformation)) < 0)
198 vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshnEntity");
201 grid->SetNumberOfPoints(size);
204 if(MEDmeshGridStructRd(
206 grid->GetParentMesh()->GetName(),
207 grid->GetComputeStep().TimeIt,
208 grid->GetComputeStep().IterationIt,
211 vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshGridStructRd");
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]));
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))
221 vtkErrorMacro("The total number of points of a Curvilinear grid should "
222 << "be the product of each axis size!");
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;
234 entity.EntityType = MED_CELL;
236 switch(grid->GetDimension())
239 entity.GeometryType = MED_POINT1;
242 entity.GeometryType = MED_SEG2;
245 entity.GeometryType = MED_QUAD4;
248 entity.GeometryType = MED_HEXA8;
251 vtkErrorMacro("Unsupported dimension for curvilinear grid : "
252 << grid->GetDimension());
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);
263 // this triggers the creation of undefined families
264 this->LoadFamilyIds(array);
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)
277 med_bool changement, transformation;
279 const char* meshName = grid->GetParentMesh()->GetName();
281 const vtkMedComputeStep& cs = grid->GetComputeStep();
283 med_int nentity = MEDmeshnEntity(
295 for(med_int geotypeit = 1; geotypeit <= nentity; geotypeit++)
297 // read cell informations
299 entity.EntityType = entityType;
301 char geometryName[MED_NAME_SIZE+1] = "";
303 // this gives us the med_geometry_type
304 if( MEDmeshEntityInfo( FileId, meshName,
310 &entity.GeometryType) < 0)
312 vtkErrorMacro("MEDmeshEntityInfo");
316 entity.GeometryName = geometryName;
319 if(entity.GeometryType == MED_POLYGON)
321 // read the number of cells of this type
322 ncell = MEDmeshnEntity(this->FileId,
331 &transformation ) - 1;
333 else if(entity.GeometryType == MED_POLYHEDRON)
335 // read the number of cells of this type
336 ncell = MEDmeshnEntity(this->FileId,
345 &transformation ) - 1;
349 ncell = MEDmeshnEntity(this->FileId,
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);
370 // this triggers the creation of undefined families
371 this->LoadFamilyIds(array);
377 // load all Information data associated with this unstructured grid.
378 void vtkMedDriver30::ReadUnstructuredGridInformation(
379 vtkMedUnstructuredGrid* grid)
383 vtkMedMesh *mesh = grid->GetParentMesh();
385 const char *meshName = mesh->GetName();
386 med_connectivity_mode connectivity;
389 med_bool transformation;
392 char profilename[MED_NAME_SIZE+1];
393 memset(profilename, '\0', MED_NAME_SIZE+1);
395 vtkMedComputeStep cs = grid->GetComputeStep();
397 // first check if we have points
398 vtkIdType npoints = MEDmeshnEntityWithProfile(
415 grid->SetNumberOfPoints(npoints);
419 if(grid->GetPreviousGrid() == NULL)
421 vtkErrorMacro("No point and no previous grid");
423 grid->SetUsePreviousCoordinates(true);
424 grid->SetNumberOfPoints(grid->GetPreviousGrid()->GetNumberOfPoints());
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);
437 // create the point vtkMedEntityArray support
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);
448 this->LoadFamilyIds(pea);
451 void vtkMedDriver30::ReadFamilyInformation(vtkMedMesh* mesh, vtkMedFamily* family)
457 char* groupNames = NULL;
458 const char* meshName = mesh->GetName();
460 ngroup = MEDnFamilyGroup(FileId, meshName, family->GetMedIterator());
462 bool has_no_group = false;
467 vtkErrorMacro("Error while reading the number of groups");
473 groupNames = new char[ngroup * MED_LNAME_SIZE + 1];
474 memset(groupNames, '\0', ngroup * MED_LNAME_SIZE + 1);
476 // special case for files written by med < 3,
477 // I have to use the 23 interface
478 if(mesh->GetParentFile()->GetVersionMajor() < 3)
480 med_int *attributenumber = NULL;
481 med_int *attributevalue = NULL;
482 char *attributedes = NULL;
484 med_int nattr = MEDnFamily23Attribute(
487 family->GetMedIterator());
491 vtkErrorMacro("MEDnFamily23Attribute");
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);
502 char familyName[MED_NAME_SIZE+1] = "";
504 if(MEDfamily23Info (this->FileId,
506 family->GetMedIterator(),
514 vtkDebugMacro("MEDfamily23Info");
517 family->SetName(familyName);
519 if(attributenumber != NULL)
520 delete[] attributenumber;
521 if(attributevalue != NULL)
522 delete[] attributevalue;
523 if(attributedes != NULL)
524 delete[] attributedes;
528 char familyName[MED_NAME_SIZE+1] = "";
529 if(MEDfamilyInfo( this->FileId,
531 family->GetMedIterator(),
537 "vtkMedDriver30::ReadInformation(vtkMedFamily* family)"
538 << " cannot read family informations.");
541 family->SetName(familyName);
544 family->SetId(familyid);
546 if( familyid == 0 ) {
547 family->SetPointOrCell( vtkMedUtilities::OnCell );
548 mesh->AppendCellFamily( family );
550 //rnv: improve algorithm to determine entity of the family:
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
556 med_int size = MEDmeshnEntity(
568 vtkErrorMacro( "vtkMedDriver30::ReadInformation(vtkMedFamily* family)"
569 <<" cannot read nb Nodes" );
573 vector<med_int> n_fams;
574 n_fams.resize( size );
576 med_int ret_val = MEDmeshEntityFamilyNumberRd( this->FileId,
583 // Remove ZERO FAMILY
584 remove( n_fams.begin(),n_fams.end(), 0 );
586 bool isOnPoints = ( ret_val >= 0) && (find( n_fams.begin(), n_fams.end(), familyid ) != n_fams.end() );
588 family->SetPointOrCell(vtkMedUtilities::OnPoint);
589 mesh->AppendPointFamily(family);
591 family->SetPointOrCell(vtkMedUtilities::OnCell);
592 mesh->AppendCellFamily(family);
596 family->AllocateNumberOfGroup(ngroup);
597 // if there where no group, set the name to the default value
600 memcpy(groupNames, vtkMedUtilities::NoGroupName,
601 strlen(vtkMedUtilities::NoGroupName));
604 for(int index = 0; index < ngroup; index++)
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(),
613 family->SetGroup(index, group);
620 vtkMedFamily* famzero = vtkMedFamily::New();
621 mesh->AppendPointFamily(famzero);
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++)
631 vtkMedGroup* group = mesh->GetOrCreateGroup(
632 vtkMedUtilities::OnPoint,
633 family->GetGroup(gid)->GetName());
634 famzero->SetGroup(gid, group);
635 mesh->AppendPointGroup(group);
640 void vtkMedDriver30::ReadLinkInformation(vtkMedLink* link)
643 char linkMeshName[MED_NAME_SIZE+1] = "";
644 if(MEDlinkInfo(this->FileId,
645 link->GetMedIterator(),
649 vtkErrorMacro("MEDlinkInfo");
652 link->SetMeshName(linkMeshName);
656 char* path = new char[size + 1];
657 memset(path, '\0', size+1);
658 if(MEDlinkRd(this->FileId, link->GetMeshName(), path) < 0)
660 vtkErrorMacro("MEDlinkRd");
661 memset(path, '\0', size+1);
669 void vtkMedDriver30::ReadFileInformation(vtkMedFile* file)
673 char comment[MED_COMMENT_SIZE+1] = "";
675 MEDfileCommentRd(this->FileId,
678 file->SetComment(comment);
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);
686 int nlink = MEDnLink(this->FileId);
687 file->AllocateNumberOfLink(nlink);
688 for(int linkid=0; linkid<nlink; linkid++)
690 vtkMedLink* link = file->GetLink(linkid);
691 link->SetMedIterator(linkid+1);
692 this->ReadLinkInformation(link);
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))
701 vtkWarningMacro("ATTENTION: The MED Reader cannot read profiles when used in parallel");
704 file->AllocateNumberOfProfile(nprof);
705 for(int i = 0; i < nprof; i++)
707 vtkMedProfile* profile = file->GetProfile(i);
708 profile->SetMedIterator(i + 1);
709 profile->SetParentFile(file);
710 this->ReadProfileInformation(profile);
713 int nloc = MEDnLocalization(this->FileId);
714 file->AllocateNumberOfLocalization(nloc);
715 for(int i = 0; i < nloc; i++)
717 vtkMedLocalization* loc = file->GetLocalization(i);
718 loc->SetMedIterator(i + 1);
719 loc->SetParentFile(file);
720 this->ReadLocalizationInformation(loc);
723 int nsupportmesh = MEDnSupportMesh(this->FileId);
724 file->AllocateNumberOfSupportMesh(nsupportmesh);
725 for(int i = 0; i < nsupportmesh; i++)
727 vtkMedMesh* supportmesh = file->GetSupportMesh(i);
728 supportmesh->SetMedIterator(i + 1);
729 supportmesh->SetParentFile(file);
730 this->ReadSupportMeshInformation(supportmesh);
733 int nmesh = MEDnMesh(this->FileId);
734 file->AllocateNumberOfMesh(nmesh);
735 for(int i = 0; i < nmesh; i++)
737 vtkMedMesh* mesh = file->GetMesh(i);
738 mesh->SetMedIterator(i + 1);
739 mesh->SetParentFile(file);
740 this->ReadMeshInformation(mesh);
743 int nfields = MEDnField(this->FileId);
744 file->AllocateNumberOfField(nfields);
745 for(int i = 0; i < nfields; i++)
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())
754 vtkMedField* newfield = vtkMedField::New();
755 int type = field->GetFirstType();
756 newfield->ExtractFieldType(field, type);
757 file->AppendField(newfield);
762 int nstruct = MEDnStructElement(this->FileId);
764 file->AllocateNumberOfStructElement(nstruct);
765 for(int i = 0; i < nstruct; i++)
767 vtkMedStructElement* structelem = file->GetStructElement(i);
768 structelem->SetMedIterator(i+1);
769 structelem->SetParentFile(file);
770 this->ReadStructElementInformation(structelem);
775 void vtkMedDriver30::ReadStructElementInformation(
776 vtkMedStructElement* structelem)
781 char modelname[MED_NAME_SIZE+1] = "";
782 med_geometry_type mgeotype;
784 char supportmeshname[MED_NAME_SIZE+1] = "";
785 med_entity_type sentitytype;
788 med_geometry_type sgeotype;
789 med_int nbofconstantattribute;
791 med_int nbofvariableattribute;
793 if(MEDstructElementInfo (this->FileId,
794 structelem->GetMedIterator(),
803 &nbofconstantattribute,
805 &nbofvariableattribute) < 0)
807 vtkErrorMacro("Error in MEDstructElementInfo");
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);
822 for(int attit = 0; attit < nbofconstantattribute; attit ++)
824 vtkMedConstantAttribute* constatt = structelem->GetConstantAttribute(attit);
825 constatt->SetMedIterator(attit+1);
826 constatt->SetParentStructElement(structelem);
827 this->ReadConstantAttributeInformation(constatt);
830 for(int attit = 0; attit < nbofvariableattribute; attit ++)
832 vtkMedVariableAttribute* varatt = structelem->GetVariableAttribute(attit);
833 varatt->SetMedIterator(attit+1);
834 varatt->SetParentStructElement(structelem);
835 this->ReadVariableAttributeInformation(varatt);
839 void vtkMedDriver30::ReadConstantAttributeInformation(
840 vtkMedConstantAttribute* constAttr)
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] = "";
852 if(MEDstructElementConstAttInfo(
854 constAttr->GetParentStructElement()->GetName(),
855 constAttr->GetMedIterator(),
863 vtkErrorMacro("MEDstructElementConstAttInfo error");
867 constAttr->SetName(constattname);
868 constAttr->SetAttributeType(constatttype);
869 constAttr->SetNumberOfComponent(nbofcomponent);
870 constAttr->SetSupportEntityType(sentitytype);
871 constAttr->SetProfileName(profilename);
872 constAttr->SetProfileSize(profilesize);
874 vtkAbstractArray* values = vtkMedUtilities::NewArray(constatttype);
877 constAttr->SetValues(values);
880 values->SetNumberOfComponents(nbofcomponent);
881 vtkIdType ntuple = 0;
882 if((strcmp(profilename, MED_NO_PROFILE) != 0) &&
883 (strcmp(profilename, "\0") != 0))
885 ntuple = profilesize;
887 else if(constAttr->GetSupportEntityType() == MED_CELL)
889 ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfCell();
893 ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfNode();
895 values->SetNumberOfTuples(ntuple);
898 vtkSmartPointer<vtkCharArray> buffer = vtkSmartPointer<vtkCharArray>::New();
899 if(constatttype != MED_ATT_NAME)
901 ptr = values->GetVoidPointer(0);
905 buffer->SetNumberOfValues(MED_NAME_SIZE*nbofcomponent*ntuple);
906 ptr = buffer->GetVoidPointer(0);
909 if(MEDstructElementConstAttRd (this->FileId,
910 constAttr->GetParentStructElement()->GetName(),
911 constAttr->GetName(), ptr) < 0)
913 vtkErrorMacro("MEDstructElementConstAttRd");
917 if(constatttype == MED_ATT_NAME)
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++)
924 memset(name, '\0', MED_NAME_SIZE+1);
925 strncpy(name, nameptr + id * MED_NAME_SIZE, MED_NAME_SIZE);
926 names->SetValue(id, name);
933 void vtkMedDriver30::ReadVariableAttributeInformation(
934 vtkMedVariableAttribute* varAttr)
939 char varattname[MED_NAME_SIZE+1] = "";
940 med_attribute_type varatttype;
941 med_int nbofcomponent;
943 if(MEDstructElementVarAttInfo (
945 varAttr->GetParentStructElement()->GetName(),
946 varAttr->GetMedIterator(),
951 vtkErrorMacro("MEDstructElementVarAttInfo");
955 varAttr->SetName(varattname);
956 varAttr->SetAttributeType(varatttype);
957 varAttr->SetNumberOfComponent(nbofcomponent);
962 void vtkMedDriver30::ReadProfileInformation(vtkMedProfile* profile)
967 char profileName[MED_NAME_SIZE+1] = "";
969 if(MEDprofileInfo(this->FileId,
970 profile->GetMedIterator(),
974 vtkErrorMacro("cannot read information on profile"
975 << profile->GetMedIterator());
977 profile->SetName(profileName);
978 profile->SetNumberOfElement(nelem);
981 void vtkMedDriver30::ReadFieldInformation(vtkMedField* field)
985 if (field->GetMedIterator() == 0)
988 int ncomp = MEDfieldnComponent(FileId, field->GetMedIterator());
992 field->SetNumberOfComponent(-1);
993 vtkErrorMacro("cannot read the number of component of field "
994 << field->GetMedIterator())
998 field->SetNumberOfComponent(ncomp);
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);
1005 //med_type_champ dataType;
1006 med_field_type dataType;
1010 char name[MED_NAME_SIZE+1] = "";
1011 char meshName[MED_NAME_SIZE+1] = "";
1012 char unit[MED_SNAME_SIZE+1] = "";
1014 if( MEDfieldInfo( FileId,
1015 field->GetMedIterator(),
1025 vtkErrorMacro("cannot read the informations on field "
1026 << field->GetMedIterator())
1030 field->SetName(name);
1031 field->SetMeshName(meshName);
1032 field->SetTimeUnit(unit);
1033 field->SetDataType(dataType);
1034 field->SetLocal(localmesh);
1036 for(int comp = 0; comp < ncomp; comp++)
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);
1042 char compName[MED_SNAME_SIZE + 1] = "";
1043 memcpy(compName, componentNames + MED_SNAME_SIZE * comp, MED_SNAME_SIZE
1045 field->GetComponentName()->SetValue(comp, compName);
1049 delete[] componentNames;
1051 med_int ninterp = MEDfieldnInterp(FileId, field->GetName());
1054 vtkErrorMacro("Error in MEDfieldnInterp");
1058 field->AllocateNumberOfInterpolation(ninterp);
1060 for(med_int interpit=0; interpit<ninterp; interpit++)
1062 vtkMedInterpolation* interp = field->GetInterpolation(interpit);
1063 interp->SetMedIterator(interpit + 1);
1064 this->ReadInterpolationInformation(interp);
1067 vtkMedFieldStep* previousStep = NULL;
1069 for(med_int csit = 0; csit < nstep; csit++)
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;
1082 void vtkMedDriver30::ReadFieldStepInformation(vtkMedFieldStep* step, bool readAllEntityInfo)
1084 vtkMedComputeStep cs;
1085 vtkMedComputeStep meshcs;
1086 vtkMedField* field = step->GetParentField();
1088 FileOpen open(this);
1090 if( MEDfieldComputingStepMeshInfo(
1093 step->GetMedIterator(),
1096 &cs.TimeOrFrequency,
1098 &meshcs.IterationIt) < 0)
1100 vtkErrorMacro("Error in MEDfieldComputingStepMeshInfo");
1104 step->SetComputeStep(cs);
1105 step->SetMeshComputeStep(meshcs);
1107 if(!readAllEntityInfo || step->GetEntityInfoLoaded())
1110 step->SetEntityInfoLoaded(1);
1112 vtkMedFile* file = field->GetParentFile();
1113 vtkMedMesh* mesh = file->GetMesh(field->GetMeshName());
1118 std::set<vtkMedEntity> tmp_entities;
1119 std::set<vtkMedEntity> entities;
1120 mesh->GatherMedEntities(tmp_entities);
1122 std::set<vtkMedEntity>::iterator tmp_entity_it = tmp_entities.begin();
1123 while(tmp_entity_it != tmp_entities.end())
1125 vtkMedEntity entity = *tmp_entity_it;
1127 entities.insert(entity);
1128 if(entity.EntityType == MED_CELL)
1130 vtkMedEntity newEntity;
1131 newEntity.EntityType = MED_NODE_ELEMENT;
1132 newEntity.GeometryType = entity.GeometryType;
1133 newEntity.GeometryName = entity.GeometryName;
1134 entities.insert(newEntity);
1138 vtkMedEntity pointEntity;
1139 pointEntity.EntityType = MED_NODE;
1140 pointEntity.GeometryType = MED_NONE;
1141 pointEntity.GeometryName = MED_NAME_BLANK;
1142 entities.insert(pointEntity);
1144 std::set<vtkMedEntity>::iterator entity_it = entities.begin();
1145 while(entity_it != entities.end())
1147 vtkMedEntity entity = *entity_it;
1150 med_int nvalues = 0;
1152 char profilename[MED_NAME_SIZE+1] = "";
1153 char localizationname[MED_NAME_SIZE+1] = "";
1155 nprofile = MEDfieldnProfile(
1158 step->GetComputeStep().TimeIt,
1159 step->GetComputeStep().IterationIt,
1161 entity.GeometryType,
1166 vtkErrorMacro("MEDfieldnProfile");
1170 bool hasprofile = (nprofile > 0);
1176 med_int profilesize;
1177 med_int nintegrationpoint;
1179 for(int pid=0; pid<nprofile; pid++)
1181 med_int medid = (hasprofile ? pid+1 : -1);
1182 nvalues = MEDfieldnValueWithProfile(
1185 step->GetComputeStep().TimeIt,
1186 step->GetComputeStep().IterationIt,
1188 entity.GeometryType,
1190 MED_COMPACT_PFLMODE,
1194 &nintegrationpoint);
1198 vtkErrorMacro("MEDfieldnValueWithProfile");
1201 else if(nvalues > 0)
1203 // I have found a profile with values, stop the loop here
1210 vtkMedFieldOverEntity* fieldOverEntity = vtkMedFieldOverEntity::New();
1211 step->AppendFieldOverEntity(fieldOverEntity);
1212 fieldOverEntity->Delete();
1214 fieldOverEntity->SetParentStep(step);
1215 fieldOverEntity->SetEntity(entity);
1217 this->ReadFieldOverEntityInformation(fieldOverEntity);
1222 void vtkMedDriver30::ReadFieldOverEntityInformation(vtkMedFieldOverEntity* fieldOverEntity)
1224 FileOpen open(this);
1226 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
1227 vtkMedField* field = step->GetParentField();
1228 vtkMedEntity entity = fieldOverEntity->GetEntity();
1230 const char* fieldName = field->GetName();
1231 const vtkMedComputeStep& cs = step->GetComputeStep();
1233 char profilename[MED_NAME_SIZE+1] = "";
1234 char localizationname[MED_NAME_SIZE+1] = "";
1236 med_int nProfiles = MEDfieldnProfile(
1242 entity.GeometryType,
1248 vtkErrorMacro("MEDfieldnProfile");
1250 else if(nProfiles == 0)
1252 fieldOverEntity->SetHasProfile(0);
1257 fieldOverEntity->SetHasProfile(1);
1259 fieldOverEntity->AllocateNumberOfFieldOnProfile(nProfiles);
1260 for(int profit = 0; profit < nProfiles; profit++)
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);
1270 void vtkMedDriver30::ReadFieldOnProfileInformation(vtkMedFieldOnProfile* fop)
1272 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
1273 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
1274 vtkMedField* field = step->GetParentField();
1276 const vtkMedComputeStep& cs = step->GetComputeStep();
1277 med_int profilesize;
1278 med_int nbofintegrationpoint;
1280 char profileName[MED_NAME_SIZE+1] = "";
1281 char localizationName[MED_NAME_SIZE+1] = "";
1283 med_int nvalue = MEDfieldnValueWithProfile(FileId,
1287 fieldOverEntity->GetEntity().EntityType,
1288 fieldOverEntity->GetEntity().GeometryType,
1289 fop->GetMedIterator(),
1294 &nbofintegrationpoint);
1298 vtkErrorMacro("Error while reading MEDfieldnValueWithProfile");
1301 fop->SetProfileName(profileName);
1302 fop->SetLocalizationName(localizationName);
1303 fop->SetNumberOfValues(nvalue);
1304 fop->SetNumberOfIntegrationPoint(nbofintegrationpoint);
1305 fop->SetProfileSize(profilesize);
1308 void vtkMedDriver30::ReadMeshInformation(vtkMedMesh* mesh)
1310 FileOpen open(this);
1314 med_mesh_type meshtype;
1316 med_sorting_type sortingtype;
1318 med_axis_type axistype;
1321 if ( (naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator())) <0 )
1323 vtkDebugMacro("Error reading mesh axis number");
1327 naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator());
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] = "";
1335 if( MEDmeshInfo( this->FileId,
1336 mesh->GetMedIterator(),
1349 vtkErrorMacro("Error reading mesh");
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);
1362 for(int axis = 0; axis < naxis; axis++)
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);
1372 char universalName[MED_LNAME_SIZE+1] = "";
1374 if(MEDmeshUniversalNameRd(this->FileId, name,
1377 vtkDebugMacro("MEDmeshUniversalNameRd < 0");
1379 mesh->SetUniversalName(universalName);
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.
1385 med_int nfam = MEDnFamily(this->FileId, name);
1387 for(int index = 0; index < nfam; index++)
1389 vtkMedFamily* family = vtkMedFamily::New();
1390 family->SetMedIterator(index + 1);
1391 this->ReadFamilyInformation(mesh, family);
1395 // this creates a family 0 if none has been read
1396 vtkMedFamily* familyZeroOnCell = mesh->GetOrCreateCellFamilyById(0);
1397 vtkMedFamily* familyZeroOnPoint = mesh->GetOrCreatePointFamilyById(0);
1399 // Load Information regarding the grid type
1400 if(meshtype == MED_STRUCTURED_MESH)
1402 // Do it for structured data
1404 if(MEDmeshGridTypeRd(FileId, name, &mtg) < 0)
1406 vtkErrorMacro("Error during structured grid Information loading.");
1409 mesh->SetStructuredGridType(mtg);
1412 vtkMedGrid* previousGrid = NULL;
1413 for(int gid=1; gid <= nstep; gid++)
1415 vtkMedComputeStep cs;
1416 if(MEDmeshComputationStepInfo(FileId,
1421 &cs.TimeOrFrequency) < 0)
1423 vtkErrorMacro("MEDmeshComputationStepInfo error");
1425 // Load Information regarding the grid type
1426 vtkMedGrid* grid = NULL;
1427 if(meshtype == MED_STRUCTURED_MESH)
1429 switch(mesh->GetStructuredGridType())
1431 case MED_CARTESIAN_GRID:
1432 grid = vtkMedCartesianGrid::New();
1434 case MED_POLAR_GRID:
1435 grid = vtkMedPolarGrid::New();
1437 case MED_CURVILINEAR_GRID:
1438 grid = vtkMedCurvilinearGrid::New();
1441 vtkErrorMacro("Unknown structured grid type " << mesh->GetStructuredGridType());
1445 else //(mesh->GetType() == MED_STRUCTURED_MESH)
1447 grid = vtkMedUnstructuredGrid::New();
1449 grid->SetParentMesh(mesh);
1450 grid->SetComputeStep(cs);
1451 this->ReadGridInformation(grid);
1452 mesh->AddGridStep(grid);
1454 grid->SetPreviousGrid(previousGrid);
1455 previousGrid = grid;
1459 void vtkMedDriver30::ReadLocalizationInformation(vtkMedLocalization* loc)
1461 FileOpen open(this);
1464 med_int spaceDimension;
1465 med_geometry_type type_geo;
1466 med_geometry_type sectiongeotype;
1467 med_int nsectionmeshcell;
1469 char name[MED_NAME_SIZE+1] = "";
1470 char interpolationName[MED_NAME_SIZE+1] = "";
1471 char sectionName[MED_NAME_SIZE+1] = "";
1473 if(MEDlocalizationInfo(
1475 loc->GetMedIterator(),
1483 §iongeotype ) < 0)
1485 vtkErrorMacro("Reading information on quadrature points definition : "
1486 << loc->GetMedIterator());
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);
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()];
1502 if(MEDlocalizationRd(FileId,
1509 vtkErrorMacro("MEDlocalizationRd : " << loc->GetName());
1512 vtkDoubleArray* lc = loc->GetPointLocalCoordinates();
1513 vtkDoubleArray *pqlc = loc->GetQuadraturePointLocalCoordinates();
1514 vtkDoubleArray *w = loc->GetWeights();
1516 lc->SetNumberOfValues(loc->GetSizeOfPointLocalCoordinates());
1517 for(int i=0; i<loc->GetSizeOfPointLocalCoordinates(); i++)
1519 lc->SetValue(i, localCoordinates[i]);
1522 pqlc->SetNumberOfValues(loc->GetSizeOfQuadraturePointLocalCoordinates());
1523 for(int i=0; i<loc->GetSizeOfQuadraturePointLocalCoordinates(); i++)
1525 pqlc->SetValue(i, pqLocalCoordinates[i]);
1528 w->SetNumberOfValues(loc->GetSizeOfWeights());
1529 for(int i=0; i<loc->GetSizeOfWeights(); i++)
1531 w->SetValue(i, weights[i]);
1535 void vtkMedDriver30::ReadInterpolationInformation(vtkMedInterpolation* interp)
1538 med_geometry_type geotype;
1540 med_int nbofbasisfunc;
1541 med_int nbofvariable;
1545 char name[MED_NAME_SIZE+1] = "";
1547 if(MEDinterpInfo (this->FileId,
1548 interp->GetMedIterator(),
1550 &geotype, &cellnode, &nbofbasisfunc,
1551 &nbofvariable, &maxdegree, &nmaxcoef) < 0)
1553 vtkErrorMacro("MEDinterpInfo");
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);
1565 for(int basisid=0; basisid < nbofbasisfunc; basisid++)
1567 vtkMedFraction* func = interp->GetBasisFunction(basisid);
1568 func->SetNumberOfVariable(nbofvariable);
1570 med_int ncoef = MEDinterpBaseFunctionCoefSize (
1574 func->SetNumberOfCoefficients(ncoef);
1576 if(ncoef <= 0 || nbofvariable <= 0)
1579 med_int *power = new med_int[nbofvariable * ncoef];
1580 med_float *coefficient = new med_float[ncoef];
1582 if(MEDinterpBaseFunctionRd (
1590 vtkErrorMacro("MEDinterpBaseFunctionRd");
1593 vtkDoubleArray* coeffs = func->GetCoefficients();
1594 for(int cid=0; cid < ncoef; cid++)
1596 coeffs->SetValue(cid, coefficient[cid]);
1598 vtkIntArray* powers = func->GetPowers();
1599 for(int pid=0; pid < ncoef*nbofvariable; pid++)
1601 powers->SetValue(pid, power[pid]);
1605 delete[] coefficient;
1609 void vtkMedDriver30::ReadSupportMeshInformation(
1610 vtkMedMesh* supportMesh)
1612 FileOpen open(this);
1614 char supportmeshname[MED_NAME_SIZE+1] = "";
1615 char description[MED_COMMENT_SIZE+1] = "";
1618 med_axis_type axistype;
1619 char axisname[3*MED_SNAME_SIZE+1] = "";
1620 char axisunit[3*MED_SNAME_SIZE+1] = "";
1622 if(MEDsupportMeshInfo (this->FileId,
1623 supportMesh->GetMedIterator(),
1632 vtkErrorMacro("MEDsupportMeshInfo");
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++)
1642 char axisname_dim[MED_SNAME_SIZE+1] = "";
1643 char axisunit_dim[MED_SNAME_SIZE+1] = "";
1645 strncpy(axisname_dim, axisname+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
1646 strncpy(axisunit_dim, axisunit+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
1648 supportMesh->GetAxisName()->SetValue(dim, axisname_dim);
1649 supportMesh->GetAxisUnit()->SetValue(dim, axisunit_dim);
1655 void vtkMedDriver30::LoadFamilyIds(vtkMedEntityArray* array)
1657 if(array->IsFamilyIdsLoaded())
1660 FileOpen open(this);
1662 vtkMedGrid* grid = array->GetParentGrid();
1664 vtkMedComputeStep cs = grid->GetComputeStep();
1666 // first, find if the family ids are implicit or explicit
1667 med_bool changement, transformation;
1669 med_int nfamid = MEDmeshnEntity(this->FileId,
1670 grid->GetParentMesh()->GetName(),
1673 array->GetEntity().EntityType,
1674 array->GetEntity().GeometryType,
1680 if(nfamid == array->GetNumberOfEntity())
1683 vtkMedIntArray* famIds = vtkMedIntArray::New();
1684 array->SetFamilyIds(famIds);
1687 famIds->SetNumberOfTuples(nfamid);
1689 if ( MEDmeshEntityFamilyNumberRd(
1691 grid->GetParentMesh()->GetName(),
1694 array->GetEntity().EntityType,
1695 array->GetEntity().GeometryType,
1696 famIds->GetPointer(0) ) < 0)
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);
1707 vtkDebugMacro("NumberOfEntity != Number of family ids");
1708 array->SetFamilyIds(NULL);
1711 array->ComputeFamilies();
1714 void vtkMedDriver30::LoadPointGlobalIds(vtkMedGrid* grid)
1716 if(grid->IsPointGlobalIdsLoaded())
1719 FileOpen open(this);
1721 vtkMedIntArray* globalIds = vtkMedIntArray::New();
1722 grid->SetPointGlobalIds(globalIds);
1723 globalIds->Delete();
1725 globalIds->SetNumberOfTuples(grid->GetNumberOfPoints());
1727 if( MEDmeshEntityNumberRd (
1729 grid->GetParentMesh()->GetName(),
1730 grid->GetComputeStep().TimeIt,
1731 grid->GetComputeStep().IterationIt,
1734 globalIds->GetPointer(0) ) < 0)
1736 grid->SetPointGlobalIds(NULL);
1740 void vtkMedDriver30::LoadCoordinates(vtkMedGrid* grid)
1742 if(grid->IsCoordinatesLoaded())
1745 vtkMedRegularGrid* rgrid = vtkMedRegularGrid::SafeDownCast(grid);
1748 this->LoadRegularGridCoordinates(rgrid);
1752 vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
1753 vtkMedCurvilinearGrid* cgrid = vtkMedCurvilinearGrid::SafeDownCast(grid);
1754 if(ugrid == NULL && cgrid == NULL)
1756 //TODO : deal with structured grids
1757 vtkWarningMacro("this kind of grid is not yet supported");
1761 if(grid->GetUsePreviousCoordinates())
1763 vtkMedGrid* previousgrid = grid->GetPreviousGrid();
1764 if(previousgrid == NULL)
1766 vtkErrorMacro("coordiantes have not changed, "
1767 << "but there is no previous grid!");
1771 this->LoadCoordinates(previousgrid);
1773 ugrid->SetCoordinates(vtkMedUnstructuredGrid::SafeDownCast(previousgrid)
1774 ->GetCoordinates());
1776 cgrid->SetCoordinates(vtkMedCurvilinearGrid::SafeDownCast(previousgrid)
1777 ->GetCoordinates());
1782 FileOpen open(this);
1784 vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
1786 ugrid->SetCoordinates(coords);
1788 cgrid->SetCoordinates(coords);
1791 vtkMedComputeStep cs = grid->GetComputeStep();
1793 coords->SetNumberOfComponents(grid->GetParentMesh()->GetSpaceDimension());
1794 coords->SetNumberOfTuples(grid->GetNumberOfPoints());
1796 if ( MEDmeshNodeCoordinateRd( this->FileId,
1797 grid->GetParentMesh()->GetName(),
1801 (med_float*) coords->GetVoidPointer(0) ) < 0)
1803 vtkErrorMacro("Load Coordinates for mesh "
1804 << grid->GetParentMesh()->GetName());
1809 void vtkMedDriver30::LoadProfile(vtkMedProfile* profile)
1811 if(!profile || profile->IsLoaded())
1814 FileOpen open(this);
1816 vtkMedIntArray* indices = vtkMedIntArray::New();
1817 profile->SetIds(indices);
1820 indices->SetNumberOfTuples(profile->GetNumberOfElement());
1822 char name[MED_NAME_SIZE+1] = "";
1824 if( MEDprofileRd(this->FileId,
1826 indices->GetPointer(0) ) < 0)
1828 vtkErrorMacro("Reading profile indices ");
1832 void vtkMedDriver30::LoadConnectivity(vtkMedEntityArray* array)
1834 if(array->IsConnectivityLoaded())
1837 FileOpen open(this);
1839 vtkMedGrid* grid = array->GetParentGrid();
1841 grid = array->GetParentGrid();
1843 const char* meshName = grid->GetParentMesh()->GetName();
1845 vtkMedIntArray* conn = vtkMedIntArray::New();
1846 array->SetConnectivityArray(conn);
1849 vtkMedComputeStep cs = grid->GetComputeStep();
1852 med_bool transformation;
1854 if(array->GetEntity().GeometryType == MED_POLYGON)
1856 // first check if we have points
1857 med_int connSize = MEDmeshnEntity(
1862 array->GetEntity().EntityType,
1865 array->GetConnectivity(),
1871 vtkErrorMacro(<< "Error while reading polygons connectivity size."
1876 conn->SetNumberOfTuples(connSize);
1878 // How many polygon in the mesh in nodal connectivity mode
1879 // For the polygons, we get the size of array index
1881 if ((indexsize = MEDmeshnEntity(this->FileId,
1885 array->GetEntity().EntityType,
1888 array->GetConnectivity(),
1890 &transformation )) < 0)
1892 vtkErrorMacro(<< "Error while reading polygons array index." << endl );
1896 vtkMedIntArray* index = vtkMedIntArray::New();
1897 array->SetFaceIndex(index);
1900 index->SetNumberOfTuples( indexsize );
1902 if ( MEDmeshPolygonRd(this->FileId,
1906 array->GetEntity().EntityType,
1907 array->GetConnectivity(),
1908 index->GetPointer(0),
1909 conn->GetPointer(0) ) < 0)
1911 vtkErrorMacro(<< "MEDmeshPolygonRd");
1915 else if(array->GetEntity().GeometryType == MED_POLYHEDRON)
1918 vtkIdType connSize = MEDmeshnEntity(this->FileId,
1920 grid->GetComputeStep().TimeIt,
1921 grid->GetComputeStep().IterationIt,
1922 array->GetEntity().EntityType,
1925 array->GetConnectivity(),
1930 vtkErrorMacro(<< "Error while reading polyhedrons connectivity size."
1935 conn->SetNumberOfTuples(connSize);
1937 vtkMedIntArray* faceIndex = vtkMedIntArray::New();
1938 array->SetFaceIndex(faceIndex);
1939 faceIndex->Delete();
1941 vtkMedIntArray* nodeIndex = vtkMedIntArray::New();
1942 array->SetNodeIndex(nodeIndex);
1943 nodeIndex->Delete();
1945 vtkIdType np = array->GetNumberOfEntity() + 1;
1946 faceIndex->SetNumberOfTuples(np);
1948 med_int nodeIndexSize;
1950 if ((nodeIndexSize = MEDmeshnEntity(this->FileId,
1952 grid->GetComputeStep().TimeIt,
1953 grid->GetComputeStep().IterationIt,
1954 array->GetEntity().EntityType,
1957 array->GetConnectivity(),
1959 &transformation )) < 0)
1961 vtkErrorMacro(<< "Error while reading polygons array index." << endl );
1965 nodeIndex->SetNumberOfTuples(nodeIndexSize);
1967 if (MEDmeshPolyhedronRd(this->FileId,
1971 array->GetEntity().EntityType,
1972 array->GetConnectivity(),
1973 faceIndex->GetPointer(0),
1974 nodeIndex->GetPointer(0),
1975 conn->GetPointer(0) ) < 0)
1977 vtkErrorMacro("Error while reading connectivity of polyhedrons");
1984 bool doReadConnectivity = true;
1985 if(array->GetConnectivity() == MED_NODAL)
1987 if(array->GetEntity().EntityType == MED_STRUCT_ELEMENT)
1989 if(array->GetStructElement() == NULL)
1991 vtkErrorMacro("Entity type = MED_STRUCT_ELEMENT, but StructElement is not set!");
1994 vtkIdType ntuple = array->GetNumberOfEntity()
1995 * array->GetStructElement()->GetConnectivitySize();
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 )
2001 for(vtkIdType cellId = 0; cellId < ntuple; cellId++)
2003 conn->SetValue(cellId, cellId+1);
2005 doReadConnectivity = false;
2010 conn->SetNumberOfTuples(array->GetNumberOfEntity()
2011 * vtkMedUtilities::GetNumberOfPoint(
2012 array->GetEntity().GeometryType));
2017 conn->SetNumberOfTuples(array->GetNumberOfEntity()
2018 * vtkMedUtilities::GetNumberOfSubEntity(
2019 array->GetEntity().GeometryType));
2022 if (this->ParallelFileId == -1) // also (array->GetFilter() == NULL)
2024 if ( (MEDmeshElementConnectivityRd(
2029 array->GetEntity().EntityType,
2030 array->GetEntity().GeometryType,
2031 array->GetConnectivity(),
2033 conn->GetPointer(0)) ) < 0)
2035 vtkErrorMacro("Error while load connectivity of cells "
2036 << array->GetEntity().GeometryType);
2041 med_filter filter = MED_FILTER_INIT;
2048 array->GetFilter()->GetFilterSizes(start, stride, count,
2049 blocksize, lastblocksize );
2051 med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
2052 array->GetEntity().GeometryType);
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,
2065 (med_size)blocksize,
2066 (med_size)lastblocksize,
2069 vtkErrorMacro("Filter creation ");
2072 if ( (MEDmeshElementConnectivityAdvancedRd(
2073 this->ParallelFileId,
2077 array->GetEntity().EntityType,
2078 array->GetEntity().GeometryType,
2079 array->GetConnectivity(),
2081 conn->GetPointer(0)) ) < 0)
2083 vtkErrorMacro("Error while load connectivity of cells "
2084 << array->GetEntity().GeometryType);
2087 if ( MEDfilterClose( &filter ) < 0)
2089 vtkErrorMacro("ERROR : filter closing ...");
2096 void vtkMedDriver30::LoadCellGlobalIds(vtkMedEntityArray* array)
2098 if(array->IsGlobalIdsLoaded())
2101 FileOpen open(this);
2103 vtkMedIntArray* globalIds = vtkMedIntArray::New();
2104 array->SetGlobalIds(globalIds);
2105 globalIds->Delete();
2107 globalIds->SetNumberOfTuples(array->GetNumberOfEntity());
2109 vtkMedGrid* grid = array->GetParentGrid();
2110 vtkMedComputeStep cs = grid->GetComputeStep();
2112 if( MEDmeshEntityNumberRd (
2114 grid->GetParentMesh()->GetName(),
2117 array->GetEntity().EntityType,
2118 array->GetEntity().GeometryType,
2119 globalIds->GetPointer(0) ) < 0)
2121 array->SetGlobalIds(NULL);
2125 void vtkMedDriver30::LoadField(vtkMedFieldOnProfile* fop, med_storage_mode mode)
2127 FileOpen open(this);
2129 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2130 vtkMedFieldStep *step = fieldOverEntity->GetParentStep();
2131 vtkMedField* field = step->GetParentField();
2132 const vtkMedComputeStep& cs = step->GetComputeStep();
2134 vtkDataArray* data = vtkMedUtilities::NewArray(field->GetDataType());
2139 if(mode == MED_COMPACT_STMODE)
2141 size = fop->GetNumberOfValues();
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,
2153 fieldOverEntity->GetEntity().EntityType,
2154 fieldOverEntity->GetEntity().GeometryType,
2155 fop->GetMedIterator(),
2160 &nbofintegrationpoint);
2163 if(fop->GetNumberOfIntegrationPoint() > 1)
2165 size *= fop->GetNumberOfIntegrationPoint();
2168 data->SetNumberOfComponents(field->GetNumberOfComponent());
2169 data->SetNumberOfTuples(size);
2170 if (this->ParallelFileId == -1)
2172 if ( MEDfieldValueWithProfileRd(
2177 fieldOverEntity->GetEntity().EntityType,
2178 fieldOverEntity->GetEntity().GeometryType,
2180 fop->GetProfileName(),
2182 MED_ALL_CONSTITUENT,
2183 (unsigned char*) data->GetVoidPointer(0) ) < 0)
2185 vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
2190 if (field->GetFieldType() == vtkMedField::CellField)
2192 med_filter filter = MED_FILTER_INIT;
2199 fop->GetFilter()->GetFilterSizes(start, stride, count,
2200 blocksize, lastblocksize );
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,
2213 (med_size)blocksize,
2214 (med_size)lastblocksize,
2217 vtkErrorMacro("Filter creation ");
2220 if ( MEDfieldValueAdvancedRd(
2221 this->ParallelFileId,
2225 fieldOverEntity->GetEntity().EntityType,
2226 fieldOverEntity->GetEntity().GeometryType,
2228 (unsigned char*) data->GetVoidPointer(0) ) < 0)
2230 vtkErrorMacro("Error on MEDfieldValueAdvancedRd");
2233 if ( MEDfilterClose( &filter ) < 0)
2235 vtkErrorMacro("ERROR : filter closing ...");
2239 {//TODO : option utilisateur pour desactiver ou non les champs avec profile en //
2240 if ( MEDfieldValueWithProfileRd(
2245 fieldOverEntity->GetEntity().EntityType,
2246 fieldOverEntity->GetEntity().GeometryType,
2248 fop->GetProfileName(),
2250 MED_ALL_CONSTITUENT,
2251 (unsigned char*) data->GetVoidPointer(0) ) < 0)
2253 vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
2259 void vtkMedDriver30::LoadVariableAttribute(vtkMedVariableAttribute* varatt,
2260 vtkMedEntityArray* array)
2262 FileOpen open(this);
2266 vtkAbstractArray* valuearray = array->GetVariableAttributeValue(varatt);
2267 // first test if this is already loaded
2268 if(valuearray != NULL && valuearray->GetNumberOfTuples() > 0)
2271 if(valuearray == NULL)
2273 valuearray = vtkMedUtilities::NewArray(varatt->GetAttributeType());
2274 array->SetVariableAttributeValues(varatt, valuearray);
2275 valuearray->Delete();
2278 valuearray->SetNumberOfComponents(varatt->GetNumberOfComponent());
2279 valuearray->SetNumberOfTuples(array->GetNumberOfEntity());
2280 valuearray->SetName(varatt->GetName());
2282 vtkSmartPointer<vtkCharArray> chararray = vtkSmartPointer<vtkCharArray>::New();
2284 if(varatt->GetAttributeType() != MED_ATT_NAME)
2286 value = valuearray->GetVoidPointer(0);
2290 chararray->SetNumberOfValues(varatt->GetNumberOfComponent() *
2291 array->GetNumberOfEntity() *
2294 value = chararray->GetVoidPointer(0);
2297 vtkMedComputeStep cs = array->GetParentGrid()->GetComputeStep();
2299 if(MEDmeshStructElementVarAttRd(
2301 array->GetParentGrid()->GetParentMesh()->GetName(),
2304 varatt->GetParentStructElement()->GetGeometryType(),
2309 if(cs.IterationIt == MED_NO_IT && cs.TimeIt == MED_NO_DT && cs.TimeOrFrequency == MED_UNDEF_DT)
2311 vtkErrorMacro("MEDmeshStructElementVarAttRd");
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)
2324 nocs_array = array->GetParentGrid()->GetParentMesh()->GetGridStep(0)->GetEntityArray(array->GetEntity());
2327 if(nocs_array == NULL || nocs_array == array)
2329 // try to force load the default compute step.
2330 if(MEDmeshStructElementVarAttRd(
2332 array->GetParentGrid()->GetParentMesh()->GetName(),
2335 varatt->GetParentStructElement()->GetGeometryType(),
2339 vtkErrorMacro("MEDmeshStructElementVarAttRd");
2345 this->LoadVariableAttribute(varatt, nocs_array);
2346 array->SetVariableAttributeValues(varatt, nocs_array->GetVariableAttributeValue(varatt));
2351 // If I am here, it means that I read the values
2352 if(varatt->GetAttributeType() == MED_ATT_NAME)
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++)
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);
2368 void vtkMedDriver30::PrintSelf(ostream& os, vtkIndent indent)
2370 this->Superclass::PrintSelf(os, indent);