src/MEDMEM/MEDMEM_EnsightUtils.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 // File      : MEDMEM_EnsightUtils.cxx
  21 // Created   : Fri Jun  6 16:08:32 2008
  22 // Author    : Edward AGAPOV (eap)
  23 //
  24 #include "MEDMEM_Field.hxx" // important to include it before MEDMEM_EnsightUtils.hxx,
  25 // in order not to redefine DBL_MIN
  26
  27 #include "MEDMEM_EnsightUtils.hxx"
  28
  29 #include "MEDMEM_EnsightMeshDriver.hxx"
  30 #include "MEDMEM_EnsightFieldDriver.hxx"
  31 #include "MEDMEM_DriverTools.hxx"
  32
  33 #ifdef WIN32
  34 #include <io.h>
  35 #else
  36 #include <unistd.h>
  37 #endif
  38 #include <fcntl.h>
  39 //#include <stdio.h>
  40
  41 using namespace std;
  42 using namespace MEDMEM;
  43
  44   // TODO:
  45 // -2) in transient model, get time values for meshes form fields on them
  46 // -1) why after reading /data/eap/ENSIGHT/SAMPLES_SRC/EnSight/crash.case only
  47 //     geo EnSIght files are written
  48 // -0) test selecting imaginary part to read
  49 // 0) appending fields to write with same names in Case file
  50 // 1) write grid
  51 // 2) if mesh description is longer that maxLen, distribute it among two lines
  52 // 3) care of mesh description at reading
  53 // 4) at writing, cut off " t=..." from mesh name and set to meshDriver->_meshName
  54 // 5) MEDMEM_EXPORT
  55 // 4.5) do something with zzzz121b.med pb
  56 // 6) stripe in getLine()
  57 // 7) care of elem numbers that were read
  58 // 8) read real and imaginary of complex field by MED driver
  59 // 9) read "mesured" and care of variable on it
  60 // 10) clear theInterMedMap and delete _InterMed's there
  61 // 11) not to fill MESH made in getMeshData()
  62 // 12) care of changing nb of '*' in names
  63 // 14) analize _TimeSet to write in shortest form
  64 // 15) add LOCALIZED where no method name specified
  65 // 17) add checks of input mesh and field like in gibi
  66 // 18) check writing interlace in all modes
  67 // 19) define creation by MEDMEM by case file and not by geo
  68 // 21) there maybe '\n' in binary file after string end
  69 // 22) optimize _ASCIIFileReader::eof() at file end
  70 // 23) What if to create CELL entity with 1-3 dim and a CELL group?
  71 // 24) not to write CELL connectivity if a group on all cells is present
  72 // 25) compact imed (sortedNodeIDs etc.)
  73 // 26) no need in converting EnSight variables to eventual types in convertReals
  74 //     since a result array is not used directly but occupies more memory
  75 // 27) try to exclude merged nodes in DTree and compare time on ChampsDarcy.med
  76
  77
  78 namespace MEDMEM
  79 {
  80   // ==============================================================================
  81   /*!
  82     Writing format.
  83    */
  84   // ==============================================================================
  85
  86   EnSightFormat theEnSightFormatForWriting = ENSIGHT_6;
  87   bool          theBinaryFormatForWriting  = false;
  88   bool          theIgnoreIncompatibility   = false;
  89
  90   void          setEnSightFormatForWriting (EnSightFormat format, bool isBinary)
  91   {
  92     theEnSightFormatForWriting = format;
  93     theBinaryFormatForWriting  = isBinary;
  94   }
  95   EnSightFormat getEnSightFormatForWriting()
  96   {
  97     return theEnSightFormatForWriting;
  98   }
  99   bool          isBinaryEnSightFormatForWriting()
 100   {
 101     return theBinaryFormatForWriting;
 102   }
 103   // ==============================================================================
 104   /*!
 105     To raise or not if MEDMEM-EnSight incompatibility encounters or suspected
 106   */
 107   // ==============================================================================
 108
 109   void          setIgnoreIncompatibility(bool toIgnore)
 110   {
 111     theIgnoreIncompatibility = toIgnore;
 112   }
 113 }
 114
 115 #ifdef WIN32
 116 #define FILE_SEPARATOR '\\'
 117 #else
 118 #define FILE_SEPARATOR '/'
 119 #endif
 120
 121 #define _ATOI( str ) atoi((str).c_str())
 122 #define _ATOF( str ) atof((str).c_str())
 123
 124 #define BUFFER_SIZE 16184  // for non-stream input
 125
 126 namespace MEDMEM_ENSIGHT
 127 {
 128   //--------------------------------------------------------------------------------
 129   /*!
 130    * \brief Registry of drivers creating _CaseFileDriver
 131    *
 132    * The problem is that if the mesh/field driver is created by the user, it is
 133    * to write Case file, but if it is created by MED driver, it should not, in
 134    * order not to overwrite the Case file written by MED driver.
 135    * To assure this feature, we register all mesh/field drivers added to a case
 136    * file by a MED driver and then ignore thier calls to _CaseFileDriver.
 137    */
 138   set< const _CaseFileDriver_User* > theCaseUsers;
 139
 140   //--------------------------------------------------------------------------------
 141   /*!
 142    * \brief Storage of data read by mesh driver to be used by field drivers
 143    *
 144    * The field driver needs info read by the mesh driver (nb of elements by type in
 145    * each part and what place the element gets in MEDMED support).
 146    * So we store needed data unless it is no more needed (i.e. no driver is alive).
 147    * The map key is "<case_file>:<mesh index>".
 148    */
 149   map< string, _InterMed* > theInterMedMap;
 150
 151   //--------------------------------------------------------------------------------
 152   /*!
 153    * \brief Add a driver to the registry and return true if was already in
 154    */
 155   bool isToIgnore(const _CaseFileDriver_User* driver);
 156   bool isToIgnore(const _CaseFileDriver_User* driver)
 157   {
 158     return ! theCaseUsers.insert( driver ).second;
 159   }
 160
 161   //--------------------------------------------------------------------------------
 162   /*!
 163    * \brief Remove a driver from the registry
 164    */
 165   void unregister(const _CaseFileDriver_User* driver);
 166   void unregister(const _CaseFileDriver_User* driver)
 167   {
 168     theCaseUsers.erase( driver );
 169     if ( theCaseUsers.empty() )
 170     {
 171       for ( map< string, _InterMed* >::iterator str_imed = theInterMedMap.begin();
 172             str_imed != theInterMedMap.end();
 173             str_imed++ )
 174         delete str_imed->second;
 175       theInterMedMap.clear();
 176     }
 177   }
 178
 179   //--------------------------------------------------------------------------------
 180   /*!
 181    * \brief Return mesh data needed by field driver
 182    */
 183   _InterMed* getMeshData( const string& key );
 184   _InterMed* getMeshData( const string& key )
 185   {
 186     // find existing data
 187     map< string, _InterMed* >::iterator kimed = theInterMedMap.find( key );
 188     if ( kimed != theInterMedMap.end() )
 189       return kimed->second;
 190
 191     // create a new data
 192     MESH * mesh = new MESH;
 193     string caseFile, meshIndex;
 194     _ASCIIFileReader::split( key, caseFile, meshIndex, ':' );
 195     ENSIGHT_MESH_RDONLY_DRIVER meshDrv(caseFile, mesh, _ATOI( meshIndex ));
 196     meshDrv.read();
 197     kimed = theInterMedMap.find( key );
 198     if ( kimed == theInterMedMap.end() )
 199       return 0;
 200     _InterMed* imed     = kimed->second;
 201     imed->_medMesh      = mesh;
 202     imed->_isOwnMedMesh = true;
 203     return imed;
 204   }
 205   // ---------------------------------------------------------------
 206   /*!
 207    * \brief Prepend "EnSight-MEDMEM compatibility problem" to the text
 208    *        of an exception
 209    */
 210   STRING compatibilityPb(const string& exceptionText)
 211   {
 212     return STRING("EnSight-MEDMEM compatibility problem:\n") << exceptionText;
 213   }
 214   // ---------------------------------------------------------------
 215   /*!
 216    * \brief To ignore incompatibility or not
 217    */
 218   bool   toIgnoreIncompatibility()
 219   {
 220     return theIgnoreIncompatibility;
 221   }
 222
 223   //--------------------------------------------------------------------------------
 224   /*!
 225    * \brief Return true if set index is empty or corresponds to an existing set
 226    */
 227   template <class TSet> bool isValidIndex(const string& index, const map<int,TSet>& aMap)
 228   {
 229     if ( index.empty() ) return true;
 230     return ( aMap.find( _ATOI( index ) ) != aMap.end() );
 231   }
 232
 233   //--------------------------------------------------------------------------------
 234   /*!
 235    * \brief Return EnSight type corresponding to med one
 236    */
 237   const TEnSightElemType& getEnSightType(medGeometryElement medType)
 238   {
 239     static TEnSightElemType theEnSightType;
 240
 241     int nbNodes = medType % 100;
 242     theEnSightType._medType = medType;
 243
 244     switch ( medType ) {
 245     case MED_NONE    :
 246     case MED_POINT1  : {
 247       theEnSightType._name = "point";
 248       theEnSightType._medIndex.resize(1,0);
 249       break ;
 250     }
 251     case MED_SEG2    : {
 252       theEnSightType._name = "bar2";
 253       int conn [2] = {0,1};
 254       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 255       break ;
 256     }
 257     case MED_SEG3    : {
 258       theEnSightType._name = "bar3";
 259       int conn [3] = {0,2,1};
 260       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 261       break ;
 262     }
 263     case MED_TRIA3   : {
 264       theEnSightType._name = "tria3";
 265       int conn [3] = {0,2,1};
 266       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 267       break ;
 268     }
 269     case MED_QUAD4   : {
 270       theEnSightType._name = "quad4";
 271       int conn [4] = {0,3,2,1};
 272       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 273       break ;
 274     }
 275     case MED_TRIA6   : {
 276       theEnSightType._name = "tria6";
 277       int conn [6] = {0,2,1,5,4,3};
 278       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 279       break ;
 280     }
 281     case MED_QUAD8   : {
 282       theEnSightType._name = "quad8";
 283       int conn [8] = {0,3,2,1,7,6,5,4};
 284       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 285       break ;
 286     }
 287     case MED_TETRA4  : {
 288       theEnSightType._name = "tetra4";
 289       int conn [4] = {0,1,3,2};
 290       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 291       break ;
 292     }
 293     case MED_PYRA5   : {
 294       theEnSightType._name = "pyramid5";
 295       int conn [5] = {0,3,2,1,4};
 296       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 297       break ;
 298     }
 299     case MED_PENTA6  : {
 300       theEnSightType._name = "penta6";
 301       int conn [6] = {0,2,1,3,5,4};
 302       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 303       break ;
 304     }
 305     case MED_HEXA8   : {
 306       theEnSightType._name = "hexa8";
 307       int conn [8] = {0,3,2,1,4,7,6,5};
 308       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 309       break ;
 310     }
 311     case MED_TETRA10 : {
 312       theEnSightType._name = "tetra10";
 313       int conn [10] = {0,2,1,3,6,5,4,7,9,8};
 314       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 315       break ;
 316     }
 317     case MED_PYRA13  : {
 318       theEnSightType._name = "pyramid13";
 319       int conn [13] = {0,3,2,1,4,8,7,6,5,9,12,11,10};
 320       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 321       break ;
 322     }
 323     case MED_PENTA15 : {
 324       theEnSightType._name = "penta15";
 325       int conn [15] = {0,2,1,3,5,4,8,7,6,11,10,12,14,13};
 326       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 327       break ;
 328     }
 329     case MED_HEXA20  : {
 330       theEnSightType._name = "hexa20";
 331       int conn [20] = {0,3,2,1,4,7,6,5,11,10,9,8,15,14,13,12,16,19,18,17};
 332       theEnSightType._medIndex.assign( conn, conn+nbNodes );
 333       break ;
 334     }
 335     case MED_POLYGON : {
 336       theEnSightType._name = "nsided";
 337       theEnSightType._medIndex.clear();
 338       break ;
 339     }
 340     case MED_POLYHEDRA : {
 341       theEnSightType._name = "nfaced";
 342       theEnSightType._medIndex.clear();
 343       break ;
 344     }
 345     default: {
 346       theEnSightType._name = "";
 347       theEnSightType._medIndex.clear();
 348     }
 349     }
 350     return theEnSightType;
 351   }
 352
 353   //--------------------------------------------------------------------------------
 354   /*!
 355    * \brief Return EnSight type having a given name
 356    */
 357   const TEnSightElemType& getEnSightType(const string& theTypeName)
 358   {
 359     string typeName = theTypeName;
 360     if ( isGhostType( typeName ))
 361       typeName = string( &typeName[2] ); // ghost type
 362
 363     static map<string, TEnSightElemType> name2Type;
 364
 365     map<string, TEnSightElemType>::iterator nameType = name2Type.find( typeName );
 366     if ( nameType != name2Type.end() )
 367       return nameType->second;
 368
 369     const list<medGeometryElement> & allMedTypes = MED_EN::meshEntities[MED_CELL];
 370     list< medGeometryElement >::const_iterator medType = allMedTypes.begin();
 371     for ( ; medType != allMedTypes.end(); ++medType )
 372     {
 373       const TEnSightElemType& enSightType = getEnSightType( *medType );
 374       if ( enSightType._name == typeName )
 375         return name2Type[ typeName ] = enSightType;
 376     }
 377     return getEnSightType(MED_ALL_ELEMENTS);
 378   }
 379
 380 // =============================================================== _CaseFileDriver
 381 /*!
 382  * \brief Case file driver constructor
 383  */
 384 //================================================================================
 385
 386 _CaseFileDriver::_CaseFileDriver(const string &              fileName,
 387                                  const _CaseFileDriver_User* creator)
 388   : _fileName( fileName ), _directory("."), _user( creator )
 389 {
 390   // Find out if the driver is blocked
 391   _blocked = isToIgnore( creator );
 392   if ( creator->getAccessMode() == MED_EN::RDONLY )
 393     _blocked = false;
 394   if ( !_blocked ) {
 395     // set directory
 396     string::size_type sepPos = _fileName.rfind( FILE_SEPARATOR );
 397     if ( sepPos != _fileName.npos ) {
 398       _directory = _fileName.substr( 0, sepPos );
 399     }
 400     _format = getEnSightFormatForWriting();
 401   }
 402 }
 403
 404 // =============================================================== _CaseFileDriver
 405 /*!
 406  * \brief Case file driver destructor
 407  */
 408 //================================================================================
 409
 410 _CaseFileDriver::~_CaseFileDriver()
 411 {
 412   if ( !_blocked )
 413   {
 414     // to write case file again by same DRIVER
 415     unregister( _user );
 416     ((_CaseFileDriver_User*)_user)->_imed = 0;
 417   }
 418 }
 419
 420 #define READ_NEXT_LINE continue
 421 #define RAISE_EXCEPTION break
 422
 423 //================================================================================
 424 /*!
 425  * \brief Read Case file
 426  */
 427 //================================================================================
 428
 429 void _CaseFileDriver::read() throw (MEDEXCEPTION)
 430 {
 431   if ( _blocked )
 432     return;
 433
 434   _ASCIIFileReader reader( _fileName );
 435
 436   STRING badFile("Invalid Case file ");
 437   badFile << _fileName << "\n";
 438
 439   list<_FileSet> fileSets;
 440   list<_TimeSet> timeSets;
 441
 442   set<string> varNames; // to detect equal variable names
 443
 444   string section = "_";
 445   int lineNb = 0;
 446   while ( !reader.eof() )
 447   {
 448     string line = reader.getLine();
 449     ++lineNb;
 450
 451     // cut off comments
 452     string comment;
 453     reader.split( line, line, comment, '#');
 454     if ( line.empty() )
 455       READ_NEXT_LINE;
 456
 457     string key, value; // parts of a line splited by column
 458     reader.split( line, key, value, ':');
 459
 460     STRING errorMsg;
 461
 462     // analyse not empty lines
 463     switch ( section[0] ) {
 464     case 'V': {
 465       // --------------------------------------------------------------------------------
 466       // VARIABLE
 467
 468       string type, s, ts, fs;
 469       reader.split(key,type,s);
 470       int newVarIndex = _variables.empty() ? 1 : _variables.rbegin()->first + 1;
 471
 472       if ( type == "scalar" || type == "vector" || type == "tensor" )
 473       {
 474         // scalar per node:        [ts] [fs] description filename
 475         // vector per node:        [ts] [fs] description filename
 476         // tensor symm per node:       [ts]  [fs]   description filename
 477         // tensor asym per node:       [ts]  [fs]   description filename
 478         // scalar per element:         [ts]  [fs]   description filename
 479         // vector per element:         [ts]  [fs]   description filename
 480         // tensor symm per element:    [ts]  [fs]   description filename
 481         // tensor asym per element:    [ts]  [fs]   description filename
 482         // scalar per measured node:   [ts]  [fs]   description filename
 483         // vector per measured node:   [ts]  [fs]   description filename
 484         list<string> parts;
 485         int nbParts = reader.split( value, parts );
 486         if ( nbParts < 2 ) {
 487           errorMsg << "invalid variable format:\n" << line;
 488           RAISE_EXCEPTION;
 489         }
 490         if ( contains( "per measured node", s.c_str() )) {
 491           //cout << "Skip not supported data type: " << key << endl;
 492           READ_NEXT_LINE;
 493         }
 494         list<string>::reverse_iterator p = parts.rbegin();
 495         _Variable& var      = _variables[newVarIndex];
 496         var._type           = key;
 497         var._fileNameOrData = *p++;
 498         var._name           = *p++;
 499         if ( nbParts == 3 ) {
 500           var._timeSetNumber = *p;
 501         }
 502         else if ( nbParts == 4 ) {
 503           var._fileSetNumber = *p++;
 504           var._timeSetNumber = *p;
 505         }
 506         varNames.insert( var._name );
 507         READ_NEXT_LINE;
 508       }
 509       else if ( type == "constant" )
 510       {
 511         // constant per case:      [ts]      description const_value(s)
 512         // constant per case file: [ts]      description cvfilename
 513         reader.split(value,s,value);
 514         if ( reader.isDigit( s )) {
 515           ts = s;
 516           reader.split(value,s,value);
 517         }
 518         _Variable& var      = _variables[newVarIndex];
 519         var._type           = key;
 520         var._name           = s;
 521         var._fileNameOrData = value;
 522         var._timeSetNumber  = ts;
 523         if ( var._name.empty() || var._fileNameOrData.empty() ) {
 524           errorMsg << "invalid variable format:\n" << line;
 525           RAISE_EXCEPTION;
 526         }
 527         varNames.insert( var._name );
 528         READ_NEXT_LINE;
 529       }
 530       else if ( type == "complex" )
 531       {
 532         // complex scalar per node:    [ts]  [fs]   description Re_fn   Im_fn freq
 533         // complex vector per node:    [ts]  [fs]   description Re_fn   Im_fn freq
 534         // complex scalar per element:    [ts] [fs] description Re_fn   Im_fn freq
 535         // complex vector per element:    [ts] [fs] description Re_fn   Im_fn freq
 536         reader.split(value,s,value);
 537         if ( reader.isDigit( s )) {
 538           ts = s;
 539           reader.split(value,s,value);
 540           if ( reader.isDigit( s )) {
 541             fs = s;
 542             reader.split(value,s,value);
 543           }
 544         }
 545         list<string> parts;
 546         int nbParts = reader.split( value, parts );
 547         if ( nbParts < 3 ) {
 548           errorMsg << "invalid variable format:\n" << line;
 549           RAISE_EXCEPTION;
 550         }
 551         // a variable contains two fields. We leave one slot in _variables empty
 552         // in order to have last key in _variables equal to number of fields.
 553         // Variable index equal to missing slot key corresponds to real part (Re_fn)
 554         _Variable& var     = _variables[++newVarIndex];
 555         var._type          = key;
 556         var._name          = s;
 557         var._timeSetNumber = ts;
 558         var._fileSetNumber = fs;
 559         var._fileNameOrData = value;
 560         varNames.insert( var._name );
 561         READ_NEXT_LINE;
 562       }
 563       break;
 564     }
 565     case 'T': {
 566       // --------------------------------------------------------------------------------
 567       // TIME
 568       //    time set:              ts [description]
 569       //    number of steps:       ns
 570       //    filename start number: fs
 571       //    filename increment:    fi
 572       //    time values:           time_1 time_2 .... time_ns
 573       // or
 574       //    time set:              ts [description]
 575       //    number of steps:       ns
 576       //    filename numbers:      fn
 577       //    time values:           time_1 time_2 .... time_ns
 578       // or
 579       //    time set:              ts [description]
 580       //    number of steps:       ns
 581       //    filename numbers file: fnfilename
 582       //    time values file:      tvfilename
 583       _TimeSet* timeSet = timeSets.empty() ? 0 : & timeSets.back();
 584       // ---------------------------------------------------------
 585       if ( key == "time set" )
 586       {
 587         int nb = _ATOI( value );
 588         if ( nb < 1 ) {
 589           errorMsg << "Invalid time set number: " << value;
 590           RAISE_EXCEPTION;
 591         }
 592         timeSets.push_back( _TimeSet() );
 593         timeSets.back()._number = nb;
 594         READ_NEXT_LINE;
 595       }
 596       // ---------------------------------------------------------
 597       if ( key == "number of steps" )
 598       {
 599         if ( !timeSet || !timeSet->_times.empty() ) {
 600           errorMsg << "Unexpected command: " << key;
 601           RAISE_EXCEPTION;
 602         }
 603         int nbSteps = _ATOI( value );
 604         if ( nbSteps < 1 ) {
 605           errorMsg << "invalid number of steps: " << value;
 606           RAISE_EXCEPTION;
 607         }
 608         timeSet->_times.resize( nbSteps );
 609         timeSet->_fileIndex.resize( nbSteps );
 610         READ_NEXT_LINE;
 611       }
 612       // ---------------------------------------------------------
 613       if ( key == "filename start number" )
 614       {
 615         if ( !timeSet || timeSet->_fileIndex.empty() ) {
 616           errorMsg << "Unexpected command: " << key;
 617           RAISE_EXCEPTION;
 618         }
 619         if ( !reader.isDigit( value )) {
 620           errorMsg << "invalid " << line;
 621           RAISE_EXCEPTION;
 622         }
 623         timeSet->_fileIndex[0] = value;
 624         READ_NEXT_LINE;
 625       }
 626       // ---------------------------------------------------------
 627       if ( key == "filename increment" ) {
 628         int incr = _ATOI( value );
 629         if ( incr == 0 ) {
 630           errorMsg << "invalid " << line;
 631           RAISE_EXCEPTION;
 632         }
 633         if ( !timeSet ||
 634              timeSet->_fileIndex.empty() ||
 635              timeSet->_fileIndex[0].empty() ) {
 636           errorMsg << "Unexpected command: " << key;
 637           RAISE_EXCEPTION;
 638         }
 639         int index = incr + _ATOI( timeSet->_fileIndex[0] );
 640         int nbSteps = timeSet->_fileIndex.size();
 641         for ( int i = 1; i < nbSteps; ++i, index += incr )
 642           timeSet->_fileIndex[i] = STRING( index );
 643         READ_NEXT_LINE;
 644       }
 645       // ---------------------------------------------------------
 646       if ( key == "time values" )
 647       {
 648         if ( !timeSet || timeSet->_times.empty() ) {
 649           errorMsg << "Unexpected command: " << key;
 650           RAISE_EXCEPTION;
 651         }
 652         list<string> times;
 653         int i, nbTimes = reader.split( value, times );
 654         list<string>::iterator t = times.begin();
 655         for ( i = 0; i < nbTimes; ++i, ++t )
 656           timeSet->_times[i] = *t;
 657         while ( nbTimes != (int)timeSet->_times.size() ) {
 658           value = reader.getLine();
 659           ++lineNb;
 660           nbTimes += reader.split( value, times );
 661           for (t = times.begin(); i < nbTimes; ++i, ++t ) {
 662             if ( ! reader.isDigit( *t, /*real=*/true ))
 663               break;
 664             timeSet->_times[i] = *t;
 665           }
 666         }
 667         if ( nbTimes != (int)timeSet->_times.size() ) {
 668           errorMsg << "incorrect number of times in time set " << timeSet->_number;
 669           RAISE_EXCEPTION;
 670         }
 671         for ( i = 1; i < nbTimes; ++i )
 672           if ( _ATOF( timeSet->_times[ i ]) <= _ATOF(timeSet->_times[ i-1 ] ))
 673             break;
 674         if ( i < nbTimes ) { // break from the previous loop occured
 675           errorMsg << "time values are not in ascending order in time set " << timeSet->_number;
 676           RAISE_EXCEPTION;
 677         }
 678         READ_NEXT_LINE;
 679       }
 680       // ---------------------------------------------------------
 681       if ( key == "filename numbers" )
 682       {
 683         if ( !timeSet || timeSet->_fileIndex.empty()) {
 684           errorMsg << "Unexpected command: " << key;
 685           RAISE_EXCEPTION;
 686         }
 687         list<string> numbers;
 688         int i, nb = reader.split( value, numbers );
 689         int nbFiles = timeSet->_fileIndex.size();
 690         timeSet->_fileIndex.insert(timeSet->_fileIndex.begin(), numbers.begin(), numbers.end() );
 691         while ( nb != nbFiles ) {
 692           value = reader.getLine();
 693           ++lineNb;
 694           i = nb;
 695           nb += reader.split( value, numbers );
 696           list<string>::iterator n = numbers.begin();
 697           for ( ; i < nb; ++i, ++n ) {
 698             if ( ! reader.isDigit( *n ))
 699               break;
 700             timeSet->_fileIndex[i] = *n;
 701           }
 702         }
 703         if ( nb != nbFiles ) {
 704           errorMsg << "incorrect number of " << key << " in time set " << timeSet->_number;
 705           RAISE_EXCEPTION;
 706         }
 707         READ_NEXT_LINE;
 708       }
 709       // ---------------------------------------------------------
 710       if ( key == "filename numbers file" ||
 711            key == "time values file" )
 712       {
 713         if ( !timeSet || timeSet->_fileIndex.empty()) {
 714           errorMsg << "Unexpected command: " << key;
 715           RAISE_EXCEPTION;
 716         }
 717         string fileName = _directory + FILE_SEPARATOR + value;
 718         if ( !_user->canOpenFile( fileName, RDONLY )) {
 719           errorMsg << "Can not open file " << fileName;
 720           RAISE_EXCEPTION;
 721         }
 722         _ASCIIFileReader file( fileName );
 723         list<string> numbers;
 724         while ( !file.eof() )
 725           numbers.push_back( file.getWord() );
 726         int nb = numbers.size();
 727         if ( nb != (int)timeSet->_times.size() ) {
 728           errorMsg << "incorrect number of values in file " << value;
 729           RAISE_EXCEPTION;
 730         }
 731         if ( key[0] == 'f' )
 732           timeSet->_fileIndex.assign( numbers.begin(), numbers.end() );
 733         else
 734           timeSet->_times.assign( numbers.begin(), numbers.end() );
 735         READ_NEXT_LINE;
 736       }
 737       break;
 738     }
 739     case 'F': {
 740       if ( section[1] == 'I' ) {
 741         // --------------------------------------------------------------------------------
 742         // FILE
 743         // file set:               fs
 744         // filename index:         fi # Note: only used when data continues in other files
 745         // number of steps:        ns
 746         if ( key == "file set" ) {
 747           int nb = _ATOI( value );
 748           if ( nb < 1 ) {
 749             errorMsg << "Invalid file set number: " << value;
 750           }
 751           else {
 752             fileSets.push_back( _FileSet() );
 753             fileSets.back()._number = nb;
 754             READ_NEXT_LINE;
 755           }
 756         }
 757         else if ( key == "filename index" ) {
 758           if ( fileSets.empty() ) {
 759             errorMsg << "'filename index' before 'file set'";
 760           }
 761           else if ( !value.empty() ) {
 762             fileSets.back()._fileIndex.push_back( value );
 763             READ_NEXT_LINE;
 764           }
 765         }
 766         else if ( key == "number of steps" ) {
 767           if ( fileSets.empty() ) {
 768             errorMsg << "'number of steps' before 'file set'";
 769           }
 770           else if ( value.empty() ) {
 771             errorMsg << "number of steps omitted: " << line;
 772           }
 773           else if ( !reader.isDigit( value )) {
 774             errorMsg << "invalid number of steps: " << value;
 775           }
 776           else {
 777             int n = _ATOI( value );
 778             if ( n < 1 ) {
 779               errorMsg << "invalid number of steps: " << value;
 780               RAISE_EXCEPTION;
 781             }
 782             fileSets.back()._nbStepsInFile.push_back( n );
 783             READ_NEXT_LINE;
 784           }
 785         }
 786       }
 787       else {
 788         // --------------------------------------------------------------------------------
 789         // FORMAT
 790         // type:  ensight
 791         if ( key != "type" ) {
 792           errorMsg << "unexpected info in section " << section << ":\n" << line;
 793           RAISE_EXCEPTION;
 794         }
 795         else {
 796           if ( value == "ensight gold" ) {
 797             _format = ENSIGHT_GOLD;
 798           }
 799           else if ( value == "ensight" ) {
 800             _format = ENSIGHT_6;
 801           }
 802           else {
 803             errorMsg << "Unsupported EnSight format: " << value;
 804             RAISE_EXCEPTION;
 805           }
 806           section = "_";
 807           READ_NEXT_LINE;
 808         }
 809       }
 810       break;
 811     }
 812     case 'G': {
 813       // --------------------------------------------------------------------------------
 814       // GEOMETRY
 815       // model:    [ts] [fs] filename [change_coords_only]
 816       // measured: [ts] [fs] filename [change_coords_only]
 817       // match:              filename
 818       // boundary:           filename
 819       // rigid_body:         filename
 820       if ( key == "measured" || key == "match" || key == "boundary" || key == "rigid_body") {
 821         //errorMsg << key << " geometry not supported";
 822         cout << "Warning: " << key << " geomerty not supported" << endl;
 823         READ_NEXT_LINE;
 824       }
 825       else if ( key == "model" ) {
 826         list<string> parts;
 827         reader.split( value, parts );
 828         list<string>::reverse_iterator s = parts.rbegin();
 829         for ( ; s != parts.rend(); ++s )
 830         {
 831           if ( *s == "change_coords_only" )
 832             _model._change_coords_only = *s;
 833           else if ( _model._fileName.empty() )
 834             _model._fileName = *s;
 835           else if ( _model._fileSetNumber.empty() )
 836             _model._fileSetNumber = *s;
 837           else
 838             _model._timeSetNumber = *s;
 839         }
 840         if ( _model._timeSetNumber.empty() && !_model._fileSetNumber.empty() )
 841           swap( _model._timeSetNumber, _model._fileSetNumber );
 842         if ( _model._fileName.empty() ) {
 843           errorMsg << "invalid model: " << value;
 844           RAISE_EXCEPTION;
 845         }
 846         READ_NEXT_LINE;
 847       }
 848       break;
 849     }
 850     case 'M': { // MATERIAL
 851       string keyWord1;
 852       reader.split(key,keyWord1,key);
 853       if ( keyWord1 == "material" || keyWord1 == "species" )
 854         READ_NEXT_LINE;
 855       break;
 856     }
 857     } // end switch (section[0])
 858
 859     if ( !errorMsg.empty() ) {
 860       throw MEDEXCEPTION(STRING("Invalid Case file ") << _fileName
 861                          << ":" << lineNb << "\n" << errorMsg );
 862     }
 863
 864     // we are here if a line was not recognized to belong to a current section
 865     if ( line == "FORMAT"   ||
 866          line == "GEOMETRY" ||
 867          line == "VARIABLE" ||
 868          line == "TIME"     ||
 869          line == "FILE"     ||
 870          line == "MATERIAL")
 871       section = line;
 872     else
 873       throw MEDEXCEPTION(STRING() << "Invalid format of Case file " << _fileName
 874                          << "\nwrong line: " << line );
 875   }
 876
 877   if ( _model._fileName.empty() )
 878     throw MEDEXCEPTION(badFile << "no supported geometry information found" );
 879
 880   // store time sets and file sets
 881   list<_FileSet>::iterator fs = fileSets.begin();
 882   for ( ; fs != fileSets.end(); ++fs ) {
 883     if ( fs->_nbStepsInFile.size() > 1 && fs->_fileIndex.empty() )
 884       throw MEDEXCEPTION(badFile << "missing file indices in a file set " << fs->_number );
 885     _fileSets[ fs->_number ] = *fs;
 886   }
 887   list<_TimeSet>::iterator ts = timeSets.begin();
 888   for ( ; ts != timeSets.end(); ++ts )
 889     _timeSets[ ts->_number ] = *ts;
 890
 891   // check validity of ts and fs
 892
 893   if ( !isValidIndex( _model._timeSetNumber, _timeSets ))
 894     throw MEDEXCEPTION(badFile << "bad time set index:" << _model._timeSetNumber );
 895   if ( !isValidIndex( _model._fileSetNumber, _fileSets ))
 896     throw MEDEXCEPTION(badFile << "bad file set index:" << _model._timeSetNumber );
 897
 898   map< int, _Variable>::iterator ivars = _variables.begin();
 899   for ( ; ivars != _variables.end(); ++ivars ) {
 900     if ( !isValidIndex( ivars->second._timeSetNumber, _timeSets ))
 901       throw MEDEXCEPTION(badFile << "bad time set index:" << ivars->second._timeSetNumber );
 902     if ( !isValidIndex( ivars->second._fileSetNumber, _fileSets ))
 903       throw MEDEXCEPTION(badFile << "bad file set index:" << ivars->second._fileSetNumber );
 904   }
 905
 906   // check uniqueness of variable names
 907   if ( varNames.size() != _variables.size() )
 908     cout << badFile <<
 909       "Warning: there are different fields with equal names, you may have problems!" << endl;
 910 //     throw MEDEXCEPTION(badFile );
 911
 912   // As variable does not refer to time set if there is one step (issue 0020113),
 913   // we try to restore the reference
 914
 915   for ( ivars = _variables.begin(); ivars != _variables.end(); ++ivars ) {
 916     _Variable & var = ivars->second;
 917     if ( var._timeSetNumber.empty() )
 918     {
 919       // try to find time set with id equal to variable order number
 920       map< int, _TimeSet >::iterator iTs = _timeSets.find( ivars->first );
 921       if ( iTs != _timeSets.end() && iTs->second._times.size() == 1 )
 922         var._timeSetNumber = STRING( iTs->second._number );
 923       else {
 924         // find any time set with 1 time value
 925         for ( iTs = _timeSets.begin(); iTs != _timeSets.end(); ++iTs )
 926           if ( iTs->second._times.size() == 1 )
 927             var._timeSetNumber = STRING( iTs->second._number );
 928       }
 929     }
 930   }
 931 }
 932
 933 //================================================================================
 934 /*!
 935  * \brief return number of time steps of a model
 936  */
 937 //================================================================================
 938
 939 int _CaseFileDriver::getNbMeshes() const
 940 {
 941   if ( _blocked || checkWasRead())
 942     return 0 ;
 943   if ( _model._timeSetNumber.empty() )
 944     return 1;
 945   int ts = _ATOI( _model._timeSetNumber );
 946   map< int, _TimeSet >::const_iterator its = _timeSets.find( ts );
 947   if ( its == _timeSets.end() )
 948     throw MEDEXCEPTION(STRING() << "Invalid format of Case file " << _fileName
 949                        << "\n Inexistent time set number of a model" );
 950   return its->second._times.size();
 951 }
 952
 953 //================================================================================
 954 /*!
 955  * \brief  Sets all data necessary for meshDriver::read()
 956  *  \param meshIndex - time step index
 957  *  \param meshDriver - driver
 958  */
 959 //================================================================================
 960
 961 void _CaseFileDriver::setDataFileName(const int                   meshIndex,
 962                                       ENSIGHT_MESH_RDONLY_DRIVER* meshDriver)
 963 {
 964   if ( _blocked || checkWasRead())
 965     return;
 966   isToIgnore( meshDriver ); // remeber
 967
 968   meshDriver->_dataFileName    = _directory + FILE_SEPARATOR + _model._fileName;
 969   meshDriver->_indexInDataFile = fixWildCardName( meshIndex,
 970                                                   _model._timeSetNumber,
 971                                                   _model._fileSetNumber,
 972                                                   meshDriver->_dataFileName,
 973                                                   meshDriver->_time);
 974   meshDriver->_isGoldFormat    = ( _format == ENSIGHT_GOLD );
 975   meshDriver->_transientMode   = ( !_model._timeSetNumber.empty() );
 976   meshDriver->_singleFileMode  = ( !_fileSets.empty() );
 977   meshDriver->_imedMapKey      = STRING(_fileName)<<":"<<meshIndex;
 978
 979   GMESH* ptrMesh = meshDriver->getMesh();
 980   ptrMesh->setName(STRING("EnSight mesh ") << meshIndex);
 981 }
 982
 983 //================================================================================
 984 /*!
 985  * \brief Return nb of ensight variables
 986  */
 987 //================================================================================
 988
 989 int  _CaseFileDriver::getNbVariables() const
 990 {
 991   return _variables.empty() ? 0 : _variables.rbegin()->first;
 992 }
 993
 994 //================================================================================
 995 /*!
 996  * \brief return number of time steps of a variable
 997  */
 998 //================================================================================
 999
1000 int _CaseFileDriver::getNbVarSteps(const int variableIndex)
1001 {
1002   if ( _blocked || checkWasRead() )
1003     return 0;
1004
1005   map< int, _Variable>::const_iterator ivar = _variables.find( variableIndex );
1006   if ( ivar == _variables.end() ) {
1007     // it can be index of real part of complex variable
1008     ivar = _variables.find( variableIndex+1 );
1009     if ( ivar == _variables.end() || !contains( "complex", ivar->second._type.c_str() )) {
1010       throw MEDEXCEPTION(STRING( "_CaseFileDriver::getNbVarSteps(): invalid variable index: ")
1011                          << variableIndex);
1012     }
1013   }
1014   const _Variable & var = ivar->second;
1015   if ( var._timeSetNumber.empty() )
1016     return 1;
1017
1018   const _TimeSet & ts = _timeSets[ _ATOI( var._timeSetNumber)];
1019   return ts._times.size();
1020 }
1021
1022 //================================================================================
1023 /*!
1024  * \brief return variable index by variable name, return 0 if none found
1025  */
1026 //================================================================================
1027
1028 int _CaseFileDriver::getVariableIndex(const string & varName) const
1029 {
1030   if ( _blocked || checkWasRead() )
1031     return 0;
1032
1033   map< int, _Variable>::const_iterator ivar = _variables.begin();
1034   for ( ; ivar != _variables.end(); ++ivar )
1035   {
1036     if ( ivar->second._name == varName ) {
1037       if ( contains( "complex", ivar->second._type.c_str() ))
1038         return ivar->first - 1; // real part of complex variable
1039       return ivar->first;
1040     }
1041   }
1042   // maybe varName is "<true_varName>_Im"
1043   size_t _ImBeg = varName.size() - 3;
1044   if ( varName[ _ImBeg + 0 ] == '_' &&
1045        varName[ _ImBeg + 1 ] == 'I' &&
1046        varName[ _ImBeg + 2 ] == 'm' )
1047   {
1048     int i = getVariableIndex( varName.substr( 0, _ImBeg ));
1049     return ( i ? i + 1 : i ); // getVariableIndex() returns index for a real part
1050   }
1051   return 0;
1052 }
1053
1054 //================================================================================
1055 /*!
1056  * \brief sets all data necessary for fieldDriver::read()
1057  */
1058 //================================================================================
1059
1060 int _CaseFileDriver::setDataFileName(const int                    varIndex,
1061                                      const int                    stepIndex,
1062                                      ENSIGHT_FIELD_RDONLY_DRIVER* fieldDriver)
1063 {
1064   const char* LOC = "_CaseFileDriver::setDataFileName(): ";
1065   if ( _blocked || checkWasRead() )
1066     return 0;
1067   isToIgnore( fieldDriver ); // remeber
1068
1069   map< int, _Variable>::iterator ivar = _variables.find( varIndex );
1070   if ( ivar == _variables.end() ) {
1071     // it can be index of real part of complex variable
1072     ivar = _variables.find( varIndex+1 );
1073     if ( ivar == _variables.end() || !contains( "complex", ivar->second._type.c_str() ))
1074       throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "invalid variable index: " << varIndex));
1075   }
1076   _Variable & var = ivar->second;
1077
1078   const bool isComplex = contains( "complex", var._type.c_str() );
1079   const bool isImaginary = ( isComplex && varIndex == ivar->first );
1080
1081   // complex scalar per element:    [ts] [fs] description Re_fn   Im_fn freq
1082   // complex scalar per node:       [ts] [fs] description Re_fn   Im_fn freq
1083   // complex vector per element:    [ts] [fs] description Re_fn   Im_fn freq
1084   // complex vector per node:       [ts] [fs] description Re_fn   Im_fn freq
1085   // constant per case file:        [ts]      description cvfilename
1086   // constant per case:             [ts]      description const_value(s)
1087   // scalar per element:            [ts] [fs] description filename
1088   // scalar per measured node:      [ts] [fs] description filename
1089   // scalar per node:               [ts] [fs] description filename
1090   // tensor asym per element:       [ts] [fs] description filename
1091   // tensor asym per node:          [ts] [fs] description filename
1092   // tensor symm per element:       [ts] [fs] description filename
1093   // tensor symm per node:          [ts] [fs] description filename
1094   // vector per element:            [ts] [fs] description filename
1095   // vector per measured node:      [ts] [fs] description filename
1096   // vector per node:               [ts] [fs] description filename
1097
1098   FIELD_* field = fieldDriver->getField();
1099
1100   // field name
1101   if ( field->getName().empty() ) {
1102     if ( isImaginary )
1103       field->setName( var._name + "_Im" );
1104     else
1105       field->setName( var._name );
1106   }
1107   // nb of components
1108   list<string> type_parts;
1109   _ASCIIFileReader::split( var._type, type_parts );
1110   string type1 = type_parts.front(); type_parts.pop_front();
1111   string type2 = type_parts.front(); type_parts.pop_front();
1112   int nbComponents = 1;
1113   if ( type1 == "vector" || type2 == "vector" )
1114     nbComponents = 3;
1115   else if ( type1 == "tensor" )
1116     nbComponents = ( type2 == "symm" ) ? 6 : 9;
1117   field->setNumberOfComponents( nbComponents );
1118
1119   // component names
1120   vector<string> compNames( nbComponents );
1121   switch ( nbComponents ) {
1122   case 1:
1123     compNames[0] = type1;
1124     break;
1125   case 3: {
1126     const char* xyz[3] = { "X comp","Y comp","Z comp" };
1127     compNames.assign(xyz, xyz+3);
1128     break;
1129   }
1130   case 6: {
1131     const char* xyz[6] = { "11 comp", "22 comp", "33 comp","12 comp", "13 comp", "23 comp" };
1132     compNames.assign(xyz, xyz+6);
1133     break;
1134   }
1135   case 9: {
1136     const char* xyz[9] = { "11 comp", "12 comp", "13 comp", "21 comp", "22 comp", "23 comp",
1137                            "31 comp", "32 comp", "33 comp" };
1138     compNames.assign(xyz, xyz+9);
1139     break;
1140   }}
1141   field->setComponentsNames( & compNames[0] );
1142
1143   // componet units
1144   vector<UNIT> units( nbComponents );
1145   vector<string> unitS( nbComponents ), descriptions( nbComponents );
1146   field->setComponentsUnits(&units[0]);
1147   field->setMEDComponentsUnits(&unitS[0]);
1148   field->setComponentsDescriptions(&descriptions[0]);
1149
1150   // support type
1151   SUPPORT* sup = const_cast<SUPPORT*>( field->getSupport());
1152   if ( !sup ) {
1153     field->setSupport( sup = new SUPPORT );
1154     sup->removeReference(); // sup belongs to field
1155   }
1156   medEntityMesh entity = ( type_parts.back() == "node" ) ? MED_NODE : MED_CELL;
1157   sup->setEntity( entity );
1158
1159   // data file name etc.
1160   list<string> fileData_parts;
1161   _ASCIIFileReader::split( var._fileNameOrData, fileData_parts );
1162   if ( isImaginary )
1163     fieldDriver->_dataFileName  = _directory + FILE_SEPARATOR + *(++fileData_parts.begin());
1164   else
1165     fieldDriver->_dataFileName  = _directory + FILE_SEPARATOR + fileData_parts.front();
1166   fieldDriver->_indexInDataFile = fixWildCardName( stepIndex,
1167                                                    var._timeSetNumber,
1168                                                    var._fileSetNumber,
1169                                                    fieldDriver->_dataFileName,
1170                                                    fieldDriver->_time);
1171   fieldDriver->_isGoldFormat    = ( _format == ENSIGHT_GOLD );
1172   fieldDriver->_transientMode   = ( !var._timeSetNumber.empty() );
1173   fieldDriver->_singleFileMode  = ( !_fileSets.empty() );
1174
1175   if ( type1 == "constant" ) {
1176     if ( type_parts.back() == "file" ) {
1177       // append constant values from cvfilename to fileData_parts
1178       string cvfilename = _directory + FILE_SEPARATOR + fileData_parts.back();
1179       _ASCIIFileReader cvfile( cvfilename );
1180       fileData_parts.pop_back();
1181       while ( !cvfile.eof() )
1182         fileData_parts.push_back( cvfile.getWord() );
1183     }
1184     if ( (int)fileData_parts.size() < stepIndex )
1185       throw MEDEXCEPTION(LOCALIZED(STRING(LOC) << "can't find value for step " << stepIndex
1186                                    << " of " << var._type << " " << var._name));
1187     list<string>::iterator value = fileData_parts.begin();
1188     advance( value, stepIndex-1 );
1189     fieldDriver->setConstantValue( *value );
1190   }
1191
1192   // mesh data
1193   int meshIndex = 1;
1194   if ( _model._fileName.find('*') != _model._fileName.npos ) {
1195     _TimeSet& ts = _timeSets[ _ATOI( _model._timeSetNumber )];
1196     vector<string>::iterator t = find( ts._times.begin(), ts._times.end(), fieldDriver->_time );
1197     if ( t != ts._times.end() )
1198       meshIndex += distance( ts._times.begin(), t );
1199   }
1200   fieldDriver->_imedMapKey = STRING(_fileName) << ":" << meshIndex;
1201
1202   // time
1203   if ( fieldDriver->_transientMode ) {
1204     field->setTime( _ATOF( fieldDriver->_time ));
1205     //field->setOrderNumber( stepIndex );
1206     field->setIterationNumber( stepIndex );
1207   }
1208   return meshIndex;
1209 }
1210
1211 //================================================================================
1212 /*!
1213  * \brief Throws if case file has not been read else return false
1214  */
1215 //================================================================================
1216
1217 bool _CaseFileDriver::checkWasRead() const throw (MEDEXCEPTION)
1218 {
1219   if ( _model._fileName.empty() )
1220     throw MEDEXCEPTION(STRING("Case file ") << _fileName << " has not been successfully read");
1221   return false;
1222 }
1223
1224 //================================================================================
1225 /*!
1226  * \brief replace '*' in file name if any and return index in file
1227  */
1228 //================================================================================
1229
1230 int  _CaseFileDriver::fixWildCardName(const int      timeStep,
1231                                       const string & ts,
1232                                       const string & fs,
1233                                       string &       fileName,
1234                                       string &       time)
1235 {
1236   int indexInFile = 0;
1237   string fileIndex;
1238
1239   STRING badFile("Invalid Case file ");
1240   badFile << _fileName << "\n";
1241
1242   if ( !fs.empty() ) { // single file mode
1243     const _FileSet & fileSet = _fileSets[ _ATOI( fs ) ];
1244     if ( fileSet._fileIndex.empty() ) { // no file continuation
1245       indexInFile = timeStep;
1246     }
1247     else {
1248       list<int>::const_iterator   nbStepsIt = fileSet._nbStepsInFile.begin();
1249       list<string>::const_iterator fIndexIt = fileSet._fileIndex.begin();
1250       int nbSteps = 0;
1251       for ( ; nbStepsIt != fileSet._nbStepsInFile.end(); ++nbStepsIt ) {
1252         if ( nbSteps + *nbStepsIt <= timeStep )
1253           break;
1254         nbSteps += *nbStepsIt;
1255       }
1256       if ( nbStepsIt == fileSet._nbStepsInFile.end() )
1257         throw MEDEXCEPTION(LOCALIZED(badFile << "Cant'f find file index for time step " <<
1258                                      timeStep << " in file set " << fs ));
1259       indexInFile = timeStep - nbSteps;
1260       fileIndex   = *fIndexIt;
1261     }
1262   }
1263
1264   string head, queue;
1265   _ASCIIFileReader::split( fileName, head, queue, '*' );
1266   int indexWidth = fileName.size() - head.size() - queue.size();
1267
1268   if ( indexWidth > 0 || !ts.empty() || timeStep > 1 ) {
1269     int tsId = ts.empty() ? 1 : _ATOI( ts );
1270     const _TimeSet& timeSet = _timeSets[ tsId ];
1271     if ( timeStep > (int)timeSet._times.size() )
1272       throw MEDEXCEPTION(LOCALIZED(badFile << "Cant'f find time for time step " <<
1273                                    timeStep << " in time set " << ts ));
1274     time = timeSet._times[ timeStep-1 ];
1275     if ( timeStep-1 < (int)timeSet._fileIndex.size()  )
1276       fileIndex = timeSet._fileIndex[ timeStep-1 ];
1277     if ( !indexInFile )
1278       indexInFile = 1;
1279   }
1280
1281   if ( indexWidth > 0 )
1282   {
1283     if ( fileIndex.empty() ) {
1284       throw MEDEXCEPTION(LOCALIZED(badFile << "Can't find file index for time step " <<
1285                                    timeStep << " in time set <" << ts <<
1286                                    "> and file set <" << fs << ">"));
1287     }
1288     if ( indexWidth ==  (int)fileIndex.size() ) {
1289       fileName = head + fileIndex + queue;
1290     }
1291     else {
1292       fileName = (STRING(head) << setw(indexWidth) << setfill('0') << fileIndex << queue);
1293     }
1294   }
1295   return indexInFile;
1296 }
1297
1298 //================================================================================
1299 /*!
1300  * \brief add a mesh to the Case file
1301  */
1302 //================================================================================
1303
1304 void _CaseFileDriver::addMesh(const ENSIGHT_MESH_WRONLY_DRIVER* meshDriver)
1305 {
1306   if ( _blocked )
1307     return;
1308
1309   _meshDrivers.push_back( const_cast<ENSIGHT_MESH_WRONLY_DRIVER*>( meshDriver ));
1310
1311   if ( _format == ENSIGHT_6 )
1312   {
1313     const GMESH* mesh = _meshDrivers.back()->getMesh();
1314     if ( mesh->getNumberOfElements(MED_CELL, MED_POLYGON) > 0 ||
1315          mesh->getNumberOfElements(MED_FACE, MED_POLYGON) > 0 ||
1316          mesh->getNumberOfElements(MED_CELL, MED_POLYHEDRA) > 0 )
1317       throw MEDEXCEPTION
1318         ( compatibilityPb(STRING("Can't write mesh <") << mesh->getName() <<
1319                           "> since Ensight6 format does not support poly elements,"
1320                           " use Ensight Gold format instead: call "
1321                           "setEnSightFormatForWriting( ENSIGHT_GOLD )"));
1322   }
1323
1324   isToIgnore( meshDriver ); // remeber
1325 }
1326
1327 //================================================================================
1328 /*!
1329  * \brief add a field to the Case file
1330  */
1331 //================================================================================
1332
1333 void _CaseFileDriver::addField(const ENSIGHT_FIELD_WRONLY_DRIVER * theFieldDriver)
1334 {
1335   if ( _blocked )
1336     return;
1337
1338   ENSIGHT_FIELD_WRONLY_DRIVER * fieldDriver =
1339     const_cast<ENSIGHT_FIELD_WRONLY_DRIVER*>( theFieldDriver );
1340
1341   FIELD_* field = fieldDriver->getField();
1342
1343   STRING problem;
1344   if ( field->getNumberOfValues() == 0 )
1345     problem << "getNumberOfValues() == 0";
1346   else if ( field->getGaussPresence() )
1347     problem << compatibilityPb("Gauss points are not supported by EnSight v8");
1348   else if ( !field->getSupport() )
1349     problem << "it has NULL support";
1350
1351   switch ( field->getNumberOfComponents() ) {
1352   case 1:
1353   case 3:
1354   case 6:
1355   case 9: break; // ok, supported
1356   case 2:
1357     if ( const GMESH* mesh = field->getSupport()->getMesh() )
1358       if ( mesh->getSpaceDimension() == 2 )
1359         break; // we add one component to both mesh and field
1360   default:
1361     problem <<
1362       compatibilityPb(STRING("it has ") << field->getNumberOfComponents()
1363                       << " components but only 1,3,6 and 9 components are supported by EnSight");
1364   }
1365   if ( !problem.empty() )
1366     throw MEDEXCEPTION(STRING("Can't write field <") << field->getName() <<
1367                        "> to EnSight: " << problem);
1368
1369   string fieldName = fieldDriver->getFieldName();
1370   if ( fieldName.empty() )
1371     fieldName = field->getName();
1372   if ( fieldName.empty() )
1373     fieldName = STRING("med_field_")<<_fieldDrivers.size();
1374   else { // replace illegal characters
1375     string::size_type pos = fieldName.find_first_of( ILLEGAL_FIELD_NAME_CHARACTERS );
1376     while ( pos != fieldName.npos ) {
1377       fieldName[ pos ] = '_';
1378       pos = fieldName.find_first_of( ILLEGAL_FIELD_NAME_CHARACTERS );
1379     }
1380   }
1381
1382   _fieldDrivers[ fieldName ].push_back( fieldDriver );
1383
1384   isToIgnore( fieldDriver ); // remeber
1385 }
1386
1387 //================================================================================
1388 /*!
1389  * \brief writing Case file
1390  */
1391 //================================================================================
1392
1393 void _CaseFileDriver::write() throw (MEDEXCEPTION)
1394 {
1395   if ( _blocked )
1396     return;
1397
1398   // Make case file data from added drivers
1399
1400   bool dataIsRead = !_model._fileName.empty();
1401   if ( !dataIsRead && _meshDrivers.empty() )
1402     throw MEDEXCEPTION("no mesh to write into Case file");
1403
1404   const int defaultNbDigits = 3; // in file index
1405
1406   string path, name, ext;
1407   _ASCIIFileReader::split( _fileName, path, name, FILE_SEPARATOR, /*fromBack=*/true);
1408   _ASCIIFileReader::split( name, name, ext, '.', /*fromBack=*/true);
1409   if ( name.empty())
1410     name = ext; // _fileName: "path/.name"
1411
1412   list<ENSIGHT_MESH_WRONLY_DRIVER*>::iterator mDrv = _meshDrivers.begin();
1413   TFieldDriversByName::iterator               fDrv = _fieldDrivers.begin();
1414
1415   map< string, string > fileToRenameTo;
1416
1417   int i, nbOldMeshes = 0, nbNewMeshes = _meshDrivers.size();
1418
1419   if ( nbNewMeshes > 0 )
1420   {
1421     if ( dataIsRead )
1422     {
1423       // A mesh is going to be added into an existing case file.
1424       // Check that number of parts is same in the existing model
1425       // and in the added ones
1426       string geoFileName     = _directory + FILE_SEPARATOR + _model._fileName;
1427       string realGeoFileName = geoFileName, time;
1428       fixWildCardName( 1,
1429                        _model._timeSetNumber,
1430                        _model._fileSetNumber,
1431                        realGeoFileName,
1432                        time);
1433       int nbParts = ENSIGHT_MESH_RDONLY_DRIVER::countParts(realGeoFileName,
1434                                                            !_model._fileSetNumber.empty());
1435       for ( ; mDrv != _meshDrivers.end(); ++mDrv )
1436       {
1437         int nbNewParts = (*mDrv)->nbPartsToWrite();
1438         if ( nbParts != nbNewParts )
1439           throw MEDEXCEPTION(compatibilityPb("Can't add a mesh ") << (*mDrv)->getMesh()->getName()
1440                              << "to " << _fileName << " as number of parts in the file "
1441                              "differs from number of parts going to be written");
1442       }
1443       if ( !_model._timeSetNumber.empty() && !_variables.empty() ) {
1444         // Make time set of the model independent of that of variables as
1445         // the number of time steps is apparently becoming different
1446         map< int, _Variable>::iterator ivar = _variables.begin();
1447         for ( ; ivar != _variables.end(); ++ivar ) {
1448           if ( ivar->second._timeSetNumber == _model._timeSetNumber ) {
1449             int tsNum    = _timeSets.rbegin()->first + 1;
1450             _TimeSet& ts = _timeSets[ tsNum ];
1451             ts._number   = tsNum;
1452             ts           = _timeSets[ _ATOI( _model._timeSetNumber )];
1453             _model._timeSetNumber = STRING(tsNum);
1454             break;
1455           }
1456         }
1457       }
1458       string::size_type digitPos = _model._fileName.find('*');
1459       bool isSeveralFiles = ( digitPos != _model._fileName.npos );
1460       int nbDigits = defaultNbDigits;
1461       if ( isSeveralFiles ) {
1462         nbDigits = 1;
1463         while ( _model._fileName[ ++digitPos ] == '*' )
1464           nbDigits++;
1465       }
1466       // Update time set and file set of the model
1467       int nbMeshes;
1468       if ( _model._timeSetNumber.empty() ) {
1469         // old model is static, create a time set
1470         nbMeshes = nbNewMeshes + 1;
1471         int tsNum    = _timeSets.empty() ? 1 : _timeSets.rbegin()->first + 1;
1472         _TimeSet& ts = _timeSets[ tsNum ];
1473         ts._number   = tsNum;
1474         _model._timeSetNumber = STRING(tsNum);
1475         // fill the new time set
1476         ts._fileIndex.resize( nbMeshes );
1477         ts._times.resize    ( nbMeshes );
1478         for ( i = 0; i < nbMeshes; ++i ) {
1479           ts._fileIndex[ i ] = (STRING() << setw(nbDigits) << setfill('0') << i);
1480           ts._times    [ i ] = STRING(i);
1481         }
1482         // not to create equal time sets
1483         map< int, _TimeSet >::iterator its, tsEnd = _timeSets.end();
1484         for ( --tsEnd, its = _timeSets.begin(); its != tsEnd; ++its ) {
1485           if ( ts == its->second ) {
1486             _model._timeSetNumber = STRING( its->first );
1487             _timeSets.erase( tsEnd );
1488             break;
1489           }
1490         }
1491       }
1492       else {
1493         // old model is transient, add times and file indices for new meshes
1494         _TimeSet& ts = _timeSets[ _ATOI( _model._timeSetNumber )];
1495         nbOldMeshes = ts._times.size();
1496         nbMeshes = nbNewMeshes + nbOldMeshes;
1497         ts._times.resize( nbMeshes );
1498         double time = 1. + _ATOF( ts._times[ nbOldMeshes-1 ] );
1499         for ( i = nbOldMeshes; i < nbMeshes; ++i, time+=1. )
1500           ts._times[ i ] = STRING( time );
1501         if ( _model._fileSetNumber.empty() ) { // multi-file mode
1502           ts._fileIndex.resize( nbMeshes, string(nbDigits,'0'));
1503           i = nbOldMeshes;
1504           int index = 1 + _ATOI( ts._fileIndex[ i-1 ] );
1505           for ( ; i < nbMeshes; ++i, ++index )
1506             ts._fileIndex[ i ] = (STRING() << setw(nbDigits) << setfill('0') << index);
1507         }
1508         else { // single-file mode
1509           _FileSet & fs = _fileSets[ _ATOI( _model._fileSetNumber )];
1510           for ( i = 0; i < nbNewMeshes; ++i )
1511             fs._nbStepsInFile.push_back( 1 );
1512           int index = 1;
1513           if ( fs._fileIndex.empty() )
1514             fs._fileIndex.push_back(string(nbDigits,'0'));
1515           else
1516             index += _ATOI( fs._fileIndex.back() );
1517           for ( i = fs._fileIndex.size(); i < nbMeshes; ++i, ++index )
1518             ts._fileIndex[ i ] = (STRING() << setw(3) << setfill('0') << index);
1519         }
1520       }
1521       // file of the old model is to be renamed
1522       if ( !isSeveralFiles ) {
1523         _model._fileName += string(nbDigits,'*');
1524         fileToRenameTo[ geoFileName ] = geoFileName + string(nbDigits,'0');
1525       }
1526     }
1527     else if ( nbNewMeshes > 1 )
1528     {
1529       // Store meshes into a new case file: create a time set
1530       int tsNum    = _timeSets.empty() ? 1 : _timeSets.rbegin()->first + 1;
1531       _TimeSet& ts = _timeSets[ tsNum ];
1532       ts._number   = tsNum;
1533       _model._timeSetNumber = STRING(tsNum);
1534       _model._fileName      = name + ".geo" + string(defaultNbDigits, '*');
1535       // fill the new time set
1536       ts._fileIndex.resize( nbNewMeshes );
1537       ts._times.resize    ( nbNewMeshes );
1538       for ( i = 0; i < nbNewMeshes; ++i ) {
1539         ts._fileIndex[ i ] = (STRING() << setw(defaultNbDigits) << setfill('0') << i);
1540         ts._times    [ i ] = STRING(i);
1541       }
1542     }
1543     else {
1544       // One mesh in a new file
1545       _model._fileName = name + ".geo";
1546     }
1547
1548     // Set data to mesh drivers
1549     i = nbOldMeshes + 1;
1550     for ( mDrv = _meshDrivers.begin(); mDrv != _meshDrivers.end(); ++mDrv, ++i )
1551     {
1552       _CaseFileDriver_User* meshDriver = (*mDrv);
1553       meshDriver->_dataFileName    = _directory + FILE_SEPARATOR + _model._fileName;
1554       meshDriver->_indexInDataFile = fixWildCardName( i,
1555                                                       _model._timeSetNumber,
1556                                                       _model._fileSetNumber,
1557                                                       meshDriver->_dataFileName,
1558                                                       meshDriver->_time);
1559       meshDriver->_isGoldFormat    = ( _format == ENSIGHT_GOLD );
1560       meshDriver->_transientMode   = ( meshDriver->_indexInDataFile > 0 );
1561       meshDriver->_singleFileMode  = ( !_fileSets.empty() );
1562       meshDriver->_imedMapKey      = STRING(_fileName)<<":"<<i;
1563     }
1564   }
1565
1566   typedef map< double, ENSIGHT_FIELD_WRONLY_DRIVER* > FDriverByDouble;
1567
1568   //bool isVarsRead = ( !_variables.empty() );
1569   for ( ; fDrv != _fieldDrivers.end(); ++fDrv )
1570   {
1571     const string &                       fieldName = fDrv->first;
1572     list< ENSIGHT_FIELD_WRONLY_DRIVER* > & drivers = fDrv->second;
1573
1574     int nbNewSteps = drivers.size();
1575
1576     // find out by which parameter fields differ and sort them by the parameter
1577     FDriverByDouble timeMap, iterMap, orderMap;
1578     list< ENSIGHT_FIELD_WRONLY_DRIVER* >::iterator drv;
1579     for ( drv = drivers.begin(); drv != drivers.end(); ++drv ) {
1580       FIELD_* field = (*drv)->getField();
1581       double time = field->getTime();
1582       int    iter = field->getIterationNumber();
1583       int    ordr = field->getOrderNumber();
1584       timeMap.insert  ( make_pair( time, *drv ));
1585       iterMap.insert  ( make_pair( iter, *drv ));
1586       orderMap.insert ( make_pair( ordr, *drv ));
1587     }
1588     FDriverByDouble * sortedDrivers;
1589     FDriverByDouble::iterator tDrv;
1590     if ( (int)timeMap.size() == nbNewSteps )
1591       sortedDrivers = & timeMap;
1592     else if ( (int)iterMap.size() == nbNewSteps )
1593       sortedDrivers = & iterMap;
1594     else if ( (int)orderMap.size() == nbNewSteps )
1595       sortedDrivers = & orderMap;
1596     else {
1597       timeMap.clear();
1598       sortedDrivers = & timeMap;
1599       for ( drv = drivers.begin(); drv != drivers.end(); ++drv ) {
1600         double time = (*drv)->getField()->getTime();
1601         if ( ! timeMap.insert( make_pair( time, *drv )).second )
1602           timeMap.insert( make_pair( timeMap.rbegin()->first + 1., *drv ));
1603       }
1604     }
1605
1606 //     if ( isVarsRead ) {
1607 //       int iVar = getVariableIndex( fieldName );
1608 //       if ( iVar > 0 ) {
1609 //         // A variable with fieldName already exists,
1610 //         // add more time steps to it
1611 //         _Variable& var = _variables[ iVar ];
1612 //         _TimeSet& ts = _timeSets[ _ATOI( var._timeSetNumber )];
1613 //         int nbOldSteps = ts._times.size();
1614 //       }
1615 //     }
1616     FIELD_* field = drivers.front()->getField();
1617     int varNum = _variables.size() + 1;
1618     _Variable& var = _variables[ varNum ];
1619     var._name = fieldName;
1620     switch ( field->getNumberOfComponents() ) {
1621     case 1: var._type = "scalar "; break;
1622     case 2: var._type = "vector "; break;// we add one component to a vector in 2d space
1623     case 3: var._type = "vector "; break;
1624     case 6: var._type = "tensor symm "; break;
1625     case 9: var._type = "tensor asym "; break;
1626     }
1627     if ( field->getSupport()->getEntity() == MED_NODE )
1628       var._type += "per node";
1629     else
1630       var._type += "per element";
1631     var._fileNameOrData = name + "." + fieldName;
1632
1633     // always create Time set to store time
1634     int nbDigits = defaultNbDigits;
1635     int tsNum = _timeSets.empty() ? 1 : _timeSets.rbegin()->first + 1;
1636     var._timeSetNumber = STRING( tsNum );
1637     _TimeSet & ts = _timeSets[ tsNum ];
1638     ts._number    = tsNum;
1639     ts._fileIndex.resize( nbNewSteps );
1640     ts._times.resize    ( nbNewSteps );
1641     tDrv = sortedDrivers->begin();
1642     for ( i = 0; tDrv != sortedDrivers->end(); ++tDrv, ++i ) {
1643       ts._times    [ i ] = (STRING( tDrv->first ));
1644       ts._fileIndex[ i ] = (STRING() << setw(nbDigits) << setfill('0') << i);
1645     }
1646     if ( nbNewSteps > 1 )
1647       var._fileNameOrData += string( nbDigits, '*' );
1648     else
1649       ts._fileIndex.clear();
1650     // not to create equal time sets
1651     map< int, _TimeSet >::iterator its, tsEnd = _timeSets.end();
1652     for ( --tsEnd, its = _timeSets.begin(); its != tsEnd; ++its ) {
1653       if ( ts == its->second ) {
1654         tsNum = its->first;
1655         var._timeSetNumber = STRING( tsNum );
1656         _timeSets.erase( tsEnd );
1657         break;
1658       }
1659     }
1660     tDrv = sortedDrivers->begin();
1661     for ( i = 1; tDrv != sortedDrivers->end(); ++tDrv, ++i ) {
1662       _CaseFileDriver_User* fieldDriver = tDrv->second;
1663       fieldDriver->_dataFileName    = _directory + FILE_SEPARATOR + var._fileNameOrData;
1664       fieldDriver->_indexInDataFile = fixWildCardName( i,
1665                                                        var._timeSetNumber,
1666                                                        var._fileSetNumber,
1667                                                        fieldDriver->_dataFileName,
1668                                                        fieldDriver->_time);
1669       fieldDriver->_isGoldFormat    = ( _format == ENSIGHT_GOLD );
1670       fieldDriver->_transientMode   = ( fieldDriver->_indexInDataFile > 0 );
1671       fieldDriver->_singleFileMode  = ( !_fileSets.empty() );
1672     }
1673
1674     // do not refer to time set if there is one step (issue 0020113)
1675     if ( nbNewSteps == 1 )
1676       var._timeSetNumber = "";
1677
1678   } // loop on _fieldDrivers
1679
1680   if ( nbOldMeshes + nbNewMeshes > 1 && !_variables.empty() )
1681   {
1682     // if there are variables on a transient model, model should not change much,
1683     // nb of entities in all parts with values must be same
1684     if ( nbOldMeshes == 0 ) {
1685       // TODO: check consistency
1686     }
1687     if ( toIgnoreIncompatibility() )
1688       cout << "Warning: EnSight file will be probably invalid " << _fileName << endl;
1689     else
1690       throw MEDEXCEPTION
1691         (compatibilityPb(STRING("EnSight file will be invalid if fields refer to "
1692                                 "the second mesh, which differ from the first one too much.\n"
1693                                 "If you are sure in data correctness, you can suppress "
1694                                 "this exception by calling setIgnoreIncompatibility(1)\n")
1695                          << _fileName));
1696   }
1697
1698   // Write a file
1699
1700   ofstream caseFile( _fileName.c_str(), ios::out );
1701   if ( !caseFile )
1702     throw MEDEXCEPTION(STRING("Can't open file for writing ") << _fileName);
1703
1704   caseFile << "# generated by MEDMEM-to-EnSight driver" << endl << endl
1705            << "FORMAT" << endl
1706            << "type: " << (_format==ENSIGHT_GOLD ? "ensight gold" : "ensight") << endl
1707            << endl
1708            << "GEOMETRY"
1709            << endl
1710            << "model:\t"
1711            << _model._timeSetNumber << " "
1712            << _model._fileSetNumber << " "
1713            << _model._fileName << " \t"
1714            << _model._change_coords_only
1715            << endl << endl;
1716
1717   // write VARIABLE section
1718   if ( !_variables.empty() )
1719   {
1720     caseFile << "VARIABLE" << endl;
1721
1722     map< int, _Variable>::iterator ivar = _variables.begin();
1723     for ( ; ivar != _variables.end(); ++ivar )
1724     {
1725       _Variable& var = ivar->second;
1726       caseFile << var._type << ": \t"
1727                << var._timeSetNumber << " "
1728                << var._fileSetNumber << " \t"
1729                << var._name << " \t"
1730                << var._fileNameOrData << endl;
1731     }
1732     caseFile << endl;
1733   }
1734   // write TIME section
1735   if ( !_timeSets.empty() )
1736   {
1737     caseFile << "TIME" << endl;
1738
1739     map< int, _TimeSet>::iterator its = _timeSets.begin();
1740     for ( ; its != _timeSets.end(); ++its )
1741     {
1742       _TimeSet & ts = its->second;
1743       caseFile << "time set:\t"        << ts._number << endl
1744                << "number of steps:\t" << ts._times.size() << endl;
1745       if ( !ts._fileIndex.empty() ) {
1746         STRING numbers( "filename numbers:" );
1747         for ( unsigned i = 0; i < ts._fileIndex.size(); ++i ) {
1748           if (int( numbers.size() + ts._fileIndex[i].size() + 2) > MAX_LINE_LENGTH ) {
1749             caseFile << numbers << endl;
1750             numbers = STRING();
1751           }
1752           numbers << " " << ts._fileIndex[i];
1753         }
1754         caseFile << numbers << endl;
1755       }
1756       STRING times( "time values:" );
1757       for ( unsigned i = 0; i < ts._times.size(); ++i ) {
1758         if (int( times.size() + ts._times[i].size() + 2 ) > MAX_LINE_LENGTH ) {
1759           caseFile << times << endl;
1760           times = STRING();
1761         }
1762         times << " " << ts._times[i];
1763       }
1764       caseFile << times << endl;
1765     }
1766   }
1767   // write FILE section
1768   if ( !_fileSets.empty() )
1769   {
1770     caseFile << "FILE" << endl;
1771
1772     map< int, _FileSet >::iterator ifs = _fileSets.begin();
1773     for ( ; ifs != _fileSets.end(); ++ifs )
1774     {
1775       _FileSet & fs = ifs->second;
1776       caseFile << "file set: " << fs._number << endl;
1777
1778       list<int>::iterator nbSteps = fs._nbStepsInFile.begin();
1779       list<string>::iterator fileIndex = fs._fileIndex.begin();
1780       for ( ; nbSteps != fs._nbStepsInFile.end(); ++nbSteps )
1781       {
1782         if ( fileIndex != fs._fileIndex.end() )
1783           caseFile << "filename index: " << *fileIndex++;
1784         caseFile << "number of steps: " << *nbSteps << endl;
1785       }
1786     }
1787   }
1788
1789   caseFile.close();
1790
1791 } // _CaseFileDriver::write()
1792
1793
1794 //================================================================================
1795 /*!
1796  * \brief _CaseFileDriver_User constructor
1797  */
1798 //================================================================================
1799
1800 _CaseFileDriver_User::_CaseFileDriver_User(const string & caseFileName,
1801                                            med_mode_acces mode)
1802   : GENDRIVER( caseFileName, mode, ENSIGHT_DRIVER ), _imed(0)
1803 {
1804 }
1805
1806 //================================================================================
1807 /*!
1808  * \brief take data from other driver
1809  */
1810 //================================================================================
1811
1812 void _CaseFileDriver_User::merge( const GENDRIVER& driver)
1813 {
1814   const _CaseFileDriver_User* other = dynamic_cast< const _CaseFileDriver_User* >( &driver );
1815   if ( other ) {
1816     _dataFileName    = other->_dataFileName   ;
1817     _isGoldFormat    = other->_isGoldFormat   ;
1818     _transientMode   = other->_transientMode  ;
1819     _singleFileMode  = other->_singleFileMode ;
1820     _indexInDataFile = other->_indexInDataFile;
1821     _time            = other->_time           ;
1822     _imed            = other->_imed           ;
1823     _imedMapKey      = other->_imedMapKey     ;
1824   }
1825 }
1826
1827 //================================================================================
1828 /*!
1829  * \brief analyse if data file is binary
1830  */
1831 //================================================================================
1832
1833 bool _CaseFileDriver_User::isBinaryDataFile(const string& dataFileName)
1834 {
1835 #ifdef WIN32
1836   int _file = ::_open (dataFileName.c_str(), _O_RDONLY|_O_BINARY);
1837 #else
1838   int _file = ::open (dataFileName.c_str(), O_RDONLY);
1839 #endif
1840   char buf[81];
1841   int nBytesRead = ::read (_file, buf, 80);
1842
1843   bool isBinary = true;
1844
1845   const char cBin[] = "C Binary";
1846   const char fBin[] = "Fortran Binary";
1847   if ( strncmp( buf, cBin, sizeof(cBin)-1) != 0 &&
1848        strncmp( buf, fBin, sizeof(fBin)-1) != 0 )
1849   {
1850     for ( int i = nBytesRead-1; i >= 0 && isBinary; --i )
1851       isBinary = ( buf[ i ] != '\n' );
1852   }
1853
1854   ::close (_file);
1855
1856   return isBinary;
1857 }
1858
1859 //================================================================================
1860 /*!
1861  * \brief return part number to write support with, zero in failure case
1862  */
1863 //================================================================================
1864
1865 int _CaseFileDriver_User::getPartNumber(const SUPPORT* support) const
1866 {
1867   bool isGroup = ( dynamic_cast<const GROUP*>( support ));
1868   bool isForField = ( dynamic_cast<const ENSIGHT_FIELD_DRIVER*>( this ));
1869   medEntityMesh entity = support->getEntity();
1870   const GMESH* mesh = support->getMesh();
1871
1872   // for supports on all entities, reserve numbers corresponding to entity
1873   bool isOnAll = support->isOnAllElements();
1874   if (!isOnAll && mesh ) {
1875     int nbMeshElem = mesh->getNumberOfElements(entity, MED_ALL_ELEMENTS);
1876     int nbSuppElem = support->getNumberOfElements(MED_ALL_ELEMENTS);
1877     isOnAll = (  nbSuppElem == nbMeshElem );
1878   }
1879   if ( !isGroup ) {
1880     if ( !isOnAll )
1881       return 0;
1882     else if ( entity == MED_NODE )
1883       return 1 + MED_CELL; // all nodes are described with all CELLs
1884     else
1885       return 1 + entity;
1886   }
1887   if ( isForField && isOnAll ) {
1888     if ( entity == MED_NODE )
1889       return 1 + MED_CELL; // all nodes are described with all CELLs
1890     else
1891       return 1 + entity;
1892   }
1893   if ( !mesh )
1894     return 0;
1895
1896   int partNum = MED_ALL_ENTITIES + 1;
1897   for ( int ent = MED_CELL; ent < MED_ALL_ENTITIES; ++ent ) {
1898     entity = (medEntityMesh) ent;
1899     int nbGroups = mesh->getNumberOfGroups(entity);
1900     if ( entity != support->getEntity() ) {
1901       partNum += nbGroups;
1902     }
1903     else {
1904       for ( int i=1; i<=nbGroups; ++i, ++partNum)
1905         if ( support == mesh->getGroup( entity, i ))
1906           return partNum;
1907     }
1908   }
1909   throw MED_EXCEPTION
1910     ( LOCALIZED( STRING("Can't find GROUP ") << support->getName() << " in its MESH"));
1911
1912   return 0;
1913 }
1914
1915 //================================================================================
1916 /*!
1917  * \brief Return true if we strore an enetuty to EnSight
1918  */
1919 //================================================================================
1920
1921 bool _CaseFileDriver_User::isToWriteEntity(const medEntityMesh entity,
1922                                            const GMESH*        mesh)
1923 {
1924   if ( entity == MED_NODE )
1925     return mesh->getNumberOfNodes() > 0;
1926
1927   if ( mesh->getNumberOfElements( entity, MED_ALL_ELEMENTS ) < 1 )
1928     return false;
1929   if ( entity == MED_CELL )
1930     return true;
1931   int meshDim = mesh->getTypes(MED_CELL)[0] / 100;
1932   if ( entity == MED_FACE )
1933     return ( meshDim == 3 || meshDim == 5 );
1934   if ( entity == MED_EDGE )
1935     return ( meshDim == 2 || meshDim == 4 );
1936
1937   return false;
1938 }
1939
1940 //================================================================================
1941 /*!
1942  * \brief Return nodes of non-nodal support, which is not on all entities
1943  */
1944 //================================================================================
1945
1946 void _CaseFileDriver_User::getSupportNodes(const SUPPORT* support, map<int, int> & nodeIds)
1947 {
1948   medEntityMesh entity = support->getEntity();
1949
1950   const medConnectivity conType     = MED_NODAL;
1951   const medGeometryElement allGeoms = MED_ALL_ELEMENTS;
1952   const int * connectivity          = 0;
1953   const int * elemConnectivity      = 0;
1954   const int * index                 = 0;
1955   const int * number                = 0;
1956   int j;
1957
1958   if ( support->isOnAllElements() )
1959   {
1960     if ( entity == MED_NODE ) { // all NODES
1961       int numberOfCell = support->getNumberOfElements(allGeoms);
1962       while ( numberOfCell ) {
1963         nodeIds.insert( nodeIds.begin(), make_pair( numberOfCell, numberOfCell ));
1964         --numberOfCell;
1965       }
1966     }
1967     else {
1968       const MESH* mesh = support->getMesh()->convertInMESH();
1969       int conLength = 0;
1970       connectivity = mesh->getConnectivity      (conType, entity, allGeoms);
1971       conLength    = mesh->getConnectivityLength(conType, entity, allGeoms);
1972       while ( conLength-- ) nodeIds[ *connectivity++ ];
1973       mesh->removeReference();
1974     }
1975     return;
1976   }
1977
1978   if ( entity == MED_NODE )
1979   {
1980     number           = support->getNumber(MED_ALL_ELEMENTS);
1981     int numberOfCell = support->getNumberOfElements(MED_ALL_ELEMENTS);
1982     for (j=0; j < numberOfCell; j++)
1983       nodeIds.insert(nodeIds.end(), make_pair( number[j], j ));
1984     return;
1985   }
1986
1987   number           = support->getNumber(MED_ALL_ELEMENTS);
1988   int numberOfCell = support->getNumberOfElements(MED_ALL_ELEMENTS);
1989   const MESH* mesh = support->getMesh()->convertInMESH();
1990   index        = mesh->getConnectivityIndex( MED_NODAL, entity);
1991   connectivity = mesh->getConnectivity( MED_NODAL, entity, MED_ALL_ELEMENTS);
1992   for ( j = 0; j < numberOfCell; ++j )
1993   {
1994     int elem = number[j];
1995     elemConnectivity   = connectivity + index[elem-1]-1;
1996     const int* connEnd = connectivity + index[elem]-1;
1997     while ( elemConnectivity < connEnd )
1998       nodeIds[ *elemConnectivity++ ];
1999   }
2000   mesh->removeReference();
2001 }
2002
2003 //================================================================================
2004 /*!
2005  * \brief method called by mesh driver at reading to store data to be used by
2006  *        field driver
2007  */
2008 //================================================================================
2009
2010 void _CaseFileDriver_User::setInterData(_InterMed* imed )
2011 {
2012   theInterMedMap[ _imedMapKey ] = imed;
2013   if ( ENSIGHT_MESH_DRIVER* mDrv = dynamic_cast<ENSIGHT_MESH_DRIVER*>( this )) {
2014     imed->_medMesh = dynamic_cast<MESH*>( mDrv->getMesh() );
2015     imed->_isOwnMedMesh = false;
2016   }
2017   else
2018     imed->_medMesh = 0;
2019 }
2020
2021 //================================================================================
2022 /*!
2023  * \brief for field driver to get mesh data
2024  */
2025 //================================================================================
2026
2027 _InterMed* _CaseFileDriver_User::getInterData()
2028 {
2029   return _imed = getMeshData( _imedMapKey );
2030 }
2031
2032 //================================================================================
2033 /*!
2034  * \brief return _SubPart by its description
2035  */
2036 //================================================================================
2037
2038 _SubPart* _CaseFileDriver_User::getSubPart(const _SubPartDesc & descriptor)
2039   throw (MEDEXCEPTION)
2040 {
2041   if ( !_imed )
2042     _imed = getMeshData( _imedMapKey );
2043
2044   map< _SubPartDesc, _SubPart >::iterator descPart =
2045     _imed->_subPartDescribed.find( descriptor );
2046
2047   if ( descPart != _imed->_subPartDescribed.end() )
2048     return & descPart->second;
2049
2050   if ( descriptor == _SubPartDesc::globalCoordDesc() )
2051     return 0; // if the mesh is structured, there are no global coordinates (EnSight 6)
2052
2053   throw MEDEXCEPTION(LOCALIZED(STRING("No mesh info for part ") << descriptor.partNumber()
2054                                << " " << descriptor.typeName()));
2055 }
2056
2057 //================================================================================
2058 /*!
2059  * \brief return _Support by its description
2060  */
2061 //================================================================================
2062
2063 _Support* _CaseFileDriver_User::getSupport(const _SupportDesc & descriptor,
2064                                            const medEntityMesh  entity)
2065   throw (MEDEXCEPTION)
2066 {
2067   const char* LOC = "_CaseFileDriver_User::getSupport(): ";
2068   if ( !_imed )
2069     _imed = getMeshData( _imedMapKey );
2070
2071   _imed->treatGroupes();
2072
2073   if ( _imed->_supportDescribed.empty() && !_imed->_subPartDescribed.empty() )
2074   {
2075     // fill _supportDescribed with _Support's corresponding to EnSight parts
2076     for (unsigned int i=0; i < _imed->groupes.size(); ++i)
2077     {
2078       _groupe& grp = _imed->groupes[i];
2079       if ( !grp.medGroup ) continue;
2080
2081       vector<int> grpIndices(1,i);
2082       if ( !grp.groupes.empty() )
2083         grpIndices.assign( grp.groupes.begin(), grp.groupes.end());
2084
2085       _SupportDesc supDescriptor;
2086       // look for a subpart for each _groupe
2087       vector<int>::iterator grIndex = grpIndices.begin(), grIndexEnd = grpIndices.end();
2088       for ( ; grIndex != grIndexEnd; ++grIndex )
2089       {
2090         map< _SubPartDesc, _SubPart >::iterator descSub;
2091         for ( descSub  = _imed->_subPartDescribed.begin();
2092               descSub != _imed->_subPartDescribed.end();
2093               ++descSub)
2094         {
2095           if ( descSub->second.myCellGroupIndex == *grIndex ) {
2096             supDescriptor.insert( descSub->first );
2097             break;
2098           }
2099         }
2100       }
2101       if ( supDescriptor.size() == grpIndices.size() ) {
2102         _Support & sup = _imed->_supportDescribed[ supDescriptor ];
2103         sup.setGroup( & grp );
2104       }
2105     }
2106   }
2107
2108   // find the support by its description
2109   map< _SupportDesc, _Support >::iterator descSup =
2110     _imed->_supportDescribed.find( descriptor );
2111
2112   // Create a new support
2113   if ( descSup == _imed->_supportDescribed.end() || !descSup->second.medSupport( entity ))
2114   {
2115     // create a _groupe composed of groups corresponding to subparts
2116     _imed->groupes.push_back(_groupe());
2117     _groupe& grp = _imed->groupes.back();
2118     grp.groupes.reserve( descriptor.size() );
2119
2120     // to detect dimension change within a support group
2121     set<int> dimensions;
2122
2123     // fill grp with sub-group indices
2124     _SupportDesc::const_iterator subPartDesc = descriptor.begin();
2125     for ( ; subPartDesc != descriptor.end(); ++subPartDesc )
2126     {
2127       const _SubPart* subPart = getSubPart(*subPartDesc);
2128       if ( !subPart )
2129         throw MEDEXCEPTION(LOCALIZED(STRING("No mesh info for ") << *subPartDesc));
2130
2131       int groupIndex =
2132         (entity == MED_NODE) ? subPart->myNodeGroupIndex : subPart->myCellGroupIndex;
2133       if ( groupIndex < 1 ) {
2134         if ( entity == MED_NODE )
2135         {
2136           // make a _groupe of nodes
2137           _imed->groupes.push_back(_groupe());
2138           groupIndex = subPart->myNodeGroupIndex = _imed->groupes.size();
2139           _groupe & groupe = _imed->groupes.back();
2140           groupe.mailles.resize( subPart->myNbNodes );
2141
2142           _maille ma( MED_POINT1, 1 );
2143           ma.sommets.resize( 1 );
2144           map< int, _noeud >::iterator n = subPart->myFirstNode;
2145           for (int i = 0; i < subPart->myNbNodes; ++i, ++n ) {
2146             ma.sommets[0] = n;
2147             ma.setOrdre( n->second.number );
2148             groupe.mailles[i] = _imed->insert(ma);
2149           }
2150         }
2151         else
2152         {
2153           throw MEDEXCEPTION(LOCALIZED(STRING("No cell info for ") << *subPartDesc));
2154         }
2155       }
2156       grp.groupes.push_back( groupIndex );
2157
2158       // find out subpart dimension
2159       if ( entity != MED_NODE )
2160         dimensions.insert( _imed->groupes[ groupIndex-1 ].maille(0).dimensionWithPoly() );
2161     }
2162     // check if MEDMEM allows creating such a support
2163     if ( dimensions.size() > 1 )
2164       throw MEDEXCEPTION
2165         (compatibilityPb(LOC) << "can't create a SUPPORT for the field from "
2166          << _dataFileName << ", since it includes different mesh entities");
2167
2168     ENSIGHT_MESH_RDONLY_DRIVER::makeGroup( grp, *_imed );
2169
2170     // add _Support
2171     descSup = _imed->_supportDescribed.insert( make_pair( descriptor, _Support() )).first;
2172     _Support & sup = descSup->second;
2173     sup.setGroup( & grp );
2174
2175   } // new SUPPORT creation
2176
2177   _Support* sup = & descSup->second;
2178
2179   // remove temporary mesh from med SUPPORT
2180   if ( _imed->_isOwnMedMesh )
2181   {
2182     if ( sup->medSupport( entity )->getMesh() == _imed->_medMesh )
2183       _imed->_medMesh->addReference(); // don't want _medMesh to die.
2184     sup->medSupport( entity )->setMesh( 0 );
2185     sup->medSupport( entity )->setMeshName( _imed->_medMesh->getName() );
2186   }
2187
2188   return sup;
2189 }
2190
2191 //================================================================================
2192 /*!
2193  * \brief check possibility to open a file with a given mode
2194  */
2195 //================================================================================
2196
2197 bool _CaseFileDriver_User::canOpenFile(const string&  fileName,
2198                                        med_mode_acces mode)
2199 {
2200   bool ok = false;
2201   if ( mode == WRONLY ) {
2202     fstream file( fileName.c_str(),
2203                   ios::app | ios_base::out ); // not to overwrite it, just to check
2204     ok = bool(file);
2205   }
2206   else {
2207     fstream file( fileName.c_str(), ios::in );
2208     ok = bool(file);
2209   }
2210   return ok;
2211 }
2212
2213 //================================================================================
2214 /*!
2215  * \brief _CaseFileDriver_User destructor
2216  */
2217 //================================================================================
2218
2219 _CaseFileDriver_User::~_CaseFileDriver_User()
2220 {
2221   if ( _imed )
2222     _imed->_nbUsers--;
2223   unregister( this );
2224 }
2225
2226 //================================================================================
2227 /*!
2228  * \brief Search substring in a string
2229  */
2230 //================================================================================
2231
2232 bool contains( const char* what, const char* inString )
2233 {
2234   size_t whatLen = strlen( what );
2235   size_t inLen   = strlen( inString );
2236   return
2237     ( search( inString, inString+inLen, what, what + whatLen) != inString+inLen );
2238 }
2239
2240 //================================================================================
2241 /*!
2242  * \brief store subPart data
2243  */
2244 //================================================================================
2245
2246 void _InterMed::addSubPart(const _SubPart& theSubPart)
2247 {
2248   if ( _needSubParts ) {
2249     _SubPart & subPart = _subPartDescribed[ theSubPart.getDescriptor() ];
2250     subPart = theSubPart;
2251     if ( subPart.myCellGroupIndex > 0 ) {
2252       _groupe & groupe = this->groupes[ subPart.myCellGroupIndex-1 ];
2253       subPart.myFirstCell = groupe.mailles.begin();
2254     }
2255   }
2256 }
2257 //================================================================================
2258 /*!
2259  * \brief delete intermediate med data
2260  */
2261 //================================================================================
2262
2263 _InterMed::~_InterMed()
2264 {
2265   if ( _isOwnMedMesh )
2266   {
2267     // remove MEDMEM groups not belonging to _medMesh
2268     for (unsigned int i=0; i < _intermediateMED::groupes.size(); ++i)
2269     {
2270       _groupe& grp = _intermediateMED::groupes[i];
2271       if ( !grp.medGroup ) continue;
2272       vector<GROUP*> groups = _medMesh->getGroups( grp.medGroup->getEntity() );
2273       if ( find( groups.begin(), groups.end(), grp.medGroup ) == groups.end() )
2274         grp.medGroup->removeReference();
2275     }
2276     if(_medMesh) _medMesh->removeReference();
2277     _medMesh=0;
2278   }
2279 }
2280 //================================================================================
2281 /*!
2282  * \brief For a node, find its index in the supporting GROUP
2283  */
2284 //================================================================================
2285
2286 int _Support::getIndex( const pair<const int,_noeud>& inode)
2287 {
2288   if ( myNodeGroup->relocMap.empty() ) // on all and not self intersecting support
2289     return abs( inode.second.number );
2290
2291   map<unsigned,int>::iterator ordreIndex = myNodeGroup->relocMap.find( abs( inode.second.number ));
2292   if ( ordreIndex == myNodeGroup->relocMap.end() )
2293     throw MEDEXCEPTION(LOCALIZED(STRING("No index found for ") << inode.second));
2294 //   map < int, int >::iterator numIndex = myNodeRelocMap.find( node->number );
2295 //   if ( numIndex == myNodeRelocMap.end() )
2296 //     throw MEDEXCEPTION(STRING("No index found for node ") << node->number);
2297
2298   return ordreIndex->second;
2299 }
2300
2301 //================================================================================
2302 /*!
2303  * \brief For a cell, return its index in the supporting GROUP
2304  */
2305 //================================================================================
2306
2307 int _Support::getIndex( const _groupe::TMaille& cell)
2308 {
2309   if ( myCellGroup->relocMap.empty() ) // on all and not self intersecting support
2310     return cell->ordre();
2311
2312   map<unsigned,int>::iterator ordreIndex = myCellGroup->relocMap.find( cell->ordre() );
2313   if ( ordreIndex == myCellGroup->relocMap.end() )
2314     throw MEDEXCEPTION(LOCALIZED(STRING("No index found for ") << *cell));
2315
2316   return ordreIndex->second;
2317 }
2318
2319 //================================================================================
2320 /*!
2321  * \brief Return med geom type for a subPart element
2322  */
2323 //================================================================================
2324
2325 // medGeometryElement _Support::getType( const pair<const int,_noeud>& node)
2326 // {
2327 //   return myNodeGroup->medGroup.getTypes[0];
2328 // }
2329
2330 // //================================================================================
2331 // /*!
2332 //  * \brief Return med geom type for a subPart element
2333 //  */
2334 // //================================================================================
2335
2336 // medGeometryElement _Support::getType( const _groupe::TMaille& cell)
2337 // {
2338 //   return cell->geometricType;
2339 // }
2340
2341 //================================================================================
2342 /*!
2343  * \brief Set med support
2344  */
2345 //================================================================================
2346
2347 void _Support::setGroup( _groupe* group )
2348 {
2349   if ( group->medGroup ) {
2350     if ( group->medGroup->getEntity() == MED_NODE )
2351       myNodeGroup = group;
2352     else
2353       myCellGroup = group;
2354   }
2355   else {
2356     throw MEDEXCEPTION(LOCALIZED("_Support::setGroup(): med GROUP is NULL"));
2357   }
2358 }
2359
2360 //================================================================================
2361 /*!
2362  * \brief Return med group correspoding to entity
2363  */
2364 //================================================================================
2365
2366 SUPPORT* _Support::medSupport( medEntityMesh entity )
2367 {
2368   if ( entity == MED_NODE )
2369     return myNodeGroup ? myNodeGroup->medGroup : 0;
2370   else
2371     return myCellGroup ? myCellGroup->medGroup : 0;
2372 }
2373
2374 //================================================================================
2375 /*!
2376  * \brief print _SubPartDesc
2377  */
2378 //================================================================================
2379
2380 std::ostream& operator << (std::ostream& os, const _SubPartDesc& desc)
2381 {
2382   if ( desc == _SubPartDesc::globalCoordDesc() )
2383     os << "'global coordinates'";
2384   else
2385     os << "<'part " << desc.partNumber() << "', '" << desc.typeName() << "'>";
2386   return os;
2387 }
2388
2389 //================================================================================
2390 /*!
2391  * \brief Constructor of ASCIIFileReader
2392  */
2393 //================================================================================
2394
2395 _ASCIIFileReader::_ASCIIFileReader(const string& fileName) throw (MEDEXCEPTION)
2396 {
2397 #ifdef WIN32
2398   _file = ::_open (fileName.c_str(), _O_RDONLY|_O_BINARY);
2399 #else
2400   _file = ::open (fileName.c_str(), O_RDONLY);
2401 #endif
2402   if (_file >= 0)
2403   {
2404     _start = new char [BUFFER_SIZE];
2405     _ptr   = _start;
2406     _eptr  = _start;
2407   }
2408   else
2409   {
2410     throw MEDEXCEPTION(STRING("Can't read from ")<<fileName);
2411   }
2412   if ( eof() )
2413     throw MEDEXCEPTION(STRING("Empty file ")<<fileName);
2414
2415   // there must be end-of-line in ASCII file
2416   char* ptr = _start + MAX_LINE_LENGTH;
2417   bool isASCII = false;
2418   while ( !isASCII && ptr >= _start )
2419     isASCII = (*ptr-- == '\n');
2420   _isWin = ( *ptr == '\r');
2421   if ( !isASCII )
2422     throw MEDEXCEPTION(STRING("Not ASCII file ")<<fileName);
2423 }
2424
2425 //================================================================================
2426 /*!
2427  * \brief Return true if the whole file has been read
2428  */
2429 //================================================================================
2430
2431 bool _ASCIIFileReader::eof()
2432 {
2433   // Check the state of the buffer;
2434   // if there is too little left, read the next portion of data
2435   int nBytesRest = _eptr - _ptr;
2436   if (nBytesRest < 2 * MAX_LINE_LENGTH)
2437   {
2438     if (nBytesRest > 0) {
2439       char* tmpBuf = new char [nBytesRest];
2440       memcpy (tmpBuf, _ptr, nBytesRest);
2441       memcpy (_start, tmpBuf, nBytesRest);
2442       delete [] tmpBuf;
2443     } else {
2444       nBytesRest = 0;
2445     }
2446     _ptr = _start;
2447     const int nBytesRead = ::read (_file,
2448                                    &_start [nBytesRest],
2449                                    BUFFER_SIZE - nBytesRest);
2450     nBytesRest += nBytesRead;
2451     _eptr = &_start [nBytesRest];
2452
2453     // skip spaces at file end
2454     if ( nBytesRest < MAX_LINE_LENGTH ) {
2455       while ( isspace( *_ptr )) _ptr++;
2456       nBytesRest = _eptr - _ptr;
2457     }
2458   }
2459   return nBytesRest < 1;
2460 }
2461
2462 //================================================================================
2463 /*!
2464  * \brief read out given data
2465  */
2466 //================================================================================
2467
2468 void _ASCIIFileReader::skip(int nbVals, int nbPerLine, int valWidth)
2469 {
2470   int nbLines = (nbVals + nbPerLine - 1) / nbPerLine;
2471   int width = nbVals * valWidth;
2472   skip( width, nbLines );
2473 }
2474
2475 //================================================================================
2476 /*!
2477  * \brief read out width chars and nbLines line-ends
2478  */
2479 //================================================================================
2480
2481 void _ASCIIFileReader::skip(int width, int nbLines)
2482 {
2483   width += nbLines; // to skip '\n'
2484   if ( _isWin )
2485     width += nbLines; // to skip '\r'
2486
2487   _ptr += width;
2488   int nBytesRest = _eptr - _ptr;
2489   while ( nBytesRest < 0 ) {
2490     width = -nBytesRest;
2491     if ( eof() ) return;
2492     _ptr += width;
2493     nBytesRest = _eptr - _ptr;
2494   }
2495 }
2496
2497 //================================================================================
2498 /*!
2499  * \brief Read a next line from the file
2500  */
2501 //================================================================================
2502
2503 char* _ASCIIFileReader::getLine() throw (MEDEXCEPTION)
2504 {
2505   if ( eof() )
2506     throw MEDEXCEPTION("Unexpected EOF");
2507
2508   // Check the buffer for the end-of-line
2509   char * ptr = _ptr;
2510   while (true)
2511   {
2512     // Check for end-of-the-buffer, the ultimate criterion for termination
2513     if (ptr >= _eptr)
2514     {
2515       //_eptr[-1] = '\0';
2516       _eptr[0] = '\0';
2517       break;
2518     }
2519     // seek the line-feed character
2520     if (ptr[0] == '\n')
2521     {
2522       if ( ptr > _start && // avoid "Invalid read" error by valgrind
2523            ptr[-1] == '\r')
2524         ptr[-1] = '\0';
2525       ptr[0] = '\0';
2526       ++ptr;
2527       break;
2528     }
2529     ++ptr;
2530   }
2531   // Output the result
2532   char * line = _ptr;
2533   _ptr = ptr;
2534
2535   return line;
2536 }
2537
2538 //================================================================================
2539 /*!
2540  * \brief
2541  */
2542 //================================================================================
2543
2544 bool _ASCIIFileReader::lookAt( const char* text )
2545 {
2546   while ( isspace(*_ptr)) ++_ptr;
2547   return ( strncmp( _ptr, text, strlen( text )) == 0 );
2548 }
2549
2550 //================================================================================
2551 /*!
2552  * \brief Read a next word from the file
2553  */
2554 //================================================================================
2555
2556 string _ASCIIFileReader::getWord()
2557 {
2558   if ( eof() )
2559     return "";
2560
2561   // skip spaces
2562   while ( isspace(*_ptr)) ++_ptr;
2563   if ( _ptr >= _eptr )
2564     return "";
2565
2566   // skip not spaces
2567   char* word = _ptr++;
2568   while ( !isspace(*_ptr)) ++_ptr;
2569
2570   return string( word, _ptr - word );
2571 }
2572
2573 //================================================================================
2574 /*!
2575  * \brief Return true if TIME_STEP_BEG follows and skips TIME_STEP_BEG line
2576  */
2577 //================================================================================
2578
2579 bool _ASCIIFileReader::isTimeStepBeginning()
2580 {
2581   if ( eof() ) throw MEDEXCEPTION(LOCALIZED("Unexpected EOF"));
2582
2583   while ( isspace(*_ptr)) ++_ptr;
2584
2585   if ( strncmp( _ptr, TIME_STEP_BEG, TIME_STEP_BEG_LEN ) != 0 )
2586     return false;
2587
2588   _ptr += TIME_STEP_BEG_LEN;
2589   while ( isspace(*_ptr)) ++_ptr;
2590   return true;
2591 }
2592
2593 //================================================================================
2594 /*!
2595  * \brief Return true if TIME_STEP_END follows and skips TIME_STEP_END line
2596  */
2597 //================================================================================
2598
2599 bool _ASCIIFileReader::isTimeStepEnd()
2600 {
2601   if ( eof() ) return true;
2602
2603   while ( isspace(*_ptr)) ++_ptr;
2604
2605   if ( strncmp( _ptr, TIME_STEP_END, TIME_STEP_END_LEN ) != 0 )
2606     return false;
2607
2608   _ptr += TIME_STEP_END_LEN;
2609   while ( isspace(*_ptr)) ++_ptr;
2610   return true;
2611 }
2612
2613
2614 //================================================================================
2615 /*!
2616  * \brief divide a string into two parts
2617  */
2618 //================================================================================
2619
2620 int _ASCIIFileReader::split(const string& str,
2621                             string &      part1,
2622                             string &      part2,
2623                             const char    separator,
2624                             const bool    fromBack)
2625 {
2626   int nbParts = 0;
2627   string parts[2];
2628   const char* ptr1 = str.c_str();
2629   const char* back = ptr1 + str.size();
2630   for (nbParts = 0; nbParts < 2; ++nbParts ) {
2631     // skip spaces before the current part
2632     while ( isspace(*ptr1)) ++ptr1;
2633     if ( !*ptr1) break;
2634     // find end of the part and beginning of the next part
2635     const char* ptr2 = ptr1;
2636     const char* nextBeg = back;
2637     if ( nbParts > 0 ) {
2638       ptr2 = back;
2639     }
2640     else if ( fromBack ) {
2641       ptr2 = back;
2642       string::size_type pos = str.rfind( separator );
2643       if ( pos != str.npos ) {
2644         ptr2 = & str[ pos ];
2645         nextBeg = ptr2 + 1;
2646         if ( separator != ' ')
2647           while ( *nextBeg && *nextBeg == separator ) ++nextBeg;
2648       }
2649     }
2650     else if ( separator == ' ' ) {
2651       while ( *ptr2 && !isspace(*ptr2)) ++ptr2;
2652       if ( *ptr2 ) nextBeg = ptr2 + 1;
2653     }
2654     else {
2655       while ( *ptr2 && *ptr2 != separator ) ++ptr2;
2656       if ( *ptr2 ) {
2657         nextBeg = ptr2 + 1;
2658         while ( *nextBeg && *nextBeg == separator ) ++nextBeg;
2659       }
2660     }
2661     //if ( !*ptr2) --ptr2;
2662     // skip spaces after the current part
2663     while ( ptr2 > ptr1 && isspace(ptr2[-1])) --ptr2;
2664     parts[ nbParts ] = string( ptr1, ptr2-ptr1 );
2665     ptr1 = nextBeg;
2666   }
2667   part1 = parts[0];
2668   part2 = parts[1];
2669   return nbParts;
2670 }
2671
2672 //================================================================================
2673 /*!
2674  * \brief divide a string into parts, return nb of parts
2675  */
2676 //================================================================================
2677
2678 int _ASCIIFileReader::split(const string&       str,
2679                             list<string> &      parts,
2680                             const char          separator,
2681                             const bool          fromBack)
2682 {
2683   parts.clear();
2684   if ( str.empty() )
2685     return 0;
2686   int nbParts = 0;
2687   const char* ptr1 = str.c_str();
2688   const char* back = ptr1 + str.size();
2689   if ( fromBack ) {
2690     swap( ptr1, back );
2691     while (1) {
2692       // skip spaces after the current part
2693       while ( isspace(ptr1[-1])) --ptr1;
2694       if ( ptr1 <= back ) break;
2695       // find beginning of the part
2696       const char* ptr2 = ptr1 - 1;
2697       if ( separator == ' ' )
2698         while ( ptr2 > back && !isspace(ptr2[-1])) --ptr2;
2699       else
2700         while ( ptr2 > back && ptr2[-1] != separator ) --ptr2;
2701       //if ( !*ptr2) --ptr2;
2702       const char* sepPtr = ptr2;
2703       // skip spaces before the current part
2704       while ( isspace(*ptr2)) ++ptr2;
2705       parts.push_back( string( ptr2, ptr1-ptr2 ));
2706       ++nbParts;
2707       ptr1 = sepPtr - 1;
2708     }
2709   }
2710   else {
2711     while (1) {
2712       // skip spaces before the current part
2713       while ( isspace(*ptr1)) ++ptr1;
2714       if ( ptr1 >= back) break;
2715       // find end of the part
2716       const char* ptr2 = ptr1 + 1;
2717       if ( separator == ' ' )
2718         while ( *ptr2 && !isspace(*ptr2)) ++ptr2;
2719       else
2720         while ( *ptr2 && *ptr2 != separator ) ++ptr2;
2721       //if ( !*ptr2) --ptr2;
2722       const char* sepPtr = ptr2;
2723       // skip spaces after the current part
2724       while ( ptr2 > ptr1 && isspace(ptr2[-1])) --ptr2;
2725       parts.push_back( string( ptr1, ptr2-ptr1 ));
2726       ++nbParts;
2727       ptr1 = sepPtr + int( sepPtr < back );
2728     }
2729   }
2730   return nbParts;
2731 }
2732
2733 //================================================================================
2734 /*!
2735  * \brief check if a string contains only digits
2736  */
2737 //================================================================================
2738
2739 bool _ASCIIFileReader::isDigit(const string& str, const bool real)
2740 {
2741   const char* s = str.c_str();
2742   if ( real ) {
2743     while ( *s ) {
2744       char c = *s++;
2745       if ( !isdigit(c) && c!='-' && c!='+' && c!='.' && c!=',' && c!='E' && c!='e')
2746         return false;
2747     }
2748   }
2749   else {
2750     while ( *s ) {
2751       if ( !isdigit( *s++ ))
2752         return false;
2753     }
2754   }
2755   return true;
2756 }
2757
2758 //================================================================================
2759 /*!
2760  * \brief Destructor of ASCIIFileReader closes its file
2761  */
2762 //================================================================================
2763
2764 _ASCIIFileReader::~_ASCIIFileReader()
2765 {
2766   if (_file >= 0)
2767   {
2768     ::close (_file);
2769     delete [] _start;
2770   }
2771 }
2772
2773 //================================================================================
2774 /*!
2775  * \brief Constructor of _BinaryFileReader opens its file
2776  */
2777 //================================================================================
2778
2779 _BinaryFileReader::_BinaryFileReader(const string& fileName) throw (MEDEXCEPTION)
2780   : _exception(STRING("Unexpected EOF ") << fileName), _mySwapBytes(false)
2781 {
2782 #ifdef WIN32
2783   _file = ::_open (fileName.c_str(), _O_RDONLY|_O_BINARY);
2784 #else
2785   _file = ::open (fileName.c_str(), O_RDONLY);
2786 #endif
2787
2788   if (_file < 0)
2789     throw MEDEXCEPTION(STRING("Can't read from ") << fileName);
2790
2791   _maxPos = ::lseek( _file, 0, SEEK_END);
2792   _pos    = ::lseek( _file, 0, SEEK_SET);
2793 }
2794
2795 //================================================================================
2796 /*!
2797  * \brief Destructor of _BinaryFileWriter closes its file
2798  */
2799 //================================================================================
2800
2801 _BinaryFileReader::~_BinaryFileReader()
2802 {
2803   if (_file >= 0)
2804     ::close (_file);
2805 }
2806
2807 //================================================================================
2808 /*!
2809  * \brief rewind the file backward
2810  */
2811 //================================================================================
2812
2813 void _BinaryFileReader::rewind()
2814 {
2815   _pos = ::lseek( _file, 0, SEEK_SET);
2816 }
2817
2818 //================================================================================
2819 /*!
2820  * \brief size of not read file portion in sizeof(int)
2821  */
2822 //================================================================================
2823
2824 int _BinaryFileReader::moreValuesAvailable() const
2825 {
2826   return (_maxPos - _pos) / sizeof(int); // believe that sizeof(int) == sizeof(float)
2827 }
2828
2829 //================================================================================
2830 /*!
2831  * \brief Return true if the whole file has been read
2832  */
2833 //================================================================================
2834
2835 bool _BinaryFileReader::eof()
2836 {
2837   return _pos >= _maxPos;
2838 }
2839
2840 //================================================================================
2841 /*!
2842  * \brief Skips given nb of bytes
2843  */
2844 //================================================================================
2845
2846 void _BinaryFileReader::skip(int size) throw (MEDEXCEPTION)
2847 {
2848   if ( _pos + size > _maxPos )
2849     throw _exception;
2850   off_t newPos = ::lseek( _file, size, SEEK_CUR);
2851   if ( newPos < _pos + size )
2852     throw _exception;
2853   _pos = newPos;
2854 }
2855
2856 //================================================================================
2857 /*!
2858  * \brief Read lines until TIME_STEP_BEG encounters
2859  */
2860 //================================================================================
2861
2862 void _BinaryFileReader::skipTimeStepBeginning() throw (MEDEXCEPTION)
2863 {
2864   bool tsbReached = false;
2865   while ( !tsbReached ) {
2866     TStrOwner line( getLine() );
2867     tsbReached = ( strncmp( line, TIME_STEP_BEG, TIME_STEP_BEG_LEN ) == 0 );
2868   }
2869 }
2870
2871 //================================================================================
2872 /*!
2873  * \brief Constructor of _BinaryFileWriter opens its file
2874  */
2875 //================================================================================
2876
2877 _BinaryFileWriter::_BinaryFileWriter(const string& fileName)  throw (MEDEXCEPTION)
2878   : _exception(STRING("Can't write into ") << fileName)
2879 {
2880 #ifdef WIN32
2881   _file = ::_open (fileName.c_str(), _O_WRONLY|_O_BINARY|_O_TRUNC);
2882 #else
2883   _file = ::open (fileName.c_str(), O_WRONLY|O_TRUNC); //length shall be truncated to 0
2884 #endif
2885
2886   if (_file < 0)
2887     throw _exception;
2888 }
2889
2890 //================================================================================
2891 /*!
2892  * \brief Destructor of _BinaryFileWriter closes its file
2893  */
2894 //================================================================================
2895
2896 _BinaryFileWriter::~_BinaryFileWriter()
2897 {
2898   if (_file >= 0)
2899     ::close (_file);
2900 }
2901
2902 //================================================================================
2903 /*!
2904  * \brief Write string
2905  */
2906 //================================================================================
2907
2908 void _BinaryFileWriter::addString(const char* str) throw (MEDEXCEPTION)
2909 {
2910   size_t len = strlen( str );
2911   if ((int) len > MAX_LINE_LENGTH )
2912     throw MEDEXCEPTION
2913       (LOCALIZED(STRING("_BinaryFileWriter::addString(), too long string (>80):\n") << str));
2914
2915   string buffer( str, len );
2916   // avoid "Syscall param write(buf) points to uninitialised byte(s)" error by valgrind
2917   buffer += string(MAX_LINE_LENGTH, ' ');
2918   buffer[ len ] = '\0';
2919   buffer[ MAX_LINE_LENGTH-1 ] = '\0'; // to have line-end within
2920
2921   add( buffer.c_str(), MAX_LINE_LENGTH );
2922 }
2923
2924
2925
2926 } // end namespace MEDMEM_ENSIGHT