]> SALOME platform Git repositories - modules/med.git/blob - src/MEDMEM/MEDMEM_MedFieldDriver.txx
Salome HOME
Fix problem of make distcheck
[modules/med.git] / src / MEDMEM / MEDMEM_MedFieldDriver.txx
1 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
5 //
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License.
10 //
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14 // Lesser General Public License for more details.
15 //
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
19 //
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 //
22 // File      : MEDMEM_MedFieldDriver22.txx
23 // Created   : Thu Jul 22 12:52:43 2010
24 // Author    : Edward AGAPOV (eap)
25 //
26 // This file contains impelementations of methods that must be included together with
27 // the template class declaration
28
29 #ifndef MED_FIELD_DRIVER_TXX
30 #define MED_FIELD_DRIVER_TXX
31
32 #include <algorithm>
33
34 #include "MEDMEM_FieldConvert.hxx"
35 #include "MEDMEM_ArrayInterface.hxx"
36 #include "MEDMEM_ArrayConvert.hxx"
37
38 #include "MEDMEM_STRING.hxx"
39 #include "MEDMEM_Exception.hxx"
40
41 #include "MEDMEM_DriversDef.hxx"
42 #include "MEDMEM_MedFieldDriver.hxx"
43 #include "MEDMEM_Unit.hxx"
44 #include "MEDMEM_Support.hxx"
45 #include "MEDMEM_Family.hxx"
46 #include "MEDMEM_Group.hxx"
47 #include "MEDMEM_GaussLocalization.hxx"
48
49 namespace MEDMEM
50 {
51
52 /*!
53   Constructor.
54 */
55 template <class T>
56 MED_FIELD_DRIVER<T>::MED_FIELD_DRIVER():
57   GENDRIVER(MED_DRIVER),
58   _ptrField((FIELD<T> *) MED_NULL),
59   _fieldName(""),_fieldNum(MED_EN::MED_INVALID),_medIdt(MED_INVALID)
60 {}
61
62 /*!
63   Constructor.
64 */
65 template <class T>
66 template <class INTERLACING_TAG>
67 MED_FIELD_DRIVER<T>::MED_FIELD_DRIVER(const string &              fileName,
68                                       FIELD<T, INTERLACING_TAG> * ptrField,
69                                       MED_EN::med_mode_acces      accessMode)
70   : GENDRIVER(fileName, accessMode, MED_DRIVER),
71     _ptrField((FIELD<T> *) ptrField),
72     _fieldName(""),_fieldNum(MED_EN::MED_INVALID),_medIdt(MED_INVALID)
73 {
74 }
75
76 /*!
77   Copy constructor.
78 */
79 template <class T>
80 MED_FIELD_DRIVER<T>::MED_FIELD_DRIVER(const MED_FIELD_DRIVER<T> & fieldDriver):
81   GENDRIVER(fieldDriver),
82   _ptrField(fieldDriver._ptrField),
83   _fieldName(fieldDriver._fieldName),
84   _fieldNum(fieldDriver._fieldNum),
85   _medIdt(fieldDriver._medIdt)
86 {
87 }
88
89 /*!
90   Destructor.
91 */
92 template <class T>
93 MED_FIELD_DRIVER<T>::~MED_FIELD_DRIVER()
94 {
95   MESSAGE_MED("MED_FIELD_DRIVER<T>::~MED_FIELD_DRIVER() has been destroyed");
96 }
97
98 /*!
99   Open the file
100 */
101 template <class T>
102 void MED_FIELD_DRIVER<T>::open() throw (MEDEXCEPTION)
103 {
104   const char * LOC = "MED_FIELD_DRIVER::open() ";
105   BEGIN_OF_MED(LOC);
106
107   // we must set fieldname before open, because we must find field number in file (if it exist !!!)
108   if ( MED_FIELD_DRIVER<T>::_fileName == "" )
109     throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
110                                      << "_fileName is |\"\"|, please set a correct fileName before calling open()"
111                                      )
112                           );
113   if ( MED_FIELD_DRIVER<T>::_status==MED_OPENED )
114     return;
115
116   int accessMode = MED_FIELD_DRIVER<T>::_accessMode;
117   if ( accessMode == MED_EN::RDWR )
118     accessMode = med_2_3::MED_ACC_RDWR;
119   MESSAGE_MED(LOC<<"_fileName.c_str : "<< MED_FIELD_DRIVER<T>::_fileName.c_str()<<",mode : "<< MED_FIELD_DRIVER<T>::_accessMode);
120   MED_FIELD_DRIVER<T>::_medIdt = med_2_3::MEDfileOpen( MED_FIELD_DRIVER<T>::_fileName.c_str(),(med_2_3::med_access_mode) accessMode);
121   MESSAGE_MED(LOC<<"_medIdt : "<< MED_FIELD_DRIVER<T>::_medIdt );
122   if (MED_FIELD_DRIVER<T>::_medIdt > 0)
123     MED_FIELD_DRIVER<T>::_status=MED_OPENED;
124   else {
125     MED_FIELD_DRIVER<T>::_status = MED_INVALID;
126     //MED_FIELD_DRIVER<T>::_medIdt = MED_INVALID;
127     throw MED_EXCEPTION (LOCALIZED( STRING(LOC)
128                                     << "Can't open |"  << MED_FIELD_DRIVER<T>::_fileName
129                                     << "|, _medIdt : " << MED_FIELD_DRIVER<T>::_medIdt
130                                     )
131                          );
132   }
133
134   END_OF_MED(LOC);
135 }
136
137 /*!
138   Close the file
139 */
140 template <class T>
141 void MED_FIELD_DRIVER<T>::close()
142 {
143   const char* LOC = "MED_FIELD_DRIVER::close()";
144   BEGIN_OF_MED(LOC);
145   med_2_3::med_int err = 0;
146   if (MED_FIELD_DRIVER<T>::_status == MED_OPENED) {
147     err=med_2_3::MEDfileClose(MED_FIELD_DRIVER<T>::_medIdt);
148     if ( err )
149       cout << LOC << "can't CLOSE file"<<MED_FIELD_DRIVER<T>::_fileName<<endl;
150     //H5close(); // If we call H5close() all the files are closed.
151     MED_FIELD_DRIVER<T>::_status = MED_CLOSED;
152     MED_FIELD_DRIVER<T>::_medIdt = MED_INVALID;
153     MESSAGE_MED(" MED_FIELD_DRIVER::close() : MEDfermer : _medIdt= " << MED_FIELD_DRIVER<T>::_medIdt );
154     MESSAGE_MED(" MED_FIELD_DRIVER::close() : MEDfermer : err    = " << err );
155   }
156   END_OF_MED(LOC);
157 }
158
159
160 /*!
161   Constructor.
162 */
163 template <class T>
164 MED_FIELD_RDONLY_DRIVER<T>::MED_FIELD_RDONLY_DRIVER():MED_FIELD_DRIVER<T>()
165 {
166   this->GENDRIVER::_accessMode = MED_EN::RDONLY;
167 }
168
169 /*!
170   Constructor.
171 */
172 template <class T>
173 template <class INTERLACING_TAG>
174 MED_FIELD_RDONLY_DRIVER<T>::MED_FIELD_RDONLY_DRIVER(const string &              fileName,
175                                                     FIELD<T, INTERLACING_TAG> * ptrField):
176   MED_FIELD_DRIVER<T>(fileName,ptrField,MED_EN::RDONLY)
177 {
178   const char* LOC = "MED_FIELD_RDONLY_DRIVER::MED_FIELD_RDONLY_DRIVER(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)";
179   BEGIN_OF_MED(LOC);
180   END_OF_MED(LOC);
181 }
182
183 /*!
184   Copy constructor.
185 */
186 template <class T>
187 MED_FIELD_RDONLY_DRIVER<T>::MED_FIELD_RDONLY_DRIVER(const MED_FIELD_RDONLY_DRIVER & fieldDriver):
188   MED_FIELD_DRIVER<T>(fieldDriver)
189 {}
190
191 /*!
192
193   Cette méthode crée le SUPPORT du champ <fieldName> pour le
194         <n°de pas de temps,n°d'itération>=<ndt,od>.
195
196   Le SUPPORT crée à pour nom  <fieldName>Support et contient
197   la liste des types géométriques sur le premier type
198   d'entité trouvé (en MEDMEM on inderdit aux champs de reposer
199   sur plusieurs types d'entité).
200   Il contient également le nombre d'entités trouvées pour chaque
201   type géométrique.
202   Par défaut l'attribut onAll du SUPPORT est positionné à true car
203   cette routine ne lit rien de ce qui concerne les entités
204   du maillage associé.
205   La méthode renvoie true si elle réussit à créer le SUPPORT
206   demandé.
207   Le nom du maillage associé ( en MEDMEM on ne
208   supporte pas encore les maillages multiples ) est renvoyé dans <meshName>.
209   Deux tableaux directements exploitables par MEDMEMnArray sont renvoyés :
210   - numberOfElementsOfTypeC : nombres d'entités cumulés de chaque type géométrique
211   avec numberOfElementsOfTypeC[0]=1 et de taille nombre de types+1
212   - numberOfGaussPoint : nombre de points de Gauss par type géométrique
213   avec numberOfGaussPoint[0]=1 et de taille  nombre de types+1.
214   <preferEntity> argument (if not default) gives entity of already present support,
215   which is used to select entity if there are several of them. See issue
216   20582: [CEA 368] MEDMEM don't work with a same field on NODES and CELLS
217 */
218
219 template <class T> bool
220 MED_FIELD_DRIVER<T>::createFieldSupportPart1(med_2_3::med_idt        id,
221                                              const string &          fieldName,
222                                              med_2_3::med_int        ndt,
223                                              med_2_3::med_int        od,
224                                              SUPPORT &               support,
225                                              string &                meshName,
226                                              vector<int> &           numberOfElementsOfTypeC,
227                                              vector<int> &           numberOfGaussPoint,
228                                              std::vector<std::string>& gaussModelName,
229                                              std::vector<std::string>& profilName,
230                                              int &                   totalNumberOfElWg,
231                                              MED_EN::medEntityMesh & fieldMedFileEntity,
232                                              MED_EN::medEntityMesh   preferEntity
233                                              ) const throw (MEDEXCEPTION)
234 {
235
236   //EF : Gérer le meshName pour le driver 2.2
237   const char * LOC="MED_FIELD_DRIVER<T>::createFieldSupportPart1(...)";
238
239   BEGIN_OF_MED(LOC);
240
241   map<int, list<MED_EN::medGeometryElement> > CellAndNodeEntities;
242   map<int, list<MED_EN::medGeometryElement> >::iterator currentEntity;
243   CellAndNodeEntities[MED_EN::MED_CELL]  = MED_EN::meshEntities[MED_EN::MED_CELL];
244   CellAndNodeEntities[MED_EN::MED_NODE] = MED_EN::meshEntities[MED_EN::MED_NODE];
245   list< MED_EN::medGeometryElement >::const_iterator currentGeometry;
246
247   MED_EN::medEntityMesh entityCurrent;
248   MED_EN::medGeometryElement geometryCurrent;
249
250   MED_EN::medEntityMesh entity;
251   MED_EN::medEntityMesh medmem_entity;
252
253   bool alreadyFoundAnEntity=false;//,alreadyFoundPdtIt = false;
254   int  numberOfElements = 0;
255   int  numberOfGeometricType = 0;
256   MED_EN::medGeometryElement geometricType[MED_N_CELL_GEO_FIXED_CON];
257   int numberOfElementsOfType[MED_N_CELL_GEO_FIXED_CON];
258   numberOfElementsOfTypeC.clear();numberOfGaussPoint.clear();
259   numberOfElementsOfTypeC.resize(MED_N_CELL_GEO_FIXED_CON+1);
260   numberOfGaussPoint.resize(MED_N_CELL_GEO_FIXED_CON+1);
261
262   med_2_3::med_int numdt=-1, numo=-1, ngauss=0, nbPdtIt=0;
263   char dtunit[MED_LNAME_SIZE+1];
264   char maa[MED_NAME_SIZE+1];
265   med_2_3::med_float   dt=-1.0;
266   med_2_3::med_bool local;
267   med_2_3::med_err     ret=1;
268   numberOfElementsOfTypeC[0] = 1;
269   numberOfGaussPoint[0] = 1;
270   totalNumberOfElWg = 0;
271   int field_dim=0;
272
273   char pflnamep3[MED_NAME_SIZE+1] = "";
274   int profilesizep3;
275   char locnamep3[MED_NAME_SIZE+1] = "";
276
277   /* Détermine le type d'entité et la liste des types géométriques associés
278      au champ <fieldName> */
279   for (currentEntity = CellAndNodeEntities.begin();
280        currentEntity != CellAndNodeEntities.end(); currentEntity++) { // loop on [CELL,NODE]
281     for (currentGeometry  = (*currentEntity).second.begin();
282          currentGeometry != (*currentEntity).second.end(); currentGeometry++) {
283
284       entityCurrent = (*currentEntity).first ;
285       geometryCurrent = (*currentGeometry) ;
286
287       medmem_entity = entityCurrent;
288
289       // That is a difference between Med File and Med Memory (NB)
290       if (geometryCurrent == MED_EN::MED_SEG2 || geometryCurrent == MED_EN::MED_SEG3)
291         medmem_entity = MED_EN::MED_EDGE;
292
293       if (geometryCurrent == MED_EN::MED_TRIA3 || geometryCurrent == MED_EN::MED_QUAD4 ||
294           geometryCurrent == MED_EN::MED_TRIA6 || geometryCurrent == MED_EN::MED_QUAD8 || 
295           geometryCurrent == MED_EN::MED_POLYGON)
296         medmem_entity = MED_EN::MED_FACE;
297
298       // med3 porting
299       numberOfElements = med_2_3::MEDfieldnValueWithProfile
300         (id,fieldName.c_str(),ndt,od,
301          (med_2_3::med_entity_type)   entityCurrent,
302          (med_2_3::med_geometry_type) geometryCurrent,
303          /*profileit=*/1,med_2_3::MED_COMPACT_PFLMODE,
304          pflnamep3,&profilesizep3,locnamep3,&ngauss);
305
306       if ( numberOfElements <= 0 )
307         numberOfElements = med_2_3::MEDfieldnValueWithProfile
308           (id,fieldName.c_str(),ndt,od,
309            (med_2_3::med_entity_type) ( entityCurrent = medmem_entity ), // try the other entity !
310            (med_2_3::med_geometry_type) geometryCurrent,
311            /*profileit=*/1,med_2_3::MED_COMPACT_PFLMODE,
312            pflnamep3,&profilesizep3,locnamep3,&ngauss);
313
314       if ( numberOfElements <= 0 )
315         continue;
316       
317       if ( preferEntity != MED_EN::MED_ALL_ENTITIES ) //issue 20582: field on NODES and CELLS
318       {
319         bool preferNode = (preferEntity  == MED_EN::MED_NODE);
320         bool  foundNode = (medmem_entity == MED_EN::MED_NODE);
321         if ( preferNode != foundNode )
322           continue;
323       }
324
325       /* Verifie que le champ n'est pas défini sur un autre type d'entité */
326       if ( alreadyFoundAnEntity )
327       {
328         if ( medmem_entity != entity )
329         {
330           throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |"  << fieldName
331                                        << "| with (ndt,or) = (" << ndt << ","
332                                        << od << ") must not be defined on different entity types" ));
333         }
334       }
335       //else
336       { 
337         entity=medmem_entity;
338         fieldMedFileEntity = entityCurrent;
339         alreadyFoundAnEntity = true;
340
341         //at this stage, the entity is only related to the geometry dimension
342         // with no relation to the mesh dimension
343         // therefore MED_CELL refer to 3D and MED_FACE to 2D
344         // the correct entity (which depends on the mesh dimension, is only 
345         // determined later
346
347         switch (entity)
348         {
349         case MED_NODE:
350           field_dim=0;
351           break;
352         case MED_EDGE:
353           field_dim=1;
354           break;
355         case MED_FACE:
356           field_dim=2;
357           break;
358         case MED_CELL:
359           field_dim=3;
360           break;
361         }
362       }
363
364       if (/* (ret != 0)  ||*/ (ngauss < 1 ) )
365         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Error in MEDfieldnValueWithProfile for  Field |" << fieldName 
366                                      << "| with (ndt,or) = ("
367                                      << ndt << "," << od << ") for (entityType,geometricType)=("
368                                      << MED_EN::entNames[entityCurrent] << ","
369                                      << MED_EN::geoNames[geometryCurrent] << ")" )); ;
370
371       numberOfElementsOfType[numberOfGeometricType] = numberOfElements;
372       numberOfElementsOfTypeC[numberOfGeometricType+1]=
373         numberOfElementsOfTypeC[numberOfGeometricType]
374         +  numberOfElementsOfType[numberOfGeometricType];
375       numberOfGaussPoint[numberOfGeometricType+1] = ngauss;
376       geometricType[numberOfGeometricType]= geometryCurrent;
377       numberOfGeometricType++;
378       totalNumberOfElWg+=numberOfElements*ngauss;
379
380       gaussModelName.push_back( locnamep3 );
381       locnamep3[0] = '\0';
382       profilName.push_back( pflnamep3 );
383       pflnamep3[0] = '\0';
384     } // End Second For
385
386   } // End Premier For
387
388   if ( numberOfGeometricType == 0 )
389     return false;
390
391   {
392     // Read mesh name and nb of timestamps in the field
393     med_2_3::med_field_type typchap3;
394     int nbCompp3=med_2_3::MEDfieldnComponentByName(id,fieldName.c_str());
395     if ( nbCompp3 <= 0 )
396       return false;
397     char *compNamesp3=new char[nbCompp3*MED_SNAME_SIZE+1];
398     char *compUnitp3=new char[nbCompp3*MED_SNAME_SIZE+1];
399     ret = med_2_3::MEDfieldInfoByName(id,fieldName.c_str(),maa,&local,&typchap3,compNamesp3,compUnitp3,dtunit,&nbPdtIt);
400     delete [] compNamesp3;
401     delete [] compUnitp3;
402     if ( ret < 0 )
403       return false;
404
405
406     //VB commented out to allow for fields on multiple meshes
407     //  if ( ! meshName.empty() )
408     //          if ( meshName != maa ) {
409     //            throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |" << fieldName << "| with (ndt,or) = ("
410     //                                         << ndt << "," << od << ") for (entityType,geometricType)=("
411     //                                         << MED_EN::entNames[entityCurrent] << ","
412     //                                         << MED_EN::geoNames[*currentGeometry] << ")"
413     //                                         << "is defined on mesh |" << maa << "| not on mesh |" << meshName ));
414     if ( meshName.empty() )
415       meshName = maa;
416     if ( meshName.empty() )
417       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Can't get mesh name for field |"  << fieldName
418                                    << "| using MEDfieldInfoByName()" ));
419     // Read the time
420     bool timeFound = false;
421     for ( med_2_3::med_int j=1; j <= nbPdtIt; j++ )
422       {
423         ret = med_2_3::MEDfieldComputingStepInfo(id,fieldName.c_str(),j,&numdt,&numo,&dt);
424         if ( ret == MED_INVALID )
425           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) << ": can't read the "
426                                         << j << "-th timestamp info of field " << fieldName ));
427         if ( numdt == ndt && numo == od )
428           {
429             MED_FIELD_DRIVER<T>::_ptrField->setTime(dt);
430             timeFound = true;
431             break;
432           }
433       }
434     if ( !timeFound )
435       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Timestamp with (ndt,or) = ("
436                                    << ndt << "," << od << ") not found for field " << fieldName));
437   }
438    
439   //retrieves the right medmem entity type from field_dim and mesh_dim
440   int mesh_dim = MED_FIELD_DRIVER<T>::getMeshDimensionFromFile(id,meshName);
441   //if (mesh_dim==2 && field_dim==2) // 0020644: [CEA 384] Problem with 1D no structured meshing
442   if (mesh_dim==field_dim)
443     entity=MED_CELL;
444
445   if ( alreadyFoundAnEntity) {
446     support.setName(fieldName+" Support");
447     support.setMeshName(meshName); // Vérifier que les différents noms de maillages lus soient identiques
448     support.setEntity(entity);
449     // REM : Le nombre <numberOfGeometricType> dans la précédente version du Driver 
450     //       était erronée pour un champ qui ne reposait pas sur toutes les entités géométriques 
451     //       du maillage mais dont le SUPPORT a été crée à partir des informations du maillage
452     //       ( méthode qui était largement utilisée pour construire un SUPPORT).
453     support.setNumberOfGeometricType(numberOfGeometricType);
454     support.setGeometricType(geometricType); // Utile uniquement si setAll == false ?
455     support.setNumberOfElements(numberOfElementsOfType);    //setNumberOfElements effectue une copie 
456     // Par défaut considère que le champ repose sur tous les type géométriques du maillage
457     // Si ce n'est pas le cas les champs geometricType et numberOfElementsOfType du SUPPORT sont corrects
458     support.setAll(true);
459     numberOfElementsOfTypeC.resize(numberOfGeometricType+1);
460     numberOfGaussPoint.resize(numberOfGeometricType+1);
461   }
462   return alreadyFoundAnEntity;
463 }
464
465 template <class T> MED_EN::medEntityMesh
466 MED_FIELD_DRIVER<T>::getMEDMEMEntityFromMEDType(medGeometryElement type,
467                                                 int mesh_dim) const
468 {
469   int elem_dim = type/100;
470   if (type==MED_POLYGON) elem_dim=2;
471   if (type==MED_POLYHEDRA) elem_dim=3;
472
473   if (elem_dim==3)
474     return MED_CELL;
475   if (elem_dim==2)
476     {
477       if (mesh_dim==2)
478         return MED_CELL;
479       else if (mesh_dim==3)
480         return MED_FACE;
481     }
482   if (elem_dim==1)
483     return MED_EDGE;
484   if(elem_dim==0)
485     return MED_NODE;
486 }
487
488
489 template <class T> int
490 MED_FIELD_DRIVER<T>::getMeshDimensionFromFile(med_2_3::med_idt id,
491                                                 const string &   meshName) const
492 {
493   const char* LOC = "MED_FIELD_DRIVER<T>::getMeshDimensionFromFile(...)";
494   BEGIN_OF_MED(LOC);
495
496   // BEGIN Issue 0020620: [CEA 378] Abort during MED mesh reading
497   // Find out if mesh is structured and read its dim if it is
498   med_2_3::med_int      meshDim;
499   med_2_3::med_int      spaceDimp3;
500   med_2_3::med_mesh_type meshType;
501   char                  meshNameInFile[MED_NAME_SIZE+1];
502   char                  meshDescription[MED_COMMENT_SIZE+1];
503   char                  dtunit[MED_LNAME_SIZE+1];
504   med_2_3::med_sorting_type      sorttypep3;
505   med_2_3::med_int nstepp3;
506   med_2_3::med_axis_type atp3;
507   // find mesh <meshName>
508   int numberOfMeshes = med_2_3::MEDnMesh(id);
509   for (int i = 1; i <= numberOfMeshes; ++i )
510     {
511 #ifdef WIN32
512       int naxis=max(3,med_2_3::MEDmeshnAxis(id,i));
513 #else
514       int naxis=std::max(3,med_2_3::MEDmeshnAxis(id,i));
515 #endif
516       char *axisname=new char[naxis*MED_SNAME_SIZE+1];
517       char *axisunit=new char[naxis*MED_SNAME_SIZE+1];
518       med_2_3::MEDmeshInfo(id, i ,meshNameInFile, &spaceDimp3, &meshDim, &meshType, meshDescription,dtunit,&sorttypep3,&nstepp3,&atp3,axisname,axisunit);
519       delete [] axisname;
520       delete [] axisunit;
521       if ( meshName == meshNameInFile )
522         {
523           if ( meshType == med_2_3::MED_STRUCTURED_MESH )
524             return meshDim;
525           break;
526         }
527     }
528   // END Issue 0020620: [CEA 378] Abort during MED mesh reading
529
530   int numberOfGeometricType=0;
531   MED_EN::medGeometryElement geometricType[MED_N_CELL_GEO_FIXED_CON];
532   med_2_3::med_int   numberOfElements=0;
533
534   /*in MED file, all entities are regarded as MED_CELL
535     (except for those related to descending connectivities),
536     whereas in MEDMEM the distinction between MED_CELL, MED_FACE and MED_EDGE exists
537     it is therefore necessary to distinguish the MED-file entity
538     that will be used for the call to MED-file
539     and the MEDMEM entity*/
540   const MED_EN::medEntityMesh entity = MED_EN::MED_CELL;
541
542   list<MED_EN::medGeometryElement>::const_iterator currentGeometry;
543
544   med_2_3::med_int dtp3 = MED_NO_DT, itp3 = MED_NO_IT;
545   med_2_3::med_float dtfp3 = MED_NO_DT;
546   med_2_3::med_bool chp3,trp3;
547   med_2_3::MEDmeshComputationStepInfo(id,meshName.c_str(),1,&dtp3,&itp3,&dtfp3);
548
549   for (currentGeometry  = (MED_EN::meshEntities[entity]).begin();
550        currentGeometry != (MED_EN::meshEntities[entity]).end(); currentGeometry++)
551     {
552       numberOfElements =
553         med_2_3::MEDmeshnEntity( id, meshName.c_str(),
554                                  dtp3,itp3,
555                                  med_2_3::MED_CELL,(med_2_3::med_geometry_type) *currentGeometry,
556                                  med_2_3::MED_CONNECTIVITY,med_2_3::MED_NODAL,
557                                  &chp3,&trp3);
558       if (numberOfElements <= 0)
559         continue;
560
561       geometricType[numberOfGeometricType] = *currentGeometry;
562
563       numberOfGeometricType++;
564     }
565
566   //Because MEDFILE and MEDMEM differ on the definition of MED_CELL
567   //it is necessary to remove the cells that do not
568   //have maximum cell dimension in the range covered by geometricType
569   int maxdim=0;
570   for (int i=0; i<numberOfGeometricType; i++)
571   {
572     const CELLMODEL& model = CELLMODEL_Map::retrieveCellModel(geometricType[i]);
573     int dim = model.getDimension();
574     if (dim>maxdim) maxdim=dim;
575   }
576
577   return maxdim;
578
579 }
580
581 /*!
582
583   Renvoie la liste <geoType> des types géométriques définis dans le maillage <meshName>
584   pour le type d'entité <entity>.
585   * < nbOfElOfType > contient le nombre d'entités de chaque type
586   * < numberOfElementsOfTypeC > contient le nombre d'entités cumulées de chaque type
587                               avec numberOfElementsOfTypeC[0]=0;
588   * < allDimensions > controls dimension of returned types of entity == MED_CELL
589 */
590 template <class T> void
591 MED_FIELD_DRIVER<T>::getMeshGeometricTypeFromFile(med_2_3::med_idt      id,
592                                                   string &              meshName,
593                                                   MED_EN::medEntityMesh entity,
594                                                   vector<MED_EN::medGeometryElement> & geoType,
595                                                   vector<int> &         nbOfElOfType,
596                                                   vector<int> &         nbOfElOfTypeC
597                                                   ) const throw(MEDEXCEPTION)
598 {
599   const char* LOC = "MED_FIELD_DRIVER<T>::getMeshGeometricTypeFromFile(...)";
600   BEGIN_OF_MED(LOC);
601
602   int numberOfGeometricType=0;
603   MED_EN::medGeometryElement geometricType[MED_N_CELL_GEO_FIXED_CON];
604   int numberOfElementsOfType [MED_N_CELL_GEO_FIXED_CON];
605   int numberOfElementsOfTypeC[MED_N_CELL_GEO_FIXED_CON+1];
606   int dimOfType[MED_N_CELL_GEO_FIXED_CON];
607   int maxdim=0;
608   med_2_3::med_int   numberOfElements=0;
609   med_2_3::med_data_type quoi;
610
611   /*in MED file, all entities are regarded as MED_CELL
612     (except for those related to descending connectivities),
613     whereas in MEDMEM the distinction between MED_CELL, MED_FACE and MED_EDGE exists
614     it is therefore necessary to distinguish the MED-file entity
615     that will be used for the call to MED-file
616     and the MEDMEM entity*/
617   MED_EN::medEntityMesh medfile_entity;
618   if (entity==MED_EN::MED_NODE)
619   {
620     medfile_entity=MED_EN::MED_NODE;
621     quoi=med_2_3::MED_COORDINATE;
622   }
623   else
624   {
625     medfile_entity=MED_EN::MED_CELL;
626     quoi=med_2_3::MED_CONNECTIVITY;
627   }
628
629   list<MED_EN::medGeometryElement>::const_iterator currentGeometry;
630   bool alreadyFoundAnEntity = false;
631   numberOfElementsOfTypeC[0]=0;
632
633   for (currentGeometry  = (MED_EN::meshEntities[entity]).begin();
634        currentGeometry != (MED_EN::meshEntities[entity]).end(); currentGeometry++)
635   {
636     med_2_3::med_int dtp3,itp3;
637     med_2_3::med_float dtfp3;
638     med_2_3::med_bool chp3,trp3;
639     med_2_3::MEDmeshComputationStepInfo(id,meshName.c_str(),1,&dtp3,&itp3,&dtfp3);
640     numberOfElements =
641       med_2_3::MEDmeshnEntity(id,meshName.c_str(),dtp3,itp3,(med_2_3::med_entity_type) medfile_entity,(med_2_3::med_geometry_type) *currentGeometry,quoi,med_2_3::MED_NODAL,&chp3,&trp3);
642
643     if (numberOfElements <= 0)
644       continue;
645
646     alreadyFoundAnEntity = true;
647     numberOfElementsOfType[numberOfGeometricType] = numberOfElements;
648     numberOfElementsOfTypeC[numberOfGeometricType+1] =
649       numberOfElementsOfTypeC[numberOfGeometricType]+numberOfElements;
650
651     MED_EN::medGeometryElement geomType = *currentGeometry;
652     geometricType[numberOfGeometricType] = geomType;
653
654     //Because MEDFILE and MEDMEM differ on the definition of MED_CELL
655     //it is necessary to remove the cells that do not
656     //have maximum cell dimension in the range covered by geometricType
657     const CELLMODEL& model = CELLMODEL_Map::retrieveCellModel( geomType );
658     const int dim = model.getDimension();
659     dimOfType[ numberOfGeometricType ] = dim;
660     if (dim>maxdim) maxdim=dim;
661
662     numberOfGeometricType++;
663   }
664
665   nbOfElOfTypeC.push_back(0);
666   for (int i=0; i<numberOfGeometricType; i++)
667   {
668     if (dimOfType[i]==maxdim || entity != MED_CELL)
669     {
670       geoType.push_back(geometricType[i]);
671       int nbelems = numberOfElementsOfType[i];
672       nbOfElOfType.push_back(nbelems);
673       nbOfElOfTypeC.push_back(nbOfElOfTypeC[nbOfElOfTypeC.size()-1]+nbelems);
674     }
675   }
676
677   //  geoType = vector<MED_EN::medGeometryElement>(geometricType,geometricType+numberOfGeometricType);
678   //  nbOfElOfType = vector<int> (numberOfElementsOfType,numberOfElementsOfType+numberOfGeometricType);
679   //  nbOfElOfTypeC = vector<int> (numberOfElementsOfTypeC,numberOfElementsOfTypeC+numberOfGeometricType+1);
680
681 //   for (int j =0 ; j<= numberOfGeometricType;++j)
682 //       cout << "nbOfElOfTypeC["<<j<<"]="<<nbOfElOfTypeC[j]<<endl;
683
684   END_OF_MED(LOC);
685 }
686
687 /*!
688 reads the MESH object in order to retrieve the list of geometric types for a given entity
689 \param[in] meshPtr pointer to MESH
690 \param[in] entity entity for which the geom types are required
691 \param[out] geoType list of geom types
692 \param[out] nbOfElOfType vector containing the number of elements per type (size : ntype)
693 \param[out] nbOfElOfTypeC accumulated version of nbOfElType (size : ntype+1)
694  */
695
696 template <class T> void
697 MED_FIELD_DRIVER<T>::getMeshGeometricTypeFromMESH( const GMESH * meshPtr,
698                                                    MED_EN::medEntityMesh  entity,
699                                                    std::vector<MED_EN::medGeometryElement> & geoType,
700                                                    std::vector<int> &nbOfElOfType,
701                                                    std::vector<int> &nbOfElOfTypeC) const throw(MEDEXCEPTION)
702 {
703   const char LOC[] = "MED_FIELD_DRIVER<T>::getMeshGeometricTypeFromMESH(...) : ";
704   BEGIN_OF_MED(LOC);
705
706   if (!meshPtr)
707     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"ptrMesh must be non null" )); ;
708
709   // Il est plus pratique de créer un support "onAll"
710   // pour calculer les tableaux du nombre d'entités cumulées
711
712   const SUPPORT *mySupportFromMesh=meshPtr->getSupportOnAll(entity);
713   geoType = vector<MED_EN::medGeometryElement>(mySupportFromMesh->getTypes(),
714                               mySupportFromMesh->getTypes()+mySupportFromMesh->getNumberOfTypes());
715   nbOfElOfType.resize(mySupportFromMesh->getNumberOfTypes());
716   nbOfElOfTypeC.resize(mySupportFromMesh->getNumberOfTypes()+1);
717   nbOfElOfTypeC[0]=0;
718
719   for (int j=1; j<=mySupportFromMesh->getNumberOfTypes(); ++j) {
720     nbOfElOfType[j-1]=mySupportFromMesh->getNumberOfElements(geoType[j-1]);
721     nbOfElOfTypeC[j]+=nbOfElOfTypeC[j-1]+nbOfElOfType[j-1];
722   }
723   END_OF_MED(LOC);
724 }
725
726 /*--------------------- RDONLY PART -------------------------------*/
727
728 template <class T> GENDRIVER * MED_FIELD_RDONLY_DRIVER<T>::copy(void) const
729 {
730   return new MED_FIELD_RDONLY_DRIVER<T>(*this);
731 }
732
733 /*!
734   In MEDMEM, FIELDs lie on support which can be plain SUPPORT, FAMILY
735   or GROUP, while in MED-file, there is no link between the FAMILY and
736   GROUP notions and the FIELDs. FIELDs lie on profiles.
737   The problem arises from the fact that the MED write driver creates
738   profiles when treating fields that lie on MEDMEM::SUPPORT,
739   MEDMEM_FAMILY or MEDMEM::GROUP. The profile is named after the
740   support name : nameOfSupport_<type_of_geometric_entity>.
741   However, the read driver is unable to link supports and profiles
742   and it recreates a new support that corresponds to the field profile.
743
744   To avoid this support recreation, pass the mesh to the FIELD's
745   constructor, and the field driver will find appropriate FAMILY or GROUP
746   in the mesh and use it for the field.
747  */
748 template <class T> void MED_FIELD_RDONLY_DRIVER<T>::read(void)
749   throw (MEDEXCEPTION)
750 {
751   const char * LOC = " MED_FIELD_RDONLY_DRIVER::read() " ;
752   BEGIN_OF_MED(LOC);
753
754   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,NoGauss>::Array       ArrayNo;
755   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,Gauss>::Array         ArrayNoWg;
756   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,NoGauss>::Array     ArrayFull;
757   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,Gauss>::Array       ArrayFullWg;
758   typedef typename MEDMEM_ArrayInterface<T,NoInterlaceByType,NoGauss>::Array ArrayByType;
759   typedef typename MEDMEM_ArrayInterface<T,NoInterlaceByType,Gauss>::Array   ArrayByTypeWg;
760
761   if (MED_FIELD_DRIVER<T>::_status!=MED_OPENED)
762     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": Method open must be called before method read.")) ;
763
764   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
765        ( MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
766     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
767                                  <<" neither <fieldName> is set in driver nor in object FIELD.")) ;
768
769   // If _fieldName is not set in driver, try to use _ptrfield->_fieldName
770   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
771        ( !MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
772     MED_FIELD_DRIVER<T>::_fieldName=MED_FIELD_DRIVER<T>::_ptrField->_name;
773
774   if ( MED_FIELD_DRIVER<T>::_fieldName.size() > MED_NAME_SIZE )
775     {
776       SCRUTE_MED(MED_FIELD_DRIVER<T>::_fieldName.size());
777       SCRUTE_MED(MED_NAME_SIZE);
778
779 //       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
780 //                                 <<" <fieldName> size in object driver FIELD is > MED_NAME_SIZE ."));
781
782       MESSAGE_MED(LOC << "Warning <fieldName> size in object driver FIELD is > MED_NAME_SIZE .");
783     }
784
785   const string & fieldName = MED_FIELD_DRIVER<T>::_fieldName;
786
787   MED_EN::medModeSwitch interlacingType = MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType();
788   bool isFullInterlace     = ( interlacingType == MED_EN::MED_FULL_INTERLACE );
789   bool isNoInterlaceByType = ( interlacingType == MED_EN::MED_NO_INTERLACE_BY_TYPE );//PAL17011
790
791   MESSAGE_MED("###### "<<LOC<<" fieldNameDRIVER : "<< fieldName << " fieldName : "<< MED_FIELD_DRIVER<T>::_ptrField->_name);
792
793 // EF :
794 //   Si un support a été donnée au champ, pour des raisons de compatibilité avec
795 //   les versions précédentes, ce support sera utilisé pour
796 //   - Obtenir le nom du maillage sur lequel on veut lire le champ
797 //     (eventuellement on pourrait l'utiliser pour selectionner un champ qui
798 //      repose sur plusieurs maillages cf HOMARD-ASTER, ce qui n'est pas géré dans MEDMEM)
799 //   -  vérifier le type d'entité (MED_NOEUD xor MED_MAILLE xor MED_FACE xor MED_ARETE ) sur lequel
800 //      il faut lire le champ qui est également retouvé.
801 //   - Si le support défini une liste d'entité ( différente de MED_ALL_ELEMENTS), celle-ci est ignorée
802 //     à la lecture et écrasé par soit :
803 //            - onall, après avoir vérifié que la liste des types géométriques utilisés par le champ
804 //                     est égale à la liste des type géométriques définis dans le maillage associé
805 //                     pour tous le même type d'entité.
806 //            - La sous liste des types géométriques utilisés (onAll quand même, cf commenataire ci-dessous )  
807 //            - les listes de profils lus s'il en existe pour une sous liste de types
808 //              géométriques
809
810 //   Si aucun support n'a été donné au champ :
811 //   - A la lecture : Un support est crée et le type d'entité unique est lu
812 //                    (cf decision gt MED qu'un champ repose sur une entité unique ?),
813 //                    l'ensemble des types géométriques est lu,
814 //                    l'ensemble des profils par type géométrique est lu
815 //                    Le nom du maillage associé est lu mais le pointeur SUPPORT-MESH non initialisé
816
817
818   char tmpFieldName[MED_NAME_SIZE+1] ;
819   int err ;
820   int    numberOfComponents          = 0;
821   char * componentName               = (char *) MED_NULL;
822   char * unitName                    = (char *) MED_NULL;
823   med_2_3::med_field_type type ;
824   med_2_3::med_idt id = MED_FIELD_DRIVER<T>::_medIdt;
825   bool needConversionToDouble = false,needConversionToInt64 = false;
826
827   // we search for the "field med number" of <fieldName>
828   // Having found <fieldName>, variables <numberOfComponents>,
829   // <componentName>, <unitname>, <type> and attribute <_fieldNum> are set.
830   if (MED_FIELD_DRIVER<T>::_fieldNum==MED_INVALID)
831     {
832       int numberOfFields = med_2_3::MEDnField(id) ;
833       if ( numberOfFields <= 0 )
834         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": There is no field found in the file !"));
835
836       for (int i=1;i<=numberOfFields;i++)
837         {
838           numberOfComponents = med_2_3::MEDfieldnComponent(id,i) ;
839
840           if ( numberOfComponents <= 0 )
841             MESSAGE_MED(LOC<<"Be careful there is no components for field n°"<<i<<"in file |"<<MED_FIELD_DRIVER<T>::_fileName<<"| !");
842
843           componentName = new char[numberOfComponents*MED_SNAME_SIZE+1] ;
844           unitName      = new char[numberOfComponents*MED_SNAME_SIZE+1] ;
845           char meshnamep3[MED_NAME_SIZE+1];
846           med_2_3::med_bool localp3;
847           char dtunitp3[MED_LNAME_SIZE+1];
848           med_2_3::med_int nstpp3;
849           err = med_2_3::MEDfieldInfo(id, i, tmpFieldName, meshnamep3,
850                                       &localp3, &type, componentName,
851                                       unitName, dtunitp3,&nstpp3) ;
852
853           MESSAGE_MED("Field "<<i<<" : #" << tmpFieldName <<"# et recherche #"<<fieldName.c_str()<<"#");
854           if ( !strcmp(tmpFieldName,fieldName.c_str()) ) {
855             MESSAGE_MED("FOUND FIELD "<< tmpFieldName <<" : "<<i);
856             MED_FIELD_DRIVER<T>::_fieldNum = i ;
857             break ;
858           }
859           // not found : release memory and search next field !
860           delete[] componentName ;
861           delete[] unitName ;
862         }
863     }
864
865   // Si aucun champ ne correspond les variables <componentName> et <unitName> ont été correctement
866   // désallouées dans la boucle de recherche
867   if (MED_FIELD_DRIVER<T>::_fieldNum==MED_INVALID)
868     throw MEDEXCEPTION(LOCALIZED( STRING(LOC) << ": Field "<<  fieldName
869                                    << " not found in file " << MED_FIELD_DRIVER<T>::_fileName) );
870
871   MESSAGE_MED ("FieldNum : "<<MED_FIELD_DRIVER<T>::_fieldNum);
872
873   if (numberOfComponents < 1) {
874     delete[] componentName; delete[] unitName;
875     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" no component found for field "
876                                  << fieldName)) ;
877   }
878
879   // Verifie que l'on essaye pas de lire un champ double dans un FIELD<int>
880   switch ( (med_2_3::med_field_type) MED_FIELD_DRIVER<T>::_ptrField->_valueType ) {
881   case  med_2_3::MED_INT :
882   case  med_2_3::MED_INT32 :
883   case  med_2_3::MED_INT64 :
884     if ( type == ( med_2_3::MED_FLOAT64 ) ) {
885       delete[] componentName; delete[] unitName;
886       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field Type in file (" << type
887                                    <<") differs from FIELD object type (" <<
888                                    MED_FIELD_DRIVER<T>::_ptrField->_valueType << ")" )) ;
889     }
890 #if defined(IRIX64) || defined(OSF1) ||defined(VPP5000) || defined(PCLINUX64)
891     if (MED_FIELD_DRIVER<T>::_ptrField->_valueType==MED_EN::MED_INT32 )
892       needConversionToInt64=true;
893 #endif
894     break;
895   case med_2_3::MED_FLOAT64 :
896     if (type != med_2_3::MED_FLOAT64)
897       needConversionToDouble=true;
898     break;
899   default:
900     break;
901   }
902
903   string meshName="";
904   const GMESH * ptrMesh = 0;
905   bool   haveSupport = false;
906   bool   haveMesh    = false;
907   MED_EN::medEntityMesh preferEntity = MED_EN::MED_ALL_ENTITIES;
908   if ( MED_FIELD_DRIVER<T>::_ptrField->getSupport() ) {
909     // Verif à faire sur la taille du meshName
910     ptrMesh = MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getMesh();
911     if ( ptrMesh) {
912       meshName =  MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getMesh()->getName() ;
913       haveMesh = true;
914     }
915     haveSupport = true;
916     preferEntity = MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getEntity();
917   }
918   else if ( MED_FIELD_DRIVER<T>::_ptrField->_mesh ) {
919     ptrMesh = MED_FIELD_DRIVER<T>::_ptrField->_mesh;
920     meshName =  ptrMesh->getName() ;
921     haveMesh = true;
922   }
923
924   // Cherche le type d'entité, le nombre d'entité  par type géométrique sur le type d'entité
925   // (MED_MAILLE ou MED_NOEUD uniquement car MEDMEMOIRE ne gère pas la connectivité descendante).
926   // et crée le support correspondant.
927   SUPPORT *   mySupport = new SUPPORT;
928   vector<int> numberOfElementsOfTypeC;
929   vector<int> numberOfGaussPoint;
930   vector<string> gaussModelName, profilName;
931   int         totalNumberOfElWg=0;
932   MED_EN::medEntityMesh fieldMedFileEntity;
933
934   bool found = this->createFieldSupportPart1(id,fieldName,
935                                              MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber,
936                                              MED_FIELD_DRIVER<T>::_ptrField->_orderNumber,
937                                              *mySupport, meshName,
938                                              numberOfElementsOfTypeC,
939                                              numberOfGaussPoint, gaussModelName, profilName,
940                                              totalNumberOfElWg, fieldMedFileEntity, preferEntity);
941
942
943   if ( !found ) {
944     mySupport->removeReference(); delete[] componentName; delete[] unitName;
945     MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ;
946     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"  Can't find any entity for field |"
947                                  << fieldName
948                                  << "| with (it,or) = ("
949                                   << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
950                                  << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), on mesh |"
951                                  << meshName << "|" ));
952   }
953
954
955   //int mesh_dim = MED_FIELD_DRIVER<T>::getMeshDimensionFromFile(id,meshName);
956
957   MED_EN::medEntityMesh entityType = mySupport->getEntity();
958   //Si un SUPPORT était donné, récupère son nom, sa description et
959   //     le pointeur du maillage associé
960   if (! haveSupport)
961     meshName = mySupport->getMeshName();
962   else {
963     // for bug 19782. Entity of support in field was set by med driver and was taken
964     // from the file without any analysis. It can differ from entity the support will
965     // have in MEDMEM.
966 //     if ( entityType != MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getEntity() ) {
967 //       delete mySupport; delete[] componentName; delete[] unitName;
968 //       MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ;
969 //       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Given entity |"
970 //                                 << MED_EN::entNames[MED_FIELD_DRIVER<T>::_ptrField->
971 //                                                     getSupport()->getEntity()]
972 //                                 << "| for field |"
973 //                                 << fieldName
974 //                                 << "| with (it,or) = ("
975 //                                 << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
976 //                                 << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), on mesh "
977 //                                 << meshName << "| differs from found entity |"
978 //                                 << MED_EN::entNames[entityType] << "|."
979 //                                 ));
980 //     }
981     // postpone name setting until profile analisys in order not to set "ON_ALL_entity"
982     // name to a partial support
983     //if ( entityType == MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getEntity() )
984     //mySupport->setName( MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getName() );
985     mySupport->setMesh( MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getMesh() );
986     mySupport->setDescription(MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getDescription());
987   }
988
989   vector< MED_EN::medGeometryElement >  MESHgeoType;
990   vector< int >  MESHnbOfElOfType;
991   vector< int >  MESHnbOfElOfTypeC;
992   if ( haveMesh )
993     this->getMeshGeometricTypeFromMESH(ptrMesh,entityType,MESHgeoType,
994                                        MESHnbOfElOfType,MESHnbOfElOfTypeC);
995   // porting MED3
996   med_2_3::med_int spceDimp3=-1,mdimp3;
997   med_2_3::med_mesh_type mtype;
998   char desccp3[MED_COMMENT_SIZE+1];
999   char dttunittp3[MED_LNAME_SIZE+1];
1000   med_2_3::med_sorting_type sttp3;
1001   med_2_3::med_int nsteppp3;
1002   med_2_3::med_axis_type axxxppp3;
1003 #ifdef WIN32
1004   int naxis=max(3,med_2_3::MEDmeshnAxisByName(id,meshName.c_str()));
1005 #else
1006   int naxis=std::max(3,med_2_3::MEDmeshnAxisByName(id,meshName.c_str()));
1007 #endif
1008   char *annp3=new char[naxis*MED_SNAME_SIZE+1];
1009   char *auup3=new char[naxis*MED_SNAME_SIZE+1];
1010   bool fileHasMesh=(med_2_3::MEDmeshInfoByName(id,meshName.c_str(),&spceDimp3,&mdimp3,&mtype,desccp3,dttunittp3,&sttp3,&nsteppp3,&axxxppp3,annp3,auup3)==0 && spceDimp3 > 0);
1011   delete [] annp3;
1012   delete [] auup3;
1013   // end porting MED3
1014   vector< MED_EN::medGeometryElement >  meshGeoType;
1015   vector< int >  meshNbOfElOfType;
1016   vector< int >  meshNbOfElOfTypeC;
1017   // Si le maillage n'est pas trouvé les tableaux renvoyés sont vides
1018   if (fileHasMesh)
1019   {
1020 //       MED_EN::medEntityMesh entityTypeLoc = entityType;
1021 //       if (entityType == MED_EN::MED_FACE || entityType == MED_EN::MED_EDGE) entityTypeLoc = MED_EN::MED_CELL;
1022
1023     this->getMeshGeometricTypeFromFile(id,meshName,entityType,meshGeoType,
1024                                        meshNbOfElOfType,meshNbOfElOfTypeC);
1025   }
1026
1027   SCRUTE_MED(meshGeoType.size());
1028   SCRUTE_MED(MESHgeoType.size());
1029   SCRUTE_MED(meshNbOfElOfTypeC.size());
1030   SCRUTE_MED(MESHnbOfElOfTypeC.size());
1031
1032   if (meshGeoType.size() != MESHgeoType.size())
1033     {
1034       for (unsigned i = 0; i<meshGeoType.size();i++)
1035         MESSAGE_MED("debug meshGeotype " << meshGeoType[i]);
1036
1037       for (unsigned i = 0; i<MESHgeoType.size();i++)
1038         MESSAGE_MED("debug MESHgeoType. " << MESHgeoType[i]);
1039     }
1040
1041   if (meshNbOfElOfTypeC.size() == MESHnbOfElOfTypeC.size())
1042     {
1043       for (unsigned i = 0; i<meshNbOfElOfTypeC.size();i++)
1044         MESSAGE_MED("debug meshNbOfElOfTypeC " << meshNbOfElOfTypeC[i]);
1045
1046       for (unsigned i = 0; i<MESHnbOfElOfTypeC.size();i++)
1047         MESSAGE_MED("debug MESHnbOfElOfTypeC " << MESHnbOfElOfTypeC[i]);
1048     }
1049
1050   if (fileHasMesh && haveSupport )
1051     if ( ( meshGeoType != MESHgeoType ) || (meshNbOfElOfTypeC != MESHnbOfElOfTypeC) )
1052       {
1053         MESSAGE_MED("Warning MedField driver 21 while getting mesh information from file for FIELD "<< fieldName
1054                 << " on entity " << MED_EN::entNames[entityType]
1055                 << " with (it,or) = ("
1056                 << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1057                 << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1058                 << " on mesh " << meshName
1059                 << " : geometric types or number of elements by type differs from MESH object !");
1060
1061 //      throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information from file for FIELD "<< fieldName
1062 //                                    << " on entity " << MED_EN::entNames[entityType]
1063 //                                    << " with (it,or) = ("
1064 //                                    << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1065 //                                    << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1066 //                                    << " on mesh " << meshName
1067 //                                    << " : geometric types or number of elements by type differs from MESH object !"
1068 //                                    )
1069 //                         );
1070       }
1071
1072   if ( !fileHasMesh && !haveSupport )
1073     throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information for FIELD "<< fieldName
1074                                   << " on entity " << MED_EN::entNames[entityType]
1075                                   << " with (it,or) = ("
1076                                   << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1077                                   << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1078                                   << " on mesh " << meshName
1079                                   << " : SUPPORT must contain a valid MESH reference or file must contain the associated MESH."
1080                                   )
1081                        );
1082
1083
1084   if (!fileHasMesh && haveSupport) {
1085     meshNbOfElOfTypeC = MESHnbOfElOfTypeC;
1086     meshGeoType       = MESHgeoType;
1087     meshNbOfElOfType  = MESHnbOfElOfType;
1088   }
1089
1090   // Test si le Support du Champ repose ou non sur toutes les entités géométriques
1091   // du maillage associé et positionne ou non l'attribut onAll du SUPPORT.
1092   // Il ne s'agit pas de la gestion des profils
1093   vector < MED_EN::medGeometryElement > v1(  mySupport->getTypes(),
1094                                              mySupport->getTypes()+mySupport->getNumberOfTypes() );
1095   vector<int> v2(numberOfElementsOfTypeC.size());
1096   transform(numberOfElementsOfTypeC.begin(),
1097             numberOfElementsOfTypeC.end(),v2.begin(), bind2nd(plus<int>(),1));
1098
1099   if ( ( meshGeoType != v1 )  || meshNbOfElOfTypeC != v2 ) {
1100     // ATTENTION : mySupport->setAll(false);
1101     // Pb : On a envie de positionner onAll à faux si le champ n'est pas défini sur tous les
1102     //      types géométriques du maillage associé.
1103     //      Mais si onAll est false et si aucun profil n'est détecté par la suite,
1104     //      l'attribut SUPPORT->_number est censé être positionné quand même ! Que faire ?
1105     // Si on veut être compatible avec la signification première de onAll,
1106     //  il faudrait créer des profils contenant toutes les entités pour chaque type géométrique
1107     //  du SUPPORT  mais d'une part c'est dommage d'un point de vue de l'encombrement mémoire
1108     //  et d'autre part, à la réécriture du fichier MED on stockera des profils 
1109     //  alors qu'il n'y en avait pas à l'origine (fichier MED différent après lecture/écriture) !
1110     // Si on laisse setAll à vrai il faut être sûr que les utilisateurs prennent les
1111     //  informations sur les types gémétrique au niveau du support et non pas du maillage.
1112     // Solution : Signification du onAll -> onAllElements des type géométriques définis
1113     // dans le SUPPORT et non du maillage associé (dans la plupart des cas si le fichier ne
1114     // contient pas de profil, le champ est défini sur toutes les entités de tous les types
1115     // géométriques définis dans le maillage).
1116   }
1117
1118
1119   // If an error occurs while reading the field, these allocated FIELD member will be deleted
1120
1121   MED_FIELD_DRIVER<T>::_ptrField->_name                   = healName( fieldName );
1122   MED_FIELD_DRIVER<T>::_ptrField->_numberOfComponents     = numberOfComponents ;
1123   //MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes        = new int   [numberOfComponents] ;
1124   //MED_FIELD_DRIVER<T>::_ptrField->_componentsNames        = new string[numberOfComponents] ;
1125   //MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits        = new UNIT  [numberOfComponents] ;
1126   //MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions = new string[numberOfComponents] ;
1127   //MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits     = new string[numberOfComponents] ;
1128   MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes.resize(numberOfComponents);
1129   MED_FIELD_DRIVER<T>::_ptrField->_componentsNames.resize(numberOfComponents);
1130   MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits.resize(numberOfComponents);
1131   MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions.resize(numberOfComponents);
1132   MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits.resize(numberOfComponents);
1133   for (int i=0; i<numberOfComponents; i++) {
1134       MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes[i]    = 1 ;
1135       MED_FIELD_DRIVER<T>::_ptrField->_componentsNames[i]    = string(componentName+i*MED_SNAME_SIZE,MED_SNAME_SIZE) ;
1136       MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits[i] = string(unitName+i*MED_SNAME_SIZE,MED_SNAME_SIZE) ;
1137       SCRUTE_MED(MED_FIELD_DRIVER<T>::_ptrField->_componentsNames[i]);
1138       SCRUTE_MED(MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits[i]);
1139   }
1140   delete[] componentName;
1141   delete[] unitName;
1142
1143   int NumberOfTypes                       = mySupport->getNumberOfTypes() ;
1144   const MED_EN::medGeometryElement *types = mySupport->getTypes() ;
1145   T * myValues = new T[totalNumberOfElWg*numberOfComponents];
1146   const int * nbOfElOfType = mySupport->getNumberOfElements() ;
1147   bool anyProfil = false;
1148   int  pflSize=0,index=0;
1149   // Le vecteur de profil est dimensionné par rapport aux nombres de types
1150   // géométriques du champ même si le champ n'a pas de profil MED FICHIER sur
1151   // tous ses types géométriques car dans MEDMEM si onAllElement 
1152   // du SUPPORT est false il faut positionner un profil pour tous les types géométriques 
1153   // du SUPPORT
1154   int profilSizeC = 0;
1155   vector < int   >                     profilSize    (NumberOfTypes,0);
1156   vector < string >                    profilNameList(NumberOfTypes);
1157   vector < vector<med_2_3::med_int>  > profilList    (NumberOfTypes);      // IPAL13481
1158   vector < vector<med_2_3::med_int>  > profilListFromFile (NumberOfTypes); // IPAL13481
1159   //char *                               profilName = new char[MED_NAME_SIZE+1];
1160
1161   MESSAGE_MED ("NumberOfTypes      : "<< NumberOfTypes);
1162   MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues=0 ;
1163  
1164   // PAL16681 (Read no interlace field from file) ->
1165   // use medModeSwitch of a field in MEDchampLire() if there is one geometric type
1166   // to exclude array conversion
1167   med_2_3::med_switch_mode modswt = med_2_3::MED_FULL_INTERLACE;
1168   // NOTE: field can be either of 3 medModeSwitch'es, MED_NO_INTERLACE_BY_TYPE added (PAL17011)
1169   if ( ( NumberOfTypes == 1 && !isFullInterlace) || isNoInterlaceByType )
1170     modswt = med_2_3::MED_NO_INTERLACE;
1171
1172   for (int typeNo=0; typeNo<NumberOfTypes; typeNo++) {
1173
1174     int numberOfValuesWc= nbOfElOfType[typeNo]*numberOfGaussPoint[typeNo+1]*numberOfComponents;
1175     //char * gaussModelName = new char[MED_NAME_SIZE+1];
1176
1177     MESSAGE_MED ("FIELD_NAME         : "<< fieldName.c_str());
1178     MESSAGE_MED ("MESH_NAME          : "<< meshName.c_str());
1179     MESSAGE_MED ("MED_ENTITE         : "<< MED_EN::entNames[entityType]);
1180     MESSAGE_MED ("MED_GEOM           : "<< MED_EN::geoNames[types[typeNo]]);
1181     MESSAGE_MED ("Iteration          : "<< MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber());
1182     MESSAGE_MED ("Order              : "<< MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber());
1183     MESSAGE_MED ("Time               : "<< MED_FIELD_DRIVER<T>::_ptrField->getTime());
1184     MESSAGE_MED ("NumberOfElements   : "<< nbOfElOfType[typeNo]);
1185     MESSAGE_MED ("NumberOfComponents : "<< numberOfComponents);
1186     MESSAGE_MED ("NumberOfGaussPts   : "<< numberOfGaussPoint[typeNo+1]);
1187     MESSAGE_MED ("NumberOfValuesWg   : "<< nbOfElOfType[typeNo]*numberOfGaussPoint[typeNo+1]);
1188     MESSAGE_MED ("NumberOfValuesWgWc : "<< numberOfValuesWc);
1189     MESSAGE_MED ("Index              : "<< index);
1190     med_2_3::med_err ret=-1;
1191
1192     med_2_3::med_int * myValuesTmp=0;
1193     unsigned char* ptrTmp=0;
1194     if (needConversionToDouble || needConversionToInt64 ) {
1195       myValuesTmp = new med_2_3::med_int[numberOfValuesWc];
1196       ptrTmp = (unsigned char*) myValuesTmp;
1197     } else
1198       ptrTmp = (unsigned char*) &myValues[index];
1199
1200     //VERIFIER LE NBRE
1201 //     med_2_3::med_entite_maillage medfile_entity;
1202 //     if (entityType==MED_NODE) 
1203 //       medfile_entity= (med_2_3::med_entite_maillage)MED_NODE;
1204 //     else 
1205 //       medfile_entity= (med_2_3::med_entite_maillage)MED_CELL;
1206     ret=med_2_3::MEDfieldValueWithProfileRd(id,fieldName.c_str(),
1207                                             MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber(),
1208                                             MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber(),
1209                                             (med_2_3::med_entity_type) fieldMedFileEntity,
1210                                             (med_2_3::med_geometry_type)types[typeNo],
1211                                             med_2_3::MED_COMPACT_PFLMODE,
1212                                             profilName[typeNo].c_str(),modswt,MED_ALL_CONSTITUENT,
1213                                             (unsigned char*) ptrTmp);
1214
1215       if (needConversionToDouble || needConversionToInt64 ) {
1216
1217       if (needConversionToInt64 ) //utiliser un trait
1218         for(int i=0;i<numberOfValuesWc;++i)
1219           myValues[index+i]=(int)(myValuesTmp[i]);
1220       else
1221         for(int i=0;i<numberOfValuesWc;++i)
1222           myValues[index+i]=myValuesTmp[i];
1223       delete[] myValuesTmp;
1224     }
1225
1226     if (ret < 0)
1227       {
1228         // The Field can't be read then we must delete all previously allocated members in FIELD
1229         //for(int j=0; j<=i;j++)
1230         //  delete[] myValues[j];
1231         delete[] myValues;
1232         //delete[] NumberOfValues ;
1233         //delete[] profilName;
1234         //delete[] gaussModelName;
1235         //delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes ;
1236         //delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsNames ;
1237         //delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits ;
1238         //delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions ;
1239         //delete[] MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits ;
1240         //MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes = NULL ;
1241         //MED_FIELD_DRIVER<T>::_ptrField->_componentsNames = NULL ;
1242         //MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits = NULL ;
1243         //MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions = NULL ;
1244         //MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits = NULL ;
1245         MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes.clear();
1246         MED_FIELD_DRIVER<T>::_ptrField->_componentsNames.clear();
1247         MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits.clear();
1248         MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions.clear();
1249         MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits.clear();
1250         MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ; // we have not found right field, so reset the field number
1251         throw MEDEXCEPTION( LOCALIZED( STRING(LOC) <<": ERROR while reading values")) ;
1252       }
1253
1254     index += numberOfValuesWc;
1255     // Le support prend en compte le nombre de valeurs lié aux profils
1256     MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues+=
1257       nbOfElOfType[typeNo];// Ne doit pas prendre en compte les points de Gauss
1258
1259     // second et troisième test lié à un bug medfichier
1260     if ( gaussModelName[typeNo][0] )
1261       {
1262         int type_geo = (int) types[typeNo];
1263         int t1       = (type_geo%100)*(type_geo/100);
1264         int ngauss   = numberOfGaussPoint[typeNo+1];
1265         int t2       = ngauss*(type_geo/100);
1266         med_2_3::med_float * refcoo = new med_2_3::med_float[t1];
1267         med_2_3::med_float * gscoo  = new med_2_3::med_float[t2];
1268         med_2_3::med_float * wg     = new med_2_3::med_float[ngauss];
1269
1270         if (med_2_3::MEDlocalizationRd(id, gaussModelName[typeNo].c_str(), modswt, refcoo, gscoo, wg ) < 0)
1271           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while reading Gauss Model |"
1272                                       << gaussModelName[typeNo] << "| for FIELD "<< fieldName
1273                                       << " on geometric type " << MED_EN::geoNames[types[typeNo]]
1274                                       )
1275                            );
1276         if (isFullInterlace ) { //serait inutile avec un driver template du type d'entrelacement
1277           GAUSS_LOCALIZATION<FullInterlace> * loc;
1278           loc = new GAUSS_LOCALIZATION<FullInterlace>(gaussModelName[typeNo],types[typeNo],ngauss, refcoo,gscoo, wg);
1279           MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]]=loc;
1280         } else {
1281           GAUSS_LOCALIZATION<NoInterlace> * loc;
1282           loc = new GAUSS_LOCALIZATION<NoInterlace>(gaussModelName[typeNo],types[typeNo],ngauss, refcoo,gscoo, wg);
1283           MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]]=loc;
1284         }
1285 //      cout << *MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]] << endl;
1286         delete [] refcoo;delete [] gscoo; delete [] wg;
1287
1288     }
1289     //delete[] gaussModelName ;
1290
1291     if ( profilName[typeNo] != MED_NO_PROFILE ) {
1292       anyProfil = true;
1293       pflSize = med_2_3::MEDprofileSizeByName(id,profilName[typeNo].c_str());
1294       if ( pflSize  <= 0)
1295         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Error while reading the size of profil |"
1296                                      << profilName[typeNo] << "|" ));
1297
1298       profilSize[typeNo]=pflSize;
1299       profilList[typeNo].resize(pflSize);
1300       profilListFromFile[typeNo].resize(pflSize);
1301       ret = med_2_3::MEDprofileRd(id,profilName[typeNo].c_str(),&profilList[typeNo][0]); // cf item 16 Effective STL // IPAL13481
1302       profilListFromFile[typeNo] = profilList[typeNo];
1303       profilNameList[typeNo]= healName(profilName[typeNo]);
1304     }
1305   }
1306
1307   //delete[] profilName;
1308
1309   //MESSAGE_MED ("Index              : "<< index);
1310   assert(index == totalNumberOfElWg*numberOfComponents);
1311   assert(MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues == mySupport->getNumberOfElements(MED_ALL_ELEMENTS));
1312
1313   if (anyProfil)
1314   {
1315     for (int typeNo=0; typeNo < NumberOfTypes; typeNo++)
1316       {
1317         MED_EN::medGeometryElement geomType = types[typeNo];
1318
1319         // Trouve l'index du type géométrique dans la liste des types géométriques du maillage
1320         // correspondant au type géométrique du champ traité
1321         vector<MED_EN::medGeometryElement>::iterator meshTypeNoIt =
1322           find(meshGeoType.begin(),meshGeoType.end(),geomType); //Gérer l'exception
1323         if (meshTypeNoIt ==  meshGeoType.end())
1324           throw MEDEXCEPTION(LOCALIZED(STRING(LOC) <<": Can't find "<< MED_EN::geoNames[geomType]
1325                                        << " on entity " << MED_EN::entNames[entityType]
1326                                        << " in geometric type list of mesh " << meshName));
1327         int meshTypeNo = meshTypeNoIt - meshGeoType.begin();
1328
1329         if (! profilList[typeNo].empty() )
1330           {
1331             //      for (int j =0 ; j< meshGeoType.size();++j)
1332             //        cout << "--MeshTypeNo : "<<meshTypeNo<<"-> meshNbOfElOfTypeC["<<j<<"]="<<meshNbOfElOfTypeC[j]<<endl;
1333             //      cout << "--typeNo--" << typeNo << endl;
1334             //      cout << "meshNbOfElOfTypeC["<<meshTypeNo<<"]=" << meshNbOfElOfTypeC[meshTypeNo] <<endl;
1335
1336             // Transformer les numéros locaux d'entités medfichier en numéro global medmémoire
1337             for (unsigned i = 0; i < profilList[typeNo].size(); i++) {
1338               // Les numéros des entités commencent à 1 dans MEDfichier comme dans MEDmémoire
1339               // meshNbOfElOfTypeC[0]=0 ...meshNbOfEltOfTypeC[meshTypeNo]=
1340               // meshNbOfElOfTypeC[meshTypeNo-1]+nbrOfElem of meshTypeNo type
1341               // rem1 : Si le meshTypeNo trouvé est 0 (premier type géométrique du maillage
1342               // il ne faut pas décaler les numéros du profils qui commencent à 1 dans MEDFICHIER
1343               // rem2 : meshNbOfElOfTypeC[NumberOfTypes] ne devrait jamais être utilisé
1344               profilList[typeNo][i]+=meshNbOfElOfTypeC[meshTypeNo];
1345             }
1346           } else {
1347           // Créer le profil <MED_ALL> pour ce type géométrique
1348           // uniquement pour renseigner le tableau skyline avec des accesseurs directs
1349           // par type géométriques
1350           // REM : Une conséquence est qu'à la réecriture le fichier contiendra des
1351           // profils sur certains types géométriques alors qu'à la lecture il n'y en avait pas !
1352           // Solution : Stocker les noms des profils et les utiliser pour savoir si il y avait ou non
1353           //            un profil
1354           int pflSize = meshNbOfElOfType[meshTypeNo];
1355           // profil    = new int[pflSize];
1356
1357           profilList[typeNo].resize(pflSize);
1358           profilSize[typeNo] = pflSize;
1359
1360           for (int j = 1; j <= pflSize; j++) {
1361             profilList[typeNo][j-1] = meshNbOfElOfTypeC[meshTypeNo] + j ; // index MEDMEM commence à 1
1362           }
1363           profilNameList[typeNo] = MED_NO_PROFILE; //Information a utiliser pour la sauvegarde : PLUTOT MED_ALL
1364         }
1365         profilSizeC += profilList[typeNo].size();
1366       }
1367
1368     MEDSKYLINEARRAY * skyLine = new MEDSKYLINEARRAY(profilList.size(), profilSizeC );
1369     vector<int> index(NumberOfTypes+1,0);
1370     index[0]=1;
1371     for( int typeNo=0; typeNo < NumberOfTypes; typeNo++ )
1372       index[typeNo+1]=index[typeNo]+profilSize[typeNo];
1373     skyLine->setIndex(&index[0]);
1374     for (unsigned i=1; i <= profilList.size() ; i++) {
1375       vector<int> aTmp(profilList[i-1].size()); // IPAL13481
1376       for (unsigned j=0; j < profilList[i-1].size(); j++)
1377         aTmp[j] = (int) profilList[i-1][j];
1378       skyLine->setI(i,&aTmp[0]);
1379       //skyLine->setI(i,&profilList[i-1][0]);
1380     }
1381
1382     MEDSKYLINEARRAY * skyLineFromFile = new MEDSKYLINEARRAY(profilListFromFile.size(), profilSizeC );
1383     skyLineFromFile->setIndex(&index[0]);
1384     for (int i=1; i <= (int)profilListFromFile.size() ; i++) {
1385       vector<int> aTmp(profilListFromFile[i-1].size()); // IPAL13481
1386       for (unsigned j=0; j < profilListFromFile[i-1].size(); j++)
1387         aTmp[j] = (int) profilListFromFile[i-1][j];
1388       skyLineFromFile->setI(i,&aTmp[0]);
1389       //skyLineFromFile->setI(i,&profilListFromFile[i-1][0]);
1390     }
1391
1392     mySupport->setAll(false);
1393     mySupport->setpartial(skyLine,true);
1394     mySupport->setpartial_fromfile(skyLineFromFile,true);
1395     mySupport->setProfilNames(profilNameList);
1396
1397     // update pointers
1398     NumberOfTypes = mySupport->getNumberOfTypes();
1399     types         = mySupport->getTypes();
1400     nbOfElOfType  = mySupport->getNumberOfElements();
1401 //    cout << "Valeurs du skyline du SUPPORT partiel crée : " << *skyLine << endl;
1402   }
1403
1404   // Créer un driver spécifique pour les modes MED_FULL_INTERLACE et MED_NO_INTERLACE
1405   // serait plus efficace.
1406   bool anyGauss = (numberOfGaussPoint != vector<int>(numberOfGaussPoint.size(),1));
1407   SCRUTE_MED(anyGauss);
1408   MEDMEM_Array_ * Values;
1409   if (anyGauss) {
1410     SCRUTE_MED(mySupport->getNumberOfElements(MED_ALL_ELEMENTS) );
1411     SCRUTE_MED(NumberOfTypes);
1412     SCRUTE_MED(numberOfElementsOfTypeC[NumberOfTypes]-1);
1413     assert(mySupport->getNumberOfElements(MED_ALL_ELEMENTS) == (numberOfElementsOfTypeC[NumberOfTypes]-1) );
1414     // PAL16681. If NumberOfTypes == 1 then myValues is what should be
1415     // in a field value, inspite of InterlacingType
1416     if ( NumberOfTypes == 1 && modswt == med_2_3::MED_NO_INTERLACE )
1417       Values = new ArrayNoWg(myValues,
1418                              numberOfComponents,
1419                              numberOfElementsOfTypeC[NumberOfTypes]-1,
1420                              NumberOfTypes,
1421                              &numberOfElementsOfTypeC[0],
1422                              &numberOfGaussPoint[0],
1423                              true,true);
1424     else if ( isNoInterlaceByType ) // PAL17011 (MEDMEM : no_interlace_by_type fields)
1425       Values = new ArrayByTypeWg(myValues,
1426                                  numberOfComponents,
1427                                  numberOfElementsOfTypeC[NumberOfTypes]-1,
1428                                  NumberOfTypes,
1429                                  &numberOfElementsOfTypeC[0],
1430                                  &numberOfGaussPoint[0],
1431                                  true,true);
1432     else
1433       Values = new ArrayFullWg(myValues,
1434                                numberOfComponents,
1435                                numberOfElementsOfTypeC[NumberOfTypes]-1,
1436                                // Up : Prend en compte les profils et
1437                                // Ne prend pas en compte le nbre de composantes et
1438                                // le nombre de points de Gauss
1439                                NumberOfTypes,
1440                                &numberOfElementsOfTypeC[0],
1441                                &numberOfGaussPoint[0],
1442                                true,true);
1443 //     cout << "Valeurs du ArrayFullWg crée : " << endl <<
1444 //       *(static_cast<ArrayFullWg*>(Values))  << endl;
1445   }
1446   else {
1447     // PAL16681. If NumberOfTypes == 1 then myValues is what should be
1448     // in a field value, inspite of InterlacingType
1449     if ( NumberOfTypes == 1 && interlacingType == MED_EN::MED_NO_INTERLACE )
1450       Values = new ArrayNo(myValues,numberOfComponents,totalNumberOfElWg,
1451                            true,true);
1452     else if ( isNoInterlaceByType ) // PAL17011 (MEDMEM : no_interlace_by_type fields)
1453       Values = new ArrayByType(myValues,numberOfComponents,totalNumberOfElWg,
1454                                NumberOfTypes, &numberOfElementsOfTypeC[0], true,true);
1455     else
1456       Values = new ArrayFull(myValues,numberOfComponents,totalNumberOfElWg,
1457                              true,true);
1458   }
1459   if (MED_FIELD_DRIVER<T>::_ptrField->_value != NULL)
1460     delete MED_FIELD_DRIVER<T>::_ptrField->_value;
1461
1462   if ( NumberOfTypes != 1 &&  // PAL16681
1463        interlacingType == MED_EN::MED_NO_INTERLACE )
1464   {
1465     // Convert MED_FULL_INTERLACE -> MED_NO_INTERLACE
1466     if (Values->getGaussPresence())
1467       MED_FIELD_DRIVER<T>::_ptrField->_value=ArrayConvert(*static_cast<ArrayFullWg*>(Values));
1468     else
1469       MED_FIELD_DRIVER<T>::_ptrField->_value=ArrayConvert(*static_cast<ArrayFull*  >(Values));
1470     delete Values;
1471   }
1472   else
1473   {
1474     MED_FIELD_DRIVER<T>::_ptrField->_value=Values;
1475   }
1476
1477   MED_FIELD_DRIVER<T>::_ptrField->_isRead = true ;
1478
1479   const GMESH* aMesh = MED_FIELD_DRIVER<T>::_ptrField->_mesh;
1480   bool isFound = false, onlyMeshProvided = ( !haveSupport && aMesh );
1481   if ( !aMesh )
1482     aMesh = ptrMesh; // try to find FAMILY OR GROUP in ptrMesh but without being so stern
1483   if ( aMesh && anyProfil)
1484   {
1485     int it = -1;
1486     for (int typeNo = 0; (typeNo < NumberOfTypes) && (it == -1); typeNo++) {
1487       if (strcmp(profilNameList[typeNo].c_str(), MED_NO_PROFILE) != 0)
1488         it = typeNo;
1489     }
1490     // IMP 0019953: link between fields and families for MED 2.2 read driver
1491     string aPN = profilNameList[it];
1492     MED_EN::medGeometryElement aPT = types[it];
1493
1494     ostringstream typestr;
1495     typestr << "_type" << aPT;
1496     string aSuff = typestr.str();
1497
1498     //- If the field profile name is toto_PFL and a family toto exists,
1499     //  the field will point to the corresponding FAMILY object.
1500     const vector<FAMILY*> aFams = aMesh->getFamilies(entityType);
1501     for (unsigned fi = 0; fi < aFams.size() && !isFound; fi++) {
1502       FAMILY* aF = aFams[fi];
1503       string aFN_suff = aF->getName() + aSuff;
1504       if (aPN == aFN_suff) {
1505         isFound = true;
1506         //family found
1507         if(MED_FIELD_DRIVER<T>::_ptrField->_support)
1508           MED_FIELD_DRIVER<T>::_ptrField->_support->removeReference();
1509         MED_FIELD_DRIVER<T>::_ptrField->_support = aF; //Prévenir l'utilisateur ?
1510         aF->addReference();
1511       }
1512     }
1513     if (!isFound) {
1514       // - If no family was found, lookup the groups and if a group toto
1515       //   exists, the field will point to the corresponding GROUP object.
1516       const vector<GROUP*> aGrps = aMesh->getGroups(entityType);
1517       for (unsigned gi = 0; gi < aGrps.size() && !isFound; gi++) {
1518         GROUP* aG = aGrps[gi];
1519         string aGN_suff = aG->getName() + aSuff;
1520         if (aPN == aGN_suff) {
1521           isFound = true;
1522           //group found
1523           if(MED_FIELD_DRIVER<T>::_ptrField->_support)
1524             MED_FIELD_DRIVER<T>::_ptrField->_support->removeReference();
1525           MED_FIELD_DRIVER<T>::_ptrField->_support = aG; //Prévenir l'utilisateur ?
1526           aG->addReference();
1527         }
1528       }
1529     }
1530     if (!isFound) {
1531       // - If no family or group was found and the
1532       //   profile name is xxx_PFL, throw an exception
1533       int pos = aPN.rfind(aSuff);
1534       if (pos + aSuff.length() - 1 == aPN.length() && onlyMeshProvided )
1535         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1536                                      << ": Can't find appropriate support (GROUP or FAMILY)"
1537                                      << " in mesh " << meshName << " for field " << fieldName
1538                                      << ", while one of its profiles " << aPN
1539                                      << " was generated from a FAMILY or a GROUP"));
1540     }
1541     else {
1542       // - Check that the found support has correct types
1543       //   and number of elements. If not, throw an exception
1544       const SUPPORT* aSupp = MED_FIELD_DRIVER<T>::_ptrField->_support;
1545       isFound = ( aSupp->getNumberOfTypes() == NumberOfTypes );
1546       if ( !isFound && onlyMeshProvided )
1547         throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << ": Invalid support (GROUP or FAMILY) found in mesh "
1548                                      << meshName << " for field " << fieldName << " by name of profile "
1549                                      << aPN << ": different number of types in found support |"
1550                                      << aSupp->getNumberOfTypes() << "| and in required |"
1551                                      << NumberOfTypes << "|"));
1552
1553       const MED_EN::medGeometryElement* aTypes = aSupp->getTypes();
1554       for (int it = 0; it < NumberOfTypes && isFound; it++)
1555       {
1556         MED_EN::medGeometryElement aType = aTypes[it];
1557         isFound = ( aType == types[it] );
1558         if ( !isFound && onlyMeshProvided )
1559           throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << ": Invalid support (GROUP or FAMILY) found in mesh "
1560                                        << meshName << " for field " << fieldName << " by name of profile "
1561                                        << aPN << ": geometric type in found support |" << aType
1562                                        << "| differs from required type |" << types[it] << "|"));
1563
1564         isFound = ( isFound && aSupp->getNumberOfElements(aType) == nbOfElOfType[it] );
1565         if ( !isFound && onlyMeshProvided )
1566           throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << ": Invalid support (GROUP or FAMILY) found in mesh "
1567                                        << meshName << " for field " << fieldName << " by name of profile "
1568                                        << aPN << ": number of elements of type " << aType
1569                                        << " in found support |" << aSupp->getNumberOfElements(aType)
1570                                        << "| differs from required |" << nbOfElOfType[it] << "|"));
1571       }
1572     }
1573   }
1574
1575   if (!isFound)
1576   {
1577     // No corresponding support (family or group)
1578     // found in the mesh, use the newly created one
1579     const MEDMEM::SUPPORT * & fieldSupport = MED_FIELD_DRIVER<T>::_ptrField->_support;
1580     if ( fieldSupport )
1581     {
1582       if ( mySupport->getEntity() == fieldSupport->getEntity() &&
1583            //not to set "ON_ALL_entity" name to a partial support
1584            mySupport->isOnAllElements() == fieldSupport->isOnAllElements())
1585         mySupport->setName( fieldSupport->getName() );
1586
1587       fieldSupport->removeReference();
1588     }
1589     fieldSupport = mySupport; //Prévenir l'utilisateur ?
1590     if ( !fieldSupport->getMesh() && aMesh )
1591       fieldSupport->setMesh( aMesh );
1592
1593     // use one support instead of keeping many equal supports
1594     if ( aMesh )
1595     {
1596       const SUPPORT* supOnAll = aMesh->getSupportOnAll( fieldSupport->getEntity() );
1597       if ( fieldSupport->deepCompare( *supOnAll ))
1598         MED_FIELD_DRIVER<T>::_ptrField->setSupport( supOnAll );
1599     }
1600   }
1601   else {
1602     mySupport->removeReference();
1603   }
1604   END_OF_MED(LOC);
1605 }
1606
1607 template <class T> void MED_FIELD_RDONLY_DRIVER<T>::write( void ) const
1608   throw (MEDEXCEPTION)
1609 {
1610   throw MEDEXCEPTION("MED_FIELD_RDONLY_DRIVER::write : Can't write with a RDONLY driver !");
1611 }
1612
1613 /*--------------------- WRONLY PART -------------------------------*/
1614
1615 /*!
1616   Constructor.
1617 */
1618 template <class T>
1619 MED_FIELD_WRONLY_DRIVER<T>::MED_FIELD_WRONLY_DRIVER():MED_FIELD_DRIVER<T>()
1620 {
1621   this->GENDRIVER::_accessMode = MED_EN::WRONLY;
1622 }
1623
1624 /*!
1625   Constructor.
1626 */
1627 template <class T>
1628 template <class INTERLACING_TAG>
1629 MED_FIELD_WRONLY_DRIVER<T>::MED_FIELD_WRONLY_DRIVER(const string &              fileName,
1630                                                     FIELD<T, INTERLACING_TAG> * ptrField):
1631   MED_FIELD_DRIVER<T>(fileName,ptrField,MED_EN::WRONLY)
1632 {
1633   const char* LOC = "MED_FIELD_WRONLY_DRIVER::MED_FIELD_WRONLY_DRIVER(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)";
1634   BEGIN_OF_MED(LOC);
1635   END_OF_MED(LOC);
1636 }
1637
1638 /*!
1639   Copy constructor.
1640 */
1641 template <class T>
1642 MED_FIELD_WRONLY_DRIVER<T>::MED_FIELD_WRONLY_DRIVER(const MED_FIELD_WRONLY_DRIVER & fieldDriver):
1643   MED_FIELD_DRIVER<T>(fieldDriver)
1644 {}
1645
1646 /*!
1647   Destructor.
1648 */
1649 template <class T>
1650 MED_FIELD_WRONLY_DRIVER<T>::~MED_FIELD_WRONLY_DRIVER() {}
1651
1652 template <class T> GENDRIVER * MED_FIELD_WRONLY_DRIVER<T>::copy(void) const
1653 {
1654   return new MED_FIELD_WRONLY_DRIVER<T>(*this);
1655 }
1656
1657 template <class T> void MED_FIELD_WRONLY_DRIVER<T>::read (void)
1658   throw (MEDEXCEPTION)
1659 {
1660   throw MEDEXCEPTION("MED_FIELD_WRONLY_DRIVER::read : Can't read with a WRONLY driver !");
1661 }
1662
1663 template <class T> void MED_FIELD_WRONLY_DRIVER<T>::write(void) const
1664   throw (MEDEXCEPTION)
1665 {
1666   const char * LOC = "MED_FIELD_WRONLY_DRIVER::write(void) const " ;
1667   BEGIN_OF_MED(LOC);
1668   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,NoGauss>::Array   ArrayNo;
1669   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,Gauss>::Array     ArrayNoWg;
1670   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,NoGauss>::Array ArrayFull;
1671   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,Gauss>::Array   ArrayFullWg;
1672
1673   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION<FullInterlace>*> locMapFull;
1674   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION<NoInterlace>*>   locMapNo;
1675   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION_*>   locMap;
1676
1677   med_2_3::med_idt id = MED_FIELD_DRIVER<T>::_medIdt;
1678
1679   if (MED_FIELD_DRIVER<T>::_status!=MED_OPENED)
1680     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": Method open must be called before method write.")) ;
1681
1682   string fieldName;
1683   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
1684        ( MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
1685     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1686                                  <<" neither <fieldName> is set in driver nor in object FIELD.")) ;
1687
1688   // If _fieldName is not set in driver, try to use _ptrfield->_fieldName
1689   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
1690        ( !MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
1691     fieldName = healName( MED_FIELD_DRIVER<T>::_ptrField->_name );
1692   else
1693     fieldName = healName( MED_FIELD_DRIVER<T>::_fieldName );
1694
1695   //if ( ! MED_FIELD_DRIVER<T>::_ptrField->_isRead )
1696   //  throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1697   //                     <<" FIELD |"<<fieldName<<"| was not read but is being written"));
1698
1699   SCRUTE_MED(fieldName);
1700   if ( fieldName.size() > MED_NAME_SIZE ) {
1701     fieldName = fieldName.substr(0,MED_NAME_SIZE);
1702     MESSAGE_MED( "Be careful <fieldName> size must not be > MED_NAME_SIZE, using fieldName : |"<< fieldName <<"|." );
1703   }
1704
1705   const SUPPORT * mySupport = MED_FIELD_DRIVER<T>::_ptrField->getSupport() ;
1706   if ( ! mySupport )
1707     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1708                                  <<" There is no SUPPORT associated with FIELD : "
1709                                  << fieldName << "."));
1710
1711   bool onAll = mySupport->isOnAllElements();
1712   const locMap & gaussModel = MED_FIELD_DRIVER<T>::_ptrField->_gaussModel;
1713
1714
1715   string meshName = healName( mySupport->getMeshName() );
1716   SCRUTE_MED(meshName);
1717   if ( meshName.size() > MED_NAME_SIZE ) {
1718     meshName = meshName.substr(0,MED_NAME_SIZE);
1719     MESSAGE_MED( "Be careful <meshName> size must not be > MED_NAME_SIZE, using meshName : |"<< meshName <<"|." );
1720   }
1721   MED_EN::medEntityMesh entityType = mySupport->getEntity();
1722
1723   // Reconstruit les listes contigues des noms de composantes et des unités
1724   // Les noms sont tronqués à MED_SNAME_SIZE
1725   int err ;
1726   int component_count=MED_FIELD_DRIVER<T>::_ptrField->getNumberOfComponents();
1727   string   component_name(component_count*MED_SNAME_SIZE,' ') ;
1728   string   component_unit(component_count*MED_SNAME_SIZE,' ') ;
1729
1730   const string * listcomponent_name=MED_FIELD_DRIVER<T>::_ptrField->getComponentsNames() ;
1731   const string * listcomponent_unit=MED_FIELD_DRIVER<T>::_ptrField->getMEDComponentsUnits() ;
1732   if ( ! listcomponent_name || ! listcomponent_unit )
1733     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) <<" Udefined components of FIELD : "
1734                                  << fieldName << "."));
1735   int length ;
1736   for (int i=0; i < component_count ; i++) {
1737     length = min(MED_SNAME_SIZE,(int)listcomponent_name[i].size());
1738     component_name.replace(i*MED_SNAME_SIZE,length,
1739                            listcomponent_name[i],0,length);
1740     length = min(MED_SNAME_SIZE,(int)listcomponent_unit[i].size());
1741     component_unit.replace(i*MED_SNAME_SIZE,length,
1742                            listcomponent_unit[i],0,length);
1743   }
1744
1745   MESSAGE_MED("using component_name=|"<<component_name<<"|");
1746   MESSAGE_MED("using component_unit=|"<<component_unit<<"|");
1747
1748   MED_EN::med_type_champ ValueType=MED_FIELD_DRIVER<T>::_ptrField->getValueType() ;
1749
1750   MESSAGE_MED("Template Type =|"<<ValueType<<"|");
1751
1752   // Vérifier si le champ existe déjà
1753   char   champName[MED_NAME_SIZE+1];
1754   char * compName, * compUnit ;
1755   med_2_3::med_field_type type ;
1756   bool Find = false ;
1757   int n = med_2_3::MEDnField(id);
1758   int nbComp = 0;
1759   for (int i=1; i<=n; i++) {
1760     nbComp   = med_2_3::MEDfieldnComponent(id,i);
1761     compName = new char[MED_SNAME_SIZE*nbComp+1];
1762     compUnit = new char[MED_SNAME_SIZE*nbComp+1];
1763     med_2_3::med_bool localmesh;
1764     med_2_3::med_int nbstpp3;
1765     char dtunit[MED_LNAME_SIZE+1];
1766     char mmmp3[MED_NAME_SIZE+1];
1767     err = med_2_3::MEDfieldInfo(id,i,champName,mmmp3,&localmesh,&type,compName,compUnit,dtunit,&nbstpp3);
1768     if (err == 0)
1769       if (!strcmp(champName,fieldName.c_str()) ) {
1770         Find = true ;
1771         break ;
1772       }
1773     delete[] compName ;
1774     delete[] compUnit ;
1775   }
1776
1777   if (Find) {
1778     // the same ?
1779     if (nbComp != component_count)
1780       throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1781                                      <<": Field exist in file, but number of component are different : "
1782                                      <<nbComp<<" in file and "
1783                                      <<component_count<<" in memory."
1784                                      )
1785                           );
1786     // component name and unit
1787     SCRUTE_MED(nbComp);
1788     MESSAGE_MED(LOC<<" Component name in file   : "<<compName);
1789     MESSAGE_MED(LOC<<" Component name in memory : "<<component_name);
1790     MESSAGE_MED(LOC<<" Component unit in file   : "<<compUnit);
1791     MESSAGE_MED(LOC<<" Component unit in memory : "<<component_unit);
1792     delete[] compName ;
1793     delete[] compUnit ;
1794
1795   } else {
1796     // Verify the field isn't yet created
1797
1798     string dataGroupName =  "/CHA/";
1799     dataGroupName        += fieldName;
1800     MESSAGE_MED(LOC << "|" << dataGroupName << "|" );
1801     med_2_3::med_idt gid =  H5Gopen(id, dataGroupName.c_str() );
1802
1803     if ( gid < 0 )
1804       { // create field
1805         string meshNameWr = meshName;
1806         meshNameWr.resize( MED_NAME_SIZE + 1, '\0' ); // avoid "Invalid read" memory error
1807         err=med_2_3::MEDfieldCr(id,fieldName.c_str(),
1808                                 (med_2_3::med_field_type)ValueType,
1809                                 component_count,
1810                                 component_name.c_str(),
1811                                 component_unit.c_str(),
1812                                 "",
1813                                 meshNameWr.c_str());
1814
1815         if ( err < 0 )
1816           throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1817                                          << ": Error MEDchampCr : "<<err
1818                                          )
1819                               );
1820       }
1821     else H5Gclose(gid);
1822   }
1823
1824
1825
1826   // On s'assure que le champ est dans le bon type d'entrelacement.
1827   // REM : Il faudrait un driver par type d'entrelacement, ce qui eviterait
1828   // de doubler l'utilisation de la taille mémoire si le champ n'est pas dans
1829   // le bon mode.
1830   FIELD<T,FullInterlace> * myField = 0;
1831   MED_EN::medModeSwitch interlacingType = MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType();
1832   bool isFullInterlace     = ( interlacingType == MED_EN::MED_FULL_INTERLACE );
1833   bool isNoInterlaceByType = ( interlacingType == MED_EN::MED_NO_INTERLACE_BY_TYPE );//PAL17011
1834   med_2_3::med_switch_mode modswt = med_2_3::MED_FULL_INTERLACE;
1835
1836   if ( isFullInterlace ) {
1837     myField = MED_FIELD_DRIVER<T>::_ptrField;
1838   }
1839   else if ( isNoInterlaceByType ) {
1840     // PAL17011, no need to convert, that is what this improvement is needed for
1841     modswt = med_2_3::MED_NO_INTERLACE;
1842   }
1843   else {
1844     myField = FieldConvert( *((FIELD<T,NoInterlace>*) MED_FIELD_DRIVER<T>::_ptrField ));
1845   }
1846
1847   // Il est necessaire de calculer le tableau
1848   // du nombre d'entités cumulées de chaque type géométrique du maillage
1849   // pour convertir les profils de la numérotation globale
1850   // à la numérotation locale au type géométrique.
1851   // Pour celà on établit ce tableau à partir de l'objet MESH si la relation SUPPORT-MESH existe.
1852   // Si le maillage existe dans le fichier MED  on essaye également de reconstituer ce tableau
1853   // pour vérifier la cohérence des informations.
1854   // Si la relation SUPPRT-MESH n'esiste pas  on constitue le tableau uniquement à partir du fichier qui
1855   // doit alors obligatoirement contenir le maillage.
1856   const int * number, *numberIndex = 0;
1857   string         profilName;
1858   vector<string> profilNameList;
1859   vector<MED_EN::medGeometryElement> meshGeoType;
1860   vector<int> meshNbOfElOfType;
1861   vector<int> meshNbOfElOfTypeC;
1862   vector<MED_EN::medGeometryElement> fileMeshGeoType;
1863   vector<int> fileMeshNbOfElOfType;
1864   vector<int> fileMeshNbOfElOfTypeC;
1865   med_2_3::med_int fileHasMesh=0;
1866
1867   if (!onAll) {
1868
1869     number = mySupport->getNumber(MED_ALL_ELEMENTS);
1870     numberIndex = mySupport->getNumberIndex();
1871     profilNameList=mySupport->getProfilNames();
1872     // porting MED3
1873     med_2_3::med_int spceDimp3,mdimp3;
1874     med_2_3::med_mesh_type mtype;
1875     char desccp3[MED_COMMENT_SIZE+1];
1876     char dttunittp3[MED_LNAME_SIZE+1];
1877     med_2_3::med_sorting_type sttp3;
1878     med_2_3::med_int nsteppp3;
1879     med_2_3::med_axis_type axxxppp3;
1880     int naxis=med_2_3::MEDmeshnAxisByName(id,meshName.c_str());
1881     char *annp3=new char[naxis*MED_SNAME_SIZE+1];
1882     char *auup3=new char[naxis*MED_SNAME_SIZE+1];
1883     fileHasMesh=(med_2_3::MEDmeshInfoByName(id,meshName.c_str(),&spceDimp3,&mdimp3,&mtype,desccp3,dttunittp3,&sttp3,&nsteppp3,&axxxppp3,annp3,auup3)==0);
1884     delete [] annp3;
1885     delete [] auup3;
1886   // end porting MED3
1887     const GMESH * meshPtr = mySupport->getMesh();
1888 //     if(!meshPtr)
1889 //       throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1890 //                                      <<": Mesh in support is null"
1891 //                                      )
1892 //                           );
1893
1894     if (fileHasMesh)
1895       this->getMeshGeometricTypeFromFile(id, meshName,
1896                                          entityType,
1897                                          fileMeshGeoType,fileMeshNbOfElOfType,fileMeshNbOfElOfTypeC);
1898
1899     if (meshPtr) {
1900       this->getMeshGeometricTypeFromMESH( meshPtr, entityType,meshGeoType,
1901                                           meshNbOfElOfType,
1902                                           meshNbOfElOfTypeC);
1903
1904       if (fileHasMesh)
1905         if ( ( fileMeshGeoType != meshGeoType ) || (fileMeshNbOfElOfTypeC != meshNbOfElOfTypeC) )
1906           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information from file for FIELD "<< fieldName
1907                                         << " on entity " << MED_EN::entNames[entityType]
1908                                         << " with (it,or) = ("
1909                                         << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1910                                         << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1911                                         << " on mesh " << meshName
1912                                         << " : geometric types or number of elements by type differs from MESH object !"
1913                                         )
1914                              );
1915
1916     }
1917
1918     if ( !fileHasMesh && meshPtr==0 )
1919       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information for FIELD "<< fieldName
1920                                     << " on entity " << MED_EN::entNames[entityType]
1921                                     << " with (it,or) = ("
1922                                     << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1923                                     << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1924                                     << " on mesh " << meshName
1925                                     << " : SUPPORT must contain a valid MESH reference or file must contain the associated MESH."
1926                                     )
1927                          );
1928
1929
1930     if (fileHasMesh && !meshPtr) {
1931       meshNbOfElOfTypeC = fileMeshNbOfElOfTypeC;
1932       meshGeoType       = fileMeshGeoType;
1933       meshNbOfElOfType  = fileMeshNbOfElOfType;
1934     }
1935   }
1936
1937   const MED_EN::medGeometryElement * types = mySupport->getTypes() ;
1938   int numberOfTypes = mySupport->getNumberOfTypes() ;
1939   int numberOfElForMED = -1;
1940   const T   * value   = NULL;
1941   int index = 1 ;
1942
1943   //converting MEDMEM type to MEDfile type 
1944   if (entityType != MED_EN::MED_NODE)
1945     entityType = MED_EN::MED_CELL;
1946
1947   // on boucle sur tout les types pour ecrire les tableaux de valeur
1948   for (int typeNo=0;typeNo<numberOfTypes;typeNo++) {
1949
1950     int numberOfElements = mySupport->getNumberOfElements(types[typeNo]) ;
1951     //UP : prend en compte les profils, pas les points de Gauss
1952
1953     //value = MED_FIELD_DRIVER<T>::_ptrField->getRow(index) ;
1954     // rem 1 : le getRow du Array est différent de celui du FIELD si le SUPPORT contient
1955     //         des profils (les indices des valeurs ne se suivent pas forcément)
1956     // rem 2 : Afin de respecter la norme MEDFICHIER, les indices contenus dans les
1957     //         profils doivent être croissant
1958     if (onAll) {
1959
1960       if ( isNoInterlaceByType ) { //PAL17011
1961         value = MED_FIELD_DRIVER<T>::_ptrField->getValueByType(typeNo+1);
1962           //((ArrayNoByType *)MED_FIELD_DRIVER<T>::_ptrField->getArray())->getValueByType(i+1);
1963       }
1964       else {
1965         value = myField->getRow(index);
1966       }
1967       profilName=MED_NO_PROFILE;
1968       numberOfElForMED = numberOfElements;
1969
1970     } else {
1971
1972       if ( isNoInterlaceByType ) { //PAL17011
1973         value = MED_FIELD_DRIVER<T>::_ptrField->getValueByType(typeNo+1);
1974       }
1975       else {
1976         value = myField->getRow(number[index-1]);
1977       }
1978       // PAL16854(Partial support on nodes) ->
1979       //profilName = (profilNameList.size()>typeNo) ? profilNameList[typeNo].substr(0,MED_NAME_SIZE) : MED_NOPFL;
1980       if (profilNameList[typeNo].size()>MED_NAME_SIZE)
1981         profilName = healName( profilNameList[typeNo].substr(0,MED_NAME_SIZE) );
1982       else
1983         profilName = healName( profilNameList[typeNo] );
1984
1985       // Rem : Si le SUPPORT n'est pas onAll mais que pour un type géométrique donné le nom
1986       // du profil associé est MED_NOPFL alors le profil n'est pas écrit dans le fichier MED.
1987       // Car en MEDMEMOIRE si le champ repose sur des éléments de deux types géométriques
1988       // différents et est défini sur tous les éléments d'un type géométrique
1989       // mais pas de l'autre, il existe tout de même des profils sur les deux types géométriques.
1990       // Ce n'est pas le cas en MEDFICHIER.
1991       vector<med_2_3::med_int> profil(&number[index-1],&(number[index-1])+numberOfElements);
1992
1993       int meshTypeNo=0;
1994       if ( entityType != MED_EN::MED_NODE ) // PAL16854(Partial support on nodes)
1995       {
1996         // Trouve l'index du type géométrique dans la liste des types géométriques du maillage
1997         // correspondant au type géométrique du champ en cours de traitement
1998         vector<MED_EN::medGeometryElement>::iterator meshTypeNoIt =
1999           find(meshGeoType.begin(),meshGeoType.end(),types[typeNo]);
2000         if ( meshTypeNoIt ==  meshGeoType.end() )
2001           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Can't find "<< MED_EN::geoNames[types[typeNo]]
2002                                         << " on entity " << MED_EN::entNames[entityType]
2003                                         << " in geometric type list of mesh " << meshName
2004                                         )
2005                              );
2006         meshTypeNo = meshTypeNoIt -  meshGeoType.begin();
2007       }
2008
2009       if ( profilName == MED_NO_PROFILE && (int)profil.size() != meshNbOfElOfType[meshTypeNo] )
2010         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating profil for FIELD "<< fieldName 
2011                                       << " on entity " << MED_EN::entNames[entityType]
2012                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
2013                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
2014                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
2015                                       << profilName << " on mesh " << meshName
2016                                       << " : There is no profileName but profilsize (" <<profil.size()
2017                                       << ") differs from number of elements in associated MESH ("
2018                                       << meshNbOfElOfType[meshTypeNo] << ")."
2019                             )
2020                          );
2021
2022       //REM : Ce n'est pas évident, mais lorsqu'il y a un profil, le nombre de valeurs
2023       //      que l'on indique à MEDchampEcr est le nombre de valeurs sans profil, d'où
2024       //      le nombre d'éléments du maillage sur le type géométrique courant.
2025       numberOfElForMED = meshNbOfElOfType[meshTypeNo];
2026
2027       for (int ind=0;ind < numberOfElements;++ind) {
2028 //      cout << "number["<<index-1<<"]="<<number[index-1]<<endl;
2029 //      cout << "profil1["<<ind<<"]="<<profil[ind]<<endl;
2030         profil[ind]-=meshNbOfElOfTypeC[meshTypeNo];
2031 //      cout << "profil2["<<ind<<"]="<<profil[ind]<<endl;
2032       }
2033
2034       if ( profil[numberOfElements-1] > numberOfElForMED )
2035         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating profil for FIELD "<< fieldName 
2036                                       << " on entity " << MED_EN::entNames[entityType]
2037                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
2038                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
2039                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
2040                                       << profilName << " on mesh " << meshName
2041                                       << " : profil["<<numberOfElements-1<<"]=" << profil[numberOfElements-1]
2042                                       << " must not be superior to field size without profil : "
2043                                       << numberOfElForMED
2044                             )
2045                          );
2046
2047       if ( med_2_3::MEDprofileWr(id,profilName.c_str(),numberOfElements,&profil[0]) < 0)
2048
2049         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing "<< numberOfElements
2050                                       << " values for MED profil "<< profilName
2051                                       )
2052                            );
2053     }
2054
2055     bool anyGauss = MED_FIELD_DRIVER<T>::_ptrField->getGaussPresence();
2056     string locName;
2057     if ( anyGauss ) {
2058       //       cout << endl << "Nombre de points de Gauss à l'écriture de " << fieldName
2059       //         << " pour le type géométrique : " << MED_EN::geoNames[types[typeNo]]
2060       //         << " : " << myField->getNumberOfGaussPoints(types[typeNo]) << endl;
2061       //       cout << *mySupport << endl;
2062
2063       const GAUSS_LOCALIZATION_ * locPtr=0;
2064       locMap::const_iterator it;
2065       if ( ( it = gaussModel.find(types[typeNo])) != gaussModel.end() )
2066         locPtr = (*it).second;
2067       else
2068         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating Gauss Model for FIELD "<< fieldName 
2069                                       << " on entity " << MED_EN::entNames[entityType]
2070                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
2071                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
2072                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
2073                                       << profilName << " on mesh " << meshName
2074                                       << " : Can't find a Gauss localisation model for this geometric type" 
2075                                       )
2076                            );
2077
2078       int ngauss = -1;
2079       int mdim=CELLMODEL_Map::retrieveCellModel(types[typeNo]).getDimension();
2080       if ( locPtr->getInterlacingType() == MED_EN::MED_FULL_INTERLACE ) {
2081         const GAUSS_LOCALIZATION<FullInterlace> & loc=*(static_cast<const GAUSS_LOCALIZATION<FullInterlace> * >(locPtr));
2082         ngauss = loc.getNbGauss();
2083         locName=healName( loc.getName() );
2084         err=med_2_3::MEDlocalizationWr(id,locName.c_str(),
2085                                        (med_2_3::med_geometry_type)loc.getType(),
2086                                        mdim,
2087                                        loc.getRefCoo().getPtr(),
2088                                        med_2_3::MED_FULL_INTERLACE,
2089                                        ngauss,
2090                                        loc.getGsCoo().getPtr(),
2091                                        &loc.getWeight()[0],
2092                                        MED_NO_INTERPOLATION, MED_NO_MESH_SUPPORT);
2093
2094       } else {
2095         const GAUSS_LOCALIZATION<NoInterlace> & loc=*(static_cast<const GAUSS_LOCALIZATION<NoInterlace> * >(locPtr));
2096         ngauss = loc.getNbGauss();
2097         locName=healName( loc.getName() );
2098         err=med_2_3::MEDlocalizationWr(id,locName.c_str(),
2099                                        (med_2_3::med_geometry_type)loc.getType(),
2100                                        mdim,
2101                                        loc.getRefCoo().getPtr(),
2102                                        med_2_3::MED_NO_INTERLACE,
2103                                        ngauss,
2104                                        loc.getGsCoo().getPtr(),
2105                                        &loc.getWeight()[0],
2106                                        MED_NO_INTERPOLATION,
2107                                        MED_NO_MESH_SUPPORT);
2108       }
2109
2110       if ( err != 0 )
2111         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing Gauss Model for FIELD "<< fieldName
2112                                       << " on entity " << MED_EN::entNames[entityType]
2113                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] 
2114                             )
2115                          );
2116
2117       //numberOfElForMED *= mySupport->getNumberOfGaussPoints(types[typeNo]); //Deplacer la méthode dans FIELD
2118       //numberOfElForMED *= ngauss;
2119     }
2120
2121     MESSAGE_MED("MED_FIELD_DRIVER<T>::_medIdt                       : "<<id);
2122     MESSAGE_MED("meshName.c_str()                : "<<meshName.c_str());
2123     MESSAGE_MED("MED_FIELD_DRIVER<T>::_ptrField->getName()            : "<<MED_FIELD_DRIVER<T>::_ptrField->getName());
2124     MESSAGE_MED("MED_FIELD_DRIVER<T>::_fieldName                      : "<<MED_FIELD_DRIVER<T>::_fieldName);
2125     MESSAGE_MED("value                           : "<<value);
2126     MESSAGE_MED("numberOfElements                : "<<numberOfElements);
2127     MESSAGE_MED("numberOfElForMED                : "<<numberOfElForMED);
2128     MESSAGE_MED("entityType                      : "<<MED_EN::entNames[entityType]);
2129     MESSAGE_MED("types[i]                        : "<<MED_EN::geoNames[types[typeNo]]);
2130     if (myField) //myField may be NULL (PAL17011)
2131       MESSAGE_MED("NumberOfGaussPoint[i]           : "<<myField->getNumberOfGaussPoints(types[typeNo]));
2132     MESSAGE_MED("MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber() : "<<MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber());
2133     MESSAGE_MED("MED_FIELD_DRIVER<T>::_ptrField->getTime()            : "<<MED_FIELD_DRIVER<T>::_ptrField->getTime());
2134     MESSAGE_MED("MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()     : "<<MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber());
2135
2136     meshName.resize( MED_NAME_SIZE+1, '\0'); // for valgrind reporting "Invalid read of size 1"
2137
2138     // Rem 1 : le nombre d'éléments passé à MEDchampEcr ne doit pas tenir compte de la taille
2139     //         des profils : c'est la taille du champ sans profil.
2140     err=med_2_3::MEDfieldValueWithProfileWr(id,fieldName.c_str(),
2141                                             MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber(),
2142                                             MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber(),
2143                                             MED_FIELD_DRIVER<T>::_ptrField->getTime(),
2144                                             (med_2_3::med_entity_type)entityType,
2145                                             (med_2_3::med_geometry_type)types[typeNo],
2146                                             med_2_3::MED_COMPACT_PFLMODE,
2147                                             profilName.c_str(),
2148                                             locName.c_str(),
2149                                             modswt,MED_ALL_CONSTITUENT,
2150                                             numberOfElForMED,
2151                                             (unsigned char*)value);
2152     if (err < MED_VALID ) {
2153       if ( !isFullInterlace )
2154         delete myField;
2155       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing "<< numberOfElements
2156                                     << " values for FIELD "<< fieldName 
2157                                     << " on entity " << MED_EN::entNames[entityType]
2158                                     << " and geometric type " << MED_EN::geoNames[types[typeNo]]
2159                                     << " with (it,or) = ("
2160                                     << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
2161                                     << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
2162                                     << profilName << " on mesh " << meshName
2163                                     )
2164                          );
2165     }
2166
2167     index += numberOfElements ; //Ne doit pas prendre en compte le nombre de points de GAUSS
2168                                 //ni les composantes.
2169
2170   }
2171   if ( !isFullInterlace ) delete myField;
2172
2173
2174   END_OF_MED(LOC);
2175 }
2176
2177 /*--------------------- RDWR PART -------------------------------*/
2178
2179 /*!
2180   Constructor.
2181 */
2182 template <class T> 
2183 MED_FIELD_RDWR_DRIVER<T>::MED_FIELD_RDWR_DRIVER():MED_FIELD_DRIVER<T>()
2184 {
2185   this->GENDRIVER::_accessMode = MED_EN::RDWR;
2186 }
2187
2188 /*!
2189   Constructor.
2190 */
2191 template <class T> 
2192 template <class INTERLACING_TAG>
2193 MED_FIELD_RDWR_DRIVER<T>::MED_FIELD_RDWR_DRIVER(const string &              fileName,
2194                                                 FIELD<T, INTERLACING_TAG> * ptrField):
2195   MED_FIELD_DRIVER<T>(fileName,ptrField,MED_EN::RDWR),
2196   MED_FIELD_RDONLY_DRIVER<T>(fileName,ptrField),
2197   MED_FIELD_WRONLY_DRIVER<T>(fileName,ptrField)
2198 {
2199   const char* LOC = "MED_FIELD_RDWR_DRIVER::MED_FIELD_RDWR_DRIVER(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)";
2200   BEGIN_OF_MED(LOC);
2201   //_accessMode = MED_RDWR ;
2202   END_OF_MED(LOC);
2203 }
2204
2205 /*!
2206   Copy constructor.
2207 */
2208 template <class T> 
2209 MED_FIELD_RDWR_DRIVER<T>::MED_FIELD_RDWR_DRIVER(const MED_FIELD_RDWR_DRIVER<T> & fieldDriver):
2210   MED_FIELD_DRIVER<T>(fieldDriver),
2211   MED_FIELD_RDONLY_DRIVER<T>(fieldDriver),
2212   MED_FIELD_WRONLY_DRIVER<T>(fieldDriver)
2213 {}
2214
2215 /*!
2216   Destructor.
2217 */
2218 template <class T> 
2219 MED_FIELD_RDWR_DRIVER<T>::~MED_FIELD_RDWR_DRIVER()
2220 {
2221 }
2222
2223 template <class T> GENDRIVER * MED_FIELD_RDWR_DRIVER<T>::copy(void) const
2224 {
2225   return new MED_FIELD_RDWR_DRIVER<T>(*this);
2226 }
2227
2228 template <class T> void MED_FIELD_RDWR_DRIVER<T>::write(void) const
2229   throw (MEDEXCEPTION)
2230 {
2231   const char* LOC = "MED_FIELD_RDWR_DRIVER::write(void)";
2232   BEGIN_OF_MED(LOC);
2233   MED_FIELD_WRONLY_DRIVER<T>::write();
2234   END_OF_MED(LOC);
2235 }
2236
2237 template <class T> void MED_FIELD_RDWR_DRIVER<T>::read (void)
2238   throw (MEDEXCEPTION)
2239 {
2240   const char* LOC = "MED_FIELD_RDWR_DRIVER::read(void)";
2241   BEGIN_OF_MED(LOC);
2242   MED_FIELD_RDONLY_DRIVER<T>::read();
2243   END_OF_MED(LOC);
2244 }
2245
2246 } // end namespace MEDMEM
2247
2248 #endif