1 // Copyright (C) 2010-2011 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"
59 // vtkCxxRevisionMacro(vtkMedDriver30, "$Revision$")
60 vtkStandardNewMacro(vtkMedDriver30)
62 vtkMedDriver30::vtkMedDriver30() : vtkMedDriver()
66 vtkMedDriver30::~vtkMedDriver30()
71 // load all Information data associated with this cartesian grid.
72 void vtkMedDriver30::ReadRegularGridInformation(vtkMedRegularGrid* grid)
76 grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
78 for(int axis=0; axis < grid->GetDimension(); axis++)
81 med_bool coordinatechangement, geotransformation;
83 if ((size = MEDmeshnEntity(
85 grid->GetParentMesh()->GetName(),
86 grid->GetComputeStep().TimeIt,
87 grid->GetComputeStep().IterationIt,
90 (med_data_type)(MED_COORDINATE_AXIS1 + axis),
92 &coordinatechangement,
93 &geotransformation)) < 0)
95 vtkErrorMacro("ERROR : number of coordinates on X axis ...");
98 grid->SetAxisSize(axis, size);
102 if(grid->GetAxisSize(0) > 1)
103 ncell *= grid->GetAxisSize(0)-1;
104 if(grid->GetAxisSize(1) > 1)
105 ncell *= grid->GetAxisSize(1)-1;
106 if(grid->GetAxisSize(2) > 1)
107 ncell *= grid->GetAxisSize(2)-1;
110 entity.EntityType = MED_CELL;
112 switch(grid->GetDimension())
115 entity.GeometryType = MED_POINT1;
118 entity.GeometryType = MED_SEG2;
121 entity.GeometryType = MED_QUAD4;
124 entity.GeometryType = MED_HEXA8;
127 vtkErrorMacro("Unsupported dimension for curvilinear grid : "
128 << grid->GetDimension());
132 vtkMedEntityArray* array = vtkMedEntityArray::New();
133 array->SetParentGrid(grid);
134 array->SetNumberOfEntity(ncell);
135 array->SetEntity(entity);
136 array->SetConnectivity(MED_NODAL);
137 grid->AppendEntityArray(array);
139 // this triggers the creation of undefined families
140 this->LoadFamilyIds(array);
144 void vtkMedDriver30::LoadRegularGridCoordinates(vtkMedRegularGrid* grid)
148 for(int axis=0; axis < grid->GetParentMesh()->GetNumberOfAxis(); axis++)
151 vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
152 grid->SetAxisCoordinate(axis, coords);
155 coords->SetNumberOfComponents(1);
156 coords->SetNumberOfTuples(grid->GetAxisSize(axis));
158 if (MEDmeshGridIndexCoordinateRd(
160 grid->GetParentMesh()->GetName(),
161 grid->GetComputeStep().TimeIt,
162 grid->GetComputeStep().IterationIt,
164 (med_float*)coords->GetVoidPointer(0)) < 0)
166 vtkErrorMacro("ERROR : read axis " << axis << " coordinates ...");
167 grid->SetAxisCoordinate(axis, NULL);
175 // load all Information data associated with this standard grid.
176 void vtkMedDriver30::ReadCurvilinearGridInformation(vtkMedCurvilinearGrid* grid)
180 grid->SetDimension(grid->GetParentMesh()->GetNumberOfAxis());
183 med_bool coordinatechangement, geotransformation;
185 if ((size = MEDmeshnEntity(
187 grid->GetParentMesh()->GetName(),
188 grid->GetComputeStep().TimeIt,
189 grid->GetComputeStep().IterationIt,
194 &coordinatechangement,
195 &geotransformation)) < 0)
197 vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshnEntity");
200 grid->SetNumberOfPoints(size);
204 if(MEDmeshGridStructRd(
206 grid->GetParentMesh()->GetName(),
207 grid->GetComputeStep().TimeIt,
208 grid->GetComputeStep().IterationIt,
211 vtkErrorMacro("ReadCurvilinearGridInformation MEDmeshGridStructRd");
214 switch(grid->GetDimension())
219 grid->SetAxisSize(0, (axissize[0] <= 0 ? 1: axissize[0]));
222 grid->SetAxisSize(0, (axissize[0] <= 0 ? 1: axissize[0]));
223 grid->SetAxisSize(1, (axissize[1] <= 0 ? 1: axissize[1]));
226 grid->SetAxisSize(0, (axissize[0] <= 0 ? 1: axissize[0]));
227 grid->SetAxisSize(1, (axissize[1] <= 0 ? 1: axissize[1]));
228 grid->SetAxisSize(2, (axissize[2] <= 0 ? 1: axissize[2]));
231 vtkErrorMacro("Unsupported dimension for curvilinear grid : "
232 << grid->GetDimension());
236 // A test to verify the number of points : total number of
237 // points should be equal to the product of each axis size
240 if (grid->GetDimension() == 1)
242 size2 = grid->GetAxisSize(0);
245 if (grid->GetDimension() == 2)
247 size2 = grid->GetAxisSize(0)*grid->GetAxisSize(1);
250 if (grid->GetDimension() == 3)
252 size2 = grid->GetAxisSize(0)*grid->GetAxisSize(1)*grid->GetAxisSize(2);
257 vtkErrorMacro("The total number of points of a Curvilinear grid should "
258 << "be the product of each axis size!");
262 if ((grid->GetDimension() >= 1) && (grid->GetAxisSize(0) > 1))
263 ncell *= grid->GetAxisSize(0)-1;
264 if ((grid->GetDimension() >= 2) && (grid->GetAxisSize(1) > 1))
265 ncell *= grid->GetAxisSize(1)-1;
266 if ((grid->GetDimension() >= 3) && (grid->GetAxisSize(2) > 1))
267 ncell *= grid->GetAxisSize(2)-1;
270 entity.EntityType = MED_CELL;
272 switch(grid->GetDimension())
275 entity.GeometryType = MED_POINT1;
278 entity.GeometryType = MED_SEG2;
281 entity.GeometryType = MED_QUAD4;
284 entity.GeometryType = MED_HEXA8;
287 vtkErrorMacro("Unsupported dimension for curvilinear grid : "
288 << grid->GetDimension());
292 vtkMedEntityArray* array = vtkMedEntityArray::New();
293 array->SetParentGrid(grid);
294 array->SetNumberOfEntity(ncell);
295 array->SetEntity(entity);
296 array->SetConnectivity(MED_NODAL);
297 grid->AppendEntityArray(array);
299 // this triggers the creation of undefined families
300 this->LoadFamilyIds(array);
303 // Description : read the number of entity of all geometry type
304 // for a given entity type and a given connectivity mode
305 void vtkMedDriver30::ReadNumberOfEntity(
306 vtkMedUnstructuredGrid* grid,
307 med_entity_type entityType,
308 med_connectivity_mode connectivity)
312 med_bool changement, transformation;
314 const char* meshName = grid->GetParentMesh()->GetName();
316 const vtkMedComputeStep& cs = grid->GetComputeStep();
318 med_int nentity = MEDmeshnEntity(
330 for(med_int geotypeit = 1; geotypeit <= nentity; geotypeit++)
332 // read cell informations
334 entity.EntityType = entityType;
336 char geometryName[MED_NAME_SIZE+1] = "";
338 // this gives us the med_geometry_type
339 if( MEDmeshEntityInfo( FileId, meshName,
345 &entity.GeometryType) < 0)
347 vtkErrorMacro("MEDmeshEntityInfo");
351 entity.GeometryName = geometryName;
354 if(entity.GeometryType == MED_POLYGON)
356 // read the number of cells of this type
357 ncell = MEDmeshnEntity(this->FileId,
366 &transformation ) - 1;
368 else if(entity.GeometryType == MED_POLYHEDRON)
370 // read the number of cells of this type
371 ncell = MEDmeshnEntity(this->FileId,
380 &transformation ) - 1;
384 ncell = MEDmeshnEntity(this->FileId,
398 vtkMedEntityArray* array = vtkMedEntityArray::New();
399 array->SetParentGrid(grid);
400 array->SetNumberOfEntity(ncell);
401 array->SetEntity(entity);
402 array->SetConnectivity(connectivity);
403 grid->AppendEntityArray(array);
405 // this triggers the creation of undefined families
406 this->LoadFamilyIds(array);
412 // load all Information data associated with this unstructured grid.
413 void vtkMedDriver30::ReadUnstructuredGridInformation(
414 vtkMedUnstructuredGrid* grid)
418 vtkMedMesh *mesh = grid->GetParentMesh();
420 const char *meshName = mesh->GetName();
421 med_connectivity_mode connectivity;
424 med_bool transformation;
427 char profilename[MED_NAME_SIZE+1];
428 memset(profilename, '\0', MED_NAME_SIZE+1);
430 vtkMedComputeStep cs = grid->GetComputeStep();
432 // first check if we have points
433 vtkIdType npoints = MEDmeshnEntityWithProfile(
450 grid->SetNumberOfPoints(npoints);
454 if(grid->GetPreviousGrid() == NULL)
456 vtkErrorMacro("No point and no previous grid");
458 grid->SetUsePreviousCoordinates(true);
459 grid->SetNumberOfPoints(grid->GetPreviousGrid()->GetNumberOfPoints());
463 this->ReadNumberOfEntity(grid, MED_CELL, MED_NODAL);
464 this->ReadNumberOfEntity(grid, MED_CELL, MED_DESCENDING);
465 this->ReadNumberOfEntity(grid, MED_DESCENDING_FACE, MED_NODAL);
466 this->ReadNumberOfEntity(grid, MED_DESCENDING_FACE, MED_DESCENDING);
467 this->ReadNumberOfEntity(grid, MED_DESCENDING_EDGE, MED_NODAL);
468 this->ReadNumberOfEntity(grid, MED_DESCENDING_EDGE, MED_DESCENDING);
469 this->ReadNumberOfEntity(grid, MED_STRUCT_ELEMENT, MED_NODAL);
470 this->ReadNumberOfEntity(grid, MED_STRUCT_ELEMENT, MED_DESCENDING);
472 // create the point vtkMedEntityArray support
474 entity.EntityType = MED_NODE;
475 entity.GeometryType = MED_POINT1;
476 vtkMedEntityArray* pea = vtkMedEntityArray::New();
477 pea->SetEntity(entity);
478 pea->SetParentGrid(grid);
479 pea->SetNumberOfEntity(grid->GetNumberOfPoints());
480 grid->AppendEntityArray(pea);
483 this->LoadFamilyIds(pea);
486 void vtkMedDriver30::ReadFamilyInformation(vtkMedMesh* mesh, vtkMedFamily* family)
492 char* groupNames = NULL;
493 const char* meshName = mesh->GetName();
495 ngroup = MEDnFamilyGroup(FileId, meshName, family->GetMedIterator());
497 bool has_no_group = false;
502 vtkErrorMacro("Error while reading the number of groups");
508 groupNames = new char[ngroup * MED_LNAME_SIZE + 1];
509 memset(groupNames, '\0', ngroup * MED_LNAME_SIZE + 1);
511 // special case for files written by med < 3,
512 // I have to use the 23 interface
513 if(mesh->GetParentFile()->GetVersionMajor() < 3)
515 med_int *attributenumber = NULL;
516 med_int *attributevalue = NULL;
517 char *attributedes = NULL;
519 med_int nattr = MEDnFamily23Attribute(
522 family->GetMedIterator());
526 vtkErrorMacro("MEDnFamily23Attribute");
531 attributenumber = new med_int[nattr];
532 attributevalue = new med_int[nattr];
533 attributedes = new char[nattr*MED_COMMENT_SIZE+1];
534 memset(attributedes, '\0', nattr*MED_COMMENT_SIZE+1);
537 char familyName[MED_NAME_SIZE+1] = "";
539 if(MEDfamily23Info (this->FileId,
541 family->GetMedIterator(),
549 vtkDebugMacro("MEDfamily23Info");
552 family->SetName(familyName);
554 if(attributenumber != NULL)
555 delete[] attributenumber;
556 if(attributevalue != NULL)
557 delete[] attributevalue;
558 if(attributedes != NULL)
559 delete[] attributedes;
563 char familyName[MED_NAME_SIZE+1] = "";
564 if(MEDfamilyInfo( this->FileId,
566 family->GetMedIterator(),
572 "vtkMedDriver30::ReadInformation(vtkMedFamily* family)"
573 << " cannot read family informations.");
576 family->SetName(familyName);
579 family->SetId(familyid);
583 family->SetPointOrCell(vtkMedUtilities::OnCell);
584 mesh->AppendCellFamily(family);
588 family->SetPointOrCell(vtkMedUtilities::OnPoint);
589 mesh->AppendPointFamily(family);
592 family->AllocateNumberOfGroup(ngroup);
593 // if there where no group, set the name to the default value
596 memcpy(groupNames, vtkMedUtilities::NoGroupName,
597 strlen(vtkMedUtilities::NoGroupName));
600 for(int index = 0; index < ngroup; index++)
602 char realGroupName[MED_LNAME_SIZE + 1];
603 memset(realGroupName, '\0', MED_LNAME_SIZE + 1);
604 memcpy(realGroupName, groupNames + index * MED_LNAME_SIZE,
605 MED_LNAME_SIZE * sizeof(char));
606 vtkMedGroup* group = mesh->GetOrCreateGroup(family->GetPointOrCell(),
609 family->SetGroup(index, group);
616 vtkMedFamily* famzero = vtkMedFamily::New();
617 mesh->AppendPointFamily(famzero);
620 famzero->SetName(family->GetName());
621 famzero->SetMedIterator(family->GetMedIterator());
622 famzero->SetId(family->GetId());
623 famzero->SetPointOrCell(vtkMedUtilities::OnPoint);
624 famzero->AllocateNumberOfGroup(family->GetNumberOfGroup());
625 for(int gid=0; gid<family->GetNumberOfGroup(); gid++)
627 vtkMedGroup* group = mesh->GetOrCreateGroup(
628 vtkMedUtilities::OnPoint,
629 family->GetGroup(gid)->GetName());
630 famzero->SetGroup(gid, group);
631 mesh->AppendPointGroup(group);
636 void vtkMedDriver30::ReadLinkInformation(vtkMedLink* link)
639 char linkMeshName[MED_NAME_SIZE+1] = "";
640 if(MEDlinkInfo(this->FileId,
641 link->GetMedIterator(),
645 vtkErrorMacro("MEDlinkInfo");
648 link->SetMeshName(linkMeshName);
652 char* path = new char[size + 1];
653 memset(path, '\0', size+1);
654 if(MEDlinkRd(this->FileId, link->GetMeshName(), path) < 0)
656 vtkErrorMacro("MEDlinkRd");
657 memset(path, '\0', size+1);
665 void vtkMedDriver30::ReadFileInformation(vtkMedFile* file)
669 char comment[MED_COMMENT_SIZE+1] = "";
671 MEDfileCommentRd(this->FileId,
674 file->SetComment(comment);
676 med_int major, minor, release;
677 MEDfileNumVersionRd(this->FileId, &major, &minor, &release);
678 file->SetVersionMajor(major);
679 file->SetVersionMinor(minor);
680 file->SetVersionRelease(release);
682 int nlink = MEDnLink(this->FileId);
683 file->AllocateNumberOfLink(nlink);
684 for(int linkid=0; linkid<nlink; linkid++)
686 vtkMedLink* link = file->GetLink(linkid);
687 link->SetMedIterator(linkid+1);
688 this->ReadLinkInformation(link);
691 int nprof = MEDnProfile(FileId);
692 // Reading id s not possible in parallel if the file contains Profiles
693 vtkMultiProcessController* controller = vtkMultiProcessController::GetGlobalController();
694 if (controller != NULL)
695 if ((nprof != 0) && (controller->GetNumberOfProcesses() > 1))
697 vtkWarningMacro("ATTENTION: The MED Reader cannot read profiles when used in parallel");
700 file->AllocateNumberOfProfile(nprof);
701 for(int i = 0; i < nprof; i++)
703 vtkMedProfile* profile = file->GetProfile(i);
704 profile->SetMedIterator(i + 1);
705 profile->SetParentFile(file);
706 this->ReadProfileInformation(profile);
709 int nloc = MEDnLocalization(this->FileId);
710 file->AllocateNumberOfLocalization(nloc);
711 for(int i = 0; i < nloc; i++)
713 vtkMedLocalization* loc = file->GetLocalization(i);
714 loc->SetMedIterator(i + 1);
715 loc->SetParentFile(file);
716 this->ReadLocalizationInformation(loc);
719 int nsupportmesh = MEDnSupportMesh(this->FileId);
720 file->AllocateNumberOfSupportMesh(nsupportmesh);
721 for(int i = 0; i < nsupportmesh; i++)
723 vtkMedMesh* supportmesh = file->GetSupportMesh(i);
724 supportmesh->SetMedIterator(i + 1);
725 supportmesh->SetParentFile(file);
726 this->ReadSupportMeshInformation(supportmesh);
729 int nmesh = MEDnMesh(this->FileId);
730 file->AllocateNumberOfMesh(nmesh);
731 for(int i = 0; i < nmesh; i++)
733 vtkMedMesh* mesh = file->GetMesh(i);
734 mesh->SetMedIterator(i + 1);
735 mesh->SetParentFile(file);
736 this->ReadMeshInformation(mesh);
739 int nfields = MEDnField(this->FileId);
740 file->AllocateNumberOfField(nfields);
741 for(int i = 0; i < nfields; i++)
743 vtkMedField* field = file->GetField(i);
744 field->SetMedIterator(i + 1);
745 field->SetParentFile(file);
746 this->ReadFieldInformation(field);
747 field->ComputeFieldType();
748 while(field->HasManyFieldTypes())
750 vtkMedField* newfield = vtkMedField::New();
751 int type = field->GetFirstType();
752 newfield->ExtractFieldType(field, type);
753 file->AppendField(newfield);
758 int nstruct = MEDnStructElement(this->FileId);
760 file->AllocateNumberOfStructElement(nstruct);
761 for(int i = 0; i < nstruct; i++)
763 vtkMedStructElement* structelem = file->GetStructElement(i);
764 structelem->SetMedIterator(i+1);
765 structelem->SetParentFile(file);
766 this->ReadStructElementInformation(structelem);
771 void vtkMedDriver30::ReadStructElementInformation(
772 vtkMedStructElement* structelem)
777 char modelname[MED_NAME_SIZE+1] = "";
778 med_geometry_type mgeotype;
780 char supportmeshname[MED_NAME_SIZE+1] = "";
781 med_entity_type sentitytype;
784 med_geometry_type sgeotype;
785 med_int nbofconstantattribute;
787 med_int nbofvariableattribute;
789 if(MEDstructElementInfo (this->FileId,
790 structelem->GetMedIterator(),
799 &nbofconstantattribute,
801 &nbofvariableattribute) < 0)
803 vtkErrorMacro("Error in MEDstructElementInfo");
806 structelem->SetName(modelname);
807 structelem->SetGeometryType(mgeotype);
808 structelem->SetModelDimension(modeldim);
809 structelem->SetSupportMeshName(supportmeshname);
810 structelem->SetSupportEntityType(sentitytype);
811 structelem->SetSupportNumberOfNode(snbofnode);
812 structelem->SetSupportNumberOfCell(snbofcell);
813 structelem->SetSupportGeometryType(sgeotype);
814 structelem->AllocateNumberOfConstantAttribute(nbofconstantattribute);
815 structelem->AllocateNumberOfVariableAttribute(nbofvariableattribute);
816 structelem->SetAnyProfile(anyprofile);
818 for(int attit = 0; attit < nbofconstantattribute; attit ++)
820 vtkMedConstantAttribute* constatt = structelem->GetConstantAttribute(attit);
821 constatt->SetMedIterator(attit+1);
822 constatt->SetParentStructElement(structelem);
823 this->ReadConstantAttributeInformation(constatt);
826 for(int attit = 0; attit < nbofvariableattribute; attit ++)
828 vtkMedVariableAttribute* varatt = structelem->GetVariableAttribute(attit);
829 varatt->SetMedIterator(attit+1);
830 varatt->SetParentStructElement(structelem);
831 this->ReadVariableAttributeInformation(varatt);
835 void vtkMedDriver30::ReadConstantAttributeInformation(
836 vtkMedConstantAttribute* constAttr)
841 char constattname[MED_NAME_SIZE+1] = "";
842 med_attribute_type constatttype;
843 med_int nbofcomponent;
844 med_entity_type sentitytype;
845 char profilename[MED_NAME_SIZE+1] = "";
848 if(MEDstructElementConstAttInfo(
850 constAttr->GetParentStructElement()->GetName(),
851 constAttr->GetMedIterator(),
859 vtkErrorMacro("MEDstructElementConstAttInfo error");
863 constAttr->SetName(constattname);
864 constAttr->SetAttributeType(constatttype);
865 constAttr->SetNumberOfComponent(nbofcomponent);
866 constAttr->SetSupportEntityType(sentitytype);
867 constAttr->SetProfileName(profilename);
868 constAttr->SetProfileSize(profilesize);
870 vtkAbstractArray* values = vtkMedUtilities::NewArray(constatttype);
873 constAttr->SetValues(values);
876 values->SetNumberOfComponents(nbofcomponent);
877 vtkIdType ntuple = 0;
878 if((strcmp(profilename, MED_NO_PROFILE) != 0) &&
879 (strcmp(profilename, "\0") != 0))
881 ntuple = profilesize;
883 else if(constAttr->GetSupportEntityType() == MED_CELL)
885 ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfCell();
889 ntuple = constAttr->GetParentStructElement()->GetSupportNumberOfNode();
891 values->SetNumberOfTuples(ntuple);
894 vtkSmartPointer<vtkCharArray> buffer = vtkSmartPointer<vtkCharArray>::New();
895 if(constatttype != MED_ATT_NAME)
897 ptr = values->GetVoidPointer(0);
901 buffer->SetNumberOfValues(MED_NAME_SIZE*nbofcomponent*ntuple);
902 ptr = buffer->GetVoidPointer(0);
905 if(MEDstructElementConstAttRd (this->FileId,
906 constAttr->GetParentStructElement()->GetName(),
907 constAttr->GetName(), ptr) < 0)
909 vtkErrorMacro("MEDstructElementConstAttRd");
913 if(constatttype == MED_ATT_NAME)
915 char name[MED_NAME_SIZE+1] = "";
916 char* nameptr = (char*) ptr;
917 vtkStringArray* names = vtkStringArray::SafeDownCast(values);
918 for(vtkIdType id = 0; id < nbofcomponent*ntuple; id++)
920 memset(name, '\0', MED_NAME_SIZE+1);
921 strncpy(name, nameptr + id * MED_NAME_SIZE, MED_NAME_SIZE);
922 names->SetValue(id, name);
929 void vtkMedDriver30::ReadVariableAttributeInformation(
930 vtkMedVariableAttribute* varAttr)
935 char varattname[MED_NAME_SIZE+1] = "";
936 med_attribute_type varatttype;
937 med_int nbofcomponent;
939 if(MEDstructElementVarAttInfo (
941 varAttr->GetParentStructElement()->GetName(),
942 varAttr->GetMedIterator(),
947 vtkErrorMacro("MEDstructElementVarAttInfo");
951 varAttr->SetName(varattname);
952 varAttr->SetAttributeType(varatttype);
953 varAttr->SetNumberOfComponent(nbofcomponent);
958 void vtkMedDriver30::ReadProfileInformation(vtkMedProfile* profile)
963 char profileName[MED_NAME_SIZE+1] = "";
965 if(MEDprofileInfo(this->FileId,
966 profile->GetMedIterator(),
970 vtkErrorMacro("cannot read information on profile"
971 << profile->GetMedIterator());
973 profile->SetName(profileName);
974 profile->SetNumberOfElement(nelem);
977 void vtkMedDriver30::ReadFieldInformation(vtkMedField* field)
981 if (field->GetMedIterator() == 0)
984 int ncomp = MEDfieldnComponent(FileId, field->GetMedIterator());
988 field->SetNumberOfComponent(-1);
989 vtkErrorMacro("cannot read the number of component of field "
990 << field->GetMedIterator())
994 field->SetNumberOfComponent(ncomp);
996 char* units = new char[MED_SNAME_SIZE * ncomp + 1];
997 char* componentNames = new char[MED_SNAME_SIZE * ncomp + 1];
998 memset(units, '\0', MED_SNAME_SIZE * ncomp + 1);
999 memset(componentNames, '\0', MED_SNAME_SIZE * ncomp + 1);
1001 //med_type_champ dataType;
1002 med_field_type dataType;
1006 char name[MED_NAME_SIZE+1] = "";
1007 char meshName[MED_NAME_SIZE+1] = "";
1008 char unit[MED_SNAME_SIZE+1] = "";
1010 if( MEDfieldInfo( FileId,
1011 field->GetMedIterator(),
1021 vtkErrorMacro("cannot read the informations on field "
1022 << field->GetMedIterator())
1026 field->SetName(name);
1027 field->SetMeshName(meshName);
1028 field->SetTimeUnit(unit);
1029 field->SetDataType(dataType);
1030 field->SetLocal(localmesh);
1032 for(int comp = 0; comp < ncomp; comp++)
1034 char unit[MED_NAME_SIZE + 1] = "";
1035 memcpy(unit, units + MED_SNAME_SIZE * comp, MED_SNAME_SIZE * sizeof(char));
1036 field->GetUnit()->SetValue(comp, unit);
1038 char compName[MED_SNAME_SIZE + 1] = "";
1039 memcpy(compName, componentNames + MED_SNAME_SIZE * comp, MED_SNAME_SIZE
1041 field->GetComponentName()->SetValue(comp, compName);
1045 delete[] componentNames;
1047 med_int ninterp = MEDfieldnInterp(FileId, field->GetName());
1050 vtkErrorMacro("Error in MEDfieldnInterp");
1054 field->AllocateNumberOfInterpolation(ninterp);
1056 for(med_int interpit=0; interpit<ninterp; interpit++)
1058 vtkMedInterpolation* interp = field->GetInterpolation(interpit);
1059 interp->SetMedIterator(interpit + 1);
1060 this->ReadInterpolationInformation(interp);
1063 vtkMedFieldStep* previousStep = NULL;
1065 for(med_int csit = 0; csit < nstep; csit++)
1067 vtkMedFieldStep* step = vtkMedFieldStep::New();
1068 step->SetMedIterator(csit + 1);
1069 step->SetParentField(field);
1070 this->ReadFieldStepInformation(step, csit == 0);
1071 field->AddFieldStep(step);
1072 step->SetPreviousStep(previousStep);
1073 previousStep = step;
1078 void vtkMedDriver30::ReadFieldStepInformation(vtkMedFieldStep* step, bool readAllEntityInfo)
1080 vtkMedComputeStep cs;
1081 vtkMedComputeStep meshcs;
1082 vtkMedField* field = step->GetParentField();
1084 FileOpen open(this);
1086 if( MEDfieldComputingStepMeshInfo(
1089 step->GetMedIterator(),
1092 &cs.TimeOrFrequency,
1094 &meshcs.IterationIt) < 0)
1096 vtkErrorMacro("Error in MEDfieldComputingStepMeshInfo");
1100 step->SetComputeStep(cs);
1101 step->SetMeshComputeStep(meshcs);
1103 if(!readAllEntityInfo || step->GetEntityInfoLoaded())
1106 step->SetEntityInfoLoaded(1);
1108 vtkMedFile* file = field->GetParentFile();
1109 vtkMedMesh* mesh = file->GetMesh(field->GetMeshName());
1114 std::set<vtkMedEntity> entities;
1115 mesh->GatherMedEntities(entities);
1117 vtkMedEntity pointEntity;
1118 pointEntity.EntityType = MED_NODE;
1119 pointEntity.GeometryType = MED_NONE;
1120 pointEntity.GeometryName = MED_NAME_BLANK;
1121 entities.insert(pointEntity);
1123 std::set<vtkMedEntity>::iterator entity_it = entities.begin();
1124 while(entity_it != entities.end())
1126 vtkMedEntity entity = *entity_it;
1129 med_int nvalues = 0;
1131 char profilename[MED_NAME_SIZE+1] = "";
1132 char localizationname[MED_NAME_SIZE+1] = "";
1134 nprofile = MEDfieldnProfile(
1137 step->GetComputeStep().TimeIt,
1138 step->GetComputeStep().IterationIt,
1140 entity.GeometryType,
1145 vtkErrorMacro("MEDfieldnProfile");
1149 bool hasprofile = (nprofile > 0);
1155 med_int profilesize;
1156 med_int nintegrationpoint;
1158 for(int pid=0; pid<nprofile; pid++)
1160 med_int medid = (hasprofile ? pid+1 : -1);
1161 nvalues = MEDfieldnValueWithProfile(
1164 step->GetComputeStep().TimeIt,
1165 step->GetComputeStep().IterationIt,
1167 entity.GeometryType,
1169 MED_COMPACT_PFLMODE,
1173 &nintegrationpoint);
1177 vtkErrorMacro("MEDfieldnValueWithProfile");
1180 else if(nvalues > 0)
1182 // I have found a profile with values, stop the loop here
1189 vtkMedFieldOverEntity* fieldOverEntity = vtkMedFieldOverEntity::New();
1190 step->AppendFieldOverEntity(fieldOverEntity);
1191 fieldOverEntity->Delete();
1193 fieldOverEntity->SetParentStep(step);
1194 fieldOverEntity->SetEntity(entity);
1196 this->ReadFieldOverEntityInformation(fieldOverEntity);
1201 void vtkMedDriver30::ReadFieldOverEntityInformation(vtkMedFieldOverEntity* fieldOverEntity)
1203 FileOpen open(this);
1205 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
1206 vtkMedField* field = step->GetParentField();
1207 vtkMedEntity entity = fieldOverEntity->GetEntity();
1209 const char* fieldName = field->GetName();
1210 const vtkMedComputeStep& cs = step->GetComputeStep();
1212 char profilename[MED_NAME_SIZE+1] = "";
1213 char localizationname[MED_NAME_SIZE+1] = "";
1215 med_int nProfiles = MEDfieldnProfile(
1221 entity.GeometryType,
1227 vtkErrorMacro("MEDfieldnProfile");
1229 else if(nProfiles == 0)
1231 fieldOverEntity->SetHasProfile(0);
1236 fieldOverEntity->SetHasProfile(1);
1238 fieldOverEntity->AllocateNumberOfFieldOnProfile(nProfiles);
1239 for(int profit = 0; profit < nProfiles; profit++)
1241 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(profit);
1242 med_int medid = (fieldOverEntity->GetHasProfile()? profit+1: -1);
1243 fop->SetMedIterator(medid);
1244 fop->SetParentFieldOverEntity(fieldOverEntity);
1245 this->ReadFieldOnProfileInformation(fop);
1249 void vtkMedDriver30::ReadFieldOnProfileInformation(vtkMedFieldOnProfile* fop)
1251 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
1252 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
1253 vtkMedField* field = step->GetParentField();
1255 const vtkMedComputeStep& cs = step->GetComputeStep();
1256 med_int profilesize;
1257 med_int nbofintegrationpoint;
1259 char profileName[MED_NAME_SIZE+1] = "";
1260 char localizationName[MED_NAME_SIZE+1] = "";
1262 med_int nvalue = MEDfieldnValueWithProfile(FileId,
1266 fieldOverEntity->GetEntity().EntityType,
1267 fieldOverEntity->GetEntity().GeometryType,
1268 fop->GetMedIterator(),
1273 &nbofintegrationpoint);
1277 vtkErrorMacro("Error while reading MEDfieldnValueWithProfile");
1280 fop->SetProfileName(profileName);
1281 fop->SetLocalizationName(localizationName);
1282 fop->SetNumberOfValues(nvalue);
1283 fop->SetNumberOfIntegrationPoint(nbofintegrationpoint);
1284 fop->SetProfileSize(profilesize);
1287 void vtkMedDriver30::ReadMeshInformation(vtkMedMesh* mesh)
1289 FileOpen open(this);
1293 med_mesh_type meshtype;
1295 med_sorting_type sortingtype;
1297 med_axis_type axistype;
1300 if ( (naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator())) <0 )
1302 vtkDebugMacro("Error reading mesh axis number");
1306 naxis=MEDmeshnAxis(this->FileId, mesh->GetMedIterator());
1308 char axisname[3*MED_SNAME_SIZE+1] = "";
1309 char axisunit[3*MED_SNAME_SIZE+1] = "";
1310 char name[MED_NAME_SIZE+1] = "";
1311 char description[MED_COMMENT_SIZE+1] = "";
1312 char timeUnit[MED_SNAME_SIZE+1] = "";
1314 if( MEDmeshInfo( this->FileId,
1315 mesh->GetMedIterator(),
1328 vtkErrorMacro("Error reading mesh");
1331 mesh->SetName(name);
1332 mesh->SetDescription(description);
1333 mesh->SetTimeUnit(timeUnit);
1334 mesh->SetSpaceDimension(sdim);
1335 mesh->SetMeshDimension(mdim);
1336 mesh->SetMeshType(meshtype);
1337 mesh->SetSortingType(sortingtype);
1338 mesh->SetAxisType(axistype);
1339 mesh->SetNumberOfAxis(naxis);
1341 for(int axis = 0; axis < naxis; axis++)
1343 char theaxisname[MED_SNAME_SIZE+1] = "";
1344 char theaxisunit[MED_SNAME_SIZE+1] = "";
1345 strncpy(theaxisname, axisname + axis*MED_SNAME_SIZE, MED_SNAME_SIZE);
1346 strncpy(theaxisunit, axisunit + axis*MED_SNAME_SIZE, MED_SNAME_SIZE);
1347 mesh->GetAxisName()->SetValue(axis, theaxisname);
1348 mesh->GetAxisUnit()->SetValue(axis, theaxisunit);
1351 char universalName[MED_LNAME_SIZE+1] = "";
1353 if(MEDmeshUniversalNameRd(this->FileId, name,
1356 vtkDebugMacro("MEDmeshUniversalNameRd < 0");
1358 mesh->SetUniversalName(universalName);
1360 // read the Information data of all families.
1361 // writing the family 0 is optional,
1362 // but I need it, so add it if it is not here.
1364 med_int nfam = MEDnFamily(this->FileId, name);
1366 for(int index = 0; index < nfam; index++)
1368 vtkMedFamily* family = vtkMedFamily::New();
1369 family->SetMedIterator(index + 1);
1370 this->ReadFamilyInformation(mesh, family);
1374 // this creates a family 0 if none has been read
1375 vtkMedFamily* familyZeroOnCell = mesh->GetOrCreateCellFamilyById(0);
1376 vtkMedFamily* familyZeroOnPoint = mesh->GetOrCreatePointFamilyById(0);
1378 // Load Information regarding the grid type
1379 if(meshtype == MED_STRUCTURED_MESH)
1381 // Do it for structured data
1383 if(MEDmeshGridTypeRd(FileId, name, &mtg) < 0)
1385 vtkErrorMacro("Error during structured grid Information loading.");
1388 mesh->SetStructuredGridType(mtg);
1391 vtkMedGrid* previousGrid = NULL;
1392 for(int gid=1; gid <= nstep; gid++)
1394 vtkMedComputeStep cs;
1395 if(MEDmeshComputationStepInfo(FileId,
1400 &cs.TimeOrFrequency) < 0)
1402 vtkErrorMacro("MEDmeshComputationStepInfo error");
1404 // Load Information regarding the grid type
1405 vtkMedGrid* grid = NULL;
1406 if(meshtype == MED_STRUCTURED_MESH)
1408 switch(mesh->GetStructuredGridType())
1410 case MED_CARTESIAN_GRID:
1411 grid = vtkMedCartesianGrid::New();
1413 case MED_POLAR_GRID:
1414 grid = vtkMedPolarGrid::New();
1416 case MED_CURVILINEAR_GRID:
1417 grid = vtkMedCurvilinearGrid::New();
1420 vtkErrorMacro("Unknown structured grid type " << mesh->GetStructuredGridType());
1424 else //(mesh->GetType() == MED_STRUCTURED_MESH)
1426 grid = vtkMedUnstructuredGrid::New();
1428 grid->SetParentMesh(mesh);
1429 grid->SetComputeStep(cs);
1430 this->ReadGridInformation(grid);
1431 mesh->AddGridStep(grid);
1433 grid->SetPreviousGrid(previousGrid);
1434 previousGrid = grid;
1438 void vtkMedDriver30::ReadLocalizationInformation(vtkMedLocalization* loc)
1440 FileOpen open(this);
1443 med_int spaceDimension;
1444 med_geometry_type type_geo;
1445 med_geometry_type sectiongeotype;
1446 med_int nsectionmeshcell;
1448 char name[MED_NAME_SIZE+1] = "";
1449 char interpolationName[MED_NAME_SIZE+1] = "";
1450 char sectionName[MED_NAME_SIZE+1] = "";
1452 if(MEDlocalizationInfo(
1454 loc->GetMedIterator(),
1462 §iongeotype ) < 0)
1464 vtkErrorMacro("Reading information on quadrature points definition : "
1465 << loc->GetMedIterator());
1469 loc->SetInterpolationName(interpolationName);
1470 loc->SetSectionName(sectionName);
1471 loc->SetNumberOfQuadraturePoint(ngp);
1472 loc->SetGeometryType(type_geo);
1473 loc->SetSpaceDimension(spaceDimension);
1474 loc->SetNumberOfCellInSection(nsectionmeshcell);
1475 loc->SetSectionGeometryType(sectiongeotype);
1477 med_float *localCoordinates = new med_float[loc->GetSizeOfPointLocalCoordinates()];
1478 med_float *pqLocalCoordinates = new med_float[loc->GetSizeOfQuadraturePointLocalCoordinates()];
1479 med_float *weights = new med_float[loc->GetSizeOfWeights()];
1481 if(MEDlocalizationRd(FileId,
1488 vtkErrorMacro("MEDlocalizationRd : " << loc->GetName());
1491 vtkDoubleArray* lc = loc->GetPointLocalCoordinates();
1492 vtkDoubleArray *pqlc = loc->GetQuadraturePointLocalCoordinates();
1493 vtkDoubleArray *w = loc->GetWeights();
1495 lc->SetNumberOfValues(loc->GetSizeOfPointLocalCoordinates());
1496 for(int i=0; i<loc->GetSizeOfPointLocalCoordinates(); i++)
1498 lc->SetValue(i, localCoordinates[i]);
1501 pqlc->SetNumberOfValues(loc->GetSizeOfQuadraturePointLocalCoordinates());
1502 for(int i=0; i<loc->GetSizeOfQuadraturePointLocalCoordinates(); i++)
1504 pqlc->SetValue(i, pqLocalCoordinates[i]);
1507 w->SetNumberOfValues(loc->GetSizeOfWeights());
1508 for(int i=0; i<loc->GetSizeOfWeights(); i++)
1510 w->SetValue(i, weights[i]);
1514 void vtkMedDriver30::ReadInterpolationInformation(vtkMedInterpolation* interp)
1517 med_geometry_type geotype;
1519 med_int nbofbasisfunc;
1520 med_int nbofvariable;
1524 char name[MED_NAME_SIZE+1] = "";
1526 if(MEDinterpInfo (this->FileId,
1527 interp->GetMedIterator(),
1529 &geotype, &cellnode, &nbofbasisfunc,
1530 &nbofvariable, &maxdegree, &nmaxcoef) < 0)
1532 vtkErrorMacro("MEDinterpInfo");
1536 interp->SetName(name);
1537 interp->SetGeometryType(geotype);
1538 interp->SetIsCellNode(cellnode);
1539 interp->SetNumberOfVariable(nbofvariable);
1540 interp->SetMaximumDegree(maxdegree);
1541 interp->SetMaximumNumberOfCoefficient(nmaxcoef);
1542 interp->AllocateNumberOfBasisFunction(nbofbasisfunc);
1544 for(int basisid=0; basisid < nbofbasisfunc; basisid++)
1546 vtkMedFraction* func = interp->GetBasisFunction(basisid);
1547 func->SetNumberOfVariable(nbofvariable);
1549 med_int ncoef = MEDinterpBaseFunctionCoefSize (
1553 func->SetNumberOfCoefficients(ncoef);
1555 if(ncoef <= 0 || nbofvariable <= 0)
1558 med_int *power = new med_int[nbofvariable * ncoef];
1559 med_float *coefficient = new med_float[ncoef];
1561 if(MEDinterpBaseFunctionRd (
1569 vtkErrorMacro("MEDinterpBaseFunctionRd");
1572 vtkDoubleArray* coeffs = func->GetCoefficients();
1573 for(int cid=0; cid < ncoef; cid++)
1575 coeffs->SetValue(cid, coefficient[cid]);
1577 vtkIntArray* powers = func->GetPowers();
1578 for(int pid=0; pid < ncoef*nbofvariable; pid++)
1580 powers->SetValue(pid, power[pid]);
1584 delete[] coefficient;
1588 void vtkMedDriver30::ReadSupportMeshInformation(
1589 vtkMedMesh* supportMesh)
1591 FileOpen open(this);
1593 char supportmeshname[MED_NAME_SIZE+1] = "";
1594 char description[MED_COMMENT_SIZE+1] = "";
1597 med_axis_type axistype;
1598 char axisname[3*MED_SNAME_SIZE+1] = "";
1599 char axisunit[3*MED_SNAME_SIZE+1] = "";
1601 if(MEDsupportMeshInfo (this->FileId,
1602 supportMesh->GetMedIterator(),
1611 vtkErrorMacro("MEDsupportMeshInfo");
1614 supportMesh->SetName(supportmeshname);
1615 supportMesh->SetDescription(description);
1616 supportMesh->SetSpaceDimension(spacedim);
1617 supportMesh->SetMeshDimension(meshdim);
1618 supportMesh->SetAxisType(axistype);
1619 for(int dim = 0; dim < 3; dim++)
1621 char axisname_dim[MED_SNAME_SIZE+1] = "";
1622 char axisunit_dim[MED_SNAME_SIZE+1] = "";
1624 strncpy(axisname_dim, axisname+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
1625 strncpy(axisunit_dim, axisunit+dim*MED_SNAME_SIZE, MED_SNAME_SIZE);
1627 supportMesh->GetAxisName()->SetValue(dim, axisname_dim);
1628 supportMesh->GetAxisUnit()->SetValue(dim, axisunit_dim);
1634 void vtkMedDriver30::LoadFamilyIds(vtkMedEntityArray* array)
1636 if(array->IsFamilyIdsLoaded())
1639 FileOpen open(this);
1641 vtkMedGrid* grid = array->GetParentGrid();
1643 vtkMedComputeStep cs = grid->GetComputeStep();
1645 // first, find if the family ids are implicit or explicit
1646 med_bool changement, transformation;
1648 med_int nfamid = MEDmeshnEntity(this->FileId,
1649 grid->GetParentMesh()->GetName(),
1652 array->GetEntity().EntityType,
1653 array->GetEntity().GeometryType,
1659 if(nfamid == array->GetNumberOfEntity())
1662 vtkMedIntArray* famIds = vtkMedIntArray::New();
1663 array->SetFamilyIds(famIds);
1666 famIds->SetNumberOfTuples(nfamid);
1668 if ( MEDmeshEntityFamilyNumberRd(
1670 grid->GetParentMesh()->GetName(),
1673 array->GetEntity().EntityType,
1674 array->GetEntity().GeometryType,
1675 famIds->GetPointer(0) ) < 0)
1677 vtkWarningMacro("Error loading the family ids of entity "
1678 << array->GetEntity().EntityType
1679 << " " << array->GetEntity().GeometryType
1680 << " on mesh " << grid->GetParentMesh()->GetName());
1681 array->SetFamilyIds(NULL);
1686 vtkDebugMacro("NumberOfEntity != Number of family ids");
1687 array->SetFamilyIds(NULL);
1690 array->ComputeFamilies();
1693 void vtkMedDriver30::LoadPointGlobalIds(vtkMedGrid* grid)
1695 if(grid->IsPointGlobalIdsLoaded())
1698 FileOpen open(this);
1700 vtkMedIntArray* globalIds = vtkMedIntArray::New();
1701 grid->SetPointGlobalIds(globalIds);
1702 globalIds->Delete();
1704 globalIds->SetNumberOfTuples(grid->GetNumberOfPoints());
1706 if( MEDmeshEntityNumberRd (
1708 grid->GetParentMesh()->GetName(),
1709 grid->GetComputeStep().TimeIt,
1710 grid->GetComputeStep().IterationIt,
1713 globalIds->GetPointer(0) ) < 0)
1715 grid->SetPointGlobalIds(NULL);
1719 void vtkMedDriver30::LoadCoordinates(vtkMedGrid* grid)
1721 if(grid->IsCoordinatesLoaded())
1724 vtkMedRegularGrid* rgrid = vtkMedRegularGrid::SafeDownCast(grid);
1727 this->LoadRegularGridCoordinates(rgrid);
1731 vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
1732 vtkMedCurvilinearGrid* cgrid = vtkMedCurvilinearGrid::SafeDownCast(grid);
1733 if(ugrid == NULL && cgrid == NULL)
1735 //TODO : deal with structured grids
1736 vtkWarningMacro("this kind of grid is not yet supported");
1740 if(grid->GetUsePreviousCoordinates())
1742 vtkMedGrid* previousgrid = grid->GetPreviousGrid();
1743 if(previousgrid == NULL)
1745 vtkErrorMacro("coordiantes have not changed, "
1746 << "but there is no previous grid!");
1750 this->LoadCoordinates(previousgrid);
1752 ugrid->SetCoordinates(vtkMedUnstructuredGrid::SafeDownCast(previousgrid)
1753 ->GetCoordinates());
1755 cgrid->SetCoordinates(vtkMedCurvilinearGrid::SafeDownCast(previousgrid)
1756 ->GetCoordinates());
1761 FileOpen open(this);
1763 vtkDataArray* coords = vtkMedUtilities::NewCoordArray();
1765 ugrid->SetCoordinates(coords);
1767 cgrid->SetCoordinates(coords);
1770 vtkMedComputeStep cs = grid->GetComputeStep();
1772 coords->SetNumberOfComponents(grid->GetParentMesh()->GetSpaceDimension());
1773 coords->SetNumberOfTuples(grid->GetNumberOfPoints());
1775 if ( MEDmeshNodeCoordinateRd( this->FileId,
1776 grid->GetParentMesh()->GetName(),
1780 (med_float*) coords->GetVoidPointer(0) ) < 0)
1782 vtkErrorMacro("Load Coordinates for mesh "
1783 << grid->GetParentMesh()->GetName());
1788 void vtkMedDriver30::LoadProfile(vtkMedProfile* profile)
1790 if(!profile || profile->IsLoaded())
1793 FileOpen open(this);
1795 vtkMedIntArray* indices = vtkMedIntArray::New();
1796 profile->SetIds(indices);
1799 indices->SetNumberOfTuples(profile->GetNumberOfElement());
1801 char name[MED_NAME_SIZE+1] = "";
1803 if( MEDprofileRd(this->FileId,
1805 indices->GetPointer(0) ) < 0)
1807 vtkErrorMacro("Reading profile indices ");
1811 void vtkMedDriver30::LoadConnectivity(vtkMedEntityArray* array)
1813 if(array->IsConnectivityLoaded())
1816 FileOpen open(this);
1818 vtkMedGrid* grid = array->GetParentGrid();
1820 grid = array->GetParentGrid();
1822 const char* meshName = grid->GetParentMesh()->GetName();
1824 vtkMedIntArray* conn = vtkMedIntArray::New();
1825 array->SetConnectivityArray(conn);
1828 vtkMedComputeStep cs = grid->GetComputeStep();
1831 med_bool transformation;
1833 if(array->GetEntity().GeometryType == MED_POLYGON)
1835 // first check if we have points
1836 med_int connSize = MEDmeshnEntity(
1841 array->GetEntity().EntityType,
1844 array->GetConnectivity(),
1850 vtkErrorMacro(<< "Error while reading polygons connectivity size."
1855 conn->SetNumberOfTuples(connSize);
1857 // How many polygon in the mesh in nodal connectivity mode
1858 // For the polygons, we get the size of array index
1860 if ((indexsize = MEDmeshnEntity(this->FileId,
1864 array->GetEntity().EntityType,
1867 array->GetConnectivity(),
1869 &transformation )) < 0)
1871 vtkErrorMacro(<< "Error while reading polygons array index." << endl );
1875 vtkMedIntArray* index = vtkMedIntArray::New();
1876 array->SetFaceIndex(index);
1879 index->SetNumberOfTuples( indexsize );
1881 if ( MEDmeshPolygonRd(this->FileId,
1885 array->GetEntity().EntityType,
1886 array->GetConnectivity(),
1887 index->GetPointer(0),
1888 conn->GetPointer(0) ) < 0)
1890 vtkErrorMacro(<< "MEDmeshPolygonRd");
1894 else if(array->GetEntity().GeometryType == MED_POLYHEDRON)
1897 vtkIdType connSize = MEDmeshnEntity(this->FileId,
1899 grid->GetComputeStep().TimeIt,
1900 grid->GetComputeStep().IterationIt,
1901 array->GetEntity().EntityType,
1904 array->GetConnectivity(),
1909 vtkErrorMacro(<< "Error while reading polyhedrons connectivity size."
1914 conn->SetNumberOfTuples(connSize);
1916 vtkMedIntArray* faceIndex = vtkMedIntArray::New();
1917 array->SetFaceIndex(faceIndex);
1918 faceIndex->Delete();
1920 vtkMedIntArray* nodeIndex = vtkMedIntArray::New();
1921 array->SetNodeIndex(nodeIndex);
1922 nodeIndex->Delete();
1924 vtkIdType np = array->GetNumberOfEntity() + 1;
1925 faceIndex->SetNumberOfTuples(np);
1927 med_int nodeIndexSize;
1929 if ((nodeIndexSize = MEDmeshnEntity(this->FileId,
1931 grid->GetComputeStep().TimeIt,
1932 grid->GetComputeStep().IterationIt,
1933 array->GetEntity().EntityType,
1936 array->GetConnectivity(),
1938 &transformation )) < 0)
1940 vtkErrorMacro(<< "Error while reading polygons array index." << endl );
1944 nodeIndex->SetNumberOfTuples(nodeIndexSize);
1946 if (MEDmeshPolyhedronRd(this->FileId,
1950 array->GetEntity().EntityType,
1951 array->GetConnectivity(),
1952 faceIndex->GetPointer(0),
1953 nodeIndex->GetPointer(0),
1954 conn->GetPointer(0) ) < 0)
1956 vtkErrorMacro("Error while reading connectivity of polyhedrons");
1963 bool doReadConnectivity = true;
1964 if(array->GetConnectivity() == MED_NODAL)
1966 if(array->GetEntity().EntityType == MED_STRUCT_ELEMENT)
1968 std::cout << " -- MED_STRUCT_ELEMENT --" << std::endl;
1969 if(array->GetStructElement() == NULL)
1971 vtkErrorMacro("Entity type = MED_STRUCT_ELEMENT, but StructElement is not set!");
1974 vtkIdType ntuple = array->GetNumberOfEntity()
1975 * array->GetStructElement()->GetConnectivitySize();
1977 conn->SetNumberOfTuples(ntuple);
1978 // particles are special : connectivity is not stored in the med file
1979 if(strcmp(array->GetStructElement()->GetName(), MED_PARTICLE_NAME) == 0 )
1981 for(vtkIdType cellId = 0; cellId < ntuple; cellId++)
1983 conn->SetValue(cellId, cellId+1);
1985 doReadConnectivity = false;
1990 conn->SetNumberOfTuples(array->GetNumberOfEntity()
1991 * vtkMedUtilities::GetNumberOfPoint(
1992 array->GetEntity().GeometryType));
1997 conn->SetNumberOfTuples(array->GetNumberOfEntity()
1998 * vtkMedUtilities::GetNumberOfSubEntity(
1999 array->GetEntity().GeometryType));
2002 if (this->ParallelFileId == -1) // also (array->GetFilter() == NULL)
2004 if ( (MEDmeshElementConnectivityRd(
2009 array->GetEntity().EntityType,
2010 array->GetEntity().GeometryType,
2011 array->GetConnectivity(),
2013 conn->GetPointer(0)) ) < 0)
2015 vtkErrorMacro("Error while load connectivity of cells "
2016 << array->GetEntity().GeometryType);
2021 med_filter filter = MED_FILTER_INIT;
2028 array->GetFilter()->GetFilterSizes(start, stride, count,
2029 blocksize, lastblocksize );
2031 med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
2032 array->GetEntity().GeometryType);
2034 if ( MEDfilterBlockOfEntityCr( this->ParallelFileId,
2035 array->GetNumberOfEntity(),
2036 1, // one is for mesh elements, more than 1 is for fields
2037 nbofconstituentpervalue,
2038 MED_ALL_CONSTITUENT,
2045 (med_size)blocksize,
2046 (med_size)lastblocksize,
2049 vtkErrorMacro("Filter creation ");
2052 if ( (MEDmeshElementConnectivityAdvancedRd(
2053 this->ParallelFileId,
2057 array->GetEntity().EntityType,
2058 array->GetEntity().GeometryType,
2059 array->GetConnectivity(),
2061 conn->GetPointer(0)) ) < 0)
2063 vtkErrorMacro("Error while load connectivity of cells "
2064 << array->GetEntity().GeometryType);
2067 if ( MEDfilterClose( &filter ) < 0)
2069 vtkErrorMacro("ERROR : filter closing ...");
2076 void vtkMedDriver30::LoadCellGlobalIds(vtkMedEntityArray* array)
2078 if(array->IsGlobalIdsLoaded())
2081 FileOpen open(this);
2083 vtkMedIntArray* globalIds = vtkMedIntArray::New();
2084 array->SetGlobalIds(globalIds);
2085 globalIds->Delete();
2087 globalIds->SetNumberOfTuples(array->GetNumberOfEntity());
2089 vtkMedGrid* grid = array->GetParentGrid();
2090 vtkMedComputeStep cs = grid->GetComputeStep();
2092 if( MEDmeshEntityNumberRd (
2094 grid->GetParentMesh()->GetName(),
2097 array->GetEntity().EntityType,
2098 array->GetEntity().GeometryType,
2099 globalIds->GetPointer(0) ) < 0)
2101 array->SetGlobalIds(NULL);
2105 void vtkMedDriver30::LoadField(vtkMedFieldOnProfile* fop, med_storage_mode mode)
2107 FileOpen open(this);
2109 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2110 vtkMedFieldStep *step = fieldOverEntity->GetParentStep();
2111 vtkMedField* field = step->GetParentField();
2112 const vtkMedComputeStep& cs = step->GetComputeStep();
2114 vtkDataArray* data = vtkMedUtilities::NewArray(field->GetDataType());
2119 if(mode == MED_COMPACT_STMODE)
2121 size = fop->GetNumberOfValues();
2125 med_int profilesize;
2126 med_int nbofintegrationpoint;
2127 char profileName[MED_NAME_SIZE+1] = "";
2128 char localizationName[MED_NAME_SIZE+1] = "";
2129 size = MEDfieldnValueWithProfile(this->FileId,
2133 fieldOverEntity->GetEntity().EntityType,
2134 fieldOverEntity->GetEntity().GeometryType,
2135 fop->GetMedIterator(),
2140 &nbofintegrationpoint);
2143 if(fop->GetNumberOfIntegrationPoint() > 1)
2145 size *= fop->GetNumberOfIntegrationPoint();
2148 data->SetNumberOfComponents(field->GetNumberOfComponent());
2149 data->SetNumberOfTuples(size);
2150 if (this->ParallelFileId == -1)
2152 if ( MEDfieldValueWithProfileRd(
2157 fieldOverEntity->GetEntity().EntityType,
2158 fieldOverEntity->GetEntity().GeometryType,
2160 fop->GetProfileName(),
2162 MED_ALL_CONSTITUENT,
2163 (unsigned char*) data->GetVoidPointer(0) ) < 0)
2165 vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
2170 if (field->GetFieldType() == vtkMedField::CellField)
2172 med_filter filter = MED_FILTER_INIT;
2179 fop->GetFilter()->GetFilterSizes(start, stride, count,
2180 blocksize, lastblocksize );
2182 if ( MEDfilterBlockOfEntityCr( this->ParallelFileId,
2183 fop->GetNumberOfValues(),
2184 1, // one is for mesh elements, more than 1 is for fields
2185 field->GetNumberOfComponent(),
2186 MED_ALL_CONSTITUENT,
2193 (med_size)blocksize,
2194 (med_size)lastblocksize,
2197 vtkErrorMacro("Filter creation ");
2200 if ( MEDfieldValueAdvancedRd(
2201 this->ParallelFileId,
2205 fieldOverEntity->GetEntity().EntityType,
2206 fieldOverEntity->GetEntity().GeometryType,
2208 (unsigned char*) data->GetVoidPointer(0) ) < 0)
2210 vtkErrorMacro("Error on MEDfieldValueAdvancedRd");
2213 if ( MEDfilterClose( &filter ) < 0)
2215 vtkErrorMacro("ERROR : filter closing ...");
2219 {//TODO : option utilisateur pour desactiver ou non les champs avec profile en //
2220 if ( MEDfieldValueWithProfileRd(
2225 fieldOverEntity->GetEntity().EntityType,
2226 fieldOverEntity->GetEntity().GeometryType,
2228 fop->GetProfileName(),
2230 MED_ALL_CONSTITUENT,
2231 (unsigned char*) data->GetVoidPointer(0) ) < 0)
2233 vtkErrorMacro("Error on MEDfieldValueWithProfileRd");
2239 void vtkMedDriver30::LoadVariableAttribute(vtkMedVariableAttribute* varatt,
2240 vtkMedEntityArray* array)
2242 FileOpen open(this);
2246 vtkAbstractArray* valuearray = array->GetVariableAttributeValue(varatt);
2247 // first test if this is already loaded
2248 if(valuearray != NULL && valuearray->GetNumberOfTuples() > 0)
2251 if(valuearray == NULL)
2253 valuearray = vtkMedUtilities::NewArray(varatt->GetAttributeType());
2254 array->SetVariableAttributeValues(varatt, valuearray);
2255 valuearray->Delete();
2258 valuearray->SetNumberOfComponents(varatt->GetNumberOfComponent());
2259 valuearray->SetNumberOfTuples(array->GetNumberOfEntity());
2260 valuearray->SetName(varatt->GetName());
2262 vtkSmartPointer<vtkCharArray> chararray = vtkSmartPointer<vtkCharArray>::New();
2264 if(varatt->GetAttributeType() != MED_ATT_NAME)
2266 value = valuearray->GetVoidPointer(0);
2270 chararray->SetNumberOfValues(varatt->GetNumberOfComponent() *
2271 array->GetNumberOfEntity() *
2274 value = chararray->GetVoidPointer(0);
2277 vtkMedComputeStep cs = array->GetParentGrid()->GetComputeStep();
2279 if(MEDmeshStructElementVarAttRd(
2281 array->GetParentGrid()->GetParentMesh()->GetName(),
2284 varatt->GetParentStructElement()->GetGeometryType(),
2289 if(cs.IterationIt == MED_NO_IT && cs.TimeIt == MED_NO_DT && cs.TimeOrFrequency == MED_UNDEF_DT)
2291 vtkErrorMacro("MEDmeshStructElementVarAttRd");
2294 // try to see if I can reuse
2295 // the variable attributes of the NO_DT, NO_IT compute step
2296 vtkMedComputeStep nocs;
2297 nocs.IterationIt = MED_NO_IT;
2298 nocs.TimeIt = MED_NO_DT;
2299 nocs.TimeOrFrequency = MED_UNDEF_DT;
2300 vtkMedEntityArray* nocs_array =
2301 array->GetParentGrid()->GetParentMesh()->GetGridStep(nocs)->GetEntityArray(array->GetEntity());
2302 if(nocs_array == NULL)
2304 nocs_array = array->GetParentGrid()->GetParentMesh()->GetGridStep(0)->GetEntityArray(array->GetEntity());
2307 if(nocs_array == NULL || nocs_array == array)
2309 // try to force load the default compute step.
2310 if(MEDmeshStructElementVarAttRd(
2312 array->GetParentGrid()->GetParentMesh()->GetName(),
2315 varatt->GetParentStructElement()->GetGeometryType(),
2319 vtkErrorMacro("MEDmeshStructElementVarAttRd");
2325 this->LoadVariableAttribute(varatt, nocs_array);
2326 array->SetVariableAttributeValues(varatt, nocs_array->GetVariableAttributeValue(varatt));
2331 // If I am here, it means that I read the values
2332 if(varatt->GetAttributeType() == MED_ATT_NAME)
2334 char current_name[MED_NAME_SIZE+1] = "";
2335 vtkStringArray* sarray = vtkStringArray::SafeDownCast(valuearray);
2336 for(vtkIdType id = 0; id < varatt->GetNumberOfComponent() *
2337 array->GetNumberOfEntity(); id++)
2339 memset(current_name, '\0', MED_NAME_SIZE+1);
2340 strncpy(current_name, ((char*)value) + id*MED_NAME_SIZE, MED_NAME_SIZE);
2341 sarray->SetValue(id, current_name);
2348 void vtkMedDriver30::PrintSelf(ostream& os, vtkIndent indent)
2350 this->Superclass::PrintSelf(os, indent);