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
 891   if (fileHasMesh && haveSupport )
 892     if ( ( meshGeoType != MESHgeoType ) || (meshNbOfElOfTypeC != MESHnbOfElOfTypeC) )
 893       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information from file for FIELD "<< fieldName
 894                                     << " on entity " << MED_EN::entNames[entityType]
 895                                     << " with (it,or) = ("
 896                                     << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
 897                                     << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
 898                                     << " on mesh " << meshName
 899                                     << " : geometric types or number of elements by type differs from MESH object !"
 900                                     )
 901                          );
 902
 903   if ( !fileHasMesh && !haveSupport )
 904     throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information for FIELD "<< fieldName
 905                                   << " on entity " << MED_EN::entNames[entityType]
 906                                   << " with (it,or) = ("
 907                                   << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
 908                                   << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
 909                                   << " on mesh " << meshName
 910                                   << " : SUPPORT must contain a valid MESH reference or file must contain the associated MESH."
 911                                   )
 912                        );
 913
 914
 915   if (!fileHasMesh && haveSupport) {
 916     meshNbOfElOfTypeC = MESHnbOfElOfTypeC;
 917     meshGeoType       = MESHgeoType;
 918     meshNbOfElOfType  = MESHnbOfElOfType;
 919   }
 920
 921
 922   // Test si le Support du Champ repose ou non sur toutes les entités géométriques
 923   // du maillage associé et positionne ou non l'attribut onAll du SUPPORT.
 924   // Il ne s'agit pas de la gestion des profils
 925   vector < MED_EN::medGeometryElement > v1(  mySupport->getTypes(),
 926                                              mySupport->getTypes()+mySupport->getNumberOfTypes() );
 927   vector<int> v2(numberOfElementsOfTypeC.size());
 928   transform(numberOfElementsOfTypeC.begin(),
 929             numberOfElementsOfTypeC.end(),v2.begin(), bind2nd(plus<int>(),1));
 930
 931   if ( ( meshGeoType != v1 )  || meshNbOfElOfTypeC != v2 ) {
 932     // ATTENTION : mySupport->setAll(false);
 933     // Pb : On a envie de positionner onAll à faux si le champ n'est pas défini sur tous les
 934     //      types géométriques du maillage associé.
 935     //      Mais si onAll est false et si aucun profil n'est détecté par la suite,
 936     //      l'attribut SUPPORT->_number est censé être positionné quand même ! Que faire ?
 937     // Si on veut être compatible avec la signification première de onAll,
 938     //  il faudrait créer des profils contenant toutes les entités pour chaque type géométrique
 939     //  du SUPPORT  mais d'une part c'est dommage d'un point de vue de l'emcombrement mémoire
 940     //  et d'autre part, à la réécriture du fichier MED on stockera des profils
 941     //  alors qu'il n'y en avait pas à l'origine (fichier MED différent après lecture/écriture) !
 942     // Si on laisse setAll à vrai il faut être sûr que les utilisateurs prennent les
 943     //  informations sur les types gémétrique au niveau du support et non pas du maillage.
 944     // Solution : Signification du onAll -> onAllElements des type géométriques définis
 945     // dans le SUPPORT et non du maillage associé (dans la plupart des cas si le fichier ne
 946     // contient pas de profil, le champ est défini sur toutes les entités de tous les types
 947     // géométriques définis dans le maillage).
 948   }
 949
 950
 951   // If an error occurs while reading the field, these allocated FIELD member will be deleted
 952
 953   MED_FIELD_DRIVER<T>::_ptrField->_name                   = fieldName;
 954   MED_FIELD_DRIVER<T>::_ptrField->_numberOfComponents     = numberOfComponents ;
 955   MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes        = new int   [numberOfComponents] ;
 956   MED_FIELD_DRIVER<T>::_ptrField->_componentsNames        = new string[numberOfComponents] ;
 957   MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits        = new UNIT  [numberOfComponents] ;
 958   MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions = new string[numberOfComponents] ;
 959   MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits     = new string[numberOfComponents] ;
 960   for (int i=0; i<numberOfComponents; i++) {
 961       MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes[i] = 1 ;
 962       MED_FIELD_DRIVER<T>::_ptrField->_componentsNames[i] = string(componentName,i*MED_TAILLE_PNOM22,MED_TAILLE_PNOM22) ;
 963       SCRUTE(MED_FIELD_DRIVER<T>::_ptrField->_componentsNames[i]);
 964       MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits[i] = string(unitName,i*MED_TAILLE_PNOM22,MED_TAILLE_PNOM22) ;
 965       SCRUTE(MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits[i]);
 966   }
 967   delete[] componentName;
 968   delete[] unitName;
 969
 970   int NumberOfTypes                       = mySupport->getNumberOfTypes() ;
 971   const MED_EN::medGeometryElement *types = mySupport->getTypes() ;
 972   T * myValues = new T[totalNumberOfElWg*numberOfComponents];
 973   const int * nbOfElOfType = mySupport->getNumberOfElements() ;
 974   bool anyProfil = false;
 975   int  pflSize=0,index=0;
 976   // Le vecteur de profil est dimensionné par rapport aux nombres de types
 977   // géométriques du champ même si le champ n'a pas de profil MED FICHIER sur
 978   // tous ses types géométriques car dans MEDMEM si onAllElement
 979   // du SUPPORT est false il faut positionner un profil pour tous les types géométriques
 980   // du SUPPORT
 981   int profilSizeC = 0;
 982   vector < int   >        profilSize    (NumberOfTypes,0);
 983   vector < vector<int>  > profilList    (NumberOfTypes);
 984   vector < string >       profilNameList(NumberOfTypes);
 985   char * profilName = new char[MED_TAILLE_NOM+1];
 986
 987   MESSAGE ("NumberOfTypes      : "<< NumberOfTypes);
 988   MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues=0 ;
 989
 990
 991   for (int typeNo=0; typeNo<NumberOfTypes; typeNo++) {
 992
 993     int numberOfValuesWc= nbOfElOfType[typeNo]*numberOfGaussPoint[typeNo+1]*numberOfComponents;
 994     char * gaussModelName = new char[MED_TAILLE_NOM+1];
 995
 996     MESSAGE ("FIELD_NAME         : "<< fieldName.c_str());
 997     MESSAGE ("MESH_NAME          : "<< meshName.c_str());
 998     MESSAGE ("MED_ENTITE         : "<< MED_EN::entNames[entityType]);
 999     MESSAGE ("MED_GEOM           : "<< MED_EN::geoNames[types[typeNo]]);
1000     MESSAGE ("Iteration          : "<< MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber());
1001     MESSAGE ("Order              : "<< MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber());
1002     MESSAGE ("NumberOfElements   : "<< nbOfElOfType[typeNo]);
1003     MESSAGE ("NumberOfComponents : "<< numberOfComponents);
1004     MESSAGE ("NumberOfGaussPts   : "<< numberOfGaussPoint[typeNo+1]);
1005     MESSAGE ("NumberOfValuesWg   : "<< nbOfElOfType[typeNo]*numberOfGaussPoint[typeNo+1]);
1006     MESSAGE ("NumberOfValuesWgWc : "<< numberOfValuesWc);
1007     MESSAGE ("Index              : "<< index);
1008     med_2_2::med_err ret=-1;
1009
1010     med_2_2::med_int * myValuesTmp=0;
1011     unsigned char* ptrTmp=0;
1012     if (needConversionToDouble || needConversionToInt64 ) {
1013       myValuesTmp = new med_2_2::med_int[numberOfValuesWc];
1014       ptrTmp = (unsigned char*) myValuesTmp;
1015     } else
1016       ptrTmp = (unsigned char*) &myValues[index];
1017
1018     //VERIFIER LE NBRE
1019     ret=med_2_2::MEDMEMchampLire(id,const_cast <char*> (meshName.c_str() ),
1020                                  const_cast <char*> (fieldName.c_str()),
1021                                 (unsigned char*) ptrTmp,
1022                                 med_2_2::MED_FULL_INTERLACE,
1023                                 MED_ALL,
1024                                 gaussModelName,
1025                                 profilName,
1026                                 med_2_2::MED_COMPACT,
1027                                 (med_2_2::med_entite_maillage) entityType,
1028                                 (med_2_2::med_geometrie_element)types[typeNo],
1029                                 MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber(),
1030                                 MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()
1031                                 );
1032
1033       if (needConversionToDouble || needConversionToInt64 ) {
1034
1035       if (needConversionToInt64 )  //utiliser un trait
1036         for(int i=0;i<numberOfValuesWc;++i)
1037           myValues[index+i]=(int)(myValuesTmp[i]);
1038       else
1039         for(int i=0;i<numberOfValuesWc;++i)
1040           myValues[index+i]=myValuesTmp[i];
1041       delete[] myValuesTmp;
1042     }
1043
1044     if (ret < 0)
1045       {
1046         // The Field can't be read then we must delete all previously allocated members in FIELD
1047         //for(int j=0; j<=i;j++)
1048         //  delete[] myValues[j];
1049         delete[] myValues;
1050         //delete[] NumberOfValues ;
1051         delete[] profilName;
1052         delete[] gaussModelName;
1053         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes ;
1054         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsNames ;
1055         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits ;
1056         delete[] MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions ;
1057         delete[] MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits ;
1058         MED_FIELD_DRIVER<T>::_ptrField->_componentsTypes = NULL ;
1059         MED_FIELD_DRIVER<T>::_ptrField->_componentsNames = NULL ;
1060         MED_FIELD_DRIVER<T>::_ptrField->_componentsUnits = NULL ;
1061         MED_FIELD_DRIVER<T>::_ptrField->_componentsDescriptions = NULL ;
1062         MED_FIELD_DRIVER<T>::_ptrField->_MEDComponentsUnits = NULL ;
1063         MED_FIELD_DRIVER<T>::_fieldNum = MED_INVALID ; // we have not found right field, so reset the field number
1064         throw MEDEXCEPTION( LOCALIZED( STRING(LOC) <<": ERROR while reading values")) ;
1065       }
1066     index += numberOfValuesWc;
1067     // Le support prend en compte le nombre de valeurs lié aux profils
1068     MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues+=
1069       nbOfElOfType[typeNo];// Ne doit pas prendre en compte les points de Gauss
1070
1071     // second et troisième test lié à un bug medfichier
1072     if ( strcmp(gaussModelName,MED_NOGAUSS) && strcmp(gaussModelName,string(MED_TAILLE_NOM,' ').c_str() )
1073          && strcmp(gaussModelName,string(16,' ').c_str() )  ) {
1074
1075         int type_geo = (int) types[typeNo];
1076         int t1       = (type_geo%100)*(type_geo/100);
1077         int ngauss   = numberOfGaussPoint[typeNo+1];
1078         int t2       = ngauss*(type_geo/100);
1079         med_2_2::med_float * refcoo = new med_2_2::med_float[t1];
1080         med_2_2::med_float * gscoo  = new med_2_2::med_float[t2];
1081         med_2_2::med_float * wg     = new med_2_2::med_float[ngauss];
1082
1083         if (MEDgaussLire(id, refcoo, gscoo, wg, (med_2_2::med_mode_switch) interlacingType,
1084                          gaussModelName ) < 0)
1085           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while reading Gauss Model |"
1086                                       << gaussModelName << "| for FIELD "<< fieldName
1087                                       << " on geometric type " << MED_EN::geoNames[types[typeNo]]
1088                                       )
1089                            );
1090         if (isFullInterlace ) { //serait inutile avec un driver template du type d'entrelacement
1091           GAUSS_LOCALIZATION<FullInterlace> * loc;
1092           loc = new GAUSS_LOCALIZATION<FullInterlace>(gaussModelName,types[typeNo],ngauss, refcoo,gscoo, wg);
1093           MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]]=loc;
1094         } else {
1095           GAUSS_LOCALIZATION<NoInterlace> * loc;
1096           loc = new GAUSS_LOCALIZATION<NoInterlace>(gaussModelName,types[typeNo],ngauss, refcoo,gscoo, wg);
1097           MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]]=loc;
1098         }
1099 //      cout << *MED_FIELD_DRIVER<T>::_ptrField->_gaussModel[types[typeNo]] << endl;
1100         delete [] refcoo;delete [] gscoo; delete [] wg;
1101     }
1102     delete[] gaussModelName ;
1103
1104     if ( strcmp(profilName,MED_NOPFL) ) {
1105       anyProfil = true;
1106       pflSize = med_2_2::MEDnValProfil(id,profilName);
1107       if ( pflSize  <= 0)
1108         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Error while reading the profil size of |"
1109                                      << profilName << "|" ));
1110
1111       profilSize[typeNo]=pflSize;
1112       profilList[typeNo].resize(pflSize);
1113       ret = med_2_2::MEDprofilLire(id,&profilList[typeNo][0],profilName); // cf item 16 Effective STL
1114       profilNameList[typeNo]=string(profilName);
1115     }
1116   }
1117
1118   delete[] profilName ;
1119
1120   //MESSAGE ("Index              : "<< index);
1121   assert(index == totalNumberOfElWg*numberOfComponents);
1122   assert(MED_FIELD_DRIVER<T>::_ptrField->_numberOfValues ==  mySupport->getNumberOfElements(MED_ALL_ELEMENTS));
1123
1124   if (anyProfil) {
1125
1126     for (int typeNo=0; typeNo < NumberOfTypes; typeNo++) {
1127
1128       // Trouve l'index du type géométrique dans la liste des types géométriques du maillage
1129       // correspondant au type géométrique du champ traité
1130       vector<MED_EN::medGeometryElement>::iterator meshTypeNoIt =
1131         find(meshGeoType.begin(),meshGeoType.end(),types[typeNo]); //Gérer l'exception
1132       if ( meshTypeNoIt ==  meshGeoType.end() )
1133         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Can't find "<< MED_EN::geoNames[types[typeNo]]
1134                                       << " on entity " << MED_EN::entNames[entityType]
1135                                       << " in geometric type list of mesh " << meshName
1136                                       )
1137                            );
1138       int meshTypeNo = meshTypeNoIt -  meshGeoType.begin();
1139
1140       if (! profilList[typeNo].empty() ) {
1141
1142 //      for (int j =0 ; j< meshGeoType.size();++j)
1143 //        cout << "--MeshTypeNo : "<<meshTypeNo<<"-> meshNbOfElOfTypeC["<<j<<"]="<<meshNbOfElOfTypeC[j]<<endl;
1144 //      cout << "--typeNo--" << typeNo << endl;
1145 //      cout << "meshNbOfElOfTypeC["<<meshTypeNo<<"]=" << meshNbOfElOfTypeC[meshTypeNo] <<endl;
1146
1147         // Transformer les numéros locaux d'entités medfichier en numéro global medmémoire
1148         for (int i = 0 ; i < profilList[typeNo].size(); i++ ) {
1149         // Les numéros des entités commencent à 1 dans MEDfichier comme dans MEDmémoire
1150         // meshNbOfElOfTypeC[0]=0 ...meshNbOfEltOfTypeC[meshTypeNo]=
1151         // meshNbOfElOfTypeC[meshTypeNo-1]+nbrOfElem of meshTypeNo type
1152         // rem1 : Si le meshTypeNo trouvé est 0 (premier type géométrique du maillage
1153         // il ne faut pas décaler les numéros du profils qui commencent à 1 dans MEDFICHIER
1154         // rem2 : meshNbOfElOfTypeC[NumberOfTypes] ne devrait jamais être utilisé
1155           profilList[typeNo][i]+=meshNbOfElOfTypeC[meshTypeNo];
1156         }
1157       } else {
1158         // Créer le profil <MED_ALL> pour ce type géométrique
1159         // uniquement pour renseigner le tableau skyline avec des accesseurs directs
1160         // par type géométriques
1161         // REM : Une conséquence est qu'à la réecriture le fichier contiendra des
1162         // profils sur certains types géométriques alors qu'à la lecture il n'y en avait pas !
1163         // Solution : Stocker les noms des profils et les utiliser pour savoir si il y avait ou non
1164         //            un profil
1165         int pflSize   = meshNbOfElOfType[meshTypeNo];
1166         // profil    = new int[pflSize];
1167
1168         profilList[typeNo].resize(pflSize);
1169         profilSize[typeNo]=pflSize;
1170
1171         for (int j = 1; j <= pflSize; j++) {
1172           profilList[typeNo][j-1]=meshNbOfElOfTypeC[meshTypeNo] + j ; // index MEDMEM commence à 1
1173         }
1174         profilNameList[typeNo] = MED_NOPFL; //Information a utiliser pour la sauvegarde : PLUTOT MED_ALL
1175       }
1176       profilSizeC+=profilList[typeNo].size();
1177     }
1178
1179     MEDSKYLINEARRAY * skyLine = new MEDSKYLINEARRAY(profilList.size(), profilSizeC );
1180     vector<int> index(NumberOfTypes+1,0);
1181     index[0]=1;
1182     for( int typeNo=0; typeNo < NumberOfTypes; typeNo++ )
1183       index[typeNo+1]=index[typeNo]+profilSize[typeNo];
1184     skyLine->setIndex(&index[0]);
1185     for (int i=1; i <= profilList.size() ; i++)
1186       skyLine->setI(i,&profilList[i-1][0]);
1187
1188     mySupport->setAll(false);
1189     mySupport->setpartial(skyLine,true);
1190     mySupport->setProfilNames(profilNameList);
1191 //    cout << "Valeurs du skyline du SUPPORT partiel crée : " << *skyLine << endl;
1192   }
1193
1194   // Créer un driver spécifique pour les modes MED_FULL_INTERLACE et MED_NO_INTERLACE
1195   // serait plus efficace.
1196   bool anyGauss = (numberOfGaussPoint != vector<int>(numberOfGaussPoint.size(),1));
1197   SCRUTE(anyGauss);
1198   MEDMEM_Array_ * Values;
1199   if (anyGauss) {
1200     SCRUTE(mySupport->getNumberOfElements(MED_ALL_ELEMENTS) );
1201     SCRUTE(NumberOfTypes);
1202     SCRUTE(numberOfElementsOfTypeC[NumberOfTypes]-1);
1203     assert(mySupport->getNumberOfElements(MED_ALL_ELEMENTS) == (numberOfElementsOfTypeC[NumberOfTypes]-1) );
1204     Values = new ArrayFullWg(myValues,
1205                              numberOfComponents,
1206                              numberOfElementsOfTypeC[NumberOfTypes]-1,
1207                              // Up : Prend en compte les profils et
1208                              // Ne prend pas en compte le nbre de composantes et
1209                              // le nombre de points de Gauss
1210                              NumberOfTypes,
1211                              &numberOfElementsOfTypeC[0],
1212                              &numberOfGaussPoint[0],
1213                              true,true);
1214 //     cout << "Valeurs du ArrayFullWg crée : " << endl <<
1215 //       *(static_cast<ArrayFullWg*>(Values))  << endl;
1216   } else
1217     Values = new ArrayFull(myValues,numberOfComponents,totalNumberOfElWg,
1218                                        true,true);
1219   if (MED_FIELD_DRIVER<T>::_ptrField->_value != NULL)
1220     delete MED_FIELD_DRIVER<T>::_ptrField->_value;
1221
1222   if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_NO_INTERLACE )
1223     {
1224       if (Values->getGaussPresence())
1225         MED_FIELD_DRIVER<T>::_ptrField->_value=ArrayConvert(*static_cast<ArrayFullWg*>(Values));
1226       else
1227         MED_FIELD_DRIVER<T>::_ptrField->_value=ArrayConvert(*static_cast<ArrayNo*    >(Values));
1228       delete Values;
1229     }
1230   else
1231     MED_FIELD_DRIVER<T>::_ptrField->_value=Values;
1232
1233   MED_FIELD_DRIVER<T>::_ptrField->_isRead = true ;
1234
1235   MED_FIELD_DRIVER<T>::_ptrField->_support=mySupport; //Prévenir l'utilisateur ?
1236
1237
1238   END_OF(LOC);
1239 }
1240
1241 template <class T> void MED_FIELD_RDONLY_DRIVER22<T>::write( void ) const
1242   throw (MEDEXCEPTION)
1243 {
1244   throw MEDEXCEPTION("MED_FIELD_RDONLY_DRIVER22::write : Can't write with a RDONLY driver !");
1245 }
1246
1247 /*--------------------- WRONLY PART -------------------------------*/
1248
1249 template <class T> GENDRIVER * MED_FIELD_WRONLY_DRIVER22<T>::copy(void) const
1250 {
1251   return new MED_FIELD_WRONLY_DRIVER22<T>(*this);
1252 }
1253
1254 template <class T> void MED_FIELD_WRONLY_DRIVER22<T>::read (void)
1255   throw (MEDEXCEPTION)
1256 {
1257   throw MEDEXCEPTION("MED_FIELD_WRONLY_DRIVER22::read : Can't read with a WRONLY driver !");
1258 }
1259
1260 template <class T> void MED_FIELD_WRONLY_DRIVER22<T>::write(void) const
1261   throw (MEDEXCEPTION)
1262 {
1263   const char * LOC = "MED_FIELD_WRONLY_DRIVER22::write(void) const " ;
1264   BEGIN_OF(LOC);
1265   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,NoGauss>::Array   ArrayNo;
1266   typedef typename MEDMEM_ArrayInterface<T,NoInterlace,Gauss>::Array     ArrayNoWg;
1267   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,NoGauss>::Array ArrayFull;
1268   typedef typename MEDMEM_ArrayInterface<T,FullInterlace,Gauss>::Array   ArrayFullWg;
1269
1270   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION<FullInterlace>*> locMapFull;
1271   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION<NoInterlace>*>   locMapNo;
1272   typedef map<MED_EN::medGeometryElement,GAUSS_LOCALIZATION_*>   locMap;
1273
1274   med_2_2::med_idt id = MED_FIELD_DRIVER22<T>::_medIdt;
1275
1276   if (MED_FIELD_DRIVER<T>::_status!=MED_OPENED)
1277     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<": Method open must be called before method read.")) ;
1278
1279   string fieldName;
1280   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
1281        ( MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
1282     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1283                                  <<" neither <fieldName> is set in driver nor in object FIELD.")) ;
1284
1285   // If _fieldName is not set in driver, try to use _ptrfield->_fieldName
1286   if ( ( MED_FIELD_DRIVER<T>::_fieldName.empty()       ) &&
1287        ( !MED_FIELD_DRIVER<T>::_ptrField->_name.empty() )    )
1288     fieldName=MED_FIELD_DRIVER<T>::_ptrField->_name;
1289   else
1290     fieldName = MED_FIELD_DRIVER<T>::_fieldName;
1291
1292   SCRUTE(fieldName);
1293   if ( fieldName.size() > MED_TAILLE_NOM ) {
1294     fieldName.substr(0,MED_TAILLE_NOM);
1295     MESSAGE( "Be careful <fieldName> size must not be > MED_TAILLE_NOM, using fieldName : |"<< fieldName <<"|." );
1296   }
1297
1298   const SUPPORT * mySupport = MED_FIELD_DRIVER<T>::_ptrField->getSupport() ;
1299   if ( ! mySupport )
1300     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)
1301                                  <<" There is no SUPPORT associated with FIELD : "
1302                                  << fieldName << "."));
1303
1304   bool onAll = mySupport->isOnAllElements();
1305   const locMap & gaussModel = MED_FIELD_DRIVER<T>::_ptrField->_gaussModel;
1306
1307
1308   string meshName = mySupport->getMeshName();
1309   SCRUTE(meshName);
1310   if ( meshName.size() > MED_TAILLE_NOM ) {
1311     meshName = meshName.substr(0,MED_TAILLE_NOM);
1312     MESSAGE( "Be careful <meshName> size must not be > MED_TAILLE_NOM, using meshName : |"<< meshName <<"|." );
1313   }
1314   MED_EN::medEntityMesh entityType = mySupport->getEntity();
1315
1316   // Reconstruit les listes contigues des noms de composantes et des unités
1317   // Les noms sont tronqués à MED_TAILLE_PNOM22
1318   int err ;
1319   int component_count=MED_FIELD_DRIVER<T>::_ptrField->getNumberOfComponents();
1320   string   component_name(component_count*MED_TAILLE_PNOM22,' ') ;
1321   string   component_unit(component_count*MED_TAILLE_PNOM22,' ') ;
1322
1323   const string * listcomponent_name=MED_FIELD_DRIVER<T>::_ptrField->getComponentsNames() ;
1324   const string * listcomponent_unit=MED_FIELD_DRIVER<T>::_ptrField->getMEDComponentsUnits() ;
1325   int length ;
1326   for (int i=0; i < component_count ; i++) {
1327     length = min(MED_TAILLE_PNOM22,(int)listcomponent_name[i].size());
1328     component_name.replace(i*MED_TAILLE_PNOM22,length,
1329                            listcomponent_name[i],0,length);
1330     length = min(MED_TAILLE_PNOM22,(int)listcomponent_unit[i].size());
1331     component_unit.replace(i*MED_TAILLE_PNOM22,length,
1332                            listcomponent_unit[i],0,length);
1333   }
1334
1335   MESSAGE("using component_name=|"<<component_name<<"|");
1336   MESSAGE("using component_unit=|"<<component_unit<<"|");
1337
1338   MED_EN::med_type_champ ValueType=MED_FIELD_DRIVER<T>::_ptrField->getValueType() ;
1339
1340   MESSAGE("Template Type =|"<<ValueType<<"|");
1341
1342   // Vérifier si le champ existe déjà
1343   char   champName[MED_TAILLE_NOM+1];
1344   char * compName, * compUnit ;
1345   med_2_2::med_type_champ type ;
1346   bool Find = false ;
1347   int n = med_2_2::MEDnChamp(id,0);
1348   int nbComp = 0;
1349   for (int i=1; i<=n; i++) {
1350     nbComp   = med_2_2::MEDnChamp(id,i);
1351     compName = new char[MED_TAILLE_PNOM22*nbComp+1];
1352     compUnit = new char[MED_TAILLE_PNOM22*nbComp+1];
1353     err = med_2_2::MEDchampInfo(id,i,champName,&type,compName,compUnit,nbComp);
1354     if (err == 0)
1355       if (!strcmp(champName,fieldName.c_str()) ) {
1356         Find = true ;
1357         break ;
1358       }
1359     delete[] compName ;
1360     delete[] compUnit ;
1361   }
1362
1363   if (Find) {
1364     // the same ?
1365     if (nbComp != component_count)
1366       throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1367                                      <<": Field exist in file, but number of component are different : "
1368                                      <<nbComp<<" in file and "
1369                                      <<component_count<<" in memory."
1370                                      )
1371                           );
1372     // component name and unit
1373     SCRUTE(nbComp);
1374     MESSAGE(LOC<<" Component name in file   : "<<compName);
1375     MESSAGE(LOC<<" Component name in memory : "<<component_name);
1376     MESSAGE(LOC<<" Component unit in file   : "<<compUnit);
1377     MESSAGE(LOC<<" Component unit in memory : "<<component_unit);
1378     delete[] compName ;
1379     delete[] compUnit ;
1380
1381   } else {
1382     // Verify the field isn't yet created
1383
1384     string dataGroupName =  "/CHA/";
1385     dataGroupName        += fieldName;
1386     MESSAGE(LOC << "|" << dataGroupName << "|" );
1387     med_2_2::med_idt gid =  H5Gopen(id, dataGroupName.c_str() );
1388
1389     if ( gid < 0 ) {
1390       // create field :
1391       err=med_2_2::MEDchampCr(id,
1392                               const_cast <char*> (fieldName.c_str()),
1393                               (med_2_2::med_type_champ) ValueType,
1394                               const_cast <char*> ( component_name.c_str() ),
1395                               const_cast <char*> ( component_unit.c_str() ),
1396                               component_count);
1397       if ( err < 0 )
1398         throw MEDEXCEPTION( LOCALIZED (STRING(LOC)
1399                                        << ": Error MEDchampCr : "<<err
1400                                        )
1401                             );
1402     }
1403     else H5Gclose(gid);
1404   }
1405
1406
1407
1408   // On s'assure que le champ est dans le bon type d'entrelacement.
1409   // REM : Il faudrait un driver par type d'entrelacement, ce qui eviterait
1410   // de doubler l'utilisation de la taille mémoire si le champ n'est pas dans
1411   // le bon mode.
1412   FIELD<T,FullInterlace> * myField = 0;
1413   if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_FULL_INTERLACE )
1414     myField = MED_FIELD_DRIVER<T>::_ptrField;
1415   else
1416     myField = FieldConvert( *( dynamic_cast< FIELD<T,NoInterlace> * > (MED_FIELD_DRIVER<T>::_ptrField )
1417                                )
1418                             );
1419
1420
1421   // Il est necessaire de calculer le tableau
1422   // du nombre d'entités cumulées de chaque type géométrique du maillage
1423   // pour convertir les profils de la numérotation globale
1424   // à la numérotation locale au type géométrique.
1425   // Pour celà on établit ce tableau à partir de l'objet MESH si la relation SUPPORT-MESH existe.
1426   // Si le maillage existe dans le fichier MED  on essaye également de reconstituer ce tableau
1427   // pour vérifier la cohérence des informations.
1428   // Si la relation SUPPRT-MESH n'esiste pas  on constitue le tableau uniquement à partir du fichier qui
1429   // doit alors obligatoirement contenir le maillage.
1430   const int * number, *numberIndex = 0;
1431   string         profilName;
1432   vector<string> profilNameList;
1433   vector<MED_EN::medGeometryElement> meshGeoType;
1434   vector<int> meshNbOfElOfType;
1435   vector<int> meshNbOfElOfTypeC;
1436   vector<MED_EN::medGeometryElement> fileMeshGeoType;
1437   vector<int> fileMeshNbOfElOfType;
1438   vector<int> fileMeshNbOfElOfTypeC;
1439   med_2_2::med_int fileHasMesh=0;
1440
1441   if (!onAll) {
1442
1443     number = mySupport->getNumber(MED_ALL_ELEMENTS);
1444     numberIndex = mySupport->getNumberIndex();
1445     profilNameList=mySupport->getProfilNames();
1446
1447     fileHasMesh = ( med_2_2::MEDdimLire(id, const_cast<char *>(meshName.c_str())) > 0);
1448     MESH * meshPtr = mySupport->getMesh();
1449
1450     if (fileHasMesh)
1451       this->getMeshGeometricTypeFromFile(id, meshName,
1452                                          entityType,
1453                                          fileMeshGeoType,fileMeshNbOfElOfType,fileMeshNbOfElOfTypeC);
1454
1455     if (meshPtr) {
1456       this->getMeshGeometricTypeFromMESH( meshPtr, entityType,meshGeoType,
1457                                           meshNbOfElOfType,
1458                                           meshNbOfElOfTypeC);
1459
1460       if (fileHasMesh)
1461         if ( ( fileMeshGeoType != meshGeoType ) || (fileMeshNbOfElOfTypeC != meshNbOfElOfTypeC) )
1462           throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information from file for FIELD "<< fieldName
1463                                         << " on entity " << MED_EN::entNames[entityType]
1464                                         << " with (it,or) = ("
1465                                         << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1466                                         << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1467                                         << " on mesh " << meshName
1468                                         << " : geometric types or number of elements by type differs from MESH object !"
1469                                         )
1470                              );
1471
1472     }
1473
1474     if ( !fileHasMesh && meshPtr==0 )
1475       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while getting mesh information for FIELD "<< fieldName
1476                                     << " on entity " << MED_EN::entNames[entityType]
1477                                     << " with (it,or) = ("
1478                                     << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1479                                     << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << ")"
1480                                     << " on mesh " << meshName
1481                                     << " : SUPPORT must contain a valid MESH reference or file must contain the associated MESH."
1482                                     )
1483                          );
1484
1485
1486     if (fileHasMesh && !meshPtr) {
1487       meshNbOfElOfTypeC = fileMeshNbOfElOfTypeC;
1488       meshGeoType       = fileMeshGeoType;
1489       meshNbOfElOfType  = fileMeshNbOfElOfType;
1490     }
1491   }
1492
1493   const MED_EN::medGeometryElement * types = mySupport->getTypes() ;
1494   int numberOfTypes = mySupport->getNumberOfTypes() ;
1495   int numberOfElForMED = -1;
1496   const T   * value   = NULL;
1497   int index = 1 ;
1498   // on boucle sur tout les types pour ecrire les tableaux de valeur
1499   for (int typeNo=0;typeNo<numberOfTypes;typeNo++) {
1500
1501     int numberOfElements = mySupport->getNumberOfElements(types[typeNo]) ;
1502     //UP : prend en compte les profils, pas les points de Gauss
1503
1504     //value = MED_FIELD_DRIVER<T>::_ptrField->getRow(index) ;
1505     // rem 1 : le getRow du Array est différent de celui du FIELD si le SUPPORT contient
1506     //         des profils (les indices des valeurs ne se suivent pas forcément)
1507     // rem 2 : Afin de respecter la norme MEDFICHIER, les indices contenus dans les
1508     //         profils doivent être croissant
1509     if (onAll) {
1510       value = myField->getRow(index);
1511       profilName=MED_NOPFL;
1512       numberOfElForMED = numberOfElements;
1513     } else {
1514       value = myField->getRow(number[index-1]);
1515       profilName = profilNameList[typeNo].substr(0,MED_TAILLE_NOM);
1516       // Rem : Si le SUPPORT n'est pas onAll mais que pour un type géométrique donné le nom
1517       // du profil associé est MED_NOPFL alors le profil n'est pas écrit dans le fichier MED.
1518       // Car en MEDMEMOIRE si le champ repose sur des éléments de deux types géométriques
1519       // différents et est défini sur tous les éléments d'un type géométrique
1520       // mais pas de l'autre, il existe tout de même des profils sur les deux types géométriques.
1521       // Ce n'est pas le cas en MEDFICHIER.
1522       vector<int> profil(&number[index-1],&(number[index-1])+numberOfElements);
1523
1524       // Trouve l'index du type géométrique dans la liste des types géométriques du maillage
1525       // correspondant au type géométrique du champ en cours de traitement
1526       vector<MED_EN::medGeometryElement>::iterator meshTypeNoIt =
1527         find(meshGeoType.begin(),meshGeoType.end(),types[typeNo]);
1528       if ( meshTypeNoIt ==  meshGeoType.end() )
1529         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Can't find "<< MED_EN::geoNames[types[typeNo]]
1530                                       << " on entity " << MED_EN::entNames[entityType]
1531                                       << " in geometric type list of mesh " << meshName
1532                                       )
1533                            );
1534
1535       int meshTypeNo = meshTypeNoIt -  meshGeoType.begin();
1536
1537       if ( profilName == MED_NOPFL && profil.size() != meshNbOfElOfType[meshTypeNo] )
1538         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating profil for FIELD "<< fieldName
1539                                       << " on entity " << MED_EN::entNames[entityType]
1540                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
1541                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1542                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1543                                       << profilName << " on mesh " << meshName
1544                                       << " : There is no profileName but profilsize (" <<profil.size()
1545                                       << ") differs from number of elements in associated MESH ("
1546                                       << meshNbOfElOfType[meshTypeNo] << ")."
1547                             )
1548                          );
1549
1550       //REM : Ce n'est pas évident, mais lorsqu'il y a un profil, le nombre de valeurs
1551       //      que l'on indique à MEDchampEcr est le nombre de valeurs sans profil, d'où
1552       //      le nombre d'éléments du maillage sur le type géométrique courant.
1553       numberOfElForMED = meshNbOfElOfType[meshTypeNo];
1554
1555       for (int ind=0;ind < numberOfElements;++ind) {
1556 //      cout << "number["<<index-1<<"]="<<number[index-1]<<endl;
1557 //      cout << "profil1["<<ind<<"]="<<profil[ind]<<endl;
1558         profil[ind]-=meshNbOfElOfTypeC[meshTypeNo];
1559 //      cout << "profil2["<<ind<<"]="<<profil[ind]<<endl;
1560       }
1561
1562       if ( profil[numberOfElements-1] > numberOfElForMED )
1563         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating profil for FIELD "<< fieldName
1564                                       << " on entity " << MED_EN::entNames[entityType]
1565                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
1566                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1567                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1568                                       << profilName << " on mesh " << meshName
1569                                       << " : profil["<<numberOfElements-1<<"]=" << profil[numberOfElements-1]
1570                                       << " must not be superior to field size without profil : "
1571                                       << numberOfElForMED
1572                             )
1573                          );
1574
1575       if ( med_2_2::MEDMEMprofilEcr(id,
1576                                  &profil[0],
1577                                  numberOfElements,
1578                                  const_cast<char *>(profilName.c_str())) < 0)
1579
1580         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing "<< numberOfElements
1581                                       << " values for MED profil "<< profilName
1582                                     )
1583                          );
1584     }
1585
1586
1587     string locName=MED_NOGAUSS;
1588     if (myField->getGaussPresence()) {
1589 //       cout << endl << "Nombre de points de Gauss à l'écriture de " << fieldName
1590 //         << " pour le type géométrique : " << MED_EN::geoNames[types[typeNo]]
1591 //         << " : " << myField->getNumberOfGaussPoints(types[typeNo]) << endl;
1592 //       cout << *mySupport << endl;
1593
1594       const GAUSS_LOCALIZATION_ * locPtr=0;
1595       locMap::const_iterator it;
1596       if ( ( it = gaussModel.find(types[typeNo])) != gaussModel.end() )
1597         locPtr = (*it).second;
1598       else
1599         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while creating Gauss Model for FIELD "<< fieldName
1600                                       << " on entity " << MED_EN::entNames[entityType]
1601                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]] << " with (it,or) = ("
1602                                       << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1603                                       << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1604                                       << profilName << " on mesh " << meshName
1605                                       << " : Can't find a Gauss localisation model for this geometric type"
1606                             )
1607                          );
1608
1609       int ngauss = -1;
1610       if ( locPtr->getInterlacingType() == MED_EN::MED_FULL_INTERLACE ) {
1611         const GAUSS_LOCALIZATION<FullInterlace> & loc=*(static_cast<const GAUSS_LOCALIZATION<FullInterlace> * >(locPtr));
1612         ngauss = loc.getNbGauss();
1613         locName=loc.getName();
1614         err=med_2_2::MEDMEMgaussEcr(id,
1615                                (med_2_2::med_geometrie_element) loc.getType(),
1616                                (med_2_2::med_float *)           loc.getRefCoo().getPtr(),
1617                                                                 med_2_2::MED_FULL_INTERLACE,
1618                                (med_2_2::med_int)               ngauss,
1619                                (med_2_2::med_float *)           loc.getGsCoo().getPtr(),
1620                                (med_2_2::med_float *)           (&loc.getWeight()[0]),
1621                                const_cast<char *>               (locName.c_str())
1622                                );
1623       } else {
1624         const GAUSS_LOCALIZATION<NoInterlace> & loc=*(static_cast<const GAUSS_LOCALIZATION<NoInterlace> * >(locPtr));
1625         ngauss = loc.getNbGauss();
1626         locName=loc.getName();
1627         err=med_2_2::MEDMEMgaussEcr(id,
1628                                (med_2_2::med_geometrie_element) loc.getType(),
1629                                (med_2_2::med_float *)           loc.getRefCoo().getPtr(),
1630                                                                 med_2_2::MED_NO_INTERLACE,
1631                                (med_2_2::med_int)               ngauss,
1632                                (med_2_2::med_float *)           loc.getGsCoo().getPtr(),
1633                                (med_2_2::med_float *)           (&loc.getWeight()[0]),
1634                                const_cast<char *>               (locName.c_str())
1635                                );
1636
1637       };
1638
1639       if ( err != 0 )
1640         throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing Gauss Model for FIELD "<< fieldName
1641                                       << " on entity " << MED_EN::entNames[entityType]
1642                                       << " and geometric type " << MED_EN::geoNames[types[typeNo]]
1643                             )
1644                          );
1645
1646         //numberOfElForMED *= mySupport->getNumberOfGaussPoints(types[typeNo]); //Deplacer la méthode dans FIELD
1647         numberOfElForMED *= ngauss;
1648     }
1649
1650     MESSAGE("MED_FIELD_DRIVER22<T>::_medIdt                       : "<<id);
1651     MESSAGE("meshName.c_str()                : "<<meshName.c_str());
1652     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getName()            : "<<MED_FIELD_DRIVER<T>::_ptrField->getName());
1653     MESSAGE("MED_FIELD_DRIVER<T>::_fieldName                      : "<<MED_FIELD_DRIVER<T>::_fieldName);
1654     MESSAGE("value                           : "<<value);
1655     MESSAGE("numberOfElements                : "<<numberOfElements);
1656     MESSAGE("numberOfElForMED                : "<<numberOfElForMED);
1657     MESSAGE("entityType                      : "<<MED_EN::entNames[entityType]);
1658     MESSAGE("types[i]                        : "<<MED_EN::geoNames[types[typeNo]]);
1659     MESSAGE("NumberOfGaussPoint[i]           : "<<myField->getNumberOfGaussPoints(types[typeNo]));
1660     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber() : "<<MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber());
1661     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getTime()            : "<<MED_FIELD_DRIVER<T>::_ptrField->getTime());
1662     MESSAGE("MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()     : "<<MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber());
1663
1664     // Rem 1 : le nombre d'éléments passé à MEDchampEcr ne doit pas tenir compte de la taille
1665     //         des profils : c'est la taille du champ sans profil.
1666     err=med_2_2::MEDchampEcr(id,
1667                              const_cast <char*> ( meshName.c_str()) ,
1668                              const_cast <char*> ( fieldName.c_str()),
1669                              (unsigned char*)value, med_2_2::MED_FULL_INTERLACE,
1670                              numberOfElForMED,
1671                              //UP : prend en compte le nombre de points de Gauss mais
1672                              //     pas le nombre de composantes
1673                              const_cast <char*> ( locName.c_str()),
1674                              MED_ALL,
1675                              const_cast <char *> (profilName.c_str()), med_2_2::MED_COMPACT,
1676                              (med_2_2::med_entite_maillage)entityType,
1677                              (med_2_2::med_geometrie_element)types[typeNo],
1678                              MED_FIELD_DRIVER<T>::_ptrField->getIterationNumber(),
1679                              "        ",                     // A FAIRE : IMPLEMENTER L'UNITE DU PAS DE TEMPS!
1680                              MED_FIELD_DRIVER<T>::_ptrField->getTime(),
1681                              MED_FIELD_DRIVER<T>::_ptrField->getOrderNumber()
1682                              );
1683
1684     if (err < MED_VALID ) {
1685       if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_NO_INTERLACE ) delete myField;
1686       throw MEDEXCEPTION(LOCALIZED( STRING(LOC) <<": Error while writing "<< numberOfElements << " values for FIELD "<< fieldName
1687                                    << " on entity " << MED_EN::entNames[entityType]
1688                                     << " and geometric type " << MED_EN::geoNames[types[typeNo]]
1689                                    << " with (it,or) = ("
1690                                    << MED_FIELD_DRIVER<T>::_ptrField->_iterationNumber << ","
1691                                    << MED_FIELD_DRIVER<T>::_ptrField->_orderNumber << "), with profilName "
1692                                     << profilName << " on mesh " << meshName
1693                                     )
1694                          );
1695     }
1696
1697     index += numberOfElements ; //Ne doit pas prendre en compte le nombre de points de GAUSS
1698                                 //ni les composantes.
1699
1700   }
1701   if ( MED_FIELD_DRIVER<T>::_ptrField->getInterlacingType() == MED_EN::MED_NO_INTERLACE ) delete myField;
1702
1703
1704   END_OF(LOC);
1705 }
1706
1707 /*--------------------- RDWR PART -------------------------------*/
1708
1709 template <class T> GENDRIVER * MED_FIELD_RDWR_DRIVER22<T>::copy(void) const
1710 {
1711   return new MED_FIELD_RDWR_DRIVER22<T>(*this);
1712 }
1713
1714 template <class T> void MED_FIELD_RDWR_DRIVER22<T>::write(void) const
1715   throw (MEDEXCEPTION)
1716 {
1717   BEGIN_OF("MED_FIELD_RDWR_DRIVER22::write(void)");
1718   MED_FIELD_WRONLY_DRIVER22<T>::write();
1719   END_OF("MED_FIELD_RDWR_DRIVER22::write(void)");
1720 }
1721
1722 template <class T> void MED_FIELD_RDWR_DRIVER22<T>::read (void)
1723   throw (MEDEXCEPTION)
1724 {
1725   BEGIN_OF("MED_FIELD_RDWR_DRIVER22::read(void)");
1726   MED_FIELD_RDONLY_DRIVER22<T>::read();
1727   END_OF("MED_FIELD_RDWR_DRIVER22::read(void)");
1728 }
1729
1730 } //End namespace MEDMEM
1731 /*-----------------------------------------------------------------*/
1732
1733 #endif /* MED_FIELD_DRIVER_HXX */