]> SALOME platform Git repositories - modules/med.git/blob - src/MEDMEM/MEDMEM_MedFieldDriver22.hxx
Salome HOME
Join modifications from branch CEAFor_V3_2_0
[modules/med.git] / src / MEDMEM / MEDMEM_MedFieldDriver22.hxx
1 // Copyright (C) 2005  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
2 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
3 // 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either 
7 // version 2.1 of the License.
8 // 
9 // This library is distributed in the hope that it will be useful 
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of 
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU 
12 // Lesser General Public License for more details.
13 //
14 // You should have received a copy of the GNU Lesser General Public  
15 // License along with this library; if not, write to the Free Software 
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
17 //
18 // See http://www.salome-platform.org/
19 //
20 #ifndef MED_FIELD_DRIVER22_HXX
21 #define MED_FIELD_DRIVER22_HXX
22
23 #include <string>
24 #include <algorithm>
25
26 #include "MEDMEM_FieldConvert.hxx"
27 #include "MEDMEM_ArrayInterface.hxx"
28 #include "MEDMEM_ArrayConvert.hxx"
29
30 #include "MEDMEM_define.hxx"
31 #include "MEDMEM_Utilities.hxx"
32 #include "MEDMEM_STRING.hxx"
33 #include "MEDMEM_Exception.hxx"
34
35 #include "MEDMEM_DriversDef.hxx"
36 #include "MEDMEM_MedFieldDriver.hxx"
37 #include "MEDMEM_Unit.hxx"
38 #include "MEDMEM_Support.hxx"
39 #include "MEDMEM_GaussLocalization.hxx"
40
41 //includes temporaires (attente release med fichier 2.3.1)
42 #include "MEDMEM_MEDMEMgaussEcr.hxx"
43 #include "MEDMEM_MEDMEMprofilEcr.hxx"
44 #include "MEDMEM_MEDMEMchampLire.hxx"
45
46 namespace MEDMEM {
47
48 /*!
49
50   Driver Med for FIELD.
51
52   Generic part : implement open and close methods.
53
54 */
55
56 template <class T> class MED_FIELD_DRIVER22 : public virtual MED_FIELD_DRIVER<T>
57 {
58 protected:
59
60   med_2_2::med_idt        _medIdt;
61
62   bool createFieldSupportPart1(med_2_2::med_idt id,
63                           const string & fieldName,
64                           med_2_2::med_int ndt,
65                           med_2_2::med_int od,
66                           SUPPORT & support,
67                           string & meshName,
68                           vector<int> & numberOfElementsOfTypeC,
69                           vector<int> & numberOfGaussPoint,
70                           int & totalNumberOfElWg
71                           ) const throw (MEDEXCEPTION);
72
73   void getMeshGeometricTypeFromFile(med_2_2::med_idt id,
74                             string & meshName,
75                             MED_EN::medEntityMesh  entite,
76                             vector<MED_EN::medGeometryElement> & geoType,
77                             vector<int> &nbOfElOfType,
78                             vector<int> &nbOfElOfTypeC) const throw(MEDEXCEPTION);
79
80   void getMeshGeometricTypeFromMESH( MESH * meshPtr,
81                                      MED_EN::medEntityMesh  entity,
82                                      vector<MED_EN::medGeometryElement> & geoType,
83                                      vector<int> &nbOfElOfType,
84                                      vector<int> &nbOfElOfTypeC) const throw(MEDEXCEPTION);
85
86 public :
87
88   /*!
89     Constructor.
90   */
91   MED_FIELD_DRIVER22():MED_FIELD_DRIVER<T>(),_medIdt(MED_INVALID)
92   {}
93   /*!
94     Constructor.
95   */
96   template <class INTERLACING_TAG>
97   MED_FIELD_DRIVER22(const string & fileName,
98                      FIELD<T, INTERLACING_TAG> * ptrField,
99                      MED_EN::med_mode_acces accessMode)
100     : MED_FIELD_DRIVER<T>(fileName,ptrField,accessMode),_medIdt(MED_INVALID)
101   {
102   }
103
104   /*!
105     Copy constructor.
106   */
107   MED_FIELD_DRIVER22(const MED_FIELD_DRIVER22 & fieldDriver):
108     MED_FIELD_DRIVER<T>(fieldDriver),
109     _medIdt(fieldDriver._medIdt)
110   {
111   }
112
113   /*!
114     Destructor.
115   */
116   virtual ~MED_FIELD_DRIVER22() {
117   }
118
119   void open() throw (MEDEXCEPTION)
120   {
121     const char * LOC = "MED_FIELD_DRIVER22::open() ";
122     BEGIN_OF(LOC);
123
124     // we must set fieldname before open, because we must find field number in file (if it exist !!!)
125     if ( MED_FIELD_DRIVER<T>::_fileName == "" )
126       throw MED_EXCEPTION ( LOCALIZED( STRING(LOC)
127                                        << "_fileName is |\"\"|, please set a correct fileName before calling open()"
128                                        )
129                             );
130
131     MESSAGE(LOC<<"_fileName.c_str : "<< MED_FIELD_DRIVER<T>::_fileName.c_str()<<",mode : "<< MED_FIELD_DRIVER<T>::_accessMode);
132     MED_FIELD_DRIVER22<T>::_medIdt = med_2_2::MEDouvrir( (const_cast <char *> (MED_FIELD_DRIVER<T>::_fileName.c_str())),(med_2_2::med_mode_acces) MED_FIELD_DRIVER<T>::_accessMode);
133     MESSAGE(LOC<<"_medIdt : "<< MED_FIELD_DRIVER22<T>::_medIdt );
134     if (MED_FIELD_DRIVER22<T>::_medIdt > 0)
135       MED_FIELD_DRIVER<T>::_status=MED_OPENED;
136     else {
137       MED_FIELD_DRIVER<T>::_status = MED_INVALID;
138       MED_FIELD_DRIVER22<T>::_medIdt = MED_INVALID;
139       throw MED_EXCEPTION (LOCALIZED( STRING(LOC)
140                                       << "Can't open |"  << MED_FIELD_DRIVER<T>::_fileName
141                                       << "|, _medIdt : " << MED_FIELD_DRIVER22<T>::_medIdt
142                                       )
143                            );
144     }
145
146     END_OF(LOC);
147   }
148
149   void close() {
150     BEGIN_OF("MED_FIELD_DRIVER22::close()");
151     med_2_2::med_int err = 0;
152     if (MED_FIELD_DRIVER<T>::_status == MED_OPENED) {
153       err=med_2_2::MEDfermer(MED_FIELD_DRIVER22<T>::_medIdt);
154       //H5close(); // If we call H5close() all the files are closed.
155       MED_FIELD_DRIVER<T>::_status = MED_CLOSED;
156       MED_FIELD_DRIVER22<T>::_medIdt = MED_INVALID;
157       MESSAGE(" MED_FIELD_DRIVER22::close() : MEDfermer : _medIdt= " << MED_FIELD_DRIVER22<T>::_medIdt );
158       MESSAGE(" MED_FIELD_DRIVER22::close() : MEDfermer : err    = " << err );
159     }
160     END_OF("MED_FIELD_DRIVER22::close()");
161   }
162 };
163
164 /*!
165
166   Driver Med for FIELD : Read only.
167
168   Implement read method.
169
170 */
171
172   template <class T> class MED_FIELD_RDONLY_DRIVER22 : public virtual MED_FIELD_DRIVER22<T>, public virtual IMED_FIELD_RDONLY_DRIVER<T>
173 {
174
175 public :
176
177   /*!
178     Constructor.
179   */
180   MED_FIELD_RDONLY_DRIVER22():MED_FIELD_DRIVER<T>() {};
181
182   /*!
183     Constructor.
184   */
185   template <class INTERLACING_TAG>
186   MED_FIELD_RDONLY_DRIVER22(const string & fileName,
187                             FIELD<T, INTERLACING_TAG> * ptrField):
188     IMED_FIELD_RDONLY_DRIVER<T>(fileName,ptrField),
189     MED_FIELD_DRIVER22<T>(fileName,ptrField,MED_EN::MED_RDONLY),
190     MED_FIELD_DRIVER<T>(fileName,ptrField,MED_EN::MED_RDONLY)
191   {
192     BEGIN_OF("MED_FIELD_RDONLY_DRIVER22::MED_FIELD_RDONLY_DRIVER22(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)");
193     END_OF("MED_FIELD_RDONLY_DRIVER22::MED_FIELD_RDONLY_DRIVER22(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)");
194   }
195
196   /*!
197     Copy constructor.
198   */
199   MED_FIELD_RDONLY_DRIVER22(const MED_FIELD_RDONLY_DRIVER22 & fieldDriver):
200     IMED_FIELD_RDONLY_DRIVER<T>(fieldDriver),
201     MED_FIELD_DRIVER22<T>(fieldDriver),
202     MED_FIELD_DRIVER<T>(fieldDriver)
203   {}
204
205   /*!
206     Destructor.
207   */
208   virtual ~MED_FIELD_RDONLY_DRIVER22() {};
209
210   // CREER UNE METHODE POUR LIRE LA LISTE DES MAILLAGES .....
211
212   /*!
213     Return a MEDEXCEPTION : it is the read-only driver.
214   */
215   void write( void ) const throw (MEDEXCEPTION) ;
216   /*!
217     Read FIELD in the specified file.
218   */
219   void read ( void ) throw (MEDEXCEPTION) ;
220
221 private:
222   GENDRIVER * copy( void ) const ;
223
224 };
225
226 /*!
227
228   Driver Med for FIELD : Write only.
229
230   Implement write method.
231
232 */
233
234 template <class T> class MED_FIELD_WRONLY_DRIVER22 : public virtual MED_FIELD_DRIVER22<T>, public virtual IMED_FIELD_WRONLY_DRIVER<T> {
235
236 public :
237
238   /*!
239     Constructor.
240   */
241   MED_FIELD_WRONLY_DRIVER22():MED_FIELD_DRIVER<T>() {}
242
243   /*!
244     Constructor.
245   */
246   template <class INTERLACING_TAG>
247   MED_FIELD_WRONLY_DRIVER22(const string & fileName,
248                             FIELD<T, INTERLACING_TAG> * ptrField):
249     IMED_FIELD_WRONLY_DRIVER<T>(fileName,ptrField),
250     MED_FIELD_DRIVER22<T>(fileName,ptrField,MED_EN::MED_WRONLY),
251     MED_FIELD_DRIVER<T>(fileName,ptrField,MED_EN::MED_WRONLY)
252   {
253     BEGIN_OF("MED_FIELD_WRONLY_DRIVER22::MED_FIELD_WRONLY_DRIVER22(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)");
254     END_OF("MED_FIELD_WRONLY_DRIVER22::MED_FIELD_WRONLY_DRIVER22(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)");
255   }
256
257   /*!
258     Copy constructor.
259   */
260   MED_FIELD_WRONLY_DRIVER22(const MED_FIELD_WRONLY_DRIVER22 & fieldDriver):
261     IMED_FIELD_WRONLY_DRIVER<T>(fieldDriver),
262     MED_FIELD_DRIVER22<T>(fieldDriver),
263     MED_FIELD_DRIVER<T>(fieldDriver)
264   {}
265
266   /*!
267     Destructor.
268   */
269   virtual ~MED_FIELD_WRONLY_DRIVER22() {};
270
271   /*!
272     Write FIELD in the specified file.
273   */
274   void write( void ) const throw (MEDEXCEPTION) ;
275   /*!
276     Return a MEDEXCEPTION : it is the write-only driver.
277   */
278   void read ( void ) throw (MEDEXCEPTION) ;
279
280 private:
281   GENDRIVER * copy( void ) const ;
282
283 };
284
285
286 /*!
287
288   Driver Med for FIELD : Read write.
289   - Use read method from MED_FIELD_RDONLY_DRIVER
290   - Use write method from MED_FIELD_WDONLY_DRIVER
291
292 */
293
294 template <class T> class MED_FIELD_RDWR_DRIVER22 : public MED_FIELD_RDONLY_DRIVER22<T>, public MED_FIELD_WRONLY_DRIVER22<T>, public IMED_FIELD_RDWR_DRIVER<T> {
295
296 public :
297
298   /*!
299     Constructor.
300   */
301   MED_FIELD_RDWR_DRIVER22():MED_FIELD_DRIVER22<T>() {}
302
303   /*!
304     Constructor.
305   */
306   template <class INTERLACING_TAG>
307   MED_FIELD_RDWR_DRIVER22(const string & fileName,
308                           FIELD<T, INTERLACING_TAG> * ptrField):
309     MED_FIELD_WRONLY_DRIVER22<T>(fileName,ptrField),
310     MED_FIELD_RDONLY_DRIVER22<T>(fileName,ptrField),
311     IMED_FIELD_RDONLY_DRIVER<T>(fileName,ptrField),
312     IMED_FIELD_WRONLY_DRIVER<T>(fileName,ptrField),
313     MED_FIELD_DRIVER<T>(fileName,ptrField,MED_EN::MED_RDWR),
314     IMED_FIELD_RDWR_DRIVER<T>(fileName,ptrField)
315   {
316     BEGIN_OF("MED_FIELD_RDWR_DRIVER22::MED_FIELD_RDWR_DRIVER22(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)");
317     //_accessMode = MED_RDWR ;
318     END_OF("MED_FIELD_RDWR_DRIVER22::MED_FIELD_RDWR_DRIVER22(const string & fileName, const FIELD<T,INTERLACING_TAG> * ptrField)");
319   }
320
321   /*!
322     Copy constructor.
323   */
324   MED_FIELD_RDWR_DRIVER22(const MED_FIELD_RDWR_DRIVER22 & fieldDriver):
325     MED_FIELD_WRONLY_DRIVER22<T>(fieldDriver),
326     MED_FIELD_RDONLY_DRIVER22<T>(fieldDriver),
327     IMED_FIELD_RDWR_DRIVER<T>(fieldDriver),
328     IMED_FIELD_RDONLY_DRIVER<T>(fieldDriver),
329     IMED_FIELD_WRONLY_DRIVER<T>(fieldDriver),
330     MED_FIELD_DRIVER<T>(fieldDriver)
331   {};
332
333   /*!
334     Destructor.
335   */
336   ~MED_FIELD_RDWR_DRIVER22() {};
337
338   /*!
339     Write FIELD in the specified file.
340   */
341   void write(void) const throw (MEDEXCEPTION) ;
342   /*!
343     Read FIELD in the specified file.
344   */
345   void read (void) throw (MEDEXCEPTION) ;
346
347 private:
348   GENDRIVER * copy( void ) const ;
349
350 };
351
352
353 /*-------------------------*/
354 /* template implementation */
355 /*-------------------------*/
356
357 /*--------------------- DRIVER PART -------------------------------*/
358
359
360 /*!
361
362   Cette méthode crée le SUPPORT du champ <fieldName> pour le
363         <n°de pas de temps,n°d'itération>=<ndt,od>.
364
365   Le SUPPORT crée à pour nom  <fieldName>Support et contient
366   la liste des types géométriques sur le premier type
367   d'entité trouvé (en MEDMEM on inderdit aux champs de reposer
368   sur plusieurs types d'entité).
369   Il contient également le nombre d'entités trouvées pour chaque
370   type géométrique.
371   Par défaut l'attribut onAll du SUPPORT est positionné à true car
372   cette routine ne lit rien de ce qui concerne les entités
373   du maillage associé.
374   La méthode renvoie true si elle réussit à créer le SUPPORT
375   demandé.
376   Le nom du maillage associé ( en MEDMEM on ne
377   supporte pas encore les maillages multiples ) est renvoyé dans <meshName>.
378   Deux tableaux directements exploitables par MEDMEMnArray sont renvoyés :
379   - numberOfElementsOfTypeC : nombres d'entités cumulés de chaque type géométrique
380        avec numberOfElementsOfTypeC[0]=1 et de taille nombre de types+1
381   - numberOfGaussPoint : nombre de points de Gauss par type géométrique
382        avec numberOfGaussPoint[0]=1 et de taille  nombre de types+1
383 */
384
385 template <class T> bool
386 MED_FIELD_DRIVER22<T>::createFieldSupportPart1(med_2_2::med_idt id,
387                                                const string & fieldName,
388                                                med_2_2::med_int ndt,
389                                                med_2_2::med_int od,
390                                                SUPPORT & support,
391                                                string & meshName,
392                                                vector<int> & numberOfElementsOfTypeC,
393                                                vector<int> & numberOfGaussPoint,
394                                                int & totalNumberOfElWg
395                                                ) const throw (MEDEXCEPTION)
396 {
397
398   //EF : Gérer le meshName pour le driver 2.2
399   const char * LOC="MED_FIELD_DRIVER<T>::createFieldSupportPart1(...)";
400
401   map<int, list<MED_EN::medGeometryElement> > CellAndNodeEntities;
402   map<int, list<MED_EN::medGeometryElement> >::iterator currentEntity;
403   CellAndNodeEntities[MED_EN::MED_CELL]  = MED_EN::meshEntities[MED_EN::MED_CELL];
404   CellAndNodeEntities[MED_EN::MED_NODE] = MED_EN::meshEntities[MED_EN::MED_NODE];
405   list< MED_EN::medGeometryElement >::const_iterator currentGeometry;
406
407   MED_EN::medEntityMesh entity;
408   bool alreadyFoundAnEntity=false,alreadyFoundPdtIt = false;
409   int  numberOfElements = 0;
410   int  numberOfGeometricType = 0;
411   MED_EN::medGeometryElement geometricType[MED_NBR_GEOMETRIE_MAILLE];
412   int numberOfElementsOfType[MED_NBR_GEOMETRIE_MAILLE];
413   numberOfElementsOfTypeC.clear();numberOfGaussPoint.clear();
414   numberOfElementsOfTypeC.resize(MED_NBR_GEOMETRIE_MAILLE+1);
415   numberOfGaussPoint.resize(MED_NBR_GEOMETRIE_MAILLE+1);
416
417   med_2_2::med_int nmaa=0, ngauss=0, numdt=-1, numo=-1, nbPdtIt=0;
418   char dtunit[MED_TAILLE_PNOM22+1];
419   char maa[MED_TAILLE_NOM+1];
420   med_2_2::med_float   dt=-1.0;
421   med_2_2::med_booleen local;
422   med_2_2::med_err     ret=1;
423   numberOfElementsOfTypeC[0] = 1;
424   numberOfGaussPoint[0] = 1;
425   totalNumberOfElWg = 0;
426
427   /* Détermine le type d'entité et la liste des types géométriques associés
428      au champ <fieldName> */
429   for (currentEntity = CellAndNodeEntities.begin();
430        currentEntity != CellAndNodeEntities.end(); currentEntity++) {
431     for (currentGeometry  = (*currentEntity).second.begin();
432          currentGeometry != (*currentEntity).second.end(); currentGeometry++) {
433
434
435       if ( (nbPdtIt = med_2_2::MEDnPasdetemps(id, const_cast <char*> ( fieldName.c_str() ),
436                                               (med_2_2::med_entite_maillage)   (*currentEntity).first,
437                                               (med_2_2::med_geometrie_element)  *currentGeometry ))  <=  0 )
438         continue;
439
440       /* Verifie que le champ n'est pas défini sur un autre type d'entité */
441       if ( alreadyFoundAnEntity ) {
442         if (entity != (*currentEntity).first )
443           throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |"  << fieldName
444                                        << "| with (ndt,or) = (" << ndt << ","
445                                        << od << ") must not be defined on nodes and cells" ));
446
447       } else { entity=(*currentEntity).first; alreadyFoundAnEntity = true; };
448
449
450       /* Cherche le champ pour le <ndt>,<ot> demandé et détermine le nombre de points de Gauss*/
451       ret = 0; alreadyFoundPdtIt = false; ngauss =0;
452       for ( med_2_2::med_int j=1; j <= nbPdtIt; j++ ) {
453
454         // Search how many <ngauss> (<fieldName>,<ndt>,<ot>) has
455         ret += med_2_2::MEDpasdetempsInfo(id, const_cast <char*> ( fieldName.c_str() ),
456                                          (med_2_2::med_entite_maillage)   (*currentEntity).first,
457                                          (med_2_2::med_geometrie_element)  *currentGeometry,
458                                          j, &ngauss,  &numdt,  &numo, dtunit, &dt,
459                                           maa, &local, &nmaa);
460
461         if ( ndt == numdt && numo == od ) {
462           alreadyFoundPdtIt = true;
463
464           if ( nmaa > 1 ) {
465             MESSAGE(LOC<<" Field |" << fieldName << "| with (ndt,or) = ("
466                     << ndt << "," << od << ") for (entityType,geometricType)=("
467                     << MED_EN::entNames[(*currentEntity).first] << ","
468                     << MED_EN::geoNames[*currentGeometry] << ")"
469                     << "is defined on multiple meshes, using dafault mesh  |" << maa << "|" );
470           }
471
472           if ( !local) {
473             MESSAGE(" Field |" << fieldName << "| with (ndt,or) = ("
474                     << ndt << "," << od << ") for (entityType,geometricType)=("
475                     << MED_EN::entNames[(*currentEntity).first] << ","
476                     << MED_EN::geoNames[*currentGeometry] << ")"
477                     << "is using a mesh on a distant file (ignored)" );
478
479 //            throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |" << fieldName << "| with (ndt,or) = ("
480 //                                         << ndt << "," << od << ") for (entityType,geometricType)=("
481 //                                         << MED_EN::entNames[(*currentEntity).first] << ","
482 //                                         << MED_EN::geoNames[*currentGeometry] << ")"
483 //                                         << "is using a mesh on a different file which is not yet supported" ));
484           }
485
486           if ( ! meshName.empty() )
487             if ( meshName != maa ) {
488               throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |" << fieldName << "| with (ndt,or) = ("
489                                            << ndt << "," << od << ") for (entityType,geometricType)=("
490                                            << MED_EN::entNames[(*currentEntity).first] << ","
491                                            << MED_EN::geoNames[*currentGeometry] << ")"
492                                            << "is defined on mesh |" << maa << "| not on mesh |" << meshName ));
493             }
494           break;
495         }
496
497       }
498
499       if ( !alreadyFoundPdtIt )
500         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |" << fieldName << "| with (ndt,or) = ("
501                                      << ndt << "," << od << ") should be defined for (entityType,geometricType)=("
502                                      << MED_EN::entNames[(*currentEntity).first] << ","
503                                      << MED_EN::geoNames[*currentGeometry] << ")" ));
504
505       if ( (ret != 0)  || (ngauss < 1 ) )
506         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Error in MEDpasdetempsInfo for  Field |" << fieldName 
507                                      << "| with (ndt,or) = ("
508                                      << ndt << "," << od << ") for (entityType,geometricType)=("
509                                      << MED_EN::entNames[(*currentEntity).first] << ","
510                                      << MED_EN::geoNames[*currentGeometry] << ")" )); ;
511
512       if ( (numberOfElements =  med_2_2::MEDnVal(id, const_cast <char*> ( fieldName.c_str() ),
513                                                 (med_2_2::med_entite_maillage)   (*currentEntity).first,
514                                                 (med_2_2::med_geometrie_element) *currentGeometry,
515                                                  numdt, numo, maa, med_2_2::MED_COMPACT))  <=  0 )
516         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Error in MEDnVal for  Field |" << fieldName
517                                      << "| with (ndt,or) = ("
518                                      << ndt << "," << od << ") for (entityType,geometricType)=("
519                                      << MED_EN::entNames[(*currentEntity).first] << ","
520                                      << MED_EN::geoNames[*currentGeometry] << ")" )); ;
521
522       numberOfElementsOfType[numberOfGeometricType] = numberOfElements/ngauss;
523       numberOfElementsOfTypeC[numberOfGeometricType+1]=
524         numberOfElementsOfTypeC[numberOfGeometricType]
525         +  numberOfElementsOfType[numberOfGeometricType];
526       numberOfGaussPoint[numberOfGeometricType+1] = ngauss;
527       geometricType[numberOfGeometricType]= *currentGeometry;
528       numberOfGeometricType++;
529       totalNumberOfElWg+=numberOfElements;
530
531     } // End Second For
532
533   } // End Premier For
534
535   if ( alreadyFoundAnEntity) {
536     support.setName(fieldName+"Support");
537     support.setMeshName(string(maa)); // Vérifier que les différents noms de maillages lus soient identiques
538     support.setEntity(entity);
539     // REM : Le nombre <numberOfGeometricType> dans la précédente version du Driver 
540     //       était erronée pour un champ qui ne reposait pas sur toutes les entités géométriques 
541     //       du maillage mais dont le SUPPORT a été crée à partir des informations du maillage
542     //       ( méthode qui était largement utilisée pour construire un SUPPORT).
543     support.setNumberOfGeometricType(numberOfGeometricType);
544     support.setGeometricType(geometricType); // Utile uniquement si setAll == false ?
545     support.setNumberOfElements(numberOfElementsOfType);    //setNumberOfElements effectue une copie 
546     // Par défaut considère que le champ repose sur tous les type géométriques du maillage
547     // Si ce n'est pas le cas les champs geometricType et numberOfElementsOfType du SUPPORT sont corrects
548     support.setAll(true);
549     numberOfElementsOfTypeC.resize(numberOfGeometricType+1);
550     numberOfGaussPoint.resize(numberOfGeometricType+1);
551
552     return alreadyFoundAnEntity;
553   } else
554     return false;
555 }
556
557 /*!
558
559   Renvoie la liste <geoType> des types géométriques définis dans le maillage <meshName>
560   pour le type d'entité <entity>.
561   * < nbOfElOfType > contient le nombre d'entités de chaque type
562   * < numberOfElementsOfTypeC > contient le nombre d'entités cumulées de chaque type
563                               avec numberOfElementsOfTypeC[0]=0;
564
565 */
566 template <class T> void
567 MED_FIELD_DRIVER22<T>::getMeshGeometricTypeFromFile(med_2_2::med_idt id,
568                                           string & meshName,
569                                           MED_EN::medEntityMesh  entity,
570                                           vector<MED_EN::medGeometryElement> & geoType,
571                                           vector<int> &nbOfElOfType,
572                                           vector<int> &nbOfElOfTypeC
573                                          ) const throw(MEDEXCEPTION)
574 {
575   const char LOC[] = "MED_FIELD_DRIVER<T>::getMeshGeometricTypeFromFile(...)";
576
577   int numberOfGeometricType=0;
578   MED_EN::medGeometryElement geometricType[MED_NBR_GEOMETRIE_MAILLE];
579   int numberOfElementsOfType [MED_NBR_GEOMETRIE_MAILLE];
580   int numberOfElementsOfTypeC[MED_NBR_GEOMETRIE_MAILLE+1];
581   med_2_2::med_int   numberOfElements=0;
582   med_2_2::med_table quoi;
583   if (entity == MED_EN::MED_CELL) quoi=med_2_2::MED_CONN;
584   else
585     if (entity == MED_EN::MED_NODE) quoi=med_2_2::MED_COOR;
586     else
587       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Support Creation from Mesh |"  << meshName
588                                    << "| on entity " << MED_EN::entNames[entity]
589                                    << "| is impossible,  must be  on MED_NODE or MED_CELL" ));
590
591   list<MED_EN::medGeometryElement>::const_iterator currentGeometry;
592   bool alreadyFoundAnEntity = false;
593   numberOfElementsOfTypeC[0]=0;
594
595   for (currentGeometry  = (MED_EN::meshEntities[entity]).begin();
596        currentGeometry != (MED_EN::meshEntities[entity]).end(); currentGeometry++) {
597
598
599     if ( (numberOfElements =
600           med_2_2::MEDnEntMaa(id,
601                               const_cast<char*> (meshName.c_str()),
602                               quoi,
603                               (med_2_2::med_entite_maillage)   entity,
604                               (med_2_2::med_geometrie_element)  *currentGeometry,
605                               med_2_2::MED_NOD) ) <= 0)
606       continue;
607
608     alreadyFoundAnEntity = true;
609     numberOfElementsOfType[numberOfGeometricType] = numberOfElements;
610     numberOfElementsOfTypeC[numberOfGeometricType+1] =
611       numberOfElementsOfTypeC[numberOfGeometricType]+numberOfElements;
612     geometricType[numberOfGeometricType] = *currentGeometry;
613     numberOfGeometricType++;
614
615   }
616
617   geoType = vector<MED_EN::medGeometryElement>(geometricType,geometricType+numberOfGeometricType);
618   nbOfElOfType = vector<int> (numberOfElementsOfType,numberOfElementsOfType+numberOfGeometricType);
619   nbOfElOfTypeC = vector<int> (numberOfElementsOfTypeC,numberOfElementsOfTypeC+numberOfGeometricType+1);
620
621 //   for (int j =0 ; j<= numberOfGeometricType;++j)
622 //       cout << "nbOfElOfTypeC["<<j<<"]="<<nbOfElOfTypeC[j]<<endl;
623
624 }
625
626 template <class T> void
627 MED_FIELD_DRIVER22<T>::getMeshGeometricTypeFromMESH( MESH * meshPtr,
628                                           MED_EN::medEntityMesh  entity,
629                                           vector<MED_EN::medGeometryElement> & geoType,
630                                           vector<int> &nbOfElOfType,
631                                           vector<int> &nbOfElOfTypeC) const throw(MEDEXCEPTION)
632 {
633   const char LOC[] = "MED_FIELD_DRIVER<T>::getMeshGeometricTypeFromMESH(...) : ";
634   BEGIN_OF(LOC);
635
636   if (!meshPtr)
637     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"ptrMesh must be non null" )); ;
638
639   // Il est plus pratique de créer un support "onAll"
640   // pour calculer les tableaux du nombre d'entités cumulées
641
642   SUPPORT mySupportFromMesh = SUPPORT(meshPtr,"Temporary Support From Associated Mesh",
643                                       entity);
644   geoType = vector<MED_EN::medGeometryElement>(mySupportFromMesh.getTypes(),
645                               mySupportFromMesh.getTypes()+mySupportFromMesh.getNumberOfTypes());
646   nbOfElOfType.resize(mySupportFromMesh.getNumberOfTypes());
647   nbOfElOfTypeC.resize(mySupportFromMesh.getNumberOfTypes()+1);
648   nbOfElOfTypeC[0]=0;
649
650   for (int j=1; j<=mySupportFromMesh.getNumberOfTypes(); ++j) {
651     nbOfElOfType[j-1]=mySupportFromMesh.getNumberOfElements(geoType[j-1]);
652     nbOfElOfTypeC[j]+=nbOfElOfTypeC[j-1]+nbOfElOfType[j-1];
653   }
654
655   END_OF(LOC);
656 }
657
658 /*--------------------- RDONLY PART -------------------------------*/
659
660 template <class T> GENDRIVER * MED_FIELD_RDONLY_DRIVER22<T>::copy(void) const
661 {
662   return new MED_FIELD_RDONLY_DRIVER22<T>(*this);
663 }
664
665 template <class T> void MED_FIELD_RDONLY_DRIVER22<T>::read(void)
666   throw (MEDEXCEPTION)
667 {
668   const char * LOC = " MED_FIELD_RDONLY_DRIVER22::read() " ;
669   BEGIN_OF(LOC);
670
671   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,NoGauss>::Array   ArrayNo;
672   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,Gauss>::Array     ArrayNoWg;
673   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,NoGauss>::Array ArrayFull;
674   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,Gauss>::Array   ArrayFullWg;
675
676   if (MED_FIELD_DRIVER<T>::_status!=MED_OPENED)
677     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": Method open must be called before method read.")) ;
678
679   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
680        ( MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
681     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
682                                  <<" neither <fieldName> is set in driver nor in object FIELD.")) ;
683
684   // If _fieldName is not set in driver, try to use _ptrfield->_fieldName
685   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
686        ( !MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
687     MED_FIELD_DRIVER<T>::_fieldName=MED_FIELD_DRIVER<T>::_ptrField->_name;
688
689   if ( MED_FIELD_DRIVER<T>::_fieldName.size() > MED_TAILLE_NOM )
690     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
691                                  <<" <fieldName> size in object driver FIELD is > MED_TAILLE_NOM ."));
692
693   const string & fieldName = MED_FIELD_DRIVER<T>::_fieldName;
694
695   MED_EN::medModeSwitch interlacingType = MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType();
696   bool isFullInterlace = ( interlacingType == MED_EN::MED_FULL_INTERLACE );
697
698   MESSAGE("###### "<<LOC<<" fieldNameDRIVER : "<< fieldName << " fieldName : "<< MED_FIELD_DRIVER<T>::_ptrField->_name);
699
700 // EF :
701 //   Si un support a été donnée au champ, pour des raisons de compatibilité avec
702 //   les versions précédentes, ce support sera utilisé pour
703 //   - Obtenir le nom du maillage sur lequel on veut lire le champ
704 //     (eventuellement on pourrait l'utiliser pour selectionner un champ qui
705 //      repose sur plusieurs maillages cf HOMARD-ASTER, ce qui n'est pas géré dans MEDMEM)
706 //   -  vérifier le type d'entité (MED_NOEUD xor MED_MAILLE xor MED_FACE xor MED_ARETE ) sur lequel
707 //      il faut lire le champ qui est également retouvé.
708 //   - Si le support défini une liste d'entité ( différente de MED_ALL_ELEMENTS), celle-ci est ignorée
709 //     à la lecture et écrasé par soit :
710 //            - onall, après avoir vérifié que la liste des types géométriques utilisés par le champ
711 //                     est égale à la liste des type géométriques définis dans le maillage associé
712 //                     pour tous le même type d'entité.
713 //            - La sous liste des types géométriques utilisés (onAll quand même, cf commenataire ci-dessous )  
714 //            - les listes de profils lus s'il en existe pour une sous liste de types
715 //              géométriques
716
717 //   Si aucun support n'a été donné au champ :
718 //   - A la lecture : Un support est crée et le type d'entité unique est lu
719 //                    (cf decision gt MED qu'un champ repose sur une entité unique ?),
720 //                    l'ensemble des types géométriques est lu,
721 //                    l'ensemble des profils par type géométrique est lu
722 //                    Le nom du maillage associé est lu mais le pointeur SUPPORT-MESH non initialisé
723
724
725   char * tmpFieldName = new char[MED_TAILLE_NOM+1] ;
726   int err ;
727   int    numberOfComponents          = 0;
728   char * componentName               = (char *) MED_NULL;
729   char * unitName                    = (char *) MED_NULL;
730   med_2_2::med_type_champ type ;
731   med_2_2::med_idt id = MED_FIELD_DRIVER22<T>::_medIdt;
732   bool needConversionToDouble = false,needConversionToInt64 = false;
733
734   // we search for the "field med number" of <fieldName>
735   // Having found <fieldName>, variables <numberOfComponents>,
736   // <componentName>, <unitname>, <type> and attribute <_fieldNum> are set.
737   if (MED_FIELD_DRIVER<T>::_fieldNum==MED_INVALID)
738     {
739       int numberOfFields = med_2_2::MEDnChamp(id,0) ;
740       if ( numberOfFields <= 0 )
741         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": There is no field found in the file !"));
742
743       for (int i=1;i<=numberOfFields;i++)
744         {
745           numberOfComponents = med_2_2::MEDnChamp(id,i) ;
746
747           if ( numberOfComponents <= 0 )
748             MESSAGE(LOC<<"Be careful there is no compound for field n°"<<i<<"in file |"<<MED_FIELD_DRIVER<T>::_fileName<<"| !");
749
750           componentName = new char[numberOfComponents*MED_TAILLE_PNOM22+1] ;
751           unitName      = new char[numberOfComponents*MED_TAILLE_PNOM22+1] ;
752
753           err = med_2_2::MEDchampInfo(id, i, tmpFieldName, &type, componentName,
754                                       unitName, numberOfComponents) ;
755
756           MESSAGE("Field "<<i<<" : #" << tmpFieldName <<"# et recherche #"<<fieldName.c_str()<<"#");
757           if ( !strcmp(tmpFieldName,fieldName.c_str()) ) {
758             MESSAGE("FOUND FIELD "<< tmpFieldName <<" : "<<i);
759             MED_FIELD_DRIVER<T>::_fieldNum = i ;
760             break ;
761           }
762           // not found : release memory and search next field !
763           delete[] componentName ;
764           delete[] unitName ;
765         }
766     }
767
768   delete[] tmpFieldName ;
769
770   // Si aucun champ ne correspond les variables <componentName> et <unitName> ont été correctement
771   // désallouées dans la boucle de recherche
772   if (MED_FIELD_DRIVER<T>::_fieldNum==MED_INVALID)
773     throw MEDEXCEPTION(LOCALIZED( STRING(LOC) << ": Field "<<  fieldName
774                                    << " not found in file " << MED_FIELD_DRIVER<T>::_fileName) );
775
776   MESSAGE ("FieldNum : "<<MED_FIELD_DRIVER<T>::_fieldNum);
777
778   if (numberOfComponents < 1) {
779     delete[] componentName; delete[] unitName;
780     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" no component found for field "
781                                  << fieldName)) ;
782   }
783
784   // Verifie que l'on essaye pas de lire un champ double dans un FIELD<int>
785   switch ( (med_2_2::med_type_champ) MED_FIELD_DRIVER<T>::_ptrField->_valueType ) {
786   case  med_2_2::MED_INT :
787   case  med_2_2::MED_INT32 :
788   case  med_2_2::MED_INT64 :
789     if ( type == ( med_2_2::MED_FLOAT64 ) ) {
790       delete[] componentName; delete[] unitName;
791       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field Type in file (" << type
792                                    <<") differs from FIELD object type (" <<
793                                    MED_FIELD_DRIVER<T>::_ptrField->_valueType << ")" )) ;
794     }
795 #if defined(IRIX64) || defined(OSF1) ||defined(VPP5000)
796     if (_ptrField->_valueType==MED_EN::MED_INT32 )
797       needConversionToInt64=true;
798 #endif
799     break;
800   case med_2_2::MED_FLOAT64 :
801     if (type != med_2_2::MED_FLOAT64)
802       needConversionToDouble=true;
803     break;
804   default:
805     break;
806   }
807
808   string meshName="";
809   MESH * ptrMesh = 0;
810   bool   haveSupport = false;
811   bool   haveMesh    = false;
812   if ( MED_FIELD_DRIVER<T>::_ptrField->getSupport() ) {
813     // Verif à faire sur la taille du meshName
814     ptrMesh = MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getMesh();
815     if ( ptrMesh) {
816       meshName =  MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getMesh()->getName() ;
817       haveMesh = true;
818     }
819     haveSupport = true;
820   }
821
822   // Cherche le type d'entité, le nombre d'entité  par type géométrique sur le type d'entité
823   // (MED_MAILLE ou MED_NOEUD uniquement car MEDMEMOIRE ne gère pas la connectivité descendante).
824   // et crée le support correspondant.
825   SUPPORT *   mySupport = new SUPPORT();
826   vector<int> numberOfElementsOfTypeC;
827   vector<int> numberOfGaussPoint;
828   int         totalNumberOfElWg=0;
829
830   bool found = createFieldSupportPart1(id,fieldName,
831                                        MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber,
832                                        MED_FIELD_DRIVER<T>::_ptrField->_orderNumber,
833                                        *mySupport, meshName,
834                                        numberOfElementsOfTypeC, numberOfGaussPoint,totalNumberOfElWg);
835
836   if ( !found ) {
837     delete mySupport; delete[] componentName; delete[] unitName;
838     MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ;
839      throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"  Can't find any entity for field |"
840                                  << fieldName
841                                  << "| with (it,or) = ("
842                                   << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
843                                  << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), on mesh "
844                                  << meshName << "|" ));
845   }
846
847
848   MED_EN::medEntityMesh entityType = mySupport->getEntity();
849   //Si un SUPPORT était donné, récupère son nom, sa description et
850   //     le pointeur du maillage associé
851   if (! haveSupport)
852     meshName = mySupport->getMeshName();
853   else {
854     if ( mySupport->getEntity() != MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getEntity() ) {
855       delete mySupport; delete[] componentName; delete[] unitName;
856       MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ;
857       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Given entity |"
858                                    << MED_EN::entNames[MED_FIELD_DRIVER<T>::_ptrField->
859                                                        getSupport()->getEntity()]
860                                    << "| for field |"
861                                    << fieldName
862                                    << "| with (it,or) = ("
863                                    << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
864                                    << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), on mesh "
865                                    << meshName << "| differs from found entity |"
866                                    << MED_EN::entNames[entityType] << "|."
867                                    ));
868     }
869     mySupport->setName( MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getName() );
870     mySupport->setMesh( MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getMesh() );
871     mySupport->setDescription(MED_FIELD_DRIVER<T>::_ptrField->getSupport()->getDescription());
872   }
873
874   vector< MED_EN::medGeometryElement >  MESHgeoType;
875   vector< int >  MESHnbOfElOfType;
876   vector< int >  MESHnbOfElOfTypeC;
877   if ( haveMesh )
878     this->getMeshGeometricTypeFromMESH(ptrMesh,entityType,MESHgeoType,
879                                        MESHnbOfElOfType,MESHnbOfElOfTypeC);
880
881   int fileHasMesh = ( med_2_2::MEDdimLire(id, const_cast<char *>(meshName.c_str())) > 0);
882   vector< MED_EN::medGeometryElement >  meshGeoType;
883   vector< int >  meshNbOfElOfType;
884   vector< int >  meshNbOfElOfTypeC;
885   // Si le maillage n'est pas trouvé les tableaux renvoyés sont vides
886   if (fileHasMesh)
887     this->getMeshGeometricTypeFromFile(id,meshName,entityType,meshGeoType,
888                                        meshNbOfElOfType,meshNbOfElOfTypeC);
889
890   SCRUTE(meshGeoType.size());
891   SCRUTE(MESHgeoType.size());
892   SCRUTE(meshNbOfElOfTypeC.size());
893   SCRUTE(MESHnbOfElOfTypeC.size());
894
895   if (meshGeoType.size() != MESHgeoType.size())
896     {
897       for (int i = 0; i<meshGeoType.size();i++)
898         MESSAGE("debug meshGeotype " << meshGeoType[i]);
899
900       for (int i = 0; i<MESHgeoType.size();i++)
901         MESSAGE("debug MESHgeoType. " << MESHgeoType[i]);
902     }
903
904   if (meshNbOfElOfTypeC.size() == MESHnbOfElOfTypeC.size())
905     {
906       for (int i = 0; i<meshNbOfElOfTypeC.size();i++)
907         MESSAGE("debug meshNbOfElOfTypeC " << meshNbOfElOfTypeC[i]);
908
909       for (int i = 0; i<MESHnbOfElOfTypeC.size();i++)
910         MESSAGE("debug MESHnbOfElOfTypeC " << MESHnbOfElOfTypeC[i]);
911     }
912
913   if (fileHasMesh && haveSupport )
914     if ( ( meshGeoType != MESHgeoType ) || (meshNbOfElOfTypeC != MESHnbOfElOfTypeC) )
915       {
916         MESSAGE("Warning MedField driver 21 while getting mesh information from file for FIELD "<< fieldName
917                 << " on entity " << MED_EN::entNames[entityType]
918                 << " with (it,or) = ("
919                 << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
920                 << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
921                 << " on mesh " << meshName
922                 << " : geometric types or number of elements by type differs from MESH object !");
923
924 //      throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information from file for FIELD "<< fieldName
925 //                                    << " on entity " << MED_EN::entNames[entityType]
926 //                                    << " with (it,or) = ("
927 //                                    << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
928 //                                    << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
929 //                                    << " on mesh " << meshName
930 //                                    << " : geometric types or number of elements by type differs from MESH object !"
931 //                                    )
932 //                         );
933       }
934
935   if ( !fileHasMesh && !haveSupport )
936     throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information for FIELD "<< fieldName
937                                   << " on entity " << MED_EN::entNames[entityType]
938                                   << " with (it,or) = ("
939                                   << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
940                                   << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
941                                   << " on mesh " << meshName
942                                   << " : SUPPORT must contain a valid MESH reference or file must contain the associated MESH."
943                                   )
944                        );
945
946
947   if (!fileHasMesh && haveSupport) {
948     meshNbOfElOfTypeC = MESHnbOfElOfTypeC;
949     meshGeoType       = MESHgeoType;
950     meshNbOfElOfType  = MESHnbOfElOfType;
951   }
952
953
954   // Test si le Support du Champ repose ou non sur toutes les entités géométriques
955   // du maillage associé et positionne ou non l'attribut onAll du SUPPORT.
956   // Il ne s'agit pas de la gestion des profils
957   vector < MED_EN::medGeometryElement > v1(  mySupport->getTypes(),
958                                              mySupport->getTypes()+mySupport->getNumberOfTypes() );
959   vector<int> v2(numberOfElementsOfTypeC.size());
960   transform(numberOfElementsOfTypeC.begin(),
961             numberOfElementsOfTypeC.end(),v2.begin(), bind2nd(plus<int>(),1));
962
963   if ( ( meshGeoType != v1 )  || meshNbOfElOfTypeC != v2 ) {
964     // ATTENTION : mySupport->setAll(false);
965     // Pb : On a envie de positionner onAll à faux si le champ n'est pas défini sur tous les
966     //      types géométriques du maillage associé.
967     //      Mais si onAll est false et si aucun profil n'est détecté par la suite,
968     //      l'attribut SUPPORT->_number est censé être positionné quand même ! Que faire ?
969     // Si on veut être compatible avec la signification première de onAll,
970     //  il faudrait créer des profils contenant toutes les entités pour chaque type géométrique
971     //  du SUPPORT  mais d'une part c'est dommage d'un point de vue de l'emcombrement mémoire
972     //  et d'autre part, à la réécriture du fichier MED on stockera des profils 
973     //  alors qu'il n'y en avait pas à l'origine (fichier MED différent après lecture/écriture) !
974     // Si on laisse setAll à vrai il faut être sûr que les utilisateurs prennent les
975     //  informations sur les types gémétrique au niveau du support et non pas du maillage.
976     // Solution : Signification du onAll -> onAllElements des type géométriques définis
977     // dans le SUPPORT et non du maillage associé (dans la plupart des cas si le fichier ne
978     // contient pas de profil, le champ est défini sur toutes les entités de tous les types
979     // géométriques définis dans le maillage).
980   }
981
982
983   // If an error occurs while reading the field, these allocated FIELD member will be deleted
984
985   MED_FIELD_DRIVER<T>::_ptrField->_name                   = fieldName;
986   MED_FIELD_DRIVER<T>::_ptrField->_numberOfComponents     = numberOfComponents ;
987   MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes        = new int   [numberOfComponents] ;
988   MED_FIELD_DRIVER<T>::_ptrField->_componentsNames        = new string[numberOfComponents] ;
989   MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits        = new UNIT  [numberOfComponents] ;
990   MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions = new string[numberOfComponents] ;
991   MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits     = new string[numberOfComponents] ;
992   for (int i=0; i<numberOfComponents; i++) {
993       MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes[i] = 1 ;
994       MED_FIELD_DRIVER<T>::_ptrField->_componentsNames[i] = string(componentName,i*MED_TAILLE_PNOM22,MED_TAILLE_PNOM22) ;
995       SCRUTE(MED_FIELD_DRIVER<T>::_ptrField->_componentsNames[i]);
996       MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits[i] = string(unitName,i*MED_TAILLE_PNOM22,MED_TAILLE_PNOM22) ;
997       SCRUTE(MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits[i]);
998   }
999   delete[] componentName;
1000   delete[] unitName;
1001
1002   int NumberOfTypes                       = mySupport->getNumberOfTypes() ;
1003   const MED_EN::medGeometryElement *types = mySupport->getTypes() ;
1004   T * myValues = new T[totalNumberOfElWg*numberOfComponents];
1005   const int * nbOfElOfType = mySupport->getNumberOfElements() ;
1006   bool anyProfil = false;
1007   int  pflSize=0,index=0;
1008   // Le vecteur de profil est dimensionné par rapport aux nombres de types
1009   // géométriques du champ même si le champ n'a pas de profil MED FICHIER sur
1010   // tous ses types géométriques car dans MEDMEM si onAllElement 
1011   // du SUPPORT est false il faut positionner un profil pour tous les types géométriques 
1012   // du SUPPORT
1013   int profilSizeC = 0;
1014   vector < int   >        profilSize    (NumberOfTypes,0);
1015   vector < vector<int>  > profilList    (NumberOfTypes);
1016   vector < string >       profilNameList(NumberOfTypes);
1017   char * profilName = new char[MED_TAILLE_NOM+1];
1018
1019   MESSAGE ("NumberOfTypes      : "<< NumberOfTypes);
1020   MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues=0 ;
1021
1022
1023   for (int typeNo=0; typeNo<NumberOfTypes; typeNo++) {
1024
1025     int numberOfValuesWc= nbOfElOfType[typeNo]*numberOfGaussPoint[typeNo+1]*numberOfComponents;
1026     char * gaussModelName = new char[MED_TAILLE_NOM+1];
1027
1028     MESSAGE ("FIELD_NAME         : "<< fieldName.c_str());
1029     MESSAGE ("MESH_NAME          : "<< meshName.c_str());
1030     MESSAGE ("MED_ENTITE         : "<< MED_EN::entNames[entityType]);
1031     MESSAGE ("MED_GEOM           : "<< MED_EN::geoNames[types[typeNo]]);
1032     MESSAGE ("Iteration          : "<< MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber());
1033     MESSAGE ("Order              : "<< MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber());
1034     MESSAGE ("NumberOfElements   : "<< nbOfElOfType[typeNo]);
1035     MESSAGE ("NumberOfComponents : "<< numberOfComponents);
1036     MESSAGE ("NumberOfGaussPts   : "<< numberOfGaussPoint[typeNo+1]);
1037     MESSAGE ("NumberOfValuesWg   : "<< nbOfElOfType[typeNo]*numberOfGaussPoint[typeNo+1]);
1038     MESSAGE ("NumberOfValuesWgWc : "<< numberOfValuesWc);
1039     MESSAGE ("Index              : "<< index);
1040     med_2_2::med_err ret=-1;
1041
1042     med_2_2::med_int * myValuesTmp=0;
1043     unsigned char* ptrTmp=0;
1044     if (needConversionToDouble || needConversionToInt64 ) {
1045       myValuesTmp = new med_2_2::med_int[numberOfValuesWc];
1046       ptrTmp = (unsigned char*) myValuesTmp;
1047     } else
1048       ptrTmp = (unsigned char*) &myValues[index];
1049
1050     //VERIFIER LE NBRE
1051     ret=med_2_2::MEDMEMchampLire(id,const_cast <char*> (meshName.c_str() ),
1052                                  const_cast <char*> (fieldName.c_str()),
1053                                 (unsigned char*) ptrTmp,
1054                                 med_2_2::MED_FULL_INTERLACE,
1055                                 MED_ALL,
1056                                 gaussModelName,
1057                                 profilName,
1058                                 med_2_2::MED_COMPACT,
1059                                 (med_2_2::med_entite_maillage) entityType,
1060                                 (med_2_2::med_geometrie_element)types[typeNo],
1061                                 MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber(),
1062                                 MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()
1063                                 );
1064
1065       if (needConversionToDouble || needConversionToInt64 ) {
1066
1067       if (needConversionToInt64 )  //utiliser un trait
1068         for(int i=0;i<numberOfValuesWc;++i)
1069           myValues[index+i]=(int)(myValuesTmp[i]);
1070       else
1071         for(int i=0;i<numberOfValuesWc;++i)
1072           myValues[index+i]=myValuesTmp[i];
1073       delete[] myValuesTmp;
1074     }
1075
1076     if (ret < 0)
1077       {
1078         // The Field can't be read then we must delete all previously allocated members in FIELD
1079         //for(int j=0; j<=i;j++)
1080         //  delete[] myValues[j];
1081         delete[] myValues;
1082         //delete[] NumberOfValues ;
1083         delete[] profilName;
1084         delete[] gaussModelName;
1085         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes ;
1086         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsNames ;
1087         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits ;
1088         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions ;
1089         delete[] MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits ;
1090         MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes = NULL ;
1091         MED_FIELD_DRIVER<T>::_ptrField->_componentsNames = NULL ;
1092         MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits = NULL ;
1093         MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions = NULL ;
1094         MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits = NULL ;
1095         MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ; // we have not found right field, so reset the field number
1096         throw MEDEXCEPTION( LOCALIZED( STRING(LOC) <<": ERROR while reading values")) ;
1097       }
1098     index += numberOfValuesWc;
1099     // Le support prend en compte le nombre de valeurs lié aux profils
1100     MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues+=
1101       nbOfElOfType[typeNo];// Ne doit pas prendre en compte les points de Gauss
1102
1103     // second et troisième test lié à un bug medfichier
1104     if ( strcmp(gaussModelName,MED_NOGAUSS) && strcmp(gaussModelName,string(MED_TAILLE_NOM,' ').c_str() )
1105          && strcmp(gaussModelName,string(16,' ').c_str() )  ) {
1106  
1107         int type_geo = (int) types[typeNo];
1108         int t1       = (type_geo%100)*(type_geo/100);
1109         int ngauss   = numberOfGaussPoint[typeNo+1];
1110         int t2       = ngauss*(type_geo/100);
1111         med_2_2::med_float * refcoo = new med_2_2::med_float[t1];
1112         med_2_2::med_float * gscoo  = new med_2_2::med_float[t2];
1113         med_2_2::med_float * wg     = new med_2_2::med_float[ngauss];
1114
1115         if (MEDgaussLire(id, refcoo, gscoo, wg, (med_2_2::med_mode_switch) interlacingType,
1116                          gaussModelName ) < 0)
1117           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while reading Gauss Model |"
1118                                       << gaussModelName << "| for FIELD "<< fieldName
1119                                       << " on geometric type " << MED_EN::geoNames[types[typeNo]]
1120                                       )
1121                            );
1122         if (isFullInterlace ) { //serait inutile avec un driver template du type d'entrelacement
1123           GAUSS_LOCALIZATION<FullInterlace> * loc;
1124           loc = new GAUSS_LOCALIZATION<FullInterlace>(gaussModelName,types[typeNo],ngauss, refcoo,gscoo, wg);
1125           MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]]=loc;
1126         } else {
1127           GAUSS_LOCALIZATION<NoInterlace> * loc;
1128           loc = new GAUSS_LOCALIZATION<NoInterlace>(gaussModelName,types[typeNo],ngauss, refcoo,gscoo, wg);
1129           MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]]=loc;
1130         }
1131 //      cout << *MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]] << endl;
1132         delete [] refcoo;delete [] gscoo; delete [] wg;
1133     }
1134     delete[] gaussModelName ;
1135
1136     if ( strcmp(profilName,MED_NOPFL) ) {
1137       anyProfil = true;
1138       pflSize = med_2_2::MEDnValProfil(id,profilName);
1139       if ( pflSize  <= 0)
1140         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Error while reading the profil size of |"
1141                                      << profilName << "|" ));
1142
1143       profilSize[typeNo]=pflSize;
1144       profilList[typeNo].resize(pflSize);
1145       ret = med_2_2::MEDprofilLire(id,&profilList[typeNo][0],profilName); // cf item 16 Effective STL
1146       profilNameList[typeNo]=string(profilName);
1147     }
1148   }
1149
1150   delete[] profilName ;
1151
1152   //MESSAGE ("Index              : "<< index);
1153   assert(index == totalNumberOfElWg*numberOfComponents);
1154   assert(MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues ==  mySupport->getNumberOfElements(MED_ALL_ELEMENTS));
1155
1156   if (anyProfil) {
1157
1158     for (int typeNo=0; typeNo < NumberOfTypes; typeNo++) {
1159
1160       // Trouve l'index du type géométrique dans la liste des types géométriques du maillage
1161       // correspondant au type géométrique du champ traité
1162       vector<MED_EN::medGeometryElement>::iterator meshTypeNoIt =
1163         find(meshGeoType.begin(),meshGeoType.end(),types[typeNo]); //Gérer l'exception
1164       if ( meshTypeNoIt ==  meshGeoType.end() )
1165         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Can't find "<< MED_EN::geoNames[types[typeNo]]
1166                                       << " on entity " << MED_EN::entNames[entityType]
1167                                       << " in geometric type list of mesh " << meshName
1168                                       )
1169                            );
1170       int meshTypeNo = meshTypeNoIt -  meshGeoType.begin();
1171
1172       if (! profilList[typeNo].empty() ) {
1173
1174 //      for (int j =0 ; j< meshGeoType.size();++j)
1175 //        cout << "--MeshTypeNo : "<<meshTypeNo<<"-> meshNbOfElOfTypeC["<<j<<"]="<<meshNbOfElOfTypeC[j]<<endl;
1176 //      cout << "--typeNo--" << typeNo << endl;
1177 //      cout << "meshNbOfElOfTypeC["<<meshTypeNo<<"]=" << meshNbOfElOfTypeC[meshTypeNo] <<endl;
1178
1179         // Transformer les numéros locaux d'entités medfichier en numéro global medmémoire
1180         for (int i = 0 ; i < profilList[typeNo].size(); i++ ) {
1181         // Les numéros des entités commencent à 1 dans MEDfichier comme dans MEDmémoire
1182         // meshNbOfElOfTypeC[0]=0 ...meshNbOfEltOfTypeC[meshTypeNo]=
1183         // meshNbOfElOfTypeC[meshTypeNo-1]+nbrOfElem of meshTypeNo type
1184         // rem1 : Si le meshTypeNo trouvé est 0 (premier type géométrique du maillage
1185         // il ne faut pas décaler les numéros du profils qui commencent à 1 dans MEDFICHIER
1186         // rem2 : meshNbOfElOfTypeC[NumberOfTypes] ne devrait jamais être utilisé
1187           profilList[typeNo][i]+=meshNbOfElOfTypeC[meshTypeNo];
1188         }
1189       } else {
1190         // Créer le profil <MED_ALL> pour ce type géométrique
1191         // uniquement pour renseigner le tableau skyline avec des accesseurs directs
1192         // par type géométriques
1193         // REM : Une conséquence est qu'à la réecriture le fichier contiendra des
1194         // profils sur certains types géométriques alors qu'à la lecture il n'y en avait pas !
1195         // Solution : Stocker les noms des profils et les utiliser pour savoir si il y avait ou non
1196         //            un profil
1197         int pflSize   = meshNbOfElOfType[meshTypeNo];
1198         // profil    = new int[pflSize];
1199
1200         profilList[typeNo].resize(pflSize);
1201         profilSize[typeNo]=pflSize;
1202
1203         for (int j = 1; j <= pflSize; j++) {
1204           profilList[typeNo][j-1]=meshNbOfElOfTypeC[meshTypeNo] + j ; // index MEDMEM commence à 1
1205         }
1206         profilNameList[typeNo] = MED_NOPFL; //Information a utiliser pour la sauvegarde : PLUTOT MED_ALL
1207       }
1208       profilSizeC+=profilList[typeNo].size();
1209     }
1210
1211     MEDSKYLINEARRAY * skyLine = new MEDSKYLINEARRAY(profilList.size(), profilSizeC );
1212     vector<int> index(NumberOfTypes+1,0);
1213     index[0]=1;
1214     for( int typeNo=0; typeNo < NumberOfTypes; typeNo++ )
1215       index[typeNo+1]=index[typeNo]+profilSize[typeNo];
1216     skyLine->setIndex(&index[0]);
1217     for (int i=1; i <= profilList.size() ; i++)
1218       skyLine->setI(i,&profilList[i-1][0]);
1219
1220     mySupport->setAll(false);
1221     mySupport->setpartial(skyLine,true);
1222     mySupport->setProfilNames(profilNameList);
1223 //    cout << "Valeurs du skyline du SUPPORT partiel crée : " << *skyLine << endl;
1224   }
1225
1226   // Créer un driver spécifique pour les modes MED_FULL_INTERLACE et MED_NO_INTERLACE
1227   // serait plus efficace.
1228   bool anyGauss = (numberOfGaussPoint != vector<int>(numberOfGaussPoint.size(),1));
1229   SCRUTE(anyGauss);
1230   MEDMEM_Array_ * Values;
1231   if (anyGauss) {
1232     SCRUTE(mySupport->getNumberOfElements(MED_ALL_ELEMENTS) );
1233     SCRUTE(NumberOfTypes);
1234     SCRUTE(numberOfElementsOfTypeC[NumberOfTypes]-1);
1235     assert(mySupport->getNumberOfElements(MED_ALL_ELEMENTS) == (numberOfElementsOfTypeC[NumberOfTypes]-1) );
1236     Values = new ArrayFullWg(myValues,
1237                              numberOfComponents,
1238                              numberOfElementsOfTypeC[NumberOfTypes]-1,
1239                              // Up : Prend en compte les profils et
1240                              // Ne prend pas en compte le nbre de composantes et
1241                              // le nombre de points de Gauss
1242                              NumberOfTypes,
1243                              &numberOfElementsOfTypeC[0],
1244                              &numberOfGaussPoint[0],
1245                              true,true);
1246 //     cout << "Valeurs du ArrayFullWg crée : " << endl <<
1247 //       *(static_cast<ArrayFullWg*>(Values))  << endl;
1248   } else
1249     Values = new ArrayFull(myValues,numberOfComponents,totalNumberOfElWg,
1250                                        true,true);
1251   if (MED_FIELD_DRIVER<T>::_ptrField->_value != NULL)
1252     delete MED_FIELD_DRIVER<T>::_ptrField->_value;
1253
1254   if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_NO_INTERLACE )
1255     {
1256       if (Values->getGaussPresence())
1257         MED_FIELD_DRIVER<T>::_ptrField->_value=ArrayConvert(*static_cast<ArrayFullWg*>(Values));
1258       else
1259         MED_FIELD_DRIVER<T>::_ptrField->_value=ArrayConvert(*static_cast<ArrayNo*    >(Values));
1260       delete Values;
1261     }
1262   else
1263     MED_FIELD_DRIVER<T>::_ptrField->_value=Values;
1264
1265   MED_FIELD_DRIVER<T>::_ptrField->_isRead = true ;
1266
1267   MED_FIELD_DRIVER<T>::_ptrField->_support=mySupport; //Prévenir l'utilisateur ?
1268
1269
1270   END_OF(LOC);
1271 }
1272
1273 template <class T> void MED_FIELD_RDONLY_DRIVER22<T>::write( void ) const
1274   throw (MEDEXCEPTION)
1275 {
1276   throw MEDEXCEPTION("MED_FIELD_RDONLY_DRIVER22::write : Can't write with a RDONLY driver !");
1277 }
1278
1279 /*--------------------- WRONLY PART -------------------------------*/
1280
1281 template <class T> GENDRIVER * MED_FIELD_WRONLY_DRIVER22<T>::copy(void) const
1282 {
1283   return new MED_FIELD_WRONLY_DRIVER22<T>(*this);
1284 }
1285
1286 template <class T> void MED_FIELD_WRONLY_DRIVER22<T>::read (void)
1287   throw (MEDEXCEPTION)
1288 {
1289   throw MEDEXCEPTION("MED_FIELD_WRONLY_DRIVER22::read : Can't read with a WRONLY driver !");
1290 }
1291
1292 template <class T> void MED_FIELD_WRONLY_DRIVER22<T>::write(void) const
1293   throw (MEDEXCEPTION)
1294 {
1295   const char * LOC = "MED_FIELD_WRONLY_DRIVER22::write(void) const " ;
1296   BEGIN_OF(LOC);
1297   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,NoGauss>::Array   ArrayNo;
1298   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,Gauss>::Array     ArrayNoWg;
1299   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,NoGauss>::Array ArrayFull;
1300   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,Gauss>::Array   ArrayFullWg;
1301
1302   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION<FullInterlace>*> locMapFull;
1303   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION<NoInterlace>*>   locMapNo;
1304   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION_*>   locMap;
1305
1306   med_2_2::med_idt id = MED_FIELD_DRIVER22<T>::_medIdt;
1307
1308   if (MED_FIELD_DRIVER<T>::_status!=MED_OPENED)
1309     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": Method open must be called before method read.")) ;
1310
1311   string fieldName;
1312   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
1313        ( MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
1314     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1315                                  <<" neither <fieldName> is set in driver nor in object FIELD.")) ;
1316
1317   // If _fieldName is not set in driver, try to use _ptrfield->_fieldName
1318   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
1319        ( !MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
1320     fieldName=MED_FIELD_DRIVER<T>::_ptrField->_name;
1321   else
1322     fieldName = MED_FIELD_DRIVER<T>::_fieldName;
1323
1324   SCRUTE(fieldName);
1325   if ( fieldName.size() > MED_TAILLE_NOM ) {
1326     fieldName.substr(0,MED_TAILLE_NOM);
1327     MESSAGE( "Be careful <fieldName> size must not be > MED_TAILLE_NOM, using fieldName : |"<< fieldName <<"|." );
1328   }
1329
1330   const SUPPORT * mySupport = MED_FIELD_DRIVER<T>::_ptrField->getSupport() ;
1331   if ( ! mySupport )
1332     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1333                                  <<" There is no SUPPORT associated with FIELD : "
1334                                  << fieldName << "."));
1335
1336   bool onAll = mySupport->isOnAllElements();
1337   const locMap & gaussModel = MED_FIELD_DRIVER<T>::_ptrField->_gaussModel;
1338
1339
1340   string meshName = mySupport->getMeshName();
1341   SCRUTE(meshName);
1342   if ( meshName.size() > MED_TAILLE_NOM ) {
1343     meshName = meshName.substr(0,MED_TAILLE_NOM);
1344     MESSAGE( "Be careful <meshName> size must not be > MED_TAILLE_NOM, using meshName : |"<< meshName <<"|." );
1345   }
1346   MED_EN::medEntityMesh entityType = mySupport->getEntity();
1347
1348   // Reconstruit les listes contigues des noms de composantes et des unités
1349   // Les noms sont tronqués à MED_TAILLE_PNOM22
1350   int err ;
1351   int component_count=MED_FIELD_DRIVER<T>::_ptrField->getNumberOfComponents();
1352   string   component_name(component_count*MED_TAILLE_PNOM22,' ') ;
1353   string   component_unit(component_count*MED_TAILLE_PNOM22,' ') ;
1354
1355   const string * listcomponent_name=MED_FIELD_DRIVER<T>::_ptrField->getComponentsNames() ;
1356   const string * listcomponent_unit=MED_FIELD_DRIVER<T>::_ptrField->getMEDComponentsUnits() ;
1357   int length ;
1358   for (int i=0; i < component_count ; i++) {
1359     length = min(MED_TAILLE_PNOM22,(int)listcomponent_name[i].size());
1360     component_name.replace(i*MED_TAILLE_PNOM22,length,
1361                            listcomponent_name[i],0,length);
1362     length = min(MED_TAILLE_PNOM22,(int)listcomponent_unit[i].size());
1363     component_unit.replace(i*MED_TAILLE_PNOM22,length,
1364                            listcomponent_unit[i],0,length);
1365   }
1366
1367   MESSAGE("using component_name=|"<<component_name<<"|");
1368   MESSAGE("using component_unit=|"<<component_unit<<"|");
1369
1370   MED_EN::med_type_champ ValueType=MED_FIELD_DRIVER<T>::_ptrField->getValueType() ;
1371
1372   MESSAGE("Template Type =|"<<ValueType<<"|");
1373
1374   // Vérifier si le champ existe déjà
1375   char   champName[MED_TAILLE_NOM+1];
1376   char * compName, * compUnit ;
1377   med_2_2::med_type_champ type ;
1378   bool Find = false ;
1379   int n = med_2_2::MEDnChamp(id,0);
1380   int nbComp = 0;
1381   for (int i=1; i<=n; i++) {
1382     nbComp   = med_2_2::MEDnChamp(id,i);
1383     compName = new char[MED_TAILLE_PNOM22*nbComp+1];
1384     compUnit = new char[MED_TAILLE_PNOM22*nbComp+1];
1385     err = med_2_2::MEDchampInfo(id,i,champName,&type,compName,compUnit,nbComp);
1386     if (err == 0)
1387       if (!strcmp(champName,fieldName.c_str()) ) {
1388         Find = true ;
1389         break ;
1390       }
1391     delete[] compName ;
1392     delete[] compUnit ;
1393   }
1394
1395   if (Find) {
1396     // the same ?
1397     if (nbComp != component_count)
1398       throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1399                                      <<": Field exist in file, but number of component are different : "
1400                                      <<nbComp<<" in file and "
1401                                      <<component_count<<" in memory."
1402                                      )
1403                           );
1404     // component name and unit
1405     SCRUTE(nbComp);
1406     MESSAGE(LOC<<" Component name in file   : "<<compName);
1407     MESSAGE(LOC<<" Component name in memory : "<<component_name);
1408     MESSAGE(LOC<<" Component unit in file   : "<<compUnit);
1409     MESSAGE(LOC<<" Component unit in memory : "<<component_unit);
1410     delete[] compName ;
1411     delete[] compUnit ;
1412
1413   } else {
1414     // Verify the field isn't yet created
1415
1416     string dataGroupName =  "/CHA/";
1417     dataGroupName        += fieldName;
1418     MESSAGE(LOC << "|" << dataGroupName << "|" );
1419     med_2_2::med_idt gid =  H5Gopen(id, dataGroupName.c_str() );
1420
1421     if ( gid < 0 ) {
1422       // create field :
1423       err=med_2_2::MEDchampCr(id,
1424                               const_cast <char*> (fieldName.c_str()),
1425                               (med_2_2::med_type_champ) ValueType,
1426                               const_cast <char*> ( component_name.c_str() ),
1427                               const_cast <char*> ( component_unit.c_str() ),
1428                               component_count);
1429       if ( err < 0 )
1430         throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1431                                        << ": Error MEDchampCr : "<<err
1432                                        )
1433                             );
1434     }
1435     else H5Gclose(gid);
1436   }
1437
1438
1439
1440   // On s'assure que le champ est dans le bon type d'entrelacement.
1441   // REM : Il faudrait un driver par type d'entrelacement, ce qui eviterait
1442   // de doubler l'utilisation de la taille mémoire si le champ n'est pas dans
1443   // le bon mode.
1444   FIELD<T,FullInterlace> * myField = 0;
1445   if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_FULL_INTERLACE )
1446     myField = MED_FIELD_DRIVER<T>::_ptrField;
1447   else
1448     myField = FieldConvert( *( dynamic_cast< FIELD<T,NoInterlace> * > (MED_FIELD_DRIVER<T>::_ptrField )
1449                                )
1450                             );
1451
1452
1453   // Il est necessaire de calculer le tableau
1454   // du nombre d'entités cumulées de chaque type géométrique du maillage
1455   // pour convertir les profils de la numérotation globale
1456   // à la numérotation locale au type géométrique.
1457   // Pour celà on établit ce tableau à partir de l'objet MESH si la relation SUPPORT-MESH existe.
1458   // Si le maillage existe dans le fichier MED  on essaye également de reconstituer ce tableau
1459   // pour vérifier la cohérence des informations.
1460   // Si la relation SUPPRT-MESH n'esiste pas  on constitue le tableau uniquement à partir du fichier qui
1461   // doit alors obligatoirement contenir le maillage.
1462   const int * number, *numberIndex = 0;
1463   string         profilName;
1464   vector<string> profilNameList;
1465   vector<MED_EN::medGeometryElement> meshGeoType;
1466   vector<int> meshNbOfElOfType;
1467   vector<int> meshNbOfElOfTypeC;
1468   vector<MED_EN::medGeometryElement> fileMeshGeoType;
1469   vector<int> fileMeshNbOfElOfType;
1470   vector<int> fileMeshNbOfElOfTypeC;
1471   med_2_2::med_int fileHasMesh=0;
1472
1473   if (!onAll) {
1474
1475     number = mySupport->getNumber(MED_ALL_ELEMENTS);
1476     numberIndex = mySupport->getNumberIndex();
1477     profilNameList=mySupport->getProfilNames();
1478
1479     fileHasMesh = ( med_2_2::MEDdimLire(id, const_cast<char *>(meshName.c_str())) > 0);
1480     MESH * meshPtr = mySupport->getMesh();
1481
1482     if (fileHasMesh)
1483       this->getMeshGeometricTypeFromFile(id, meshName,
1484                                          entityType,
1485                                          fileMeshGeoType,fileMeshNbOfElOfType,fileMeshNbOfElOfTypeC);
1486
1487     if (meshPtr) {
1488       this->getMeshGeometricTypeFromMESH( meshPtr, entityType,meshGeoType,
1489                                           meshNbOfElOfType,
1490                                           meshNbOfElOfTypeC);
1491
1492       if (fileHasMesh)
1493         if ( ( fileMeshGeoType != meshGeoType ) || (fileMeshNbOfElOfTypeC != meshNbOfElOfTypeC) )
1494           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information from file for FIELD "<< fieldName
1495                                         << " on entity " << MED_EN::entNames[entityType]
1496                                         << " with (it,or) = ("
1497                                         << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1498                                         << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1499                                         << " on mesh " << meshName
1500                                         << " : geometric types or number of elements by type differs from MESH object !"
1501                                         )
1502                              );
1503
1504     }
1505
1506     if ( !fileHasMesh && meshPtr==0 )
1507       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information for FIELD "<< fieldName
1508                                     << " on entity " << MED_EN::entNames[entityType]
1509                                     << " with (it,or) = ("
1510                                     << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1511                                     << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1512                                     << " on mesh " << meshName
1513                                     << " : SUPPORT must contain a valid MESH reference or file must contain the associated MESH."
1514                                     )
1515                          );
1516
1517
1518     if (fileHasMesh && !meshPtr) {
1519       meshNbOfElOfTypeC = fileMeshNbOfElOfTypeC;
1520       meshGeoType       = fileMeshGeoType;
1521       meshNbOfElOfType  = fileMeshNbOfElOfType;
1522     }
1523   }
1524
1525   const MED_EN::medGeometryElement * types = mySupport->getTypes() ;
1526   int numberOfTypes = mySupport->getNumberOfTypes() ;
1527   int numberOfElForMED = -1;
1528   const T   * value   = NULL;
1529   int index = 1 ;
1530   // on boucle sur tout les types pour ecrire les tableaux de valeur
1531   for (int typeNo=0;typeNo<numberOfTypes;typeNo++) {
1532
1533     int numberOfElements = mySupport->getNumberOfElements(types[typeNo]) ;
1534     //UP : prend en compte les profils, pas les points de Gauss
1535
1536     //value = MED_FIELD_DRIVER<T>::_ptrField->getRow(index) ;
1537     // rem 1 : le getRow du Array est différent de celui du FIELD si le SUPPORT contient
1538     //         des profils (les indices des valeurs ne se suivent pas forcément)
1539     // rem 2 : Afin de respecter la norme MEDFICHIER, les indices contenus dans les
1540     //         profils doivent être croissant
1541     if (onAll) {
1542       value = myField->getRow(index);
1543       profilName=MED_NOPFL;
1544       numberOfElForMED = numberOfElements;
1545     } else {
1546       value = myField->getRow(number[index-1]);
1547       profilName = profilNameList[typeNo].substr(0,MED_TAILLE_NOM);
1548       // Rem : Si le SUPPORT n'est pas onAll mais que pour un type géométrique donné le nom
1549       // du profil associé est MED_NOPFL alors le profil n'est pas écrit dans le fichier MED.
1550       // Car en MEDMEMOIRE si le champ repose sur des éléments de deux types géométriques
1551       // différents et est défini sur tous les éléments d'un type géométrique
1552       // mais pas de l'autre, il existe tout de même des profils sur les deux types géométriques.
1553       // Ce n'est pas le cas en MEDFICHIER.
1554       vector<int> profil(&number[index-1],&(number[index-1])+numberOfElements);
1555
1556       // Trouve l'index du type géométrique dans la liste des types géométriques du maillage
1557       // correspondant au type géométrique du champ en cours de traitement
1558       vector<MED_EN::medGeometryElement>::iterator meshTypeNoIt =
1559         find(meshGeoType.begin(),meshGeoType.end(),types[typeNo]);
1560       if ( meshTypeNoIt ==  meshGeoType.end() )
1561         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Can't find "<< MED_EN::geoNames[types[typeNo]]
1562                                       << " on entity " << MED_EN::entNames[entityType]
1563                                       << " in geometric type list of mesh " << meshName
1564                                       )
1565                            );
1566
1567       int meshTypeNo = meshTypeNoIt -  meshGeoType.begin();
1568
1569       if ( profilName == MED_NOPFL && profil.size() != meshNbOfElOfType[meshTypeNo] )
1570         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating profil for FIELD "<< fieldName 
1571                                       << " on entity " << MED_EN::entNames[entityType]
1572                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
1573                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1574                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1575                                       << profilName << " on mesh " << meshName
1576                                       << " : There is no profileName but profilsize (" <<profil.size()
1577                                       << ") differs from number of elements in associated MESH ("
1578                                       << meshNbOfElOfType[meshTypeNo] << ")."
1579                             )
1580                          );
1581
1582       //REM : Ce n'est pas évident, mais lorsqu'il y a un profil, le nombre de valeurs
1583       //      que l'on indique à MEDchampEcr est le nombre de valeurs sans profil, d'où
1584       //      le nombre d'éléments du maillage sur le type géométrique courant.
1585       numberOfElForMED = meshNbOfElOfType[meshTypeNo];
1586
1587       for (int ind=0;ind < numberOfElements;++ind) {
1588 //      cout << "number["<<index-1<<"]="<<number[index-1]<<endl;
1589 //      cout << "profil1["<<ind<<"]="<<profil[ind]<<endl;
1590         profil[ind]-=meshNbOfElOfTypeC[meshTypeNo];
1591 //      cout << "profil2["<<ind<<"]="<<profil[ind]<<endl;
1592       }
1593
1594       if ( profil[numberOfElements-1] > numberOfElForMED )
1595         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating profil for FIELD "<< fieldName 
1596                                       << " on entity " << MED_EN::entNames[entityType]
1597                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
1598                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1599                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1600                                       << profilName << " on mesh " << meshName
1601                                       << " : profil["<<numberOfElements-1<<"]=" << profil[numberOfElements-1]
1602                                       << " must not be superior to field size without profil : "
1603                                       << numberOfElForMED
1604                             )
1605                          );
1606
1607       if ( med_2_2::MEDMEMprofilEcr(id,
1608                                  &profil[0],
1609                                  numberOfElements,
1610                                  const_cast<char *>(profilName.c_str())) < 0)
1611
1612         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing "<< numberOfElements
1613                                       << " values for MED profil "<< profilName
1614                                     )
1615                          );
1616     }
1617
1618
1619     string locName=MED_NOGAUSS;
1620     if (myField->getGaussPresence()) {
1621 //       cout << endl << "Nombre de points de Gauss à l'écriture de " << fieldName
1622 //         << " pour le type géométrique : " << MED_EN::geoNames[types[typeNo]]
1623 //         << " : " << myField->getNumberOfGaussPoints(types[typeNo]) << endl;
1624 //       cout << *mySupport << endl;
1625
1626       const GAUSS_LOCALIZATION_ * locPtr=0;
1627       locMap::const_iterator it;
1628       if ( ( it = gaussModel.find(types[typeNo])) != gaussModel.end() )
1629         locPtr = (*it).second;
1630       else
1631         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating Gauss Model for FIELD "<< fieldName 
1632                                       << " on entity " << MED_EN::entNames[entityType]
1633                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
1634                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1635                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1636                                       << profilName << " on mesh " << meshName
1637                                       << " : Can't find a Gauss localisation model for this geometric type" 
1638                             )
1639                          );
1640
1641       int ngauss = -1;
1642       if ( locPtr->getInterlacingType() == MED_EN::MED_FULL_INTERLACE ) {
1643         const GAUSS_LOCALIZATION<FullInterlace> & loc=*(static_cast<const GAUSS_LOCALIZATION<FullInterlace> * >(locPtr));
1644         ngauss = loc.getNbGauss();
1645         locName=loc.getName();
1646         err=med_2_2::MEDMEMgaussEcr(id,
1647                                (med_2_2::med_geometrie_element) loc.getType(),
1648                                (med_2_2::med_float *)           loc.getRefCoo().getPtr(),
1649                                                                 med_2_2::MED_FULL_INTERLACE,
1650                                (med_2_2::med_int)               ngauss,
1651                                (med_2_2::med_float *)           loc.getGsCoo().getPtr(),
1652                                (med_2_2::med_float *)           (&loc.getWeight()[0]),
1653                                const_cast<char *>               (locName.c_str())
1654                                );
1655       } else {
1656         const GAUSS_LOCALIZATION<NoInterlace> & loc=*(static_cast<const GAUSS_LOCALIZATION<NoInterlace> * >(locPtr));
1657         ngauss = loc.getNbGauss();
1658         locName=loc.getName();
1659         err=med_2_2::MEDMEMgaussEcr(id,
1660                                (med_2_2::med_geometrie_element) loc.getType(),
1661                                (med_2_2::med_float *)           loc.getRefCoo().getPtr(),
1662                                                                 med_2_2::MED_NO_INTERLACE,
1663                                (med_2_2::med_int)               ngauss,
1664                                (med_2_2::med_float *)           loc.getGsCoo().getPtr(),
1665                                (med_2_2::med_float *)           (&loc.getWeight()[0]),
1666                                const_cast<char *>               (locName.c_str())
1667                                );
1668
1669       };
1670
1671       if ( err != 0 )
1672         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing Gauss Model for FIELD "<< fieldName
1673                                       << " on entity " << MED_EN::entNames[entityType]
1674                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] 
1675                             )
1676                          );
1677
1678         //numberOfElForMED *= mySupport->getNumberOfGaussPoints(types[typeNo]); //Deplacer la méthode dans FIELD
1679         numberOfElForMED *= ngauss;
1680     }
1681
1682     MESSAGE("MED_FIELD_DRIVER22<T>::_medIdt                       : "<<id);
1683     MESSAGE("meshName.c_str()                : "<<meshName.c_str());
1684     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getName()            : "<<MED_FIELD_DRIVER<T>::_ptrField->getName());
1685     MESSAGE("MED_FIELD_DRIVER<T>::_fieldName                      : "<<MED_FIELD_DRIVER<T>::_fieldName);
1686     MESSAGE("value                           : "<<value);
1687     MESSAGE("numberOfElements                : "<<numberOfElements);
1688     MESSAGE("numberOfElForMED                : "<<numberOfElForMED);
1689     MESSAGE("entityType                      : "<<MED_EN::entNames[entityType]);
1690     MESSAGE("types[i]                        : "<<MED_EN::geoNames[types[typeNo]]);
1691     MESSAGE("NumberOfGaussPoint[i]           : "<<myField->getNumberOfGaussPoints(types[typeNo]));
1692     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber() : "<<MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber());
1693     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getTime()            : "<<MED_FIELD_DRIVER<T>::_ptrField->getTime());
1694     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()     : "<<MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber());
1695
1696     // Rem 1 : le nombre d'éléments passé à MEDchampEcr ne doit pas tenir compte de la taille
1697     //         des profils : c'est la taille du champ sans profil.
1698     err=med_2_2::MEDchampEcr(id,
1699                              const_cast <char*> ( meshName.c_str()) ,
1700                              const_cast <char*> ( fieldName.c_str()),
1701                              (unsigned char*)value, med_2_2::MED_FULL_INTERLACE,
1702                              numberOfElForMED,
1703                              //UP : prend en compte le nombre de points de Gauss mais 
1704                              //     pas le nombre de composantes
1705                              const_cast <char*> ( locName.c_str()),
1706                              MED_ALL,
1707                              const_cast <char *> (profilName.c_str()), med_2_2::MED_COMPACT,
1708                              (med_2_2::med_entite_maillage)entityType,
1709                              (med_2_2::med_geometrie_element)types[typeNo],
1710                              MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber(),
1711                              "        ",                     // A FAIRE : IMPLEMENTER L'UNITE DU PAS DE TEMPS!
1712                              MED_FIELD_DRIVER<T>::_ptrField->getTime(),
1713                              MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()
1714                              );
1715
1716     if (err < MED_VALID ) {
1717       if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_NO_INTERLACE ) delete myField;
1718       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing "<< numberOfElements << " values for FIELD "<< fieldName 
1719                                    << " on entity " << MED_EN::entNames[entityType]
1720                                     << " and geometric type " << MED_EN::geoNames[types[typeNo]]
1721                                    << " with (it,or) = ("
1722                                    << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1723                                    << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1724                                     << profilName << " on mesh " << meshName
1725                                     )
1726                          );
1727     }
1728
1729     index += numberOfElements ; //Ne doit pas prendre en compte le nombre de points de GAUSS
1730                                 //ni les composantes.
1731
1732   }
1733   if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_NO_INTERLACE ) delete myField;
1734
1735
1736   END_OF(LOC);
1737 }
1738
1739 /*--------------------- RDWR PART -------------------------------*/
1740
1741 template <class T> GENDRIVER * MED_FIELD_RDWR_DRIVER22<T>::copy(void) const
1742 {
1743   return new MED_FIELD_RDWR_DRIVER22<T>(*this);
1744 }
1745
1746 template <class T> void MED_FIELD_RDWR_DRIVER22<T>::write(void) const
1747   throw (MEDEXCEPTION)
1748 {
1749   BEGIN_OF("MED_FIELD_RDWR_DRIVER22::write(void)");
1750   MED_FIELD_WRONLY_DRIVER22<T>::write();
1751   END_OF("MED_FIELD_RDWR_DRIVER22::write(void)");
1752 }
1753
1754 template <class T> void MED_FIELD_RDWR_DRIVER22<T>::read (void)
1755   throw (MEDEXCEPTION)
1756 {
1757   BEGIN_OF("MED_FIELD_RDWR_DRIVER22::read(void)");
1758   MED_FIELD_RDONLY_DRIVER22<T>::read();
1759   END_OF("MED_FIELD_RDWR_DRIVER22::read(void)");
1760 }
1761
1762 } //End namespace MEDMEM
1763 /*-----------------------------------------------------------------*/
1764
1765 #endif /* MED_FIELD_DRIVER_HXX */