src/MEDMEM/MEDMEM_GibiMeshDriver.cxx

   1 //  Copyright (C) 2007-2008  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 //  Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 //  CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 //  This library is free software; you can redistribute it and/or
   7 //  modify it under the terms of the GNU Lesser General Public
   8 //  License as published by the Free Software Foundation; either
   9 //  version 2.1 of the License.
  10 //
  11 //  This library is distributed in the hope that it will be useful,
  12 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 //  Lesser General Public License for more details.
  15 //
  16 //  You should have received a copy of the GNU Lesser General Public
  17 //  License along with this library; if not, write to the Free Software
  18 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 //  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22 #include <algorithm>
  23 #include <queue>
  24
  25 #include "MEDMEM_GibiMeshDriver.hxx"
  26
  27 #include "MEDMEM_DriversDef.hxx"
  28
  29 #include "MEDMEM_Med.hxx"
  30 #include "MEDMEM_Family.hxx"
  31 #include "MEDMEM_Field.hxx"
  32 #include "MEDMEM_Group.hxx"
  33 #include "MEDMEM_Coordinate.hxx"
  34 #include "MEDMEM_Connectivity.hxx"
  35 #include "MEDMEM_Mesh.hxx"
  36 #include "MEDMEM_CellModel.hxx"
  37 #include "MEDMEM_define.hxx"
  38 #include "MEDMEM_DriverTools.hxx"
  39
  40 class MEDMEM::FIELD_;
  41
  42 #include <stdio.h>
  43 #include <fcntl.h>
  44 #ifdef WNT
  45 #include <io.h>
  46 #else
  47 #include <unistd.h>
  48 #endif
  49
  50 #include <float.h>
  51
  52 /////
  53 using namespace std;
  54 using namespace MED_EN;
  55 using namespace MEDMEM;
  56 /////
  57
  58 // read or not non-named fields
  59 //#define GIBI_READ_ONLY_NAMED_FIELD
  60
  61 // to throw an exception when try to write a name longer than 8 or non-unique encounters
  62 #define THROW_ON_BAD_NAME
  63
  64 // to see full dump of RESULTATS STRUCTURE INTERMEDIAIRES
  65 #ifdef _DEBUG_
  66 // #undef MESSAGE_MED
  67 // #define MESSAGE_MED(txt) std::cout << txt << endl;
  68 // #undef INFOS_MED
  69 // #define INFOS_MED(txt) std::cout << txt << endl;
  70 #endif
  71
  72 // Every memory allocation made in the MedDriver members function are desallocated in the Mesh destructor
  73
  74 /////
  75 //const size_t GIBI_MESH_DRIVER::nb_geometrie_gibi;
  76
  77 const medGeometryElement GIBI_MESH_DRIVER::geomGIBItoMED[nb_geometrie_gibi] =
  78   {   /*1 */ MED_POINT1 ,/*2 */ MED_SEG2   ,/*3 */ MED_SEG3   ,/*4 */ MED_TRIA3  ,/*5 */ MED_NONE   ,
  79       /*6 */ MED_TRIA6  ,/*7 */ MED_NONE   ,/*8 */ MED_QUAD4  ,/*9 */ MED_NONE   ,/*10*/ MED_QUAD8  ,
  80       /*11*/ MED_NONE   ,/*12*/ MED_NONE   ,/*13*/ MED_NONE   ,/*14*/ MED_HEXA8  ,/*15*/ MED_HEXA20 ,
  81       /*16*/ MED_PENTA6 ,/*17*/ MED_PENTA15,/*18*/ MED_NONE   ,/*19*/ MED_NONE   ,/*20*/ MED_NONE   ,
  82       /*21*/ MED_NONE   ,/*22*/ MED_NONE   ,/*23*/ MED_TETRA4 ,/*24*/ MED_TETRA10,/*25*/ MED_PYRA5  ,
  83       /*26*/ MED_PYRA13 ,/*27*/ MED_NONE   ,/*28*/ MED_NONE   ,/*29*/ MED_NONE   ,/*30*/ MED_NONE   ,
  84       /*31*/ MED_NONE   ,/*32*/ MED_NONE   ,/*33*/ MED_NONE   ,/*34*/ MED_NONE   ,/*35*/ MED_NONE   ,
  85       /*36*/ MED_NONE   ,/*37*/ MED_NONE   ,/*38*/ MED_NONE   ,/*39*/ MED_NONE   ,/*40*/ MED_NONE   ,
  86       /*41*/ MED_NONE   ,/*42*/ MED_NONE   ,/*43*/ MED_NONE   ,/*44*/ MED_NONE   ,/*45*/ MED_NONE   ,
  87       /*46*/ MED_NONE   ,/*47*/ MED_NONE   };
  88
  89 //=======================================================================
  90 //function : gibi2medGeom
  91 //purpose  :
  92 //=======================================================================
  93
  94 medGeometryElement GIBI_MESH_DRIVER::gibi2medGeom( size_t gibiTypeNb )
  95 {
  96   if ( gibiTypeNb < 1 || gibiTypeNb > 47 )
  97     return MED_NONE;
  98
  99   return geomGIBItoMED[ gibiTypeNb - 1 ];
 100 }
 101
 102 //=======================================================================
 103 //function : med2gibiGeom
 104 //purpose  :
 105 //=======================================================================
 106
 107 int GIBI_MESH_DRIVER::med2gibiGeom( medGeometryElement medGeomType )
 108 {
 109   for ( int i = 0; i < nb_geometrie_gibi; i++ )
 110     if ( geomGIBItoMED[ i ] == medGeomType )
 111       return i + 1;
 112
 113   return -1;
 114 }
 115
 116 //=======================================================================
 117 //function : getGroupId
 118 //purpose  :
 119 //=======================================================================
 120
 121 static int getGroupId(const vector<int>& support_ids, _intermediateMED*  medi)
 122 {
 123   int group_id = 0;
 124   vector<int>::const_iterator sb = support_ids.begin(), se = support_ids.end();
 125   if (support_ids.size() == 1 || // one or equal support ids
 126       *std::max_element( sb, se ) == *std::min_element( sb, se ))
 127   {
 128     group_id = support_ids[0] - 1;
 129   }
 130   else
 131   {
 132     // try to find an existing group with the same sub-groups
 133     set<int> sup_set;
 134     sup_set.insert( sb, se );
 135
 136     for ( group_id = 0; group_id < (int)medi->groupes.size(); ++group_id )
 137     {
 138       if (sup_set.size() == medi->groupes[ group_id ].groupes.size() &&
 139           std::equal (sup_set.begin(), sup_set.end(),
 140                       medi->groupes[ group_id ].groupes.begin()))
 141         break;
 142     }
 143     if ( group_id == medi->groupes.size() )
 144     {
 145       // no such a group, add a new one
 146       medi->groupes.push_back( _groupe() );
 147       _groupe& new_grp = medi->groupes.back();
 148       //new_grp.nom = string( group_id % 10 + 1, 'G' );
 149       new_grp.groupes.reserve( sup_set.size() );
 150       for ( set<int>::iterator it = sup_set.begin(); it != sup_set.end(); it++ ) {
 151         new_grp.groupes.push_back( *it );
 152         //new_grp.nom += "_" + medi->groupes[ *it - 1 ].nom;
 153       }
 154     }
 155   }
 156   return group_id;
 157 }
 158
 159 //=======================================================================
 160 //function : isNamedObject
 161 //purpose  :
 162 //=======================================================================
 163
 164 //#ifdef GIBI_READ_ONLY_NAMED_FIELD
 165 static inline bool isNamedObject( int obj_index, const vector<int>& indices_objets_nommes )
 166 {
 167   return ( std::find( indices_objets_nommes.begin(), indices_objets_nommes.end(), obj_index)
 168            != indices_objets_nommes.end() );
 169 }
 170 //#endif
 171
 172 //=======================================================================
 173 //function : setFieldNames
 174 //purpose  : set field names
 175 //=======================================================================
 176
 177 static void setFieldNames(const vector< _fieldBase* >& fields,
 178                           const vector<string>&        objets_nommes,
 179                           const vector<int>&           indices_objets_nommes)
 180 {
 181   // set field names
 182   int i;
 183   set<string> fnames;
 184   for ( i = 0; i < indices_objets_nommes.size(); ++i ) {
 185     int fieldIndex = indices_objets_nommes[ i ];
 186     fnames.insert( objets_nommes[ i ]);
 187     if ( fields[ fieldIndex - 1 ] )
 188       fields[ fieldIndex - 1 ]->_name = objets_nommes[ i ];
 189   }
 190   int noNameIndex = 0;
 191   for ( int i = 0; i < fields.size(); ++i ) {
 192     if ( !fields[ i ] ) {
 193       if ( !isNamedObject( i+1, indices_objets_nommes ))
 194         ++noNameIndex;
 195     }
 196     else if ( fields[ i ]->_name.empty() ) {
 197       do {
 198         fields[ i ]->_name = STRING("F_") << ++noNameIndex;
 199       } while ( !fnames.insert( fields[ i ]->_name ).second );
 200     }
 201   }
 202 }
 203
 204 //=======================================================================
 205 //function : read
 206 //purpose  :
 207 //=======================================================================
 208
 209 #define GIBI_EQUAL(var_str, stat_str) \
 210   (strncmp (var_str, stat_str, strlen(stat_str)) == 0)
 211 #define DUMP_LINE_NB " on line " << _lineNb
 212
 213 bool GIBI_MESH_RDONLY_DRIVER::readFile (_intermediateMED* medi, bool readFields )
 214 {
 215   const char* LOC = "GIBI_MESH_RDONLY_DRIVER::readFile() : ";
 216   BEGIN_OF_MED(LOC);
 217
 218   // LECTURE DES DONNEES DS FICHIER GIBI
 219
 220   enum Readable_Piles {
 221     PILE_SOUS_MAILLAGE=1,
 222     PILE_NODES_FIELD  =2,
 223     PILE_NOEUDS       =32,
 224     PILE_COORDONNEES  =33,
 225     PILE_FIELD        =39,
 226     PILE_LAST_READABLE=39
 227     };
 228   Readable_Piles readable_Piles [] = {
 229     PILE_SOUS_MAILLAGE,
 230     PILE_NODES_FIELD,
 231     PILE_NOEUDS,
 232     PILE_COORDONNEES,
 233     PILE_FIELD,
 234     PILE_LAST_READABLE
 235     };
 236   // other known piles:
 237   // PILE NUMERO  26 - Integers
 238   // PILE NUMERO  25 - Floats
 239   // PILE NUMERO  27 - Strings
 240   // PILE NUMERO  10 - Tables
 241   char* ligne; // pour lire une ligne
 242   const char* enregistrement_type=" ENREGISTREMENT DE TYPE";
 243   vector<int> numero_noeuds; // tableau de travail (indices)
 244   set<int> donePiles; // already read piles
 245   unsigned space_dimension = 0;
 246
 247   while ( getNextLine(ligne, false)) // boucle externe de recherche de "ENREGISTREMENT DE TYPE"
 248   {
 249     if ( !GIBI_EQUAL( ligne, enregistrement_type ))
 250       continue; // "ENREGISTREMENT DE TYPE" non trouvé -> on lit la ligne suivante
 251
 252     // lecture du numéro d'enregistrement
 253     int numero_enregistrement = atoi( ligne + strlen(enregistrement_type) + 1 );
 254
 255     enum { ENREG_TYPE_2=2, ENREG_TYPE_4=4}; // énumération des types d'enregistrement traités
 256     int numero_pile, nb_objets_nommes, nb_objets, nb_indices;
 257     vector<int> indices_objets_nommes;
 258     vector<string> objets_nommes;
 259
 260     if (numero_enregistrement == ENREG_TYPE_4)
 261     {
 262       getNextLine(ligne);
 263       const char* s = " NIVEAU  15 NIVEAU ERREUR   0 DIMENSION";
 264       space_dimension = atoi( ligne + strlen( s ) + 1 );
 265       if ( !GIBI_EQUAL( ligne, " NIVEAU" ) || space_dimension < 1 ) {
 266         INFOS_MED( " Could not read file: syntax error in type 4 record" << DUMP_LINE_NB );
 267         return false;
 268       }
 269     }
 270     else if (numero_enregistrement == ENREG_TYPE_2 )
 271     {
 272       if ( space_dimension == 0 ) {
 273         INFOS_MED( "Missing ENREGISTREMENT DE TYPE   4" << DUMP_LINE_NB );
 274         return false;
 275       }
 276       // FORMAT(' PILE NUMERO',I4,'NBRE OBJETS NOMMES',I8,'NBRE OBJETS',I8)
 277       getNextLine(ligne);
 278       const char *s1 = " PILE NUMERO", *s2 = "NBRE OBJETS NOMMES", *s3 = "NBRE OBJETS";
 279       if ( ! GIBI_EQUAL( ligne, s1 ) ) {
 280         INFOS_MED( " Could not read file: error in type 2 record. " << ligne << DUMP_LINE_NB );
 281         return false;
 282       }
 283       ligne = ligne + strlen(s1);
 284       numero_pile = atoi( ligne );
 285       ligne = ligne + 4 + strlen(s2);
 286       nb_objets_nommes = atoi( ligne );
 287       ligne = ligne + 8 + strlen(s3);
 288       nb_objets = atoi( ligne );
 289       if ( nb_objets_nommes<0 || nb_objets<0  ) {
 290         INFOS_MED(" Could not read file: " << nb_objets << " " <<nb_objets_nommes << DUMP_LINE_NB );
 291         return false;
 292       }
 293       if ( !donePiles.insert( numero_pile ).second ) // piles may repeat
 294         continue;
 295
 296       if ( numero_pile > PILE_LAST_READABLE )
 297         break; // stop file reading
 298
 299       // skip not readable piles
 300       int i = -1;
 301       while ( readable_Piles[ ++i ] != PILE_LAST_READABLE )
 302         if ( readable_Piles[ i ] == numero_pile )
 303           break;
 304       if ( readable_Piles[ i ] != numero_pile )
 305         continue;
 306
 307       // lecture des objets nommés et de leurs indices
 308       objets_nommes.resize(nb_objets_nommes);
 309       indices_objets_nommes.resize(nb_objets_nommes);
 310       for ( initNameReading( nb_objets_nommes ); more(); next() ) {
 311         objets_nommes[ index() ] = getName();
 312       }
 313       for ( initIntReading( nb_objets_nommes ); more(); next() )
 314         indices_objets_nommes[ index() ] = getInt();
 315
 316       // boucle interne : lecture de la pile
 317
 318       MESSAGE_MED(PREFIX_MED << "---- Traitement pile " << numero_pile);
 319
 320       // -----------------------------------
 321       //                        MESH GROUPS
 322       // -----------------------------------
 323
 324       if (numero_pile == PILE_SOUS_MAILLAGE ) // PILE NUMERO   1
 325       {
 326         map<int,int> strangeGroupType;
 327         medi->groupes.reserve(nb_objets*2); // fields may add some groups
 328         map< int , int > nbElemsByGeomType;
 329
 330         for (int objet=0; objet!=nb_objets; ++objet) // pour chaque groupe
 331         {
 332           initIntReading( 5 );
 333           unsigned type_geom_castem = getInt(); next();
 334           unsigned nb_sous_maillage = getInt(); next();
 335           unsigned nb_reference     = getInt(); next();
 336           unsigned nb_noeud         = getInt(); next();
 337           unsigned nb_elements      = getInt();
 338
 339           // le cas type_geom_castem=0 correspond aux maillages composites
 340           if (type_geom_castem<0) {
 341             INFOS_MED(" Error while reading file, bad geometric type:" <<
 342                   type_geom_castem << DUMP_LINE_NB );
 343             return false;
 344           }
 345
 346           medi->groupes.push_back(_groupe());
 347           _groupe & groupe = medi->groupes.back();
 348
 349           // si le groupe se compose de sous-maillages (ie groupe composite)
 350           if (type_geom_castem==0 && nb_sous_maillage>0)
 351           {
 352             // lecture des indices des sous-maillages, stockage.
 353             // les mailles correspondant a ces sous_maillages seront inserees a la fin du case
 354             groupe.groupes.resize( nb_sous_maillage );
 355             for ( initIntReading( nb_sous_maillage ); more(); next() ) {
 356               groupe.groupes[ index() ] = getInt();
 357             }
 358             if ( readFields )
 359               std::sort( groupe.groupes.begin(), groupe.groupes.end() );
 360           }
 361           // lecture des references (non utilisé pour MED)
 362           for ( i = 0; i < (int)nb_reference; i += 10 ) {// FORMAT(10I8)
 363             getNextLine(ligne);
 364           }
 365           // lecture des couleurs (non utilisé pour MED)
 366           for ( i = 0; i < (int)nb_elements; i += 10 ) {
 367             getNextLine(ligne);
 368           }
 369           // not a composit group
 370           if (type_geom_castem>0 && nb_sous_maillage==0)
 371           {
 372             medGeometryElement medType = gibi2medGeom(type_geom_castem);
 373
 374             initIntReading( nb_elements * nb_noeud );
 375             if ( medType == MED_NONE ) { // look for group end
 376               while ( more() )
 377                 next();
 378               strangeGroupType.insert( make_pair( objet, type_geom_castem ));
 379               continue;
 380             }
 381             if ( nbElemsByGeomType.find( medType ) == nbElemsByGeomType.end())
 382               nbElemsByGeomType[ medType ] = 0;
 383             int & order = nbElemsByGeomType[ medType ];
 384
 385             pair<set<_maille>::iterator,bool> p;
 386             pair<map<int,_noeud>::iterator,bool> p_no;
 387             _noeud no;
 388             no.coord.resize(space_dimension);
 389             _maille ma( medType, nb_noeud );
 390             ma.sommets.resize(nb_noeud);
 391             groupe.mailles.resize( nb_elements );
 392
 393             // lecture pour chaque maille des sommets et insertions
 394             for ( i = 0; i < nb_elements; ++i )
 395             {
 396               for (unsigned n = 0; n < nb_noeud; ++n, next() )
 397               {
 398                 if ( !more() ) {
 399                   INFOS_MED( " Error while reading elem nodes ");
 400                   return false;
 401                 }
 402                 no.number = getInt();
 403                 p_no=medi->points.insert(make_pair(no.number, no));
 404                 ma.sommets[n]=p_no.first;
 405               }
 406               ma.setOrdre( ++order );
 407               // on stocke dans le groupe un iterateur sur la maille
 408               groupe.mailles[i] = medi->insert(ma);
 409             }
 410           }
 411         } // loop on groups
 412
 413         // set group names
 414         for (i=0; i!=nb_objets_nommes; ++i) {
 415           int grpID = indices_objets_nommes[i];
 416           _groupe & grp = medi->groupes[ grpID-1 ];
 417           if ( !grp.nom.empty() ) // a group has several names
 418           { // create a group with subgroup grp and named grp.nom
 419             medi->groupes.push_back(_groupe());
 420             medi->groupes.back().groupes.push_back( grpID );
 421             medi->groupes.back().nom = grp.nom;
 422           }
 423           grp.nom=objets_nommes[i];
 424           map<int,int>::iterator it = strangeGroupType.find( grpID - 1 );
 425           if ( it != strangeGroupType.end() ) {
 426 #ifdef _DEBUG
 427             cout << "Skip " << grp.nom << " of not supported CASTEM type: " << it->second << endl;
 428 #endif
 429             //INFOS_MED( "Skip " << grp.nom << " of not supported CASTEM type: " << it->second );
 430           }
 431         }
 432
 433       }// Fin case PILE_SOUS_MAILLAGE
 434
 435       // ---------------------------------
 436       //                            NODES
 437       // ---------------------------------
 438
 439       else if ( numero_pile == PILE_NOEUDS ) // PILE NUMERO  32
 440       {
 441         getNextLine( ligne );
 442         std::vector<int> place_noeuds;
 443         nb_indices = atoi ( ligne );
 444         if (nb_indices != nb_objets)
 445         {
 446           INFOS_MED("Erreur de lecture dans enregistrement de pile " <<
 447                 PILE_NOEUDS << DUMP_LINE_NB );
 448           return false;
 449         }
 450
 451         place_noeuds.resize(nb_objets);
 452         for ( initIntReading( nb_objets ); more(); next() )
 453           place_noeuds[ index() ] = getInt();
 454         int max=(* std::max_element(place_noeuds.begin(),place_noeuds.end()));
 455
 456         // numero_noeuds contient pour chacun des max noeuds qu'on va lire dans le case PILE_COORDONNEES
 457         // son indice dans la connectivite du maillage. Cet indice correspond egalement a la cle du map
 458         // medi->points ou l'on stocke les noeuds.
 459         numero_noeuds.resize(max,-1);
 460         for (unsigned i=0; i!=place_noeuds.size(); ++i)
 461           numero_noeuds[place_noeuds[i]-1]=i+1;
 462       }
 463
 464       // ---------------------------------------
 465       //                            COORDINATES
 466       // ---------------------------------------
 467
 468       else if ( numero_pile == PILE_COORDONNEES )// PILE NUMERO  33
 469       {
 470         getNextLine( ligne );
 471         unsigned nb_reels = atoi( ligne );
 472         // PROVISOIRE : certains fichier gibi n`ont
 473         if (nb_reels < numero_noeuds.size()*(space_dimension)) {
 474           INFOS_MED("Erreur de lecture dans enregistrement de pile " <<
 475                 PILE_COORDONNEES << DUMP_LINE_NB );
 476           return false;
 477         }
 478         initDoubleReading( nb_reels );
 479         map< int, _noeud >::iterator pIt;
 480         for (unsigned i=0; i!=numero_noeuds.size(); ++i)
 481         {
 482           // si le noeud est utilisé dans le maillage,
 483           //on lit ses coordonnées et on les stocke dans la structure
 484           if (( numero_noeuds[i] != -1 ) &&
 485               (( pIt = medi->points.find(numero_noeuds[i])) != medi->points.end()))
 486           {
 487             for (unsigned j=0; j!=space_dimension; ++j, next())
 488               pIt->second.coord[j] = getDouble();
 489             next(); // on ne conserve pas la densite
 490           }
 491           else // sinon, on passe au noeud suivant
 492           {
 493             for (unsigned j=0; j!=space_dimension+1; ++j)
 494               next();
 495           }
 496         }
 497       }
 498
 499       // ---------------------------------------
 500       //                            NODE FIELDS
 501       // ---------------------------------------
 502
 503       else if ( numero_pile == PILE_NODES_FIELD && readFields ) // PILE NUMERO   2
 504       {
 505         vector< _fieldBase* > fields( nb_objets );
 506         for (int objet=0; objet!=nb_objets; ++objet) // pour chaque field
 507         {
 508           bool ignoreField = false;
 509 #ifdef GIBI_READ_ONLY_NAMED_FIELD
 510           ignoreField = !isNamedObject( objet+1, indices_objets_nommes );
 511           if ( ignoreField )
 512             INFOS_MED("Skip non-named field " << objet+1 << DUMP_LINE_NB);
 513 #endif
 514
 515           // EXAMPLE ( with no values )
 516
 517           // (1)       4       7       2       1
 518           // (2)     -88       0       3     -89       0       1     -90       0       2     -91
 519           // (2)       0       1
 520           // (3) FX   FY   FZ   FZ   FX   FY   FLX
 521           // (4)       0       0       0       0       0       0       0
 522           // (5)           cree  par  muc pri
 523           // (6)
 524           // (7)       2
 525
 526           // (1): nb subcomponents, nb components(total), IFOUR, nb attributes
 527           initIntReading( 4 );
 528           int i_sub, nb_sub         = getInt(); next();
 529           int i_comp, total_nb_comp = getInt(); next();
 530           next(); // ignore IFOUR
 531           int nb_attr = getInt();
 532           if ( nb_sub < 0 || total_nb_comp < 0 || nb_attr < 0 ) {
 533             INFOS_MED("Error of field reading: wrong nb of components "
 534                   << nb_sub << " " << total_nb_comp << DUMP_LINE_NB);
 535             return false;
 536           }
 537           // (2) loop on subcomponents of a field, for each read
 538           // (a) support, (b) number of values and (c) number of components
 539           vector<int>     support_ids( nb_sub );
 540           vector<int>     nb_values  ( nb_sub );
 541           vector<int>     nb_comps   ( nb_sub );
 542           int total_nb_values = 0;
 543           initIntReading( nb_sub * 3 );
 544           for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
 545           {
 546             support_ids[ i_sub ] = -getInt(); next(); // (a) reference to support
 547             if ( support_ids[ i_sub ] < 1 || support_ids[ i_sub ] > medi->groupes.size() ) {
 548               INFOS_MED("Error of field reading: wrong mesh reference "<<
 549                     support_ids[ i_sub ] << DUMP_LINE_NB );
 550               return false;
 551             }
 552             nb_values[ i_sub ] = getInt(); next();    // (b) nb points
 553             total_nb_values += nb_values[ i_sub ];
 554             if ( nb_values[ i_sub ] < 0 ) {
 555               INFOS_MED(" Wrong nb of points: " << nb_values[ i_sub ]  << DUMP_LINE_NB );
 556               return false;
 557             }
 558             nb_comps[ i_sub ] = getInt(); next();     // (c) nb of components in i_sub
 559           }
 560           // create a field if there are values
 561           _field<double>* fdouble = 0;
 562           if ( total_nb_values > 0 && !ignoreField )
 563           {
 564             fdouble = new _field<double>( MED_REEL64, nb_sub, total_nb_comp );
 565             medi->fields.push_back( fields[ objet ] = fdouble );
 566           }
 567           // (3) component names
 568           initNameReading( total_nb_comp, 4 );
 569           for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
 570           {
 571             // store support id and nb components of a sub
 572             if ( fdouble )
 573               fdouble->_sub[ i_sub ].setData( nb_comps[ i_sub ], support_ids[ i_sub ] );
 574             for ( i_comp = 0; i_comp < nb_comps[ i_sub ]; ++i_comp, next() )
 575             {
 576               ASSERT_MED( more() );
 577               // store component name
 578               if ( fdouble )
 579                 fdouble->_sub[ i_sub ].compName( i_comp ) = getName();
 580             }
 581           }
 582           // (4) nb harmonics ( ignored )
 583           for ( initIntReading( total_nb_comp ); more(); next() )
 584             ;
 585           // (5) TYPE ( ignored )
 586           getNextLine( ligne );
 587           // (6) TITRE ( ignored )
 588           getNextLine( ligne );
 589           // (7) attributes ( ignored )
 590           for ( initIntReading( nb_attr ); more(); next() )
 591             ;
 592
 593           for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
 594           {
 595             // loop on components: read values
 596             initDoubleReading( nb_values[ i_sub ] * nb_comps[ i_sub ] );
 597             for ( i_comp = 0; i_comp < nb_comps[ i_sub ]; ++i_comp )
 598             {
 599               vector<double>* vals = 0;
 600               if ( fdouble ) vals = & fdouble->addComponent( nb_values[ i_sub ] );
 601               for ( int i = 0; more() && i < nb_values[ i_sub ]; next(), ++i ) {
 602                 if ( vals ) (*vals)[ i ] = getDouble();
 603               }
 604             }
 605           } // loop on subcomponents of a field
 606
 607           // set id of a group including all subs supports but only
 608           // if all subs have the same components
 609           if ( fdouble && fdouble->hasSameComponentsBySupport() )
 610             fdouble->_group_id = getGroupId( support_ids, medi );
 611
 612         } // end loop on field objects
 613
 614         // set field names
 615         setFieldNames( fields, objets_nommes, indices_objets_nommes );
 616
 617       }  // Fin numero_pile == PILE_NODES_FIELD
 618
 619       // -------------------------------------------------
 620       //                                           FIELDS
 621       // -------------------------------------------------
 622
 623       else if ( numero_pile == PILE_FIELD && readFields ) // PILE NUMERO  39
 624       {
 625         // REAL EXAMPLE
 626
 627         // (1)        1       2       6      16
 628         // (2)                                                         CARACTERISTIQUES
 629         // (3)      -15  317773       4       0       0       0      -2       0       3
 630         // (4)             317581
 631         // (5)  0
 632         // (6)   317767  317761  317755  317815
 633         // (7)  YOUN     NU       H        SIGY
 634         // (8)  REAL*8            REAL*8            REAL*8            REAL*8
 635         // (9)        1       1       0       0
 636         // (10)  2.00000000000000E+05
 637         // (11)       1       1       0       0
 638         // (12)  3.30000000000000E-01
 639         // (13)       1       1       0       0
 640         // (14)  1.00000000000000E+04
 641         // (15)       6     706       0       0
 642         // (16)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
 643         // (17)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
 644         // (18)  ...
 645         vector< _fieldBase* > fields( nb_objets, (_fieldBase*)0 );
 646         for (int objet=0; objet!=nb_objets; ++objet) // pour chaque field
 647         {
 648           bool ignoreField = false;
 649 #ifdef GIBI_READ_ONLY_NAMED_FIELD
 650           ignoreField = !isNamedObject( objet+1, indices_objets_nommes );
 651           if ( ignoreField )
 652             INFOS_MED("Skip non-named field " << objet+1 << DUMP_LINE_NB);
 653 #endif
 654           initIntReading( 4 );
 655           int i_sub, nb_sub = getInt(); // (1) <nb_sub> 2 6 <title length>
 656           if ( nb_sub < 1 ) {
 657             INFOS_MED("Error of field reading: wrong nb of subcomponents " <<
 658                   nb_sub << DUMP_LINE_NB );
 659             return false;
 660           }
 661           next(); next(); next(); // skip (1) <nb_sub> 2 6
 662           int title_length = getInt(); // <title length>
 663           string description;
 664           if ( title_length ) {
 665             getNextLine( ligne ); // (2) title
 666             const int len = 72; // line length
 667             // title is right justified
 668             description = string(ligne + len - title_length, title_length);
 669           }
 670           // look for a line starting with '-' : <reference to support>
 671           do {
 672             initIntReading( nb_sub * 9 );
 673           } while ( getInt() >= 0 );
 674
 675           int total_nb_comp = 0;
 676           vector<int> support_ids( nb_sub ), nb_comp( nb_sub );
 677           for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
 678           {                                           // (3)
 679             support_ids[ i_sub ] = -getInt(); next(); // <reference to support>
 680             next();                                   // ignore <address>
 681             nb_comp    [ i_sub ] =  getInt(); next(); // <nb of components in the sub>
 682             for ( i = 0; i < 6; ++i )                 // ignore 6 ints, in example 0 0 0 -2 0 3
 683               next();
 684             if ( support_ids[ i_sub ] < 1 || support_ids[ i_sub ] > medi->groupes.size() ) {
 685               INFOS_MED("Error of field reading: wrong mesh reference "<<
 686                     support_ids[ i_sub ] << DUMP_LINE_NB );
 687               return false;
 688             }
 689             if ( nb_comp[ i_sub ] < 0 /*1*/ ) { // [SALOME platform 0019886]
 690               INFOS_MED("Error of field reading: wrong nb of components " <<
 691                     nb_comp[ i_sub ] << DUMP_LINE_NB );
 692               return false;
 693             }
 694             total_nb_comp += nb_comp[ i_sub ];
 695           }
 696           for ( initNameReading( nb_sub, 17 ); more(); next() )
 697             ; // (4) dummy strings
 698           for ( initNameReading( nb_sub ); more(); next() )
 699             ; // (5) dummy strings
 700
 701           // loop on subcomponents of a field, each of which refers to
 702           // a certain support and has its own number of components;
 703           // read component values
 704           _field<double>* fdouble = 0;
 705           _field<int>*    fint    = 0;
 706           _fieldBase *    fbase   = 0;
 707           for ( i_sub = 0; i_sub < nb_sub; ++ i_sub )
 708           {
 709             vector<string> comp_names( nb_comp[ i_sub ]), comp_type( nb_comp[ i_sub ]);
 710             for ( initIntReading( nb_comp[ i_sub ] ); more(); next() )
 711               ;  // (6) nb_comp addresses of MELVAL structure
 712
 713             // (7) component names
 714             for ( initNameReading( nb_comp[ i_sub ] ); more(); next() )
 715               comp_names[ index() ] = getName();
 716
 717             // (8) component type
 718             for ( initNameReading( nb_comp[ i_sub ], 17 ); more(); next() ) { // 17 is name width
 719               comp_type[ index() ] = getName();
 720               // component types must be the same
 721               if ( index() > 0 && comp_type[ index() ] != comp_type[ index() - 1] ) {
 722                 INFOS_MED( "Error of field reading: diff component types <"
 723                       << comp_type[ index() ] << "> != <" << comp_type[ index() - 1 ]
 724                        << ">" << DUMP_LINE_NB );
 725                 return false;
 726               }
 727             }
 728             // now type is known, create a field, one for all subs
 729             bool isReal = (nb_comp[i_sub] > 0) ? (comp_type[0] == "REAL*8") : true;
 730             if ( !ignoreField && !fbase && total_nb_comp ) {
 731               if ( !isReal ) {
 732                 fbase = fint = new _field<int>( MED_INT32, nb_sub, total_nb_comp );
 733                 INFOS_MED( "Warning: read NOT REAL field, type <" << comp_type[0] << ">"
 734                       << DUMP_LINE_NB);
 735               }
 736               else
 737                 fbase = fdouble = new _field<double>( MED_REEL64, nb_sub, total_nb_comp );
 738               medi->fields.push_back( fields[ objet ] = fbase ); // medi->fields is a std::list
 739             }
 740             // store support id and nb components of a sub
 741             if ( fbase ) {
 742               fbase->_sub[ i_sub ].setData( nb_comp[ i_sub ], support_ids[ i_sub ]);
 743               fbase->_description = description;
 744             }
 745             // loop on components: read values
 746             for ( int i_comp = 0; i_comp < nb_comp[ i_sub ]; ++i_comp )
 747             {
 748               // (9) nb of values
 749               initIntReading( 2 );
 750               int nb_val_by_elem = getInt(); next();
 751               int nb_values      = getInt();
 752               fbase->_sub[ i_sub ]._nb_gauss[ i_comp ] = nb_val_by_elem;
 753               // (10) values
 754               nb_values *= nb_val_by_elem;
 755               if ( fbase ) {
 756                 if ( isReal ) {
 757                   vector<double> & vals = fdouble->addComponent( nb_values );
 758                   for ( initDoubleReading( nb_values ); more(); next()) {
 759                     vals[ index() ] = getDouble();
 760                   }
 761                 }
 762                 else {
 763                   vector<int> & vals = fint->addComponent( nb_values );
 764                   for ( initIntReading( nb_values ); more(); next() ) {
 765                     vals[ index() ] = getInt();
 766                   }
 767                 }
 768                 // store component name
 769                 fbase->_sub[ i_sub ].compName( i_comp ) = comp_names[ i_comp ];
 770               }
 771               else {
 772                 for ( isReal ? initDoubleReading( nb_values ) : initIntReading( nb_values );
 773                       more();
 774                       next() ) ;
 775               }
 776             }
 777           } // loop on subcomponents of a field
 778
 779           // set id of a group including all sub supports but only
 780           // if all subs have the same nb of components
 781           if ( fbase && fbase->hasSameComponentsBySupport() )
 782             fbase->_group_id = getGroupId( support_ids, medi );
 783
 784         } // end loop on field objects
 785
 786         // set field names
 787         setFieldNames( fields, objets_nommes, indices_objets_nommes );
 788
 789       } // numero_pile == PILE_FIELD && readFields
 790
 791       else if ( numero_pile >= PILE_LAST_READABLE )
 792         break; // stop file reading
 793
 794     } // Fin case ENREG_TYPE_2
 795   } //  fin de la boucle while de lecture externe
 796
 797   // check if all needed piles present
 798   if ( donePiles.find( PILE_SOUS_MAILLAGE ) != donePiles.end() )
 799   {
 800     if (donePiles.find( PILE_NOEUDS ) == donePiles.end() ) {
 801       INFOS_MED( " Missing pile " << PILE_NOEUDS );
 802       return false;
 803     }
 804     if (donePiles.find( PILE_COORDONNEES ) == donePiles.end()) {
 805       INFOS_MED( " Missing pile " << PILE_COORDONNEES );
 806       return false;
 807     }
 808   }
 809
 810   END_OF_MED(LOC);
 811   return true;
 812 }
 813
 814 GIBI_MESH_DRIVER::GIBI_MESH_DRIVER():
 815        GENDRIVER(GIBI_DRIVER),
 816        _ptrMesh(( MESH *) NULL),
 817        // A VOIR _medIdt(MED_INVALID),
 818        _meshName("")
 819 {
 820   MESSAGE_MED("GIBI_MESH_DRIVER()");
 821 }
 822
 823 GIBI_MESH_DRIVER::GIBI_MESH_DRIVER(const string & fileName,
 824                                    MESH * ptrMesh,
 825                                    MED_EN::med_mode_acces accessMode):
 826   GENDRIVER(fileName, accessMode, GIBI_DRIVER),
 827   _ptrMesh(ptrMesh)
 828   // A VOIR _medIdt(MED_INVALID),
 829 {
 830   MESSAGE_MED( "GIBI_MESH_DRIVER(" << fileName <<","<<accessMode );
 831 //   _meshName=fileName.substr(0,fileName.rfind("."));
 832     // mesh name construction from fileName
 833     const string ext=".sauve"; // expected extension
 834     string::size_type pos=fileName.find(ext,0);
 835     string::size_type pos1=fileName.rfind('/');
 836     _meshName = string(fileName,pos1+1,pos-pos1-1); //get rid of directory & extension
 837     SCRUTE_MED(_meshName);
 838 }
 839
 840 GIBI_MESH_DRIVER::GIBI_MESH_DRIVER(const GIBI_MESH_DRIVER & driver):
 841   GENDRIVER(driver),
 842   _ptrMesh(driver._ptrMesh),
 843   // A VOIR _medIdt(MED_INVALID),
 844   _meshName(driver._meshName)
 845 {
 846   MESSAGE_MED("GIBI_MESH_DRIVER(const GIBI_MESH_DRIVER & driver)");
 847 }
 848
 849 GIBI_MESH_DRIVER::~GIBI_MESH_DRIVER()
 850 {
 851   MESSAGE_MED("~GIBI_MESH_DRIVER()");
 852 }
 853 void    GIBI_MESH_DRIVER::setMeshName(const string & meshName) { _meshName = meshName; };
 854 string  GIBI_MESH_DRIVER::getMeshName() const { return _meshName; };
 855
 856
 857 //---------------------------------- RDONLY PART -------------------------------------------------------------
 858
 859 GIBI_MESH_RDONLY_DRIVER::GIBI_MESH_RDONLY_DRIVER():
 860        GIBI_MESH_DRIVER(),
 861        _File (-1),_start(0L),_ptr  (0L),_eptr (0L)
 862 {
 863 }
 864 GIBI_MESH_RDONLY_DRIVER::GIBI_MESH_RDONLY_DRIVER(const string & fileName,MESH * ptrMesh):
 865        GIBI_MESH_DRIVER(fileName,ptrMesh,RDONLY),
 866        _File (-1),_start(0L),_ptr  (0L),_eptr (0L)
 867 {
 868     MESSAGE_MED("GIBI_MESH_RDONLY_DRIVER::GIBI_MESH_RDONLY_DRIVER"
 869             "(const string & fileName, MESH * ptrMesh) has been created, "
 870             << fileName << ", " << RDONLY);
 871 }
 872 GIBI_MESH_RDONLY_DRIVER::GIBI_MESH_RDONLY_DRIVER(const GIBI_MESH_RDONLY_DRIVER & driver):
 873 GIBI_MESH_DRIVER(driver)
 874 {
 875 }
 876 GIBI_MESH_RDONLY_DRIVER::~GIBI_MESH_RDONLY_DRIVER()
 877 {
 878   const char* LOC = "~GIBI_MESH_RDONLY_DRIVER()";
 879   BEGIN_OF_MED(LOC);
 880   if (_File >= 0)
 881   {
 882     ::close (_File);
 883     if (_start != 0L)
 884       delete [] _start;
 885   }
 886   MESSAGE_MED("GIBI_MESH_RDONLY_DRIVER::~GIBI_MESH_RDONLY_DRIVER() has been destroyed");
 887 }
 888 GENDRIVER * GIBI_MESH_RDONLY_DRIVER::copy(void) const
 889 {
 890   return new GIBI_MESH_RDONLY_DRIVER(*this);
 891 }
 892
 893 //=======================================================================
 894 //function : open
 895 //purpose  :
 896 //=======================================================================
 897
 898 const int GIBI_MaxOutputLen   = 150;
 899 const int GIBI_BufferSize     = 16184; // for non-stream input
 900
 901 void GIBI_MESH_RDONLY_DRIVER::open()
 902   //     throw (MEDEXCEPTION)
 903 {
 904   if( _status == MED_OPENED )
 905     return;
 906
 907   const char * LOC = "GIBI_MESH_RDONLY_DRIVER::open()" ;
 908   BEGIN_OF_MED(LOC);
 909
 910 //   MED_EN::med_mode_acces aMode = getAccessMode();
 911 //   if ( aMode != MED_EN::MED_LECT && aMode != MED_EN::MED_REMP )
 912 //     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << " Bad file mode access ! " << aMode));
 913
 914 #ifdef WNT
 915   _File = ::_open (_fileName.c_str(), _O_RDONLY|_O_BINARY);
 916 #else
 917   _File = ::open (_fileName.c_str(), O_RDONLY);
 918 #endif
 919   if (_File >= 0)
 920   {
 921     _start = new char [GIBI_BufferSize];
 922     _ptr   = _start;
 923     _eptr  = _start;
 924     _status = MED_OPENED;
 925     _lineNb = 0;
 926   }
 927   else
 928   {
 929     _status = MED_CLOSED;
 930     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Could not open file "<<_fileName
 931                                  << " fd: " << _File));
 932   }
 933   END_OF_MED(LOC);
 934 }
 935
 936 //=======================================================================
 937 //function : close
 938 //purpose  :
 939 //=======================================================================
 940
 941 void GIBI_MESH_RDONLY_DRIVER::close()
 942 {
 943   const char* LOC = "GIBI_MESH_DRIVER::close() ";
 944   BEGIN_OF_MED(LOC);
 945   if ( _status == MED_OPENED) {
 946     if (_File >= 0) {
 947       ::close (_File);
 948       if (_start != 0L)
 949         delete [] _start;
 950       _File = -1;
 951     }
 952     _status = MED_CLOSED;
 953   }
 954   END_OF_MED(LOC);
 955 }
 956
 957 //=======================================================================
 958 //function : getLine
 959 //purpose  :
 960 //=======================================================================
 961
 962 bool GIBI_MESH_RDONLY_DRIVER::getLine(char* & aLine)
 963 {
 964   bool aResult = true;
 965   // Check the state of the buffer;
 966   // if there is too little left, read the next portion of data
 967   int nBytesRest = _eptr - _ptr;
 968   if (nBytesRest < GIBI_MaxOutputLen)
 969   {
 970     if (nBytesRest > 0) {
 971       //memcpy (_start, _ptr, nBytesRest);
 972       char* tmpBuf = new char [nBytesRest];
 973       memcpy (tmpBuf, _ptr, nBytesRest);
 974       memcpy (_start, tmpBuf, nBytesRest);
 975       delete [] tmpBuf;
 976     } else
 977       nBytesRest = 0;
 978     _ptr = _start;
 979     const int nBytesRead = ::read (_File,
 980                                    &_start [nBytesRest],
 981                                    GIBI_BufferSize - nBytesRest);
 982     nBytesRest += nBytesRead;
 983     _eptr = &_start [nBytesRest];
 984   }
 985   // Check the buffer for the end-of-line
 986   char * ptr = _ptr;
 987   while (~0)
 988   {
 989     // Check for end-of-the-buffer, the ultimate criterion for termination
 990     if (ptr >= _eptr)
 991     {
 992       if (nBytesRest <= 0)
 993         aResult = false;
 994       else
 995         _eptr[-1] = '\0';
 996       break;
 997     }
 998     // seek the line-feed character
 999     if (ptr[0] == '\n')
1000     {
1001       if (ptr[-1] == '\r')
1002         ptr[-1] = '\0';
1003       ptr[0] = '\0';
1004       ++ptr;
1005       break;
1006     }
1007     ++ptr;
1008   }
1009   // Output the result
1010   aLine = _ptr;
1011   _ptr = ptr;
1012   _lineNb++;
1013
1014   return aResult;
1015 }
1016
1017 //=======================================================================
1018 //function : init
1019 //purpose  :
1020 //=======================================================================
1021
1022 void GIBI_MESH_RDONLY_DRIVER::init( int nbToRead, int nbPosInLine, int width, int shift )
1023 {
1024   _nbToRead = nbToRead;
1025   _nbPosInLine = nbPosInLine;
1026   _width = width;
1027   _shift = shift;
1028   _iPos = _iRead = 0;
1029   if ( _nbToRead ) {
1030     getNextLine( _curPos );
1031     _curPos = _curPos + _shift;
1032   }
1033   else
1034     _curPos = 0;
1035 }
1036
1037 //=======================================================================
1038 //function : next
1039 //purpose  : line getting
1040 //=======================================================================
1041
1042 void GIBI_MESH_RDONLY_DRIVER::next()
1043 {
1044   if ( !more() ) throw MEDEXCEPTION(LOCALIZED("!more()"));
1045   ++_iRead;
1046   ++_iPos;
1047   if ( _iRead < _nbToRead ) {
1048     if ( _iPos >= _nbPosInLine ) {
1049       getNextLine( _curPos );
1050       _curPos = _curPos + _shift;
1051       _iPos = 0;
1052     }
1053     else
1054       _curPos = _curPos + _width + _shift;
1055   }
1056   else
1057     _curPos = 0;
1058 }
1059
1060 //=======================================================================
1061 //function : getName
1062 //purpose  : names reading
1063 //=======================================================================
1064
1065 string GIBI_MESH_RDONLY_DRIVER::getName() const
1066 {
1067   int len = _width;
1068   while (( _curPos[len-1] == ' ' || _curPos[len-1] == 0) && len > 0 )
1069     len--;
1070   return string( _curPos, len );
1071 }
1072
1073 //=======================================================================
1074 //function : read
1075 //purpose  :
1076 //=======================================================================
1077
1078 void GIBI_MESH_RDONLY_DRIVER::read(void) throw (MEDEXCEPTION)
1079 {
1080   const char * LOC = "_GIBI_RDONLY_DRIVER::read() : " ;
1081   BEGIN_OF_MED(LOC);
1082
1083   if (_status!=MED_OPENED)
1084     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "The _idt of file " << _fileName << " is : "
1085                                  <<  " (the file is not opened)." )) ;
1086   if ( ! _ptrMesh->isEmpty() )
1087     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"Mesh object not empty : can't fill it!"));
1088
1089   _intermediateMED medi;
1090   try {
1091     if ( readFile( &medi, false )) {
1092       // impression résultats
1093       MESSAGE_MED(LOC << "GIBI_MESH_RDONLY_DRIVER::read : RESULTATS STRUCTURE INTERMEDIAIRES : ");
1094       MESSAGE_MED(LOC <<  medi );
1095
1096       fillMesh( &medi );
1097     }
1098   }
1099   catch (MEDEXCEPTION &ex)
1100   {
1101     INFOS_MED( ex.what() );
1102   }
1103   END_OF_MED(LOC);
1104 }
1105
1106 //=======================================================================
1107 //function : getReverseVector
1108 //purpose  :
1109 //=======================================================================
1110
1111 static void getReverseVector (const medGeometryElement type,
1112                               vector<pair<int,int> > & swapVec )
1113 {
1114   const char* LOC = "void getReverseVector()";
1115   BEGIN_OF_MED(LOC);
1116   swapVec.clear();
1117
1118   switch ( type ) {
1119   case MED_TETRA4:
1120     swapVec.resize(1);
1121     swapVec[0] = make_pair( 1, 2 );
1122     break;
1123   case MED_PYRA5:
1124     swapVec.resize(1);
1125     swapVec[0] = make_pair( 1, 3 );
1126     break;
1127   case MED_PENTA6:
1128     swapVec.resize(2);
1129     swapVec[0] = make_pair( 1, 2 );
1130     swapVec[1] = make_pair( 4, 5 );
1131     break;
1132   case MED_HEXA8:
1133     swapVec.resize(2);
1134     swapVec[0] = make_pair( 1, 3 );
1135     swapVec[1] = make_pair( 5, 7 );
1136     break;
1137   case MED_TETRA10:
1138     swapVec.resize(3);
1139     swapVec[0] = make_pair( 1, 2 );
1140     swapVec[1] = make_pair( 4, 5 );
1141     swapVec[2] = make_pair( 8, 9 );
1142     break;
1143   case MED_PYRA13:
1144     swapVec.resize(4);
1145     swapVec[0] = make_pair( 1, 3 );
1146     swapVec[1] = make_pair( 5, 8 );
1147     swapVec[2] = make_pair( 6, 7 );
1148     swapVec[3] = make_pair( 10, 12 );
1149     break;
1150   case MED_PENTA15:
1151     swapVec.resize(4);
1152     swapVec[0] = make_pair( 1, 2 );
1153     swapVec[1] = make_pair( 4, 5 );
1154     swapVec[2] = make_pair( 6, 8 );
1155     swapVec[3] = make_pair( 9, 11 );
1156     break;
1157   case MED_HEXA20:
1158     swapVec.resize(7);
1159     swapVec[0] = make_pair( 1, 3 );
1160     swapVec[1] = make_pair( 5, 7 );
1161     swapVec[2] = make_pair( 8, 11 );
1162     swapVec[3] = make_pair( 9, 10 );
1163     swapVec[4] = make_pair( 12, 15 );
1164     swapVec[5] = make_pair( 13, 14 );
1165     swapVec[6] = make_pair( 17, 19 );
1166     break;
1167 //   case MED_SEG3: no need to reverse edges
1168 //     swapVec.resize(1);
1169 //     swapVec[0] = make_pair( 1, 2 );
1170 //     break;
1171   case MED_TRIA6:
1172     swapVec.resize(2);
1173     swapVec[0] = make_pair( 1, 2 );
1174     swapVec[1] = make_pair( 3, 5 );
1175     break;
1176   case MED_QUAD8:
1177     swapVec.resize(3);
1178     swapVec[0] = make_pair( 1, 3 );
1179     swapVec[1] = make_pair( 4, 7 );
1180     swapVec[2] = make_pair( 5, 6 );
1181     break;
1182   default:;
1183   }
1184   END_OF_MED(LOC);
1185 }
1186
1187 //=======================================================================
1188 //function : reverse
1189 //purpose  : inverse element orientation using vector of indices to swap
1190 //=======================================================================
1191
1192 static void reverse(const _maille & aMaille, const vector<pair<int,int> > & swapVec )
1193 {
1194   _maille* ma = (_maille*) & aMaille;
1195   for ( int i = 0; i < swapVec.size(); ++i ) {
1196     std::swap( ma->sommets[ swapVec[i].first ],
1197                ma->sommets[ swapVec[i].second ]);
1198   }
1199   if ( swapVec.empty() )
1200     ma->reverse = true;
1201   else
1202     ma->reverse = false;
1203 }
1204
1205 //=======================================================================
1206 //function : getGibi2MedConnectivity
1207 //purpose  : return array of indices to transform GIBI connectivity to MED one
1208 //=======================================================================
1209
1210 static const int * getGibi2MedConnectivity( const medGeometryElement type )
1211 {
1212   static vector<const int*> conn;
1213   static int hexa20 [] = {0,6,4,2, 12,18,16,14, 7,5,3,1, 19,17,15,13, 8,11,10,9};
1214   static int penta15[] = {0,2,4, 9,11,13, 1,3,5, 10,12,14, 6,7,3};
1215   static int pyra13 [] = {0,2,4,6, 12, 1,3,5,7, 8,9,10,11};
1216   static int tetra10[] = {0,2,4, 9, 1,3,5, 6,7,8};
1217   static int quad8  [] = {0,2,4,6, 1,3,5,7};
1218   static int tria6  [] = {0,2,4, 1,3,5};
1219   static int seg3   [] = {0,2,1};
1220   if ( conn.empty() ) {
1221     conn.resize( MED_HEXA20 + 1, 0 );
1222     conn[ MED_HEXA20 ] = hexa20;
1223     conn[ MED_PENTA15] = penta15;
1224     conn[ MED_PYRA13 ] = pyra13;
1225     conn[ MED_TETRA10] = tetra10;
1226     conn[ MED_SEG3   ] = seg3;
1227     conn[ MED_TRIA6  ] = tria6;
1228     conn[ MED_QUAD8  ] = quad8;
1229   }
1230   return conn[ type ];
1231 }
1232
1233 //=======================================================================
1234 //function : fixConnectivity
1235 //purpose  : GIBI connectivity -> MED one
1236 //=======================================================================
1237
1238 static inline void fixConnectivity(const _maille & aMaille )
1239 {
1240   if ( const int * conn = getGibi2MedConnectivity( aMaille.geometricType )) {
1241     _maille* ma = (_maille*) & aMaille;
1242     //cout << "###### BEFORE fixConnectivity() " << *ma << endl;
1243     vector< _maille::TNoeud > newSommets( ma->sommets.size() );
1244     for ( int i = 0; i < newSommets.size(); ++i )
1245       newSommets[ i ] = ma->sommets[ conn[ i ]];
1246     ma->sommets = newSommets;
1247     //cout << "###### AFTER fixConnectivity() " << *ma << endl;
1248   }
1249 }
1250
1251 //=======================================================================
1252 //function : orientElements
1253 //purpose  :
1254 //=======================================================================
1255
1256 static void orientElements( _intermediateMED& medi )
1257 {
1258   MESSAGE_MED("orientElements()");
1259
1260   set<_maille>::const_iterator elemIt, elemEnd;
1261
1262   int type = -100;
1263   vector< pair<int,int> > swapVec;
1264
1265   bool isQuadratic = false;
1266   _maillageByDimIterator maIt(medi);
1267   while ( !isQuadratic && maIt.nextType() )
1268     isQuadratic = getGibi2MedConnectivity( maIt.type() );
1269
1270   if ( medi.points.begin()->second.coord.size() == 2 ) { // space dimension
1271
1272     // fix connectivity of quadratic edges
1273     if ( isQuadratic ) {
1274       if ( const set<_maille > * edges = _maillageByDimIterator( medi, 1 ).nextType() ) {
1275         elemIt = edges->begin(), elemEnd = edges->end();
1276         for ( ; elemIt != elemEnd; ++elemIt )
1277           fixConnectivity( *elemIt );
1278       }
1279     }
1280     // --------------------------
1281     // Orient 2D faces clockwise
1282     // --------------------------
1283
1284     _maillageByDimIterator faceMailIt( medi, 2 );
1285     while ( const set<_maille > * faces = faceMailIt.nextType() )
1286     {
1287       elemIt = faces->begin(), elemEnd = faces->end();
1288       for ( ; elemIt != elemEnd; elemIt++ )
1289       {
1290         // fix connectivity of quadratic faces
1291         if ( isQuadratic )
1292           fixConnectivity( *elemIt );
1293
1294         // look for index of the most left node
1295         int iLeft = 0, iNode, nbNodes = elemIt->sommets.size();
1296         if ( nbNodes > 4 ) // quadratic face
1297           nbNodes /= 2;
1298         double minX = elemIt->sommets[0]->second.coord[0];
1299         for ( iNode = 1; iNode < nbNodes; ++iNode )
1300         {
1301           if ( minX > elemIt->sommets[ iNode ]->second.coord[ 0 ]) {
1302             minX = elemIt->sommets[ iNode ]->second.coord[ 0 ];
1303             iLeft = iNode;
1304           }
1305         }
1306         // indeces of the nodes neighboring the most left one
1307         int iPrev = ( iLeft - 1 < 0 ) ? nbNodes - 1 : iLeft - 1;
1308         int iNext = ( iLeft + 1 == nbNodes ) ? 0 : iLeft + 1;
1309         // find components of prev-left and left-next vectors
1310         double xP = elemIt->sommets[ iPrev ]->second.coord[ 0 ];
1311         double yP = elemIt->sommets[ iPrev ]->second.coord[ 1 ];
1312         double xN = elemIt->sommets[ iNext ]->second.coord[ 0 ];
1313         double yN = elemIt->sommets[ iNext ]->second.coord[ 1 ];
1314         double xL = elemIt->sommets[ iLeft ]->second.coord[ 0 ];
1315         double yL = elemIt->sommets[ iLeft ]->second.coord[ 1 ];
1316         double xPL = xL - xP, yPL = yL - yP; // components of prev-left vector
1317         double xLN = xN - xL, yLN = yN - yL; // components of left-next vector
1318         // normalise y of the vectors
1319         double modPL = sqrt ( xPL * xPL + yPL * yPL );
1320         double modLN = sqrt ( xLN * xLN + yLN * yLN );
1321         if ( modLN > DBL_MIN && modPL > DBL_MIN )
1322         {
1323           yPL /= modPL;
1324           yLN /= modLN;
1325           // summury direction of neighboring links must be positive
1326           bool clockwise = ( yPL + yLN > 0 );
1327           if ( !clockwise ) {
1328             if ( elemIt->geometricType != type ) {
1329               type = elemIt->geometricType;
1330               getReverseVector( type, swapVec );
1331             }
1332             reverse( *elemIt, swapVec );
1333           }
1334         }
1335       }
1336     }
1337   }
1338   else {
1339
1340     // --------------------------------------
1341     // orient equally all connected 3D faces
1342     // --------------------------------------
1343     // connectivity of quadratic faces will be fixed by fixConnectivity()
1344     // in the next loop on elements
1345
1346     // fill map of links and their faces
1347     set<const _maille*> faces;
1348     map<const _maille*, _groupe*> fgm;
1349     map<_link, list<const _maille*> > linkFacesMap;
1350     map<_link, list<const _maille*> >::iterator lfIt, lfIt2;
1351
1352     medi.treatGroupes(); // erase groupes that wont be converted
1353     for (unsigned int i=0; i!=medi.groupes.size(); ++i)
1354     {
1355       _groupe& grp = medi.groupes[i];
1356       _groupe::TMailleIter maIt=grp.mailles.begin();
1357       if ( maIt==grp.mailles.end() || (*maIt)->dimensionWithPoly() != 2 )
1358         continue;
1359       for(; maIt!=grp.mailles.end(); ++maIt) {
1360         if ( faces.insert( &(**maIt )).second ) {
1361           for ( int j = 0; j < (int)(*maIt)->sommets.size(); ++j )
1362             linkFacesMap[ (*maIt)->link( j ) ].push_back( &(**maIt) );
1363           fgm.insert( make_pair( &(**maIt), &grp ));
1364         }
1365       }
1366     }
1367     // dump linkFacesMap
1368 //     for ( lfIt = linkFacesMap.begin(); lfIt!=linkFacesMap.end(); lfIt++) {
1369 //       cout<< "LINK: " << lfIt->first.first << "-" << lfIt->first.second << endl;
1370 //       list<const _maille*> & fList = lfIt->second;
1371 //       list<const _maille*>::iterator fIt = fList.begin();
1372 //       for ( ; fIt != fList.end(); fIt++ )
1373 //         cout << "\t" << **fIt << fgm[*fIt]->nom << endl;
1374 //     }
1375
1376     // Each oriented link must appear in one face only, else a face is reversed.
1377
1378     queue<const _maille*> faceQueue; // the queue contains well oriented faces
1379     // whose neighbors orientation is to be checked
1380
1381     bool manifold = true;
1382     while ( !linkFacesMap.empty() )
1383     {
1384       if ( faceQueue.empty() ) {
1385         ASSERT_MED( !linkFacesMap.begin()->second.empty() );
1386         faceQueue.push( linkFacesMap.begin()->second.front() );
1387       }
1388       while ( !faceQueue.empty() )
1389       {
1390         const _maille* face = faceQueue.front();
1391         faceQueue.pop();
1392
1393         // loop on links of <face>
1394         for ( int i = 0; i < (int)face->sommets.size(); ++i ) {
1395           _link link = face->link( i );
1396           // find the neighbor faces
1397           lfIt = linkFacesMap.find( link );
1398           int nbFaceByLink = 0;
1399           list< const _maille* > ml;
1400           if ( lfIt != linkFacesMap.end() )
1401           {
1402             list<const _maille*> & fList = lfIt->second;
1403             list<const _maille*>::iterator fIt = fList.begin();
1404             ASSERT_MED( fIt != fList.end() );
1405             for ( ; fIt != fList.end(); fIt++, nbFaceByLink++ ) {
1406               ml.push_back( *fIt );
1407               if ( *fIt != face ) // wrongly oriented neighbor face
1408               {
1409                 const _maille* badFace = *fIt;
1410                 // reverse and remove badFace from linkFacesMap
1411                 for ( int j = 0; j < (int)badFace->sommets.size(); ++j ) {
1412                   _link badlink = badFace->link( j );
1413                   if ( badlink == link ) continue;
1414                   lfIt2 = linkFacesMap.find( badlink );
1415                   if ( lfIt2 != linkFacesMap.end() ) {
1416                     list<const _maille*> & ff = lfIt2->second;
1417                     ff.erase( find( ff.begin(), ff.end(), badFace ));
1418                         if ( ff.empty() )
1419                           linkFacesMap.erase( lfIt2 );
1420                   }
1421                 }
1422                 badFace->reverse = true; // reverse
1423                 //INFOS_MED( "REVERSE " << *badFace );
1424                 faceQueue.push( badFace );
1425               }
1426             }
1427             linkFacesMap.erase( lfIt );
1428           }
1429           // add good neighbors to the queue
1430           _link revLink( link.second, link.first );
1431           lfIt = linkFacesMap.find( revLink );
1432           if ( lfIt != linkFacesMap.end() )
1433           {
1434             list<const _maille*> & fList = lfIt->second;
1435             list<const _maille*>::iterator fIt = fList.begin();
1436             for ( ; fIt != fList.end(); fIt++, nbFaceByLink++ ) {
1437               ml.push_back( *fIt );
1438               if ( *fIt != face )
1439                 faceQueue.push( *fIt );
1440             }
1441             linkFacesMap.erase( lfIt );
1442           }
1443           if ( nbFaceByLink > 2 ) {
1444             if ( manifold ) {
1445               list<const _maille*>::iterator i = ml.begin();
1446               INFOS_MED(nbFaceByLink << " faces by 1 link:");
1447               for( ; i!= ml.end(); i++ ) {
1448                 INFOS_MED("in object " << fgm[ *i ]->nom);
1449                 INFOS_MED( **i );
1450               }
1451             }
1452             manifold = false;
1453           }
1454         } // loop on links of the being checked face
1455       } // loop on the face queue
1456     } // while ( !linkFacesMap.empty() )
1457
1458     if ( !manifold )
1459       INFOS_MED(" -> Non manifold mesh, faces orientation may be incorrect");
1460
1461
1462     // ---------------------------------------------------
1463     // fix connectivity of quadratic elements
1464     // ---------------------------------------------------
1465
1466     for ( int dim = 1; dim <= 3; ++dim )
1467     {
1468       _maillageByDimIterator mailIt( medi, dim );
1469       while ( const set<_maille > * elems = mailIt.nextType() )
1470       {
1471         elemIt = elems->begin(), elemEnd = elems->end();
1472         for ( ; elemIt != elemEnd; elemIt++ )
1473         {
1474           // GIBI connectivity -> MED one
1475           if ( isQuadratic )
1476             fixConnectivity( *elemIt );
1477
1478           // reverse quadratic faces
1479           if ( elemIt->reverse ) {
1480             if ( elemIt->geometricType != type ) {
1481               type = elemIt->geometricType;
1482               getReverseVector( type, swapVec );
1483             }
1484             reverse ( *elemIt, swapVec );
1485           }
1486         }
1487       }
1488     }
1489
1490     // ---------------------------------------------------
1491     // Orient volumes according to MED conventions:
1492     // normal of a bottom (first) face should be downward
1493     // ---------------------------------------------------
1494
1495     _maillageByDimIterator mailIt( medi, 3 );
1496     while ( const set<_maille > * elems = mailIt.nextType() )
1497     {
1498       elemIt = elems->begin(), elemEnd = elems->end();
1499       int nbBottomNodes = 0;
1500       switch ( elemIt->geometricType ) {
1501       case MED_TETRA4:
1502       case MED_TETRA10:
1503       case MED_PENTA6:
1504       case MED_PENTA15:
1505         nbBottomNodes = 3; break;
1506       case MED_PYRA5:
1507       case MED_PYRA13:
1508       case MED_HEXA8:
1509       case MED_HEXA20:
1510         nbBottomNodes = 4; break;
1511       default: continue;
1512       }
1513       // find a normal to the bottom face
1514       const _noeud* n[4] = {
1515         &elemIt->sommets[0]->second, // 3 bottom nodes
1516         &elemIt->sommets[1]->second,
1517         &elemIt->sommets[2]->second,
1518         &elemIt->sommets[nbBottomNodes]->second };// a top node
1519       double vec01 [3] = { // vector n[0]-n[1]
1520         n[1]->coord[0] - n[0]->coord[0],
1521         n[1]->coord[1] - n[0]->coord[1],
1522         n[1]->coord[2] - n[0]->coord[2], };
1523       double vec02 [3] = { // vector n[0]-n[2]
1524         n[2]->coord[0] - n[0]->coord[0],
1525         n[2]->coord[1] - n[0]->coord[1],
1526         n[2]->coord[2] - n[0]->coord[2] };
1527       double normal [3] = { // vec01 ^ vec02
1528         vec01[1] * vec02[2] - vec01[2] * vec02[1],
1529         vec01[2] * vec02[0] - vec01[0] * vec02[2],
1530         vec01[0] * vec02[1] - vec01[1] * vec02[0] };
1531       // check if the 102 angle is convex
1532       if ( nbBottomNodes > 3 ) {
1533         const _noeud* n3 = &elemIt->sommets[nbBottomNodes-1]->second;// last bottom node
1534         double vec03 [3] = { // vector n[0]-n3
1535           n3->coord[0] - n[0]->coord[0],
1536           n3->coord[1] - n[0]->coord[1],
1537           n3->coord[2] - n[0]->coord[2], };
1538         if ( fabs( normal[0]+normal[1]+normal[2] ) <= DBL_MIN ) { // vec01 || vec02
1539           normal[0] = vec01[1] * vec03[2] - vec01[2] * vec03[1]; // vec01 ^ vec03
1540           normal[1] = vec01[2] * vec03[0] - vec01[0] * vec03[2];
1541           normal[2] = vec01[0] * vec03[1] - vec01[1] * vec03[0];
1542         }
1543         else {
1544           double vec [3] = { // normal ^ vec01
1545             normal[1] * vec01[2] - normal[2] * vec01[1],
1546             normal[2] * vec01[0] - normal[0] * vec01[2],
1547             normal[0] * vec01[1] - normal[1] * vec01[0] };
1548           double dot2 = vec[0]*vec03[0] + vec[1]*vec03[1] + vec[2]*vec03[2]; // vec*vec03
1549           if ( dot2 < 0 ) { // concave -> reverse normal
1550             normal[0] *= -1;
1551             normal[1] *= -1;
1552             normal[2] *= -1;
1553           }
1554         }
1555       }
1556       // direction from top to bottom
1557       vector<double> tbDir(3);
1558       tbDir[0] = n[0]->coord[0] - n[3]->coord[0];
1559       tbDir[1] = n[0]->coord[1] - n[3]->coord[1];
1560       tbDir[2] = n[0]->coord[2] - n[3]->coord[2];
1561       // compare 2 directions: normal and top-bottom
1562       double dot = normal[0]*tbDir[0] + normal[1]*tbDir[1] + normal[2]*tbDir[2];
1563       if ( dot < 0. ) // need reverse
1564       {
1565         getReverseVector( elemIt->geometricType, swapVec );
1566         for ( ; elemIt != elemEnd; elemIt++ )
1567           reverse( *elemIt, swapVec );
1568       }
1569     } // loop on volumes of one geometry
1570
1571   } // space dimension == 3
1572 }
1573
1574 //=======================================================================
1575 //function : fillMesh
1576 //purpose  : load data from medi to mesh
1577 //=======================================================================
1578
1579 void GIBI_MESH_RDONLY_DRIVER::fillMesh(_intermediateMED* _ptrMedi)
1580 {
1581   const char* LOC = "GIBI_MESH_RDONLY_DRIVER::fillMesh(_intermediateMED* _ptrMedi) : ";
1582   BEGIN_OF_MED(LOC);
1583
1584   _ptrMesh->_name = _meshName;
1585
1586   if (_ptrMedi)
1587   {
1588     if (!_maillageByDimIterator(*_ptrMedi).nextType() ||
1589         _ptrMedi->groupes.empty() ||
1590         _ptrMedi->points.empty()) {
1591       INFOS_MED(" Error while reading file: the data read are not completed " ) ;
1592       return;
1593     }
1594     // fix element orientation
1595     orientElements( *_ptrMedi );
1596
1597     _ptrMesh->_spaceDimension = _ptrMedi->points.begin()->second.coord.size();
1598     _ptrMesh->_meshDimension  = _ptrMedi->getMeshDimension();
1599     _ptrMesh->_numberOfNodes  = _ptrMedi->points.size();
1600     _ptrMesh->_isAGrid        = 0;
1601     _ptrMesh->_coordinate     = _ptrMedi->getCoordinate();
1602
1603     //Construction des groupes
1604     _ptrMedi->getGroups(_ptrMesh->_groupCell,
1605                         _ptrMesh->_groupFace,
1606                         _ptrMesh->_groupEdge,
1607                         _ptrMesh->_groupNode, _ptrMesh);
1608
1609     _ptrMesh->_connectivity = _ptrMedi->getConnectivity();
1610
1611     // calcul de la connectivite d-1 complete, avec renumerotation des groupes
1612     //if (_ptrMesh->_spaceDimension==3)
1613     //    _ptrMesh->_connectivity->updateGroup(_ptrMesh->_groupFace) ;
1614     //else if (_ptrMesh->_spaceDimension==2)
1615     //    _ptrMesh->_connectivity->updateGroup(_ptrMesh->_groupEdge) ;
1616
1617     // Creation des familles à partir des groupes
1618     // NC : Cet appel pourra être différé quand la gestion de la cohérence famille/groupes sera assurée
1619     _ptrMesh->createFamilies();
1620     // TAKE CARE OF ELEMENTS ORDER IN GROUPS AFTER THEIR SPLITING INTO FAMILIES !!!!
1621     // _ptrMesh->createFamilies() breaks the order
1622 //     _ptrMedi->getFamilies(_ptrMesh->_familyCell,
1623 //                           _ptrMesh->_familyFace,
1624 //                           _ptrMesh->_familyEdge,
1625 //                           _ptrMesh->_familyNode, _ptrMesh);
1626
1627     // add attributes to families
1628     set<string> famNames;
1629     for (medEntityMesh entity=MED_CELL; entity<MED_ALL_ENTITIES; ++entity)
1630     {
1631       int i, nb = _ptrMesh->getNumberOfFamilies(entity);
1632       for ( i = 1; i <= nb; ++i ) {
1633         FAMILY* f = const_cast<FAMILY*>( _ptrMesh->getFamily( entity, i ));
1634         f->setNumberOfAttributes( 1 );
1635         int* attIDs = new int[1];
1636         attIDs[0] = 1;
1637         f->setAttributesIdentifiers( attIDs );
1638         int* attVals = new int[1];
1639         attVals[0] = 1;
1640         f->setAttributesValues( attVals );
1641         string* attDescr = new string[1];
1642         attDescr[0] = "med_family";
1643         f->setAttributesDescriptions( attDescr );
1644         // limit a name length
1645         if ( f->getName().length() > 31 ) {
1646           ostringstream name;
1647           name << "FAM" << f->getIdentifier();
1648           f->setName( name.str());
1649         }
1650         // check if family is on the whole mesh entity
1651         if (_ptrMesh->getNumberOfElements( entity, MED_ALL_ELEMENTS ) ==
1652             f->getNumberOfElements( MED_ALL_ELEMENTS ))
1653           f->setAll( true );
1654       }
1655       // setAll() for groups
1656       nb = _ptrMesh->getNumberOfGroups(entity);
1657       for ( i = 1; i <= nb; ++i ) {
1658         GROUP * g = const_cast<GROUP*>( _ptrMesh->getGroup( entity, i ));
1659         if (_ptrMesh->getNumberOfElements( entity, MED_ALL_ELEMENTS ) ==
1660             g->getNumberOfElements( MED_ALL_ELEMENTS ))
1661           g->setAll( true );
1662       }
1663     }
1664   }
1665   END_OF_MED(LOC);
1666 }
1667
1668 void GIBI_MESH_RDONLY_DRIVER::write( void ) const
1669   throw (MEDEXCEPTION)
1670 {
1671   throw MEDEXCEPTION("GIBI_MESH_RDONLY_DRIVER::write : Can't write with a RDONLY driver !");
1672 }
1673
1674
1675 /*--------------------- WRONLY PART -------------------------------*/
1676
1677 GIBI_MESH_WRONLY_DRIVER::GIBI_MESH_WRONLY_DRIVER():GIBI_MESH_DRIVER()
1678 {
1679 }
1680 GIBI_MESH_WRONLY_DRIVER::GIBI_MESH_WRONLY_DRIVER(const string & fileName,
1681                                                  MESH * ptrMesh):
1682   GIBI_MESH_DRIVER(fileName,ptrMesh,WRONLY)
1683 {
1684   MESSAGE_MED("GIBI_MESH_WRONLY_DRIVER::GIBI_MESH_WRONLY_DRIVER(const string & fileName, MESH * ptrMesh) has been created");
1685 }
1686 GIBI_MESH_WRONLY_DRIVER::GIBI_MESH_WRONLY_DRIVER(const GIBI_MESH_WRONLY_DRIVER & driver):
1687   GIBI_MESH_DRIVER(driver)
1688 {
1689 }
1690 GIBI_MESH_WRONLY_DRIVER::~GIBI_MESH_WRONLY_DRIVER()
1691 {
1692   //MESSAGE_MED("GIBI_MESH_WRONLY_DRIVER::GIBI_MESH_WRONLY_DRIVER(const string & fileName, MESH * ptrMesh) has been destroyed");
1693 }
1694 GENDRIVER * GIBI_MESH_WRONLY_DRIVER::copy(void) const
1695 {
1696   return new GIBI_MESH_WRONLY_DRIVER(*this);
1697 }
1698 void GIBI_MESH_WRONLY_DRIVER::read (void)
1699   throw (MEDEXCEPTION)
1700 {
1701   throw MEDEXCEPTION("GIBI_MESH_WRONLY_DRIVER::read : Can't read with a WRONLY driver !");
1702 }
1703
1704 //=======================================================================
1705 //function : open
1706 //purpose  :
1707 //=======================================================================
1708
1709 void GIBI_MESH_WRONLY_DRIVER::open()
1710   //     throw (MEDEXCEPTION)
1711 {
1712   if( _status == MED_OPENED )
1713     return;
1714
1715   const char * LOC = "GIBI_MESH_DRIVER::open()" ;
1716   BEGIN_OF_MED(LOC);
1717
1718   MED_EN::med_mode_acces aMode = getAccessMode();
1719   switch (aMode) {
1720   case MED_EN::RDWR:
1721   case MED_EN::WRONLY: // should never append !!
1722     _gibi.open(_fileName.c_str(), ios::out);
1723     break;
1724   default:
1725     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "Bad file mode access ! " << aMode));
1726   }
1727   //change for windows compilation
1728   if ( !_gibi ||
1729 #ifdef WNT
1730       !_gibi.is_open()
1731 #else
1732       !_gibi.rdbuf()->is_open()
1733 #endif
1734       )
1735   {
1736     _status = MED_CLOSED;
1737     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Could not open file "<<_fileName));
1738   }
1739   else
1740   {
1741     _status = MED_OPENED;
1742   }
1743   END_OF_MED(LOC);
1744 }
1745
1746 //=======================================================================
1747 //function : close
1748 //purpose  :
1749 //=======================================================================
1750
1751 void GIBI_MESH_WRONLY_DRIVER::close()
1752   //  throw (MEDEXCEPTION)
1753 {
1754   const char* LOC = "GIBI_MESH_DRIVER::close() ";
1755   BEGIN_OF_MED(LOC);
1756     if ( _status == MED_OPENED)
1757     {
1758         _gibi.close();
1759         _status = MED_CLOSED;
1760     }
1761   END_OF_MED(LOC);
1762 }
1763
1764 //=======================================================================
1765 //function : write
1766 //purpose  :
1767 //=======================================================================
1768
1769 void GIBI_MESH_WRONLY_DRIVER::write(void) const
1770   throw (MEDEXCEPTION)
1771 {
1772   const char* LOC = "void GIBI_MESH_WRONLY_DRIVER::write(void) const : ";
1773   BEGIN_OF_MED(LOC);
1774
1775   // we are going to modify the _gibi field
1776   GIBI_MESH_WRONLY_DRIVER * me = const_cast<GIBI_MESH_WRONLY_DRIVER *>(this);
1777 //  try {
1778     me->writeSupportsAndMesh();
1779     me->writeLastRecord();
1780 //   }
1781 //   catch (MEDEXCEPTION &ex)
1782 //   {
1783 //     INFOS_MED( ex.what() );
1784 //   }
1785
1786   END_OF_MED(LOC);
1787 }
1788
1789 //=======================================================================
1790 //function : getName
1791 //purpose  : return cleaned up support name
1792 //=======================================================================
1793
1794 static string cleanName( const string& theName )
1795 {
1796   string name = theName;
1797   if ( !name.empty() ) {
1798     // find a name string end
1799     int i, len = name.length();
1800     for ( i = 0; i < len; ++i ) {
1801       if ( name[i] == 0 )
1802         break;
1803     }
1804     // cut off trailing white spaces
1805     while ( i > 0 && name[i-1] == ' ' )
1806       i--;
1807     if ( i != len ) {
1808       name = name.substr( 0, i );
1809       len = i;
1810     }
1811   }
1812   return name;
1813 }
1814
1815 //=======================================================================
1816 //function : addSupport
1817 //purpose  :
1818 //=======================================================================
1819
1820 bool GIBI_MESH_WRONLY_DRIVER::addSupport( const SUPPORT * support )
1821 {
1822   if ( !support )
1823     return false;
1824   map<const SUPPORT*,supportData>::iterator su = _supports.find( support );
1825   if ( su != _supports.end() )
1826     return ( su->second.getNumberOfTypes() > 0 );
1827
1828   if ( support->getMesh() != _ptrMesh )
1829     throw MEDEXCEPTION(LOCALIZED(STRING("cant write support of other mesh" )));
1830
1831   // get sub-supports and define a support type name
1832   string supType;
1833   list<const SUPPORT*> sList;
1834   const GROUP* group = dynamic_cast< const GROUP* >(support);
1835   if ( group )
1836   {
1837     if ( group->getNumberOfTypes() > 0 || group->isOnAllElements() )
1838       sList.push_back( group );
1839     else {
1840       int iFam, nbFam = group->getNumberOfFamilies();
1841       for ( iFam = 1; iFam <= nbFam; ++iFam )
1842         sList.push_back( group->getFamily( iFam ));
1843     }
1844     supType = "group";
1845   }
1846   else
1847   {
1848     sList.push_back( support );
1849     supType = dynamic_cast< const FAMILY* >(support) ? "family" : "support";
1850   }
1851
1852   supportData & data = _supports[ support ];
1853   data._cleanName = cleanName( support->getName() );
1854
1855   // check if it is a writtable support, i.e.
1856   // nodal connectivity for a support entity exists
1857   medEntityMesh entity = support->getEntity();
1858   if ( entity != MED_NODE && !_ptrMesh->existConnectivity( MED_NODAL, entity )) {
1859     INFOS_MED("Do not save " << supType << " of entity " << entity
1860           << " named <" << data._cleanName << "> nodal connectivity not defined");
1861     return false;
1862   }
1863
1864   // fill supportData
1865   list<const SUPPORT*>::iterator sIt = sList.begin();
1866   for ( ; sIt != sList.end(); sIt++ )
1867   {
1868     bool onAll = (*sIt)->isOnAllElements();
1869     int nbTypes = 0;
1870     if ( !onAll )
1871       nbTypes = (*sIt)->getNumberOfTypes();
1872     else
1873       nbTypes = _ptrMesh->getNumberOfTypes( entity );
1874     if ( nbTypes == 0 )
1875       continue;
1876     const medGeometryElement* types = 0;
1877     if ( !onAll )
1878       types = (*sIt)->getTypes();
1879     else if ( entity != MED_NODE )
1880       types = _ptrMesh->getTypes( entity );
1881     for ( int iType = 0; iType < nbTypes; ++iType )
1882     {
1883       if ( types && types[ iType ] > MED_HEXA20 )
1884         continue; // poly
1885       medGeometryElement geomType = types ? types[ iType ] : MED_ALL_ELEMENTS;
1886       const int * ptrElemIDs = 0;
1887       int elemID1 = 0, nbElems = 0;
1888       if ( onAll ) {
1889         nbElems = _ptrMesh->getNumberOfElements( entity, geomType );
1890         elemID1 = (entity == MED_NODE) ? 1 : _ptrMesh->getGlobalNumberingIndex (entity)[ iType ];
1891       }
1892       else {
1893         nbElems = (*sIt)->getNumberOfElements( geomType );
1894         ptrElemIDs = (*sIt)->getNumber( geomType );
1895       }
1896       if ( geomType == 0 || ( entity == MED_NODE ))
1897         geomType = MED_POINT1;
1898
1899       data.addTypeData( geomType, nbElems, ptrElemIDs, elemID1 );
1900     }
1901   }
1902
1903   if ( data.getNumberOfTypes() == 0 ) {
1904     INFOS_MED("Do not save " << supType << " of entity " << entity
1905           << " named <" << data._cleanName << "> no geometric types");
1906     return false;
1907   }
1908
1909   return true;
1910 }
1911
1912 //=======================================================================
1913 //function : getSupportIndex
1914 //purpose  :
1915 //=======================================================================
1916
1917 int GIBI_MESH_WRONLY_DRIVER::getSubMeshIdAndSize(const SUPPORT *        support,
1918                                                  list<pair<int,int> > & idsAndSizes) const
1919 {
1920   idsAndSizes.clear();
1921   map<const SUPPORT*,supportData>::const_iterator su = _supports.find( support );
1922   if ( su == _supports.end() )
1923     return 0;
1924
1925   supportData * data = const_cast<supportData *>( & su->second );
1926   int id = data->_id;
1927   if ( data->getNumberObjects() > data->getNumberOfTypes() )
1928     id++;
1929   supportData::typeIterator tIt = data->_types.begin();
1930   for ( ; tIt != data->_types.end(); ++tIt )
1931   {
1932     int size = 0;
1933     list< typeData >& td = tIt->second;
1934     list< typeData >::iterator tdIt = td.begin();
1935     for ( ; tdIt != td.end(); ++tdIt )
1936       size += tdIt->_nbElems;
1937     idsAndSizes.push_back( make_pair( id++, size ));
1938   }
1939   return idsAndSizes.size();
1940 }
1941
1942 // ============================================================
1943 // the class writes endl to the file as soon as <limit> fields
1944 // have been written after the last endl
1945 // ============================================================
1946
1947 class TFieldCounter
1948 {
1949   fstream& _file;
1950   int _count, _limit;
1951  public:
1952   TFieldCounter(fstream& f, int limit=0): _file(f), _limit(limit) { init(); }
1953   void init(int limit=0) // init, is done by stop() as well
1954   { if (limit) _limit = limit; _count = 0; }
1955   void operator++(int) // next
1956   { if ( ++_count == _limit ) { _file << endl; init(); }}
1957   void stop() // init() and write endl if there was no endl after the last written field
1958   { if ( _count ) _file << endl; init(); }
1959 };
1960
1961 //=======================================================================
1962 //function : writeElements
1963 //purpose  : ptrElemIDs and elemID1 provide two alternative ways of giving
1964 //           elements to write.
1965 //           If elemSet != 0 then an element is
1966 //           ( addElemInSet ? <written and added to elemSet> : <ignored if id is in elemSet>)
1967 //=======================================================================
1968
1969 void GIBI_MESH_WRONLY_DRIVER::writeElements (medGeometryElement geomType,
1970                                              list< typeData >&  typeDataList,
1971                                              const int *        nodalConnect,
1972                                              const int *        nodalConnectIndex)
1973 {
1974   // ITYPEL : type de l'鬩ment 1=point, 2=segment ?eux noeuds...
1975   // NBSOUS : nombre de sous parties dans cet objet,
1976   //          une sous partie par type d'鬩ments le composant.
1977   // NBREF : nombre de sous r馩rences. Une r馩rence est par exemple le contour
1978   // NBNOEL : nombre de noeuds par 鬩ment
1979   // NBEL : nombre d'鬩ments
1980
1981   int castemType = GIBI_MESH_DRIVER::med2gibiGeom( geomType );
1982   const char* zeroI8 = "       0"; // FORMAT(I8)
1983   int nbElemNodes = geomType % 100;
1984
1985   // indices to transform MED connectivity to GIBI one
1986   vector< int > toGibiConn;
1987   toGibiConn.reserve( nbElemNodes );
1988   if ( const int * toMedConn = getGibi2MedConnectivity( geomType )) {
1989     toGibiConn.resize( nbElemNodes );
1990     for ( int i = 0; i < nbElemNodes; ++i )
1991       toGibiConn[ toMedConn[ i ]] = i;
1992   }
1993   else {
1994     while ( toGibiConn.size() < nbElemNodes )
1995       toGibiConn.push_back( toGibiConn.size() );
1996   }
1997
1998   // count total nb of elements
1999   int nbElements = 0;
2000   list< typeData >::iterator td = typeDataList.begin();
2001   for ( ; td != typeDataList.end(); td++ )
2002     nbElements += td->_nbElems;
2003
2004   _gibi << setw(8) << castemType <<  // ITYPE
2005     zeroI8 <<                       // NBSOUS
2006       zeroI8 <<                     // NBREF
2007         setw(8) << nbElemNodes <<   // NBNOEL
2008           setw(8) << nbElements <<  // NBEL
2009             endl;
2010
2011   MESSAGE_MED("writeElements(): geomType=" << geomType << " nbElements= " << nbElements)
2012
2013   // L 'enregistrement donnant le num? de la couleur des 鬩ments.
2014   // * 8000 FORMAT(10I8)
2015   TFieldCounter fcount( _gibi, 10 );
2016   int iElem = 0;
2017   for ( ; iElem < nbElements; ++iElem, fcount++ )
2018     _gibi << zeroI8;
2019   fcount.stop();
2020
2021   // Tableau des connectivit鳮 Description du premier 鬩ment puis du deuxi譥...
2022   // ATTENTION il ne s'agit pas de la num?tation vraie,
2023   // il faut la faire passer par le filtre du dernier tableau de la pile num? 32.
2024   //int nbSkipped = 0;
2025
2026   for ( td = typeDataList.begin(); td != typeDataList.end(); td++ )
2027   {
2028     for ( int i = 0; i < td->_nbElems; i++ )
2029     {
2030       iElem = td->_ptrElemIDs ? td->_ptrElemIDs[ i ] : td->_elemID1 + i;
2031       if ( geomType == MED_POINT1 )
2032       {
2033         _gibi << setw(8) << iElem;
2034         fcount++;
2035       }
2036       else
2037       {
2038         int nodeId = nodalConnectIndex[ iElem - 1 ] - 1;
2039         for ( int iNode = 0; iNode < nbElemNodes; ++iNode, fcount++ ) {
2040           _gibi << setw(8) << nodalConnect[ nodeId + toGibiConn[ iNode ]];
2041         }
2042       }
2043     }
2044   }
2045
2046   fcount.stop();
2047 }
2048
2049 //=======================================================================
2050 //function : addName
2051 //purpose  : make name uppercase and shorter than 9, add it to nameNbMap,
2052 //           raise if not unique
2053 //=======================================================================
2054
2055 void GIBI_MESH_WRONLY_DRIVER::addName(map<string,int>& nameMap,
2056                                       string&          theName,
2057                                       int              index,
2058                                       string           prefix)
2059 {
2060   string name = cleanName( theName );
2061   if ( !name.empty() ) {
2062     int len = name.length();
2063 #ifdef THROW_ON_BAD_NAME
2064     if ( len > 8 )
2065       throw MEDEXCEPTION(STRING("Can't write name longer than 8: ") << name );
2066
2067     for ( int i = 0; i < len; ++i )
2068       name[i] = toupper( name[i] );
2069     if ( ! nameMap.insert( make_pair( name, index )).second )
2070       throw MEDEXCEPTION(STRING("Can't write not unique name: ") << name );
2071 #else
2072     bool ok = ( len <= 8 && len > 0 );
2073     if ( ok ) {
2074       for ( int i = 0; i < len; ++i )
2075         name[i] = toupper( name[i] );
2076       ok = nameMap.insert( make_pair( name, index )).second;
2077     }
2078     if ( !ok ) {
2079       char *str=new char[ prefix.size() + 13 ];
2080       int j = 1;
2081       do {
2082         sprintf( str, "%s_%d", prefix.c_str(), nameMap.size()+j );
2083         ok = nameMap.insert( make_pair( str, index )).second;
2084         j++;
2085       } while ( !ok );
2086       INFOS_MED( "Save <" << name << "> as <" << str << ">");
2087       delete [] str;
2088     }
2089 #endif
2090   }
2091 }
2092
2093 //=======================================================================
2094 //function : writeNames
2095 //purpose  :
2096 //=======================================================================
2097
2098 void GIBI_MESH_WRONLY_DRIVER::writeNames( map<string,int>& nameNbMap )
2099 {
2100   // La pile num? 1 est celle des objets de type maillage.
2101   // La ligne suivante donne le nom des objets maillages sauv鳮
2102   // * 8001       FORMAT(8(1X,A8))
2103   if ( !nameNbMap.empty() )
2104   {
2105     TFieldCounter fcount( _gibi, 8 );
2106     _gibi << left;
2107     map<string,int>::iterator nameNbIt = nameNbMap.begin();
2108     for ( ; nameNbIt != nameNbMap.end(); nameNbIt++, fcount++ ) {
2109       _gibi << " " << setw(8) << nameNbIt->first;
2110     }
2111     fcount.stop();
2112     _gibi << right;
2113     // La ligne suivante donne les num?s d'ordre, dans la pile,
2114     // des objets nomm?cit?pr飩demment.
2115     // *  8000 FORMAT(10I8)
2116     nameNbIt = nameNbMap.begin();
2117     for ( fcount.init(10); nameNbIt != nameNbMap.end(); nameNbIt++, fcount++ )
2118       _gibi << setw(8) << nameNbIt->second;
2119     fcount.stop();
2120   }
2121 }
2122
2123 //=======================================================================
2124 //function : writeSupportsAndMesh
2125 //purpose  :
2126 //=======================================================================
2127
2128 void GIBI_MESH_WRONLY_DRIVER::writeSupportsAndMesh()
2129 {
2130   const char * LOC = "void GIBI_MESH_WRONLY_DRIVER::writeSupportsAndMesh() ";
2131   BEGIN_OF_MED(LOC);
2132
2133   if (_status!=MED_OPENED)
2134     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "file " << _fileName<<  " is not opened." ));
2135   if (!_ptrMesh)
2136     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "can't write a NULL mesh" ));
2137   if (!_ptrMesh->getConnectivityptr())
2138     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "can't write a mesh with NULL connectivity" ));
2139
2140   // fill _supports with families and groups
2141   medEntityMesh entity;
2142   for (entity=MED_CELL; entity<MED_ALL_ENTITIES; ++entity)
2143   {
2144     int i, nb = _ptrMesh->getNumberOfGroups(entity);
2145     for ( i = 1; i <= nb; ++i )
2146       addSupport( _ptrMesh->getGroup( entity, i ));
2147 //     nb = _ptrMesh->getNumberOfFamilies(entity);
2148 //     for ( i = 1; i <= nb; ++i )
2149 //       addSupport( _ptrMesh->getFamily( entity, i ));
2150   }
2151
2152   // --------------------------------------------------------------------
2153   // Count total nb of objects: an object per an element type in support
2154   // plus an object per an element type not used in _supports.
2155   // Collect object names
2156   // --------------------------------------------------------------------
2157
2158   vector<int> nbSuppElemsByType(MED_HEXA20+1,0);
2159   map<string,int> nameNbMap;
2160   map<const SUPPORT*,supportData>::iterator supIt = _supports.begin();
2161   int i, nb_objects = 0;
2162   for ( ; supIt != _supports.end(); supIt++ )
2163   {
2164     supportData & data = supIt->second;
2165     int nbSupObj = data.getNumberObjects();
2166     if ( nbSupObj == 0 )
2167       continue;
2168     data._id = nb_objects + 1;
2169     nb_objects += nbSupObj;
2170
2171     addName( nameNbMap, data._cleanName, data._id, "C" );
2172     MESSAGE_MED( "obj " << data._id << " " << data._cleanName);
2173
2174     // count elements: take into account supports on all elements and families only
2175     const SUPPORT* support = supIt->first;
2176     if ( support->isOnAllElements() || dynamic_cast< const FAMILY* >( support ))
2177     {
2178       supportData::typeIterator tIt = data._types.begin();
2179       for ( ; tIt != data._types.end(); ++tIt )
2180         if ( support->isOnAllElements() )
2181         {
2182           nbSuppElemsByType[ tIt->first ] = INT_MAX / 100;
2183         }
2184         else
2185         {
2186           list< typeData >& td = tIt->second;
2187           list< typeData >::iterator tdIt = td.begin();
2188           for ( ; tdIt != td.end(); ++tdIt )
2189             nbSuppElemsByType[ tIt->first] += tdIt->_nbElems;
2190         }
2191     }
2192   }
2193
2194   // count types of mesh elements that are not all in _supports
2195   int iType, nbTypes;
2196   entity = _ptrMesh->getConnectivityptr()->getEntity();
2197   for ( ; entity < MED_NODE; entity++ )
2198   {
2199     nbTypes = _ptrMesh->getNumberOfTypes( entity );
2200     if ( nbTypes == 0 || !_ptrMesh->existConnectivity( MED_NODAL, entity ))
2201       continue;
2202     const medGeometryElement* types = _ptrMesh->getTypes( entity );
2203     for ( iType = 0; iType < nbTypes; ++iType )
2204     {
2205       int nbElemInSups = nbSuppElemsByType[ types[ iType ]];
2206       int nbElemInMesh = _ptrMesh->getNumberOfElements(entity, types[ iType ]);
2207       if ( nbElemInSups < nbElemInMesh ) {
2208         nb_objects++;
2209         nbSuppElemsByType[ types[ iType ]] = -1; // to keep written elements of _supports
2210       }
2211     }
2212   }
2213
2214   // ------------
2215   // Write file
2216   // ------------
2217
2218   // Premier paquet dont le nombre de lignes ne varie pas.
2219   // On y trouve des indications g鮩rales.
2220   const int dim = _ptrMesh->getSpaceDimension();
2221   _gibi << " ENREGISTREMENT DE TYPE   4" << endl;
2222   _gibi << " NIVEAU  15 NIVEAU ERREUR   0 DIMENSION   " << dim <<endl;
2223   _gibi << " DENSITE  .00000E+00" << endl;
2224   _gibi << " ENREGISTREMENT DE TYPE   7" << endl;
2225   _gibi << " NOMBRE INFO CASTEM2000   8" <<endl;
2226   _gibi << " IFOUR  -1 NIFOUR   0 IFOMOD  -1 IECHO   1 IIMPI   0 IOSPI   0 ISOTYP   1" << endl;
2227   _gibi << " NSDPGE     0" << endl;
2228
2229   // Deuxi譥 paquet qui d馩nit toutes les piles
2230   // (une pile par type d'objet et certaines piles en plus).
2231   // Un enregistrement de type 2 pr鶩ent de l'飲iture d'une nouvelle pile,
2232   // celui de type 5 pr鶩ent de la fin.
2233   // * 800   FORMAT (' ENREGISTREMENT DE TYPE', I4)
2234   _gibi << " ENREGISTREMENT DE TYPE   2" << endl;
2235   // * 801     FORMAT(' PILE NUMERO',I4,'NBRE OBJETS NOMMES',I8,'NBRE OBJETS',I8)
2236   _gibi << " PILE NUMERO   1NBRE OBJETS NOMMES" << setw(8) << nameNbMap.size() <<
2237     "NBRE OBJETS" << setw(8) << nb_objects <<endl;
2238
2239   writeNames( nameNbMap );
2240
2241   // Passage ?a description des objets les uns apr賠les autres.
2242   // Le premier enregistrement de chaque objet est compos頤e 5 nombres repr鳥ntant :
2243   // ITYPEL : type de l'鬩ment 1=point, 2=segment ?eux noeuds...
2244   // NBSOUS : nombre de sous parties dans cet objet,
2245   //          une sous partie par type d'鬩ments le composant.
2246   // NBREF : nombre de sous r馩rences. Une r馩rence est par exemple le contour
2247   // NBNOEL : nombre de noeuds par 鬩ment
2248   // NBEL : nombre d'鬩ments
2249   // Si ITYPEL=0 alors NBSOUS diff鲥nt de z?. Dans ce cas on lira la liste des positions,
2250   // dans la pile des objets, des sous parties le composant.
2251   // Si NBSOUS=0, NBNOEL et NBEL sont diff鲥nts de z?, on trouve, au besoin,
2252   // la liste des r馩rences , les num?s des couleurs puis les connectivit鳮
2253
2254   TFieldCounter fcount( _gibi, 10 );
2255   const char* zeroI8 = "       0"; // FORMAT(I8)
2256   for ( supIt = _supports.begin(); supIt != _supports.end(); supIt++ )
2257   {
2258     supportData & data = supIt->second;
2259     int nbSupObj = data.getNumberObjects();
2260     if ( nbSupObj == 0 )
2261       continue;
2262     MESSAGE_MED("support " << data._id << "<" << data._cleanName << ">");
2263
2264     // write a compound object
2265     int nbTypes = data.getNumberOfTypes();
2266     if ( nbSupObj > nbTypes )
2267     {
2268       _gibi << zeroI8 << setw(8) << nbTypes << zeroI8 << zeroI8 << zeroI8 << endl;
2269       for ( int i_sub = 1; i_sub <= nbTypes; ++i_sub, fcount++ )
2270         _gibi << setw(8) << ( data._id + i_sub );
2271       fcount.stop();
2272     }
2273
2274     // write components
2275     entity = supIt->first->getEntity();
2276     const int * nodalConnect = 0, * nodalConnectIndex = 0;
2277     if ( entity != MED_NODE ) {
2278       nodalConnect = _ptrMesh->getConnectivity (MED_FULL_INTERLACE,
2279                                                 MED_NODAL,
2280                                                 entity,
2281                                                 MED_ALL_ELEMENTS);
2282       nodalConnectIndex = _ptrMesh->getConnectivityIndex (MED_NODAL,entity);
2283     }
2284     supportData::typeIterator tIt = data._types.begin();
2285     for ( ; tIt != data._types.end(); ++tIt )
2286     {
2287       writeElements (tIt->first,
2288                      tIt->second,
2289                      nodalConnect,
2290                      nodalConnectIndex);
2291     }
2292   }  // loop on _supports
2293
2294   // Write elements that are not in _supports
2295
2296   supportData data;
2297   entity = _ptrMesh->getConnectivityptr()->getEntity();
2298   for ( ; entity < MED_NODE; entity++ )
2299   {
2300     int nbTypes = _ptrMesh->getNumberOfTypes( entity );
2301     if ( nbTypes == 0 || !_ptrMesh->existConnectivity( MED_NODAL, entity ))
2302       continue;
2303     const medGeometryElement* types = _ptrMesh->getTypes( entity );
2304     const int * nbIndex = _ptrMesh->getGlobalNumberingIndex (entity);
2305     const int * nodalConnect = 0, * nodalConnectIndex = 0;
2306     nodalConnect = _ptrMesh->getConnectivity (MED_FULL_INTERLACE,
2307                                               MED_NODAL,
2308                                               entity,
2309                                               MED_ALL_ELEMENTS);
2310     nodalConnectIndex = _ptrMesh->getConnectivityIndex (MED_NODAL,entity);
2311
2312     for ( int iType = 1; iType <= nbTypes; ++iType )
2313     {
2314       int nbElements = nbIndex[ iType ] - nbIndex[ iType - 1 ];
2315       medGeometryElement geomType = types[ iType - 1 ];
2316       if ( nbSuppElemsByType[ geomType ] >= nbElements )
2317         continue; // all elements are written with _supports
2318
2319       int elemId1 = nbIndex[ iType - 1 ];
2320       data.addTypeData( geomType, nbElements, 0, elemId1 );
2321
2322       writeElements (geomType,
2323                      data._types[ geomType ],
2324                      nodalConnect,
2325                      nodalConnectIndex);
2326     }
2327   }
2328
2329   // D颵t de la pile 32 (celle des points)
2330
2331   int nbNodes = _ptrMesh->getNumberOfNodes();
2332   _gibi << " ENREGISTREMENT DE TYPE   2" << endl;
2333   _gibi << " PILE NUMERO  32NBRE OBJETS NOMMES       0" <<
2334     "NBRE OBJETS" << setw(8) << nbNodes << endl;
2335   // Liste des noms de points
2336   // * 8001       FORMAT(8(1X,A8))
2337   // No named nodes
2338   // suit le nombre de noeuds
2339   _gibi << setw(8) << nbNodes << endl;
2340   // Le tableau suivant donne le filtre pour avoir le vrai num? des noeuds
2341   // appartenant aux 鬩ments d飲its. Par exemple, si un 鬩ment, d飲it
2342   // dans la pile 1, fait r馩rence ?n num? de noeud 駡l ? il faut le
2343   // mettre 駡l ?2
2344   // * 8000 FORMAT(10I8)
2345   for ( i = 0; i < nbNodes; ++i, fcount++ )
2346     _gibi << setw(8) << i + 1;
2347   fcount.stop();
2348
2349   // D颵t de pile 33 (celle des configurations (coordonn?))
2350   _gibi << " ENREGISTREMENT DE TYPE   2" << endl;
2351   _gibi << " PILE NUMERO  33NBRE OBJETS NOMMES       0NBRE OBJETS       1" << endl;
2352   // Suit le nombre de points dont on donne les coordonn?
2353   int nbValues = nbNodes * ( dim + 1 );
2354   _gibi << setw(8) << nbValues << endl;
2355   // Les coordonn? sont donn? par noeuds. D'abord le premier puis le deuxi譥...
2356   // Pour chaque noeuds, on donne les 2 ou 3 coordonn? plus la densit頣ourante
2357   // au moment de sa cr顴ion.
2358   // * 8003   FORMAT(1P,3E22.14)
2359   _gibi.precision(14);
2360   _gibi.setf( ios_base::scientific, ios_base::floatfield );
2361   _gibi.setf( ios_base::uppercase );
2362   const double * coords = _ptrMesh->getCoordinates(MED_FULL_INTERLACE);
2363   int j = 0;
2364   const double precision = 1.e-99; // PAL12077
2365   for ( fcount.init(3),i = 0; i < nbNodes; ++i, j += dim )
2366   {
2367     for ( int iCoord = 0; iCoord < dim; ++iCoord, fcount++ ) {
2368       double coo = coords[ j + iCoord ];
2369       bool  zero = ( -precision < coo && coo < precision );
2370       _gibi << setw(22) << ( zero ? 0.0 : coo );
2371     }
2372     _gibi << setw(22) << 0.0; // densite
2373     fcount++;
2374   }
2375   fcount.stop();
2376
2377   END_OF_MED(LOC);
2378 }
2379
2380 //=======================================================================
2381 //function : writeLastRecord
2382 //purpose  :
2383 //=======================================================================
2384
2385 void GIBI_MESH_WRONLY_DRIVER::writeLastRecord()
2386 {
2387   _gibi << " ENREGISTREMENT DE TYPE   5" << endl;
2388   _gibi << "LABEL AUTOMATIQUE :   1" << endl;
2389 }
2390
2391 /*--------------------- RDWR PART -------------------------------*/
2392
2393 GIBI_MESH_RDWR_DRIVER::GIBI_MESH_RDWR_DRIVER():GIBI_MESH_DRIVER()
2394 {
2395 }
2396 GIBI_MESH_RDWR_DRIVER::GIBI_MESH_RDWR_DRIVER(const string & fileName,
2397                                              MESH * ptrMesh):
2398        GIBI_MESH_DRIVER(fileName,ptrMesh,RDWR)
2399 {
2400   MESSAGE_MED("GIBI_MESH_RDWR_DRIVER::GIBI_MESH_RDWR_DRIVER(const string & fileName, MESH * ptrMesh) has been created");
2401 }
2402 GIBI_MESH_RDWR_DRIVER::GIBI_MESH_RDWR_DRIVER(const GIBI_MESH_RDWR_DRIVER & driver):
2403   GIBI_MESH_RDONLY_DRIVER::GIBI_MESH_DRIVER(driver)
2404 {
2405   MESSAGE_MED("GIBI_MESH_RDWR_DRIVER::GIBI_MESH_RDWR_DRIVER(driver) has been created");
2406 }
2407 GIBI_MESH_RDWR_DRIVER::~GIBI_MESH_RDWR_DRIVER() {
2408   MESSAGE_MED("GIBI_MESH_RDWR_DRIVER::GIBI_MESH_RDWR_DRIVER(const string & fileName, MESH * ptrMesh) has been destroyed");
2409 }
2410 GENDRIVER * GIBI_MESH_RDWR_DRIVER::copy(void) const
2411 {
2412   const char* LOC = "GIBI_MESH_RDWR_DRIVER::copy()";
2413   BEGIN_OF_MED(LOC);
2414   GENDRIVER * driver = new GIBI_MESH_RDWR_DRIVER(*this);
2415   END_OF_MED(LOC);
2416   return driver;
2417 }
2418 void GIBI_MESH_RDWR_DRIVER::write(void) const
2419   throw (MEDEXCEPTION)
2420 {
2421   GIBI_MESH_RDWR_DRIVER * me = const_cast<GIBI_MESH_RDWR_DRIVER *>(this);
2422   me->GIBI_MESH_WRONLY_DRIVER::open();
2423   me->GIBI_MESH_WRONLY_DRIVER::write();
2424   me->GIBI_MESH_WRONLY_DRIVER::close();
2425 }
2426 void GIBI_MESH_RDWR_DRIVER::read (void)
2427   throw (MEDEXCEPTION)
2428 {
2429   const char* LOC = "GIBI_MESH_RDWR_DRIVER::read()";
2430   BEGIN_OF_MED(LOC);
2431   GIBI_MESH_RDONLY_DRIVER::open();
2432   GIBI_MESH_RDONLY_DRIVER::read();
2433   GIBI_MESH_RDONLY_DRIVER::close();
2434   END_OF_MED(LOC);
2435 }
2436 void GIBI_MESH_RDWR_DRIVER::open()
2437   // throw (MEDEXCEPTION)
2438 {
2439 }
2440 void GIBI_MESH_RDWR_DRIVER::close()
2441   // throw (MEDEXCEPTION)
2442 {
2443 }
2444
2445 //============================== ====================================================
2446 //============================== FIELD Reading Driver ==============================
2447 //============================== ====================================================
2448
2449 GIBI_MED_RDONLY_DRIVER::GIBI_MED_RDONLY_DRIVER():GIBI_MESH_RDONLY_DRIVER()
2450 {
2451 }
2452 GIBI_MED_RDONLY_DRIVER::GIBI_MED_RDONLY_DRIVER(const string & fileName, MED * ptrMed):
2453        GIBI_MESH_RDONLY_DRIVER(fileName,NULL), _med( ptrMed )
2454 {
2455   MESSAGE_MED("GIBI_MED_RDONLY_DRIVER(const string & fileName, MED * ptrMed) has been created");
2456   _fileName = fileName;
2457   _accessMode = RDONLY;
2458 }
2459 GIBI_MED_RDONLY_DRIVER::GIBI_MED_RDONLY_DRIVER(const GIBI_MED_RDONLY_DRIVER & driver)
2460 {
2461 }
2462  GIBI_MED_RDONLY_DRIVER::~GIBI_MED_RDONLY_DRIVER()
2463 {
2464 }
2465 GENDRIVER * GIBI_MED_RDONLY_DRIVER::copy ( void ) const
2466 {
2467   return new GIBI_MED_RDONLY_DRIVER(*this);
2468 }
2469
2470 //=======================================================================
2471 //function : read
2472 //purpose  :
2473 //=======================================================================
2474
2475 void GIBI_MED_RDONLY_DRIVER::read ( void ) throw (MEDEXCEPTION)
2476 {
2477   const char * LOC = "GIBI_MED_RDONLY_DRIVER::read() : " ;
2478   BEGIN_OF_MED(LOC);
2479
2480   if (_status!=MED_OPENED)
2481     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "file " << _fileName<<" is not opened." ));
2482
2483   _ptrMesh = new MESH();
2484
2485   _intermediateMED medi;
2486   try {
2487     if ( !readFile( &medi, true ) )
2488       return;
2489
2490     //MESSAGE_MED(LOC <<  medi );
2491     fillMesh( &medi );
2492     MESSAGE_MED(LOC << "GIBI_MED_RDONLY_DRIVER::read : RESULTATS STRUCTURE INTERMEDIAIRES : ");
2493     MESSAGE_MED(LOC <<  medi );
2494
2495     list< FIELD_* > fields;
2496     medi.getFields( fields );
2497     MESSAGE_MED( "nb fields: " << fields.size() );
2498
2499     if ( _ptrMesh->getName().empty() )
2500       _ptrMesh->setName( "MESH" );
2501
2502     _med->addMesh( _ptrMesh );
2503
2504     list< FIELD_* >::iterator it = fields.begin();
2505     for ( ; it != fields.end(); it++ ) {
2506       int nbComponents = (*it)->getNumberOfComponents();
2507       if(nbComponents>0) {
2508         UNIT* compoUnits = new UNIT[nbComponents];
2509         string* MEDcompoUnits = new string[nbComponents];
2510         for(int l = 0; l<nbComponents; l++) {
2511           compoUnits[l] = UNIT("", "");
2512           MEDcompoUnits[l] = "";
2513         }
2514         (*it)->setComponentsUnits(compoUnits);
2515         (*it)->setMEDComponentsUnits(MEDcompoUnits);
2516       }
2517       _med->addField( *it );
2518     }
2519   }
2520   catch (MEDEXCEPTION &ex)
2521   {
2522     INFOS_MED( ex.what() );
2523   }
2524
2525   END_OF_MED(LOC);
2526 }
2527
2528 //============================== ====================================================
2529 //============================== FIELD Writting Driver ==============================
2530 //============================== ====================================================
2531
2532 GIBI_MED_WRONLY_DRIVER::GIBI_MED_WRONLY_DRIVER():GIBI_MESH_WRONLY_DRIVER()
2533 {
2534 }
2535 GIBI_MED_WRONLY_DRIVER::GIBI_MED_WRONLY_DRIVER(const string & fileName,
2536                                                    MED * ptrMed,
2537                                                    MESH * ptrMesh)
2538      :GIBI_MESH_WRONLY_DRIVER(fileName,ptrMesh), _med( ptrMed )
2539 {
2540   const char * LOC =
2541     "GIBI_MED_WRONLY_DRIVER(const string & fileName, MED * ptrMed, MESH * ptrMesh)" ;
2542   BEGIN_OF_MED(LOC);
2543
2544   _fileName = fileName;
2545   _accessMode = WRONLY;
2546   _ptrMesh = ptrMesh;
2547   if ( !_med || !_ptrMesh )
2548     throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << " Bad params " << ptrMed << " " << ptrMesh ));
2549 }
2550 GIBI_MED_WRONLY_DRIVER::GIBI_MED_WRONLY_DRIVER(const GIBI_MED_WRONLY_DRIVER & driver)
2551 {
2552 }
2553 GIBI_MED_WRONLY_DRIVER::~GIBI_MED_WRONLY_DRIVER()
2554 {
2555 }
2556 GENDRIVER * GIBI_MED_WRONLY_DRIVER::copy ( void ) const
2557 {
2558   return new GIBI_MED_WRONLY_DRIVER(*this);
2559 }
2560
2561 //=======================================================================
2562 //function : writeDataSection
2563 //purpose  :
2564 //=======================================================================
2565
2566 template< class T, class INTERLACING_TAG>
2567 static void writeDataSection (fstream&                    file,
2568                               FIELD<T, INTERLACING_TAG> * field,
2569                               const int                   id1,
2570                               const int                   id2) throw (MEDEXCEPTION)
2571 {
2572   const char * LOC="writeDataSection (.....) :";
2573   BEGIN_OF_MED(LOC);
2574
2575   int nbGauss, nbComp = field->getNumberOfComponents();
2576
2577   typedef typename MEDMEM_ArrayInterface<T,INTERLACING_TAG,NoGauss>::Array ArrayNoGauss;
2578   typedef typename MEDMEM_ArrayInterface<T,INTERLACING_TAG,Gauss>::Array   ArrayGauss;
2579
2580   MEDMEM_Array_ * array        = field->getArray();
2581   ArrayNoGauss  * arrayNoGauss = 0;
2582   ArrayGauss    * arrayGauss   = 0;
2583
2584   if ( !array )
2585     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" Field |"<< field->getName()
2586                                  << "| not allocated"));
2587
2588   if ( array->getGaussPresence() ) {
2589     arrayGauss = field->getArrayGauss();
2590     nbGauss    = arrayGauss->getNbGauss( id1 );
2591   }
2592   else {
2593     arrayNoGauss = field->getArrayNoGauss();
2594     nbGauss      = 1;
2595   }
2596
2597   TFieldCounter fcount( file, 3 ); // 3 values on a line
2598
2599 #ifdef CASTEM_FULL_INTERLACE
2600   const int gauss_step = field->getInterlacingType() == MED_EN::MED_FULL_INTERLACE ? nbComp : 1;
2601 #endif
2602
2603   for ( int iComp = 1; iComp <= nbComp; ++iComp )
2604   {
2605     file << setw(8) << nbGauss       // nb scalar values by element
2606          << setw(8) << ( id2 - id1 ) // total nb of scalar values
2607          << setw(8) << 0
2608          << setw(8) << 0
2609          << endl;
2610     // * 8003   FORMAT(1P,3E22.14)
2611     if ( arrayNoGauss ) {
2612       for (int id = id1; id < id2; id++, fcount++ )
2613         file << setw(22) << arrayNoGauss->getIJ( id, iComp );
2614     }
2615     else {
2616 #ifdef CASTEM_FULL_INTERLACE
2617       for (int id = id1; id < id2; id++ ) {
2618         const T* val = & arrayGauss->getIJK( id, iComp, 1 );
2619         const T* valEnd = val + nbGauss * gauss_step;
2620         for ( ; val < valEnd; val += gauss_step, fcount++ )
2621           file << setw(22) << *val;
2622       }
2623 #else
2624       for ( int iGauss = 1; iGauss <= nbGauss; ++iGauss ) {
2625         for (int id = id1; id < id2; id++, fcount++  )
2626           file << setw(22) << arrayGauss->getIJK( id, iComp, iGauss );
2627       }
2628 #endif
2629     }
2630     fcount.stop();
2631   }
2632   END_OF_MED(LOC);
2633 }
2634
2635 //=======================================================================
2636 //function : write
2637 //purpose  :
2638 //=======================================================================
2639
2640 void GIBI_MED_WRONLY_DRIVER::write( void ) const throw (MEDEXCEPTION)
2641 {
2642   const char* LOC = "void GIBI_MED_WRONLY_DRIVER::write(void) const : ";
2643   BEGIN_OF_MED(LOC);
2644
2645   // we are going to modify the _gibi field
2646   GIBI_MED_WRONLY_DRIVER * me = const_cast<GIBI_MED_WRONLY_DRIVER *>(this);
2647
2648   // get all fields on _ptrMesh and add their support to be written
2649   list<FIELD_*> fields;
2650   int iField, nbFields = _med->getNumberOfFields();
2651   int nb_obj = 0;
2652   list<int> nb_sub_list;
2653   map<string,int> nameNbMap;
2654
2655   list<pair<int,int> >           subIdSizeList; // pair( <submesh id>, <submesh size> );
2656   list<pair<int,int> >::iterator idsize;
2657
2658   string *names=new string[ nbFields ];
2659   _med->getFieldNames( names );
2660   for ( iField = 0; iField < nbFields; ++iField )
2661   {
2662     int nb_sub = 0;
2663     deque<DT_IT_> dtit = _med->getFieldIteration( names[ iField ]);
2664     deque<DT_IT_>::iterator fIt = dtit.begin();
2665     for ( ; fIt != dtit.end(); fIt++ )
2666     {
2667       FIELD_ * f = _med->getField( names[ iField ], fIt->dt, fIt->it );
2668       if ( f->getValueType() != MED_EN::MED_REEL64 )
2669       {
2670         MESSAGE_MED("GIBI_MED_WRONLY_DRIVER::write( FIELD< int > ) not implemented");
2671         continue;
2672       }
2673       const SUPPORT * sup = f->getSupport();
2674       if ( me->addSupport( sup ) ) {
2675         fields.push_back( f );
2676         nb_sub += getSubMeshIdAndSize( sup, subIdSizeList );
2677       }
2678     }
2679     if ( nb_sub ) {
2680       addName( nameNbMap, names[ iField ], ++nb_obj, "F" );
2681       nb_sub_list.push_back( nb_sub );
2682     }
2683   }
2684
2685   // write mesh
2686
2687   //try {
2688     me->writeSupportsAndMesh();
2689 //   }
2690 //   catch (MEDEXCEPTION &ex)
2691 //   {
2692 //     INFOS_MED( ex.what() );
2693 //     return;
2694 //   }
2695
2696   // write fields
2697
2698   if ( !fields.empty() ) {
2699
2700     fstream & gibi = me->_gibi;
2701
2702     TFieldCounter fcount( gibi, 10 );
2703
2704     gibi << " ENREGISTREMENT DE TYPE   2" << endl;
2705     gibi << " PILE NUMERO  39NBRE OBJETS NOMMES" << setw(8) << nameNbMap.size()
2706       << "NBRE OBJETS" << setw(8) << nb_obj << endl;
2707
2708     me->writeNames( nameNbMap );
2709
2710     list<FIELD_*>::iterator itF = fields.begin();
2711     list<int>::iterator itNbSub = nb_sub_list.begin();
2712     int nb_sub = 0, cur_nb_sub = 0;
2713     for ( ; itF != fields.end(); itF++ )
2714     {
2715       if ( cur_nb_sub == nb_sub && itNbSub != nb_sub_list.end() ) {
2716         // start the next field writting
2717         nb_sub = *(itNbSub++);
2718         string description = (*itF)->getDescription();
2719         if ( description.size() > 72 )
2720           description = description.substr(0,72);
2721         gibi << setw(8) << nb_sub
2722              << setw(8) << -1
2723              << setw(8) << 6
2724              << setw(8) << 72 /*description.size()*/ << endl; // PAL19100
2725         if ( !description.empty() )
2726           gibi << setw(72) << description << endl;
2727         else
2728           gibi << setw(72) << "Field" << endl;
2729         gibi << setw(72) << " " << endl;
2730
2731         // Sub Components section
2732         list<FIELD_*>::iterator itF2 = itF;
2733         vector<int> vals( 9, 0 );
2734         vals[8] = 2;
2735         fcount.init(10);
2736         cur_nb_sub = 0;
2737         while ( itF2 != fields.end() && cur_nb_sub < nb_sub )
2738         {
2739           FIELD_* f = *itF2++;
2740           vals[2] = f->getNumberOfComponents();
2741           getSubMeshIdAndSize( f->getSupport(), subIdSizeList );
2742           for ( idsize = subIdSizeList.begin(); idsize != subIdSizeList.end(); idsize++ )
2743           {
2744             ++cur_nb_sub;
2745             vals[0] = -idsize->first; // support id
2746             for ( int i = 0; i < (int)vals.size(); ++i, fcount++ )
2747               gibi << setw(8) << vals[ i ];
2748           }
2749         }
2750         fcount.stop();
2751         cur_nb_sub = 0;
2752         // dummy strings
2753         int i_sub;
2754         for ( fcount.init(4), i_sub = 0; i_sub < nb_sub; ++i_sub, fcount++ )
2755           gibi << "                  ";
2756         fcount.stop();
2757         for ( fcount.init(8), i_sub = 0; i_sub < nb_sub; ++i_sub, fcount++ )
2758           gibi << "         ";
2759         fcount.stop();
2760       }
2761
2762       FIELD_* f = *itF;
2763       int id1 = 1;
2764       int iComp, nbComp = f->getNumberOfComponents();
2765       // loop on sub-components
2766       getSubMeshIdAndSize( f->getSupport(), subIdSizeList );
2767       for ( idsize = subIdSizeList.begin(); idsize != subIdSizeList.end(); idsize++ )
2768       {
2769         cur_nb_sub++;
2770         // component addresses
2771         for ( fcount.init(10), iComp = 0; iComp < nbComp; ++iComp, fcount++ )
2772           gibi << setw(8) << 777; // a good number
2773         fcount.stop();
2774         // component names
2775         gibi << left;
2776         for ( fcount.init(8), iComp = 0; iComp < nbComp; ++iComp, fcount++ ) {
2777           string compName = f->getComponentName( iComp + 1 );
2778           if ( compName.size() > 8 )
2779             compName = compName.substr(0, 8);
2780           gibi << " "  << setw(8) << compName;
2781         }
2782         fcount.stop();
2783         // component types
2784         for ( fcount.init(4), iComp = 0; iComp < nbComp; ++iComp, fcount++ )
2785           gibi << " "  << setw(17) << "REAL*8";
2786         fcount.stop();
2787         gibi << right;
2788
2789         // Data section
2790         int id2 = id1 + idsize->second;
2791
2792         // PAL11040
2793         //if  (f->getGaussPresence() )
2794         //throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "GibiDriver don't support Field with Gauss point"));
2795
2796         if ( f->getInterlacingType() == MED_NO_INTERLACE )
2797           writeDataSection( gibi, static_cast<FIELD<double,NoInterlace  > * >(f), id1, id2 );
2798         else if ( f->getInterlacingType() == MED_FULL_INTERLACE )
2799           writeDataSection( gibi, static_cast<FIELD<double,FullInterlace> * >(f), id1, id2 );
2800         else
2801           writeDataSection( gibi, static_cast<FIELD<double,NoInterlaceByType> * >(f), id1, id2 );
2802
2803         id1 = id2;
2804       }
2805     } // loop on fields
2806   }
2807   me->writeLastRecord();
2808   delete [] names;
2809   END_OF_MED(LOC);
2810 }