1 // Copyright (C) 2007-2012 CEA/DEN, EDF R&D, OPEN CASCADE
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.
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.
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
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 // File : MEDMEM_EnsightUtils.cxx
21 // Created : Fri Jun 6 16:08:32 2008
22 // Author : Edward AGAPOV (eap)
24 #include "MEDMEM_Field.hxx" // important to include it before MEDMEM_EnsightUtils.hxx,
25 // in order not to redefine DBL_MIN
27 #include "MEDMEM_EnsightUtils.hxx"
29 #include "MEDMEM_EnsightMeshDriver.hxx"
30 #include "MEDMEM_EnsightFieldDriver.hxx"
31 #include "MEDMEM_DriverTools.hxx"
42 using namespace MEDMEM;
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
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
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
80 // ==============================================================================
84 // ==============================================================================
86 EnSightFormat theEnSightFormatForWriting = ENSIGHT_6;
87 bool theBinaryFormatForWriting = false;
88 bool theIgnoreIncompatibility = false;
90 void setEnSightFormatForWriting (EnSightFormat format, bool isBinary)
92 theEnSightFormatForWriting = format;
93 theBinaryFormatForWriting = isBinary;
95 EnSightFormat getEnSightFormatForWriting()
97 return theEnSightFormatForWriting;
99 bool isBinaryEnSightFormatForWriting()
101 return theBinaryFormatForWriting;
103 // ==============================================================================
105 To raise or not if MEDMEM-EnSight incompatibility encounters or suspected
107 // ==============================================================================
109 void setIgnoreIncompatibility(bool toIgnore)
111 theIgnoreIncompatibility = toIgnore;
116 #define FILE_SEPARATOR '\\'
118 #define FILE_SEPARATOR '/'
121 #define _ATOI( str ) atoi((str).c_str())
122 #define _ATOF( str ) atof((str).c_str())
124 #define BUFFER_SIZE 16184 // for non-stream input
126 namespace MEDMEM_ENSIGHT
128 //--------------------------------------------------------------------------------
130 * \brief Registry of drivers creating _CaseFileDriver
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.
138 set< const _CaseFileDriver_User* > theCaseUsers;
140 //--------------------------------------------------------------------------------
142 * \brief Storage of data read by mesh driver to be used by field drivers
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>".
149 map< string, _InterMed* > theInterMedMap;
151 //--------------------------------------------------------------------------------
153 * \brief Add a driver to the registry and return true if was already in
155 bool isToIgnore(const _CaseFileDriver_User* driver);
156 bool isToIgnore(const _CaseFileDriver_User* driver)
158 return ! theCaseUsers.insert( driver ).second;
161 //--------------------------------------------------------------------------------
163 * \brief Remove a driver from the registry
165 void unregister(const _CaseFileDriver_User* driver);
166 void unregister(const _CaseFileDriver_User* driver)
168 theCaseUsers.erase( driver );
169 if ( theCaseUsers.empty() )
171 for ( map< string, _InterMed* >::iterator str_imed = theInterMedMap.begin();
172 str_imed != theInterMedMap.end();
174 delete str_imed->second;
175 theInterMedMap.clear();
179 //--------------------------------------------------------------------------------
181 * \brief Return mesh data needed by field driver
183 _InterMed* getMeshData( const string& key );
184 _InterMed* getMeshData( const string& key )
186 // find existing data
187 map< string, _InterMed* >::iterator kimed = theInterMedMap.find( key );
188 if ( kimed != theInterMedMap.end() )
189 return kimed->second;
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 ));
197 kimed = theInterMedMap.find( key );
198 if ( kimed == theInterMedMap.end() )
200 _InterMed* imed = kimed->second;
201 imed->_medMesh = mesh;
202 imed->_isOwnMedMesh = true;
205 // ---------------------------------------------------------------
207 * \brief Prepend "EnSight-MEDMEM compatibility problem" to the text
210 STRING compatibilityPb(const string& exceptionText)
212 return STRING("EnSight-MEDMEM compatibility problem:\n") << exceptionText;
214 // ---------------------------------------------------------------
216 * \brief To ignore incompatibility or not
218 bool toIgnoreIncompatibility()
220 return theIgnoreIncompatibility;
223 //--------------------------------------------------------------------------------
225 * \brief Return true if set index is empty or corresponds to an existing set
227 template <class TSet> bool isValidIndex(const string& index, const map<int,TSet>& aMap)
229 if ( index.empty() ) return true;
230 return ( aMap.find( _ATOI( index ) ) != aMap.end() );
233 //--------------------------------------------------------------------------------
235 * \brief Return EnSight type corresponding to med one
237 const TEnSightElemType& getEnSightType(medGeometryElement medType)
239 static TEnSightElemType theEnSightType;
241 int nbNodes = medType % 100;
242 theEnSightType._medType = medType;
247 theEnSightType._name = "point";
248 theEnSightType._medIndex.resize(1,0);
252 theEnSightType._name = "bar2";
253 int conn [2] = {0,1};
254 theEnSightType._medIndex.assign( conn, conn+nbNodes );
258 theEnSightType._name = "bar3";
259 int conn [3] = {0,2,1};
260 theEnSightType._medIndex.assign( conn, conn+nbNodes );
264 theEnSightType._name = "tria3";
265 int conn [3] = {0,2,1};
266 theEnSightType._medIndex.assign( conn, conn+nbNodes );
270 theEnSightType._name = "quad4";
271 int conn [4] = {0,3,2,1};
272 theEnSightType._medIndex.assign( conn, conn+nbNodes );
276 theEnSightType._name = "tria6";
277 int conn [6] = {0,2,1,5,4,3};
278 theEnSightType._medIndex.assign( conn, conn+nbNodes );
282 theEnSightType._name = "quad8";
283 int conn [8] = {0,3,2,1,7,6,5,4};
284 theEnSightType._medIndex.assign( conn, conn+nbNodes );
288 theEnSightType._name = "tetra4";
289 int conn [4] = {0,1,3,2};
290 theEnSightType._medIndex.assign( conn, conn+nbNodes );
294 theEnSightType._name = "pyramid5";
295 int conn [5] = {0,3,2,1,4};
296 theEnSightType._medIndex.assign( conn, conn+nbNodes );
300 theEnSightType._name = "penta6";
301 int conn [6] = {0,2,1,3,5,4};
302 theEnSightType._medIndex.assign( conn, conn+nbNodes );
306 theEnSightType._name = "hexa8";
307 int conn [8] = {0,3,2,1,4,7,6,5};
308 theEnSightType._medIndex.assign( conn, conn+nbNodes );
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 );
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 );
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 );
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 );
336 theEnSightType._name = "nsided";
337 theEnSightType._medIndex.clear();
340 case MED_POLYHEDRA : {
341 theEnSightType._name = "nfaced";
342 theEnSightType._medIndex.clear();
346 theEnSightType._name = "";
347 theEnSightType._medIndex.clear();
350 return theEnSightType;
353 //--------------------------------------------------------------------------------
355 * \brief Return EnSight type having a given name
357 const TEnSightElemType& getEnSightType(const string& theTypeName)
359 string typeName = theTypeName;
360 if ( isGhostType( typeName ))
361 typeName = string( &typeName[2] ); // ghost type
363 static map<string, TEnSightElemType> name2Type;
365 map<string, TEnSightElemType>::iterator nameType = name2Type.find( typeName );
366 if ( nameType != name2Type.end() )
367 return nameType->second;
369 const list<medGeometryElement> & allMedTypes = MED_EN::meshEntities[MED_CELL];
370 list< medGeometryElement >::const_iterator medType = allMedTypes.begin();
371 for ( ; medType != allMedTypes.end(); ++medType )
373 const TEnSightElemType& enSightType = getEnSightType( *medType );
374 if ( enSightType._name == typeName )
375 return name2Type[ typeName ] = enSightType;
377 return getEnSightType(MED_ALL_ELEMENTS);
380 // =============================================================== _CaseFileDriver
382 * \brief Case file driver constructor
384 //================================================================================
386 _CaseFileDriver::_CaseFileDriver(const string & fileName,
387 const _CaseFileDriver_User* creator)
388 : _fileName( fileName ), _directory("."), _user( creator )
390 // Find out if the driver is blocked
391 _blocked = isToIgnore( creator );
392 if ( creator->getAccessMode() == MED_EN::RDONLY )
396 string::size_type sepPos = _fileName.rfind( FILE_SEPARATOR );
397 if ( sepPos != _fileName.npos ) {
398 _directory = _fileName.substr( 0, sepPos );
400 _format = getEnSightFormatForWriting();
404 // =============================================================== _CaseFileDriver
406 * \brief Case file driver destructor
408 //================================================================================
410 _CaseFileDriver::~_CaseFileDriver()
414 // to write case file again by same DRIVER
416 ((_CaseFileDriver_User*)_user)->_imed = 0;
420 #define READ_NEXT_LINE continue
421 #define RAISE_EXCEPTION break
423 //================================================================================
425 * \brief Read Case file
427 //================================================================================
429 void _CaseFileDriver::read() throw (MEDEXCEPTION)
434 _ASCIIFileReader reader( _fileName );
436 STRING badFile("Invalid Case file ");
437 badFile << _fileName << "\n";
439 list<_FileSet> fileSets;
440 list<_TimeSet> timeSets;
442 set<string> varNames; // to detect equal variable names
444 string section = "_";
446 while ( !reader.eof() )
448 string line = reader.getLine();
453 reader.split( line, line, comment, '#');
457 string key, value; // parts of a line splited by column
458 reader.split( line, key, value, ':');
462 // analyse not empty lines
463 switch ( section[0] ) {
465 // --------------------------------------------------------------------------------
468 string type, s, ts, fs;
469 reader.split(key,type,s);
470 int newVarIndex = _variables.empty() ? 1 : _variables.rbegin()->first + 1;
472 if ( type == "scalar" || type == "vector" || type == "tensor" )
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
485 int nbParts = reader.split( value, parts );
487 errorMsg << "invalid variable format:\n" << line;
490 if ( contains( "per measured node", s.c_str() )) {
491 //cout << "Skip not supported data type: " << key << endl;
494 list<string>::reverse_iterator p = parts.rbegin();
495 _Variable& var = _variables[newVarIndex];
497 var._fileNameOrData = *p++;
499 if ( nbParts == 3 ) {
500 var._timeSetNumber = *p;
502 else if ( nbParts == 4 ) {
503 var._fileSetNumber = *p++;
504 var._timeSetNumber = *p;
506 varNames.insert( var._name );
509 else if ( type == "constant" )
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 )) {
516 reader.split(value,s,value);
518 _Variable& var = _variables[newVarIndex];
521 var._fileNameOrData = value;
522 var._timeSetNumber = ts;
523 if ( var._name.empty() || var._fileNameOrData.empty() ) {
524 errorMsg << "invalid variable format:\n" << line;
527 varNames.insert( var._name );
530 else if ( type == "complex" )
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 )) {
539 reader.split(value,s,value);
540 if ( reader.isDigit( s )) {
542 reader.split(value,s,value);
546 int nbParts = reader.split( value, parts );
548 errorMsg << "invalid variable format:\n" << line;
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];
557 var._timeSetNumber = ts;
558 var._fileSetNumber = fs;
559 var._fileNameOrData = value;
560 varNames.insert( var._name );
566 // --------------------------------------------------------------------------------
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
574 // time set: ts [description]
575 // number of steps: ns
576 // filename numbers: fn
577 // time values: time_1 time_2 .... time_ns
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" )
587 int nb = _ATOI( value );
589 errorMsg << "Invalid time set number: " << value;
592 timeSets.push_back( _TimeSet() );
593 timeSets.back()._number = nb;
596 // ---------------------------------------------------------
597 if ( key == "number of steps" )
599 if ( !timeSet || !timeSet->_times.empty() ) {
600 errorMsg << "Unexpected command: " << key;
603 int nbSteps = _ATOI( value );
605 errorMsg << "invalid number of steps: " << value;
608 timeSet->_times.resize( nbSteps );
609 timeSet->_fileIndex.resize( nbSteps );
612 // ---------------------------------------------------------
613 if ( key == "filename start number" )
615 if ( !timeSet || timeSet->_fileIndex.empty() ) {
616 errorMsg << "Unexpected command: " << key;
619 if ( !reader.isDigit( value )) {
620 errorMsg << "invalid " << line;
623 timeSet->_fileIndex[0] = value;
626 // ---------------------------------------------------------
627 if ( key == "filename increment" ) {
628 int incr = _ATOI( value );
630 errorMsg << "invalid " << line;
634 timeSet->_fileIndex.empty() ||
635 timeSet->_fileIndex[0].empty() ) {
636 errorMsg << "Unexpected command: " << key;
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 );
645 // ---------------------------------------------------------
646 if ( key == "time values" )
648 if ( !timeSet || timeSet->_times.empty() ) {
649 errorMsg << "Unexpected command: " << key;
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();
660 nbTimes += reader.split( value, times );
661 for (t = times.begin(); i < nbTimes; ++i, ++t ) {
662 if ( ! reader.isDigit( *t, /*real=*/true ))
664 timeSet->_times[i] = *t;
667 if ( nbTimes != (int)timeSet->_times.size() ) {
668 errorMsg << "incorrect number of times in time set " << timeSet->_number;
671 for ( i = 1; i < nbTimes; ++i )
672 if ( _ATOF( timeSet->_times[ i ]) <= _ATOF(timeSet->_times[ i-1 ] ))
674 if ( i < nbTimes ) { // break from the previous loop occured
675 errorMsg << "time values are not in ascending order in time set " << timeSet->_number;
680 // ---------------------------------------------------------
681 if ( key == "filename numbers" )
683 if ( !timeSet || timeSet->_fileIndex.empty()) {
684 errorMsg << "Unexpected command: " << key;
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();
695 nb += reader.split( value, numbers );
696 list<string>::iterator n = numbers.begin();
697 for ( ; i < nb; ++i, ++n ) {
698 if ( ! reader.isDigit( *n ))
700 timeSet->_fileIndex[i] = *n;
703 if ( nb != nbFiles ) {
704 errorMsg << "incorrect number of " << key << " in time set " << timeSet->_number;
709 // ---------------------------------------------------------
710 if ( key == "filename numbers file" ||
711 key == "time values file" )
713 if ( !timeSet || timeSet->_fileIndex.empty()) {
714 errorMsg << "Unexpected command: " << key;
717 string fileName = _directory + FILE_SEPARATOR + value;
718 if ( !_user->canOpenFile( fileName, RDONLY )) {
719 errorMsg << "Can not open file " << fileName;
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;
732 timeSet->_fileIndex.assign( numbers.begin(), numbers.end() );
734 timeSet->_times.assign( numbers.begin(), numbers.end() );
740 if ( section[1] == 'I' ) {
741 // --------------------------------------------------------------------------------
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 );
749 errorMsg << "Invalid file set number: " << value;
752 fileSets.push_back( _FileSet() );
753 fileSets.back()._number = nb;
757 else if ( key == "filename index" ) {
758 if ( fileSets.empty() ) {
759 errorMsg << "'filename index' before 'file set'";
761 else if ( !value.empty() ) {
762 fileSets.back()._fileIndex.push_back( value );
766 else if ( key == "number of steps" ) {
767 if ( fileSets.empty() ) {
768 errorMsg << "'number of steps' before 'file set'";
770 else if ( value.empty() ) {
771 errorMsg << "number of steps omitted: " << line;
773 else if ( !reader.isDigit( value )) {
774 errorMsg << "invalid number of steps: " << value;
777 int n = _ATOI( value );
779 errorMsg << "invalid number of steps: " << value;
782 fileSets.back()._nbStepsInFile.push_back( n );
788 // --------------------------------------------------------------------------------
791 if ( key != "type" ) {
792 errorMsg << "unexpected info in section " << section << ":\n" << line;
796 if ( value == "ensight gold" ) {
797 _format = ENSIGHT_GOLD;
799 else if ( value == "ensight" ) {
803 errorMsg << "Unsupported EnSight format: " << value;
813 // --------------------------------------------------------------------------------
815 // model: [ts] [fs] filename [change_coords_only]
816 // measured: [ts] [fs] filename [change_coords_only]
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;
825 else if ( key == "model" ) {
827 reader.split( value, parts );
828 list<string>::reverse_iterator s = parts.rbegin();
829 for ( ; s != parts.rend(); ++s )
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;
838 _model._timeSetNumber = *s;
840 if ( _model._timeSetNumber.empty() && !_model._fileSetNumber.empty() )
841 swap( _model._timeSetNumber, _model._fileSetNumber );
842 if ( _model._fileName.empty() ) {
843 errorMsg << "invalid model: " << value;
850 case 'M': { // MATERIAL
852 reader.split(key,keyWord1,key);
853 if ( keyWord1 == "material" || keyWord1 == "species" )
857 } // end switch (section[0])
859 if ( !errorMsg.empty() ) {
860 throw MEDEXCEPTION(STRING("Invalid Case file ") << _fileName
861 << ":" << lineNb << "\n" << errorMsg );
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" ||
873 throw MEDEXCEPTION(STRING() << "Invalid format of Case file " << _fileName
874 << "\nwrong line: " << line );
877 if ( _model._fileName.empty() )
878 throw MEDEXCEPTION(badFile << "no supported geometry information found" );
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;
887 list<_TimeSet>::iterator ts = timeSets.begin();
888 for ( ; ts != timeSets.end(); ++ts )
889 _timeSets[ ts->_number ] = *ts;
891 // check validity of ts and fs
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 );
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 );
906 // check uniqueness of variable names
907 if ( varNames.size() != _variables.size() )
909 "Warning: there are different fields with equal names, you may have problems!" << endl;
910 // throw MEDEXCEPTION(badFile );
912 // As variable does not refer to time set if there is one step (issue 0020113),
913 // we try to restore the reference
915 for ( ivars = _variables.begin(); ivars != _variables.end(); ++ivars ) {
916 _Variable & var = ivars->second;
917 if ( var._timeSetNumber.empty() )
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 );
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 );
933 //================================================================================
935 * \brief return number of time steps of a model
937 //================================================================================
939 int _CaseFileDriver::getNbMeshes() const
941 if ( _blocked || checkWasRead())
943 if ( _model._timeSetNumber.empty() )
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();
953 //================================================================================
955 * \brief Sets all data necessary for meshDriver::read()
956 * \param meshIndex - time step index
957 * \param meshDriver - driver
959 //================================================================================
961 void _CaseFileDriver::setDataFileName(const int meshIndex,
962 ENSIGHT_MESH_RDONLY_DRIVER* meshDriver)
964 if ( _blocked || checkWasRead())
966 isToIgnore( meshDriver ); // remeber
968 meshDriver->_dataFileName = _directory + FILE_SEPARATOR + _model._fileName;
969 meshDriver->_indexInDataFile = fixWildCardName( meshIndex,
970 _model._timeSetNumber,
971 _model._fileSetNumber,
972 meshDriver->_dataFileName,
974 meshDriver->_isGoldFormat = ( _format == ENSIGHT_GOLD );
975 meshDriver->_transientMode = ( !_model._timeSetNumber.empty() );
976 meshDriver->_singleFileMode = ( !_fileSets.empty() );
977 meshDriver->_imedMapKey = STRING(_fileName)<<":"<<meshIndex;
979 GMESH* ptrMesh = meshDriver->getMesh();
980 ptrMesh->setName(STRING("EnSight mesh ") << meshIndex);
983 //================================================================================
985 * \brief Return nb of ensight variables
987 //================================================================================
989 int _CaseFileDriver::getNbVariables() const
991 return _variables.empty() ? 0 : _variables.rbegin()->first;
994 //================================================================================
996 * \brief return number of time steps of a variable
998 //================================================================================
1000 int _CaseFileDriver::getNbVarSteps(const int variableIndex)
1002 if ( _blocked || checkWasRead() )
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: ")
1014 const _Variable & var = ivar->second;
1015 if ( var._timeSetNumber.empty() )
1018 const _TimeSet & ts = _timeSets[ _ATOI( var._timeSetNumber)];
1019 return ts._times.size();
1022 //================================================================================
1024 * \brief return variable index by variable name, return 0 if none found
1026 //================================================================================
1028 int _CaseFileDriver::getVariableIndex(const string & varName) const
1030 if ( _blocked || checkWasRead() )
1033 map< int, _Variable>::const_iterator ivar = _variables.begin();
1034 for ( ; ivar != _variables.end(); ++ivar )
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
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' )
1048 int i = getVariableIndex( varName.substr( 0, _ImBeg ));
1049 return ( i ? i + 1 : i ); // getVariableIndex() returns index for a real part
1054 //================================================================================
1056 * \brief sets all data necessary for fieldDriver::read()
1058 //================================================================================
1060 int _CaseFileDriver::setDataFileName(const int varIndex,
1061 const int stepIndex,
1062 ENSIGHT_FIELD_RDONLY_DRIVER* fieldDriver)
1064 const char* LOC = "_CaseFileDriver::setDataFileName(): ";
1065 if ( _blocked || checkWasRead() )
1067 isToIgnore( fieldDriver ); // remeber
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));
1076 _Variable & var = ivar->second;
1078 const bool isComplex = contains( "complex", var._type.c_str() );
1079 const bool isImaginary = ( isComplex && varIndex == ivar->first );
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
1098 FIELD_* field = fieldDriver->getField();
1101 if ( field->getName().empty() ) {
1103 field->setName( var._name + "_Im" );
1105 field->setName( var._name );
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" )
1115 else if ( type1 == "tensor" )
1116 nbComponents = ( type2 == "symm" ) ? 6 : 9;
1117 field->setNumberOfComponents( nbComponents );
1120 vector<string> compNames( nbComponents );
1121 switch ( nbComponents ) {
1123 compNames[0] = type1;
1126 const char* xyz[3] = { "X comp","Y comp","Z comp" };
1127 compNames.assign(xyz, xyz+3);
1131 const char* xyz[6] = { "11 comp", "22 comp", "33 comp","12 comp", "13 comp", "23 comp" };
1132 compNames.assign(xyz, xyz+6);
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);
1141 field->setComponentsNames( & compNames[0] );
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]);
1151 SUPPORT* sup = const_cast<SUPPORT*>( field->getSupport());
1153 field->setSupport( sup = new SUPPORT );
1154 sup->removeReference(); // sup belongs to field
1156 medEntityMesh entity = ( type_parts.back() == "node" ) ? MED_NODE : MED_CELL;
1157 sup->setEntity( entity );
1159 // data file name etc.
1160 list<string> fileData_parts;
1161 _ASCIIFileReader::split( var._fileNameOrData, fileData_parts );
1163 fieldDriver->_dataFileName = _directory + FILE_SEPARATOR + *(++fileData_parts.begin());
1165 fieldDriver->_dataFileName = _directory + FILE_SEPARATOR + fileData_parts.front();
1166 fieldDriver->_indexInDataFile = fixWildCardName( stepIndex,
1169 fieldDriver->_dataFileName,
1170 fieldDriver->_time);
1171 fieldDriver->_isGoldFormat = ( _format == ENSIGHT_GOLD );
1172 fieldDriver->_transientMode = ( !var._timeSetNumber.empty() );
1173 fieldDriver->_singleFileMode = ( !_fileSets.empty() );
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() );
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 );
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 );
1200 fieldDriver->_imedMapKey = STRING(_fileName) << ":" << meshIndex;
1203 if ( fieldDriver->_transientMode ) {
1204 field->setTime( _ATOF( fieldDriver->_time ));
1205 //field->setOrderNumber( stepIndex );
1206 field->setIterationNumber( stepIndex );
1211 //================================================================================
1213 * \brief Throws if case file has not been read else return false
1215 //================================================================================
1217 bool _CaseFileDriver::checkWasRead() const throw (MEDEXCEPTION)
1219 if ( _model._fileName.empty() )
1220 throw MEDEXCEPTION(STRING("Case file ") << _fileName << " has not been successfully read");
1224 //================================================================================
1226 * \brief replace '*' in file name if any and return index in file
1228 //================================================================================
1230 int _CaseFileDriver::fixWildCardName(const int timeStep,
1236 int indexInFile = 0;
1239 STRING badFile("Invalid Case file ");
1240 badFile << _fileName << "\n";
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;
1248 list<int>::const_iterator nbStepsIt = fileSet._nbStepsInFile.begin();
1249 list<string>::const_iterator fIndexIt = fileSet._fileIndex.begin();
1251 for ( ; nbStepsIt != fileSet._nbStepsInFile.end(); ++nbStepsIt ) {
1252 if ( nbSteps + *nbStepsIt <= timeStep )
1254 nbSteps += *nbStepsIt;
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;
1265 _ASCIIFileReader::split( fileName, head, queue, '*' );
1266 int indexWidth = fileName.size() - head.size() - queue.size();
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 ];
1281 if ( indexWidth > 0 )
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 << ">"));
1288 if ( indexWidth == (int)fileIndex.size() ) {
1289 fileName = head + fileIndex + queue;
1292 fileName = (STRING(head) << setw(indexWidth) << setfill('0') << fileIndex << queue);
1298 //================================================================================
1300 * \brief add a mesh to the Case file
1302 //================================================================================
1304 void _CaseFileDriver::addMesh(const ENSIGHT_MESH_WRONLY_DRIVER* meshDriver)
1309 _meshDrivers.push_back( const_cast<ENSIGHT_MESH_WRONLY_DRIVER*>( meshDriver ));
1311 if ( _format == ENSIGHT_6 )
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 )
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 )"));
1324 isToIgnore( meshDriver ); // remeber
1327 //================================================================================
1329 * \brief add a field to the Case file
1331 //================================================================================
1333 void _CaseFileDriver::addField(const ENSIGHT_FIELD_WRONLY_DRIVER * theFieldDriver)
1338 ENSIGHT_FIELD_WRONLY_DRIVER * fieldDriver =
1339 const_cast<ENSIGHT_FIELD_WRONLY_DRIVER*>( theFieldDriver );
1341 FIELD_* field = fieldDriver->getField();
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";
1351 switch ( field->getNumberOfComponents() ) {
1355 case 9: break; // ok, supported
1357 if ( const GMESH* mesh = field->getSupport()->getMesh() )
1358 if ( mesh->getSpaceDimension() == 2 )
1359 break; // we add one component to both mesh and field
1362 compatibilityPb(STRING("it has ") << field->getNumberOfComponents()
1363 << " components but only 1,3,6 and 9 components are supported by EnSight");
1365 if ( !problem.empty() )
1366 throw MEDEXCEPTION(STRING("Can't write field <") << field->getName() <<
1367 "> to EnSight: " << problem);
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 );
1382 _fieldDrivers[ fieldName ].push_back( fieldDriver );
1384 isToIgnore( fieldDriver ); // remeber
1387 //================================================================================
1389 * \brief writing Case file
1391 //================================================================================
1393 void _CaseFileDriver::write() throw (MEDEXCEPTION)
1398 // Make case file data from added drivers
1400 bool dataIsRead = !_model._fileName.empty();
1401 if ( !dataIsRead && _meshDrivers.empty() )
1402 throw MEDEXCEPTION("no mesh to write into Case file");
1404 const int defaultNbDigits = 3; // in file index
1406 string path, name, ext;
1407 _ASCIIFileReader::split( _fileName, path, name, FILE_SEPARATOR, /*fromBack=*/true);
1408 _ASCIIFileReader::split( name, name, ext, '.', /*fromBack=*/true);
1410 name = ext; // _fileName: "path/.name"
1412 list<ENSIGHT_MESH_WRONLY_DRIVER*>::iterator mDrv = _meshDrivers.begin();
1413 TFieldDriversByName::iterator fDrv = _fieldDrivers.begin();
1415 map< string, string > fileToRenameTo;
1417 int i, nbOldMeshes = 0, nbNewMeshes = _meshDrivers.size();
1419 if ( nbNewMeshes > 0 )
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;
1429 _model._timeSetNumber,
1430 _model._fileSetNumber,
1433 int nbParts = ENSIGHT_MESH_RDONLY_DRIVER::countParts(realGeoFileName,
1434 !_model._fileSetNumber.empty());
1435 for ( ; mDrv != _meshDrivers.end(); ++mDrv )
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");
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 ];
1452 ts = _timeSets[ _ATOI( _model._timeSetNumber )];
1453 _model._timeSetNumber = STRING(tsNum);
1458 string::size_type digitPos = _model._fileName.find('*');
1459 bool isSeveralFiles = ( digitPos != _model._fileName.npos );
1460 int nbDigits = defaultNbDigits;
1461 if ( isSeveralFiles ) {
1463 while ( _model._fileName[ ++digitPos ] == '*' )
1466 // Update time set and file set of the model
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 ];
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);
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 );
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'));
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);
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 );
1513 if ( fs._fileIndex.empty() )
1514 fs._fileIndex.push_back(string(nbDigits,'0'));
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);
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');
1527 else if ( nbNewMeshes > 1 )
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 ];
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);
1544 // One mesh in a new file
1545 _model._fileName = name + ".geo";
1548 // Set data to mesh drivers
1549 i = nbOldMeshes + 1;
1550 for ( mDrv = _meshDrivers.begin(); mDrv != _meshDrivers.end(); ++mDrv, ++i )
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,
1559 meshDriver->_isGoldFormat = ( _format == ENSIGHT_GOLD );
1560 meshDriver->_transientMode = ( meshDriver->_indexInDataFile > 0 );
1561 meshDriver->_singleFileMode = ( !_fileSets.empty() );
1562 meshDriver->_imedMapKey = STRING(_fileName)<<":"<<i;
1566 typedef map< double, ENSIGHT_FIELD_WRONLY_DRIVER* > FDriverByDouble;
1568 //bool isVarsRead = ( !_variables.empty() );
1569 for ( ; fDrv != _fieldDrivers.end(); ++fDrv )
1571 const string & fieldName = fDrv->first;
1572 list< ENSIGHT_FIELD_WRONLY_DRIVER* > & drivers = fDrv->second;
1574 int nbNewSteps = drivers.size();
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 ));
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;
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 ));
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();
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;
1627 if ( field->getSupport()->getEntity() == MED_NODE )
1628 var._type += "per node";
1630 var._type += "per element";
1631 var._fileNameOrData = name + "." + fieldName;
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 ];
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);
1646 if ( nbNewSteps > 1 )
1647 var._fileNameOrData += string( nbDigits, '*' );
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 ) {
1655 var._timeSetNumber = STRING( tsNum );
1656 _timeSets.erase( tsEnd );
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,
1667 fieldDriver->_dataFileName,
1668 fieldDriver->_time);
1669 fieldDriver->_isGoldFormat = ( _format == ENSIGHT_GOLD );
1670 fieldDriver->_transientMode = ( fieldDriver->_indexInDataFile > 0 );
1671 fieldDriver->_singleFileMode = ( !_fileSets.empty() );
1674 // do not refer to time set if there is one step (issue 0020113)
1675 if ( nbNewSteps == 1 )
1676 var._timeSetNumber = "";
1678 } // loop on _fieldDrivers
1680 if ( nbOldMeshes + nbNewMeshes > 1 && !_variables.empty() )
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
1687 if ( toIgnoreIncompatibility() )
1688 cout << "Warning: EnSight file will be probably invalid " << _fileName << endl;
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")
1700 ofstream caseFile( _fileName.c_str(), ios::out );
1702 throw MEDEXCEPTION(STRING("Can't open file for writing ") << _fileName);
1704 caseFile << "# generated by MEDMEM-to-EnSight driver" << endl << endl
1706 << "type: " << (_format==ENSIGHT_GOLD ? "ensight gold" : "ensight") << endl
1711 << _model._timeSetNumber << " "
1712 << _model._fileSetNumber << " "
1713 << _model._fileName << " \t"
1714 << _model._change_coords_only
1717 // write VARIABLE section
1718 if ( !_variables.empty() )
1720 caseFile << "VARIABLE" << endl;
1722 map< int, _Variable>::iterator ivar = _variables.begin();
1723 for ( ; ivar != _variables.end(); ++ivar )
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;
1734 // write TIME section
1735 if ( !_timeSets.empty() )
1737 caseFile << "TIME" << endl;
1739 map< int, _TimeSet>::iterator its = _timeSets.begin();
1740 for ( ; its != _timeSets.end(); ++its )
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;
1752 numbers << " " << ts._fileIndex[i];
1754 caseFile << numbers << endl;
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;
1762 times << " " << ts._times[i];
1764 caseFile << times << endl;
1767 // write FILE section
1768 if ( !_fileSets.empty() )
1770 caseFile << "FILE" << endl;
1772 map< int, _FileSet >::iterator ifs = _fileSets.begin();
1773 for ( ; ifs != _fileSets.end(); ++ifs )
1775 _FileSet & fs = ifs->second;
1776 caseFile << "file set: " << fs._number << endl;
1778 list<int>::iterator nbSteps = fs._nbStepsInFile.begin();
1779 list<string>::iterator fileIndex = fs._fileIndex.begin();
1780 for ( ; nbSteps != fs._nbStepsInFile.end(); ++nbSteps )
1782 if ( fileIndex != fs._fileIndex.end() )
1783 caseFile << "filename index: " << *fileIndex++;
1784 caseFile << "number of steps: " << *nbSteps << endl;
1791 } // _CaseFileDriver::write()
1794 //================================================================================
1796 * \brief _CaseFileDriver_User constructor
1798 //================================================================================
1800 _CaseFileDriver_User::_CaseFileDriver_User(const string & caseFileName,
1801 med_mode_acces mode)
1802 : GENDRIVER( caseFileName, mode, ENSIGHT_DRIVER ), _imed(0)
1806 //================================================================================
1808 * \brief take data from other driver
1810 //================================================================================
1812 void _CaseFileDriver_User::merge( const GENDRIVER& driver)
1814 const _CaseFileDriver_User* other = dynamic_cast< const _CaseFileDriver_User* >( &driver );
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 ;
1827 //================================================================================
1829 * \brief analyse if data file is binary
1831 //================================================================================
1833 bool _CaseFileDriver_User::isBinaryDataFile(const string& dataFileName)
1836 int _file = ::_open (dataFileName.c_str(), _O_RDONLY|_O_BINARY);
1838 int _file = ::open (dataFileName.c_str(), O_RDONLY);
1841 int nBytesRead = ::read (_file, buf, 80);
1843 bool isBinary = true;
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 )
1850 for ( int i = nBytesRead-1; i >= 0 && isBinary; --i )
1851 isBinary = ( buf[ i ] != '\n' );
1859 //================================================================================
1861 * \brief return part number to write support with, zero in failure case
1863 //================================================================================
1865 int _CaseFileDriver_User::getPartNumber(const SUPPORT* support) const
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();
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 );
1882 else if ( entity == MED_NODE )
1883 return 1 + MED_CELL; // all nodes are described with all CELLs
1887 if ( isForField && isOnAll ) {
1888 if ( entity == MED_NODE )
1889 return 1 + MED_CELL; // all nodes are described with all CELLs
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;
1904 for ( int i=1; i<=nbGroups; ++i, ++partNum)
1905 if ( support == mesh->getGroup( entity, i ))
1910 ( LOCALIZED( STRING("Can't find GROUP ") << support->getName() << " in its MESH"));
1915 //================================================================================
1917 * \brief Return true if we strore an enetuty to EnSight
1919 //================================================================================
1921 bool _CaseFileDriver_User::isToWriteEntity(const medEntityMesh entity,
1924 if ( entity == MED_NODE )
1925 return mesh->getNumberOfNodes() > 0;
1927 if ( mesh->getNumberOfElements( entity, MED_ALL_ELEMENTS ) < 1 )
1929 if ( entity == MED_CELL )
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 );
1940 //================================================================================
1942 * \brief Return nodes of non-nodal support, which is not on all entities
1944 //================================================================================
1946 void _CaseFileDriver_User::getSupportNodes(const SUPPORT* support, map<int, int> & nodeIds)
1948 medEntityMesh entity = support->getEntity();
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;
1958 if ( support->isOnAllElements() )
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 ));
1968 const MESH* mesh = support->getMesh()->convertInMESH();
1970 connectivity = mesh->getConnectivity (conType, entity, allGeoms);
1971 conLength = mesh->getConnectivityLength(conType, entity, allGeoms);
1972 while ( conLength-- ) nodeIds[ *connectivity++ ];
1973 mesh->removeReference();
1978 if ( entity == MED_NODE )
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 ));
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 )
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++ ];
2000 mesh->removeReference();
2003 //================================================================================
2005 * \brief method called by mesh driver at reading to store data to be used by
2008 //================================================================================
2010 void _CaseFileDriver_User::setInterData(_InterMed* imed )
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;
2021 //================================================================================
2023 * \brief for field driver to get mesh data
2025 //================================================================================
2027 _InterMed* _CaseFileDriver_User::getInterData()
2029 return _imed = getMeshData( _imedMapKey );
2032 //================================================================================
2034 * \brief return _SubPart by its description
2036 //================================================================================
2038 _SubPart* _CaseFileDriver_User::getSubPart(const _SubPartDesc & descriptor)
2039 throw (MEDEXCEPTION)
2042 _imed = getMeshData( _imedMapKey );
2044 map< _SubPartDesc, _SubPart >::iterator descPart =
2045 _imed->_subPartDescribed.find( descriptor );
2047 if ( descPart != _imed->_subPartDescribed.end() )
2048 return & descPart->second;
2050 if ( descriptor == _SubPartDesc::globalCoordDesc() )
2051 return 0; // if the mesh is structured, there are no global coordinates (EnSight 6)
2053 throw MEDEXCEPTION(LOCALIZED(STRING("No mesh info for part ") << descriptor.partNumber()
2054 << " " << descriptor.typeName()));
2057 //================================================================================
2059 * \brief return _Support by its description
2061 //================================================================================
2063 _Support* _CaseFileDriver_User::getSupport(const _SupportDesc & descriptor,
2064 const medEntityMesh entity)
2065 throw (MEDEXCEPTION)
2067 const char* LOC = "_CaseFileDriver_User::getSupport(): ";
2069 _imed = getMeshData( _imedMapKey );
2071 _imed->treatGroupes();
2073 if ( _imed->_supportDescribed.empty() && !_imed->_subPartDescribed.empty() )
2075 // fill _supportDescribed with _Support's corresponding to EnSight parts
2076 for (unsigned int i=0; i < _imed->groupes.size(); ++i)
2078 _groupe& grp = _imed->groupes[i];
2079 if ( !grp.medGroup ) continue;
2081 vector<int> grpIndices(1,i);
2082 if ( !grp.groupes.empty() )
2083 grpIndices.assign( grp.groupes.begin(), grp.groupes.end());
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 )
2090 map< _SubPartDesc, _SubPart >::iterator descSub;
2091 for ( descSub = _imed->_subPartDescribed.begin();
2092 descSub != _imed->_subPartDescribed.end();
2095 if ( descSub->second.myCellGroupIndex == *grIndex ) {
2096 supDescriptor.insert( descSub->first );
2101 if ( supDescriptor.size() == grpIndices.size() ) {
2102 _Support & sup = _imed->_supportDescribed[ supDescriptor ];
2103 sup.setGroup( & grp );
2108 // find the support by its description
2109 map< _SupportDesc, _Support >::iterator descSup =
2110 _imed->_supportDescribed.find( descriptor );
2112 // Create a new support
2113 if ( descSup == _imed->_supportDescribed.end() || !descSup->second.medSupport( entity ))
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() );
2120 // to detect dimension change within a support group
2121 set<int> dimensions;
2123 // fill grp with sub-group indices
2124 _SupportDesc::const_iterator subPartDesc = descriptor.begin();
2125 for ( ; subPartDesc != descriptor.end(); ++subPartDesc )
2127 const _SubPart* subPart = getSubPart(*subPartDesc);
2129 throw MEDEXCEPTION(LOCALIZED(STRING("No mesh info for ") << *subPartDesc));
2132 (entity == MED_NODE) ? subPart->myNodeGroupIndex : subPart->myCellGroupIndex;
2133 if ( groupIndex < 1 ) {
2134 if ( entity == MED_NODE )
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 );
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 ) {
2147 ma.setOrdre( n->second.number );
2148 groupe.mailles[i] = _imed->insert(ma);
2153 throw MEDEXCEPTION(LOCALIZED(STRING("No cell info for ") << *subPartDesc));
2156 grp.groupes.push_back( groupIndex );
2158 // find out subpart dimension
2159 if ( entity != MED_NODE )
2160 dimensions.insert( _imed->groupes[ groupIndex-1 ].maille(0).dimensionWithPoly() );
2162 // check if MEDMEM allows creating such a support
2163 if ( dimensions.size() > 1 )
2165 (compatibilityPb(LOC) << "can't create a SUPPORT for the field from "
2166 << _dataFileName << ", since it includes different mesh entities");
2168 ENSIGHT_MESH_RDONLY_DRIVER::makeGroup( grp, *_imed );
2171 descSup = _imed->_supportDescribed.insert( make_pair( descriptor, _Support() )).first;
2172 _Support & sup = descSup->second;
2173 sup.setGroup( & grp );
2175 } // new SUPPORT creation
2177 _Support* sup = & descSup->second;
2179 // remove temporary mesh from med SUPPORT
2180 if ( _imed->_isOwnMedMesh )
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() );
2191 //================================================================================
2193 * \brief check possibility to open a file with a given mode
2195 //================================================================================
2197 bool _CaseFileDriver_User::canOpenFile(const string& fileName,
2198 med_mode_acces mode)
2201 if ( mode == WRONLY ) {
2202 fstream file( fileName.c_str(),
2203 ios::app | ios_base::out ); // not to overwrite it, just to check
2207 fstream file( fileName.c_str(), ios::in );
2213 //================================================================================
2215 * \brief _CaseFileDriver_User destructor
2217 //================================================================================
2219 _CaseFileDriver_User::~_CaseFileDriver_User()
2226 //================================================================================
2228 * \brief Search substring in a string
2230 //================================================================================
2232 bool contains( const char* what, const char* inString )
2234 size_t whatLen = strlen( what );
2235 size_t inLen = strlen( inString );
2237 ( search( inString, inString+inLen, what, what + whatLen) != inString+inLen );
2240 //================================================================================
2242 * \brief store subPart data
2244 //================================================================================
2246 void _InterMed::addSubPart(const _SubPart& theSubPart)
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();
2257 //================================================================================
2259 * \brief delete intermediate med data
2261 //================================================================================
2263 _InterMed::~_InterMed()
2265 if ( _isOwnMedMesh )
2267 // remove MEDMEM groups not belonging to _medMesh
2268 for (unsigned int i=0; i < _intermediateMED::groupes.size(); ++i)
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();
2276 if(_medMesh) _medMesh->removeReference();
2280 //================================================================================
2282 * \brief For a node, find its index in the supporting GROUP
2284 //================================================================================
2286 int _Support::getIndex( const pair<const int,_noeud>& inode)
2288 if ( myNodeGroup->relocMap.empty() ) // on all and not self intersecting support
2289 return abs( inode.second.number );
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);
2298 return ordreIndex->second;
2301 //================================================================================
2303 * \brief For a cell, return its index in the supporting GROUP
2305 //================================================================================
2307 int _Support::getIndex( const _groupe::TMaille& cell)
2309 if ( myCellGroup->relocMap.empty() ) // on all and not self intersecting support
2310 return cell->ordre();
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));
2316 return ordreIndex->second;
2319 //================================================================================
2321 * \brief Return med geom type for a subPart element
2323 //================================================================================
2325 // medGeometryElement _Support::getType( const pair<const int,_noeud>& node)
2327 // return myNodeGroup->medGroup.getTypes[0];
2330 // //================================================================================
2332 // * \brief Return med geom type for a subPart element
2334 // //================================================================================
2336 // medGeometryElement _Support::getType( const _groupe::TMaille& cell)
2338 // return cell->geometricType;
2341 //================================================================================
2343 * \brief Set med support
2345 //================================================================================
2347 void _Support::setGroup( _groupe* group )
2349 if ( group->medGroup ) {
2350 if ( group->medGroup->getEntity() == MED_NODE )
2351 myNodeGroup = group;
2353 myCellGroup = group;
2356 throw MEDEXCEPTION(LOCALIZED("_Support::setGroup(): med GROUP is NULL"));
2360 //================================================================================
2362 * \brief Return med group correspoding to entity
2364 //================================================================================
2366 SUPPORT* _Support::medSupport( medEntityMesh entity )
2368 if ( entity == MED_NODE )
2369 return myNodeGroup ? myNodeGroup->medGroup : 0;
2371 return myCellGroup ? myCellGroup->medGroup : 0;
2374 //================================================================================
2376 * \brief print _SubPartDesc
2378 //================================================================================
2380 std::ostream& operator << (std::ostream& os, const _SubPartDesc& desc)
2382 if ( desc == _SubPartDesc::globalCoordDesc() )
2383 os << "'global coordinates'";
2385 os << "<'part " << desc.partNumber() << "', '" << desc.typeName() << "'>";
2389 //================================================================================
2391 * \brief Constructor of ASCIIFileReader
2393 //================================================================================
2395 _ASCIIFileReader::_ASCIIFileReader(const string& fileName) throw (MEDEXCEPTION)
2398 _file = ::_open (fileName.c_str(), _O_RDONLY|_O_BINARY);
2400 _file = ::open (fileName.c_str(), O_RDONLY);
2404 _start = new char [BUFFER_SIZE];
2410 throw MEDEXCEPTION(STRING("Can't read from ")<<fileName);
2413 throw MEDEXCEPTION(STRING("Empty file ")<<fileName);
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');
2422 throw MEDEXCEPTION(STRING("Not ASCII file ")<<fileName);
2425 //================================================================================
2427 * \brief Return true if the whole file has been read
2429 //================================================================================
2431 bool _ASCIIFileReader::eof()
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)
2438 if (nBytesRest > 0) {
2439 char* tmpBuf = new char [nBytesRest];
2440 memcpy (tmpBuf, _ptr, nBytesRest);
2441 memcpy (_start, tmpBuf, nBytesRest);
2447 const int nBytesRead = ::read (_file,
2448 &_start [nBytesRest],
2449 BUFFER_SIZE - nBytesRest);
2450 nBytesRest += nBytesRead;
2451 _eptr = &_start [nBytesRest];
2453 // skip spaces at file end
2454 if ( nBytesRest < MAX_LINE_LENGTH ) {
2455 while ( isspace( *_ptr )) _ptr++;
2456 nBytesRest = _eptr - _ptr;
2459 return nBytesRest < 1;
2462 //================================================================================
2464 * \brief read out given data
2466 //================================================================================
2468 void _ASCIIFileReader::skip(int nbVals, int nbPerLine, int valWidth)
2470 int nbLines = (nbVals + nbPerLine - 1) / nbPerLine;
2471 int width = nbVals * valWidth;
2472 skip( width, nbLines );
2475 //================================================================================
2477 * \brief read out width chars and nbLines line-ends
2479 //================================================================================
2481 void _ASCIIFileReader::skip(int width, int nbLines)
2483 width += nbLines; // to skip '\n'
2485 width += nbLines; // to skip '\r'
2488 int nBytesRest = _eptr - _ptr;
2489 while ( nBytesRest < 0 ) {
2490 width = -nBytesRest;
2491 if ( eof() ) return;
2493 nBytesRest = _eptr - _ptr;
2497 //================================================================================
2499 * \brief Read a next line from the file
2501 //================================================================================
2503 char* _ASCIIFileReader::getLine() throw (MEDEXCEPTION)
2506 throw MEDEXCEPTION("Unexpected EOF");
2508 // Check the buffer for the end-of-line
2512 // Check for end-of-the-buffer, the ultimate criterion for termination
2519 // seek the line-feed character
2522 if ( ptr > _start && // avoid "Invalid read" error by valgrind
2531 // Output the result
2538 //================================================================================
2542 //================================================================================
2544 bool _ASCIIFileReader::lookAt( const char* text )
2546 while ( isspace(*_ptr)) ++_ptr;
2547 return ( strncmp( _ptr, text, strlen( text )) == 0 );
2550 //================================================================================
2552 * \brief Read a next word from the file
2554 //================================================================================
2556 string _ASCIIFileReader::getWord()
2562 while ( isspace(*_ptr)) ++_ptr;
2563 if ( _ptr >= _eptr )
2567 char* word = _ptr++;
2568 while ( !isspace(*_ptr)) ++_ptr;
2570 return string( word, _ptr - word );
2573 //================================================================================
2575 * \brief Return true if TIME_STEP_BEG follows and skips TIME_STEP_BEG line
2577 //================================================================================
2579 bool _ASCIIFileReader::isTimeStepBeginning()
2581 if ( eof() ) throw MEDEXCEPTION(LOCALIZED("Unexpected EOF"));
2583 while ( isspace(*_ptr)) ++_ptr;
2585 if ( strncmp( _ptr, TIME_STEP_BEG, TIME_STEP_BEG_LEN ) != 0 )
2588 _ptr += TIME_STEP_BEG_LEN;
2589 while ( isspace(*_ptr)) ++_ptr;
2593 //================================================================================
2595 * \brief Return true if TIME_STEP_END follows and skips TIME_STEP_END line
2597 //================================================================================
2599 bool _ASCIIFileReader::isTimeStepEnd()
2601 if ( eof() ) return true;
2603 while ( isspace(*_ptr)) ++_ptr;
2605 if ( strncmp( _ptr, TIME_STEP_END, TIME_STEP_END_LEN ) != 0 )
2608 _ptr += TIME_STEP_END_LEN;
2609 while ( isspace(*_ptr)) ++_ptr;
2614 //================================================================================
2616 * \brief divide a string into two parts
2618 //================================================================================
2620 int _ASCIIFileReader::split(const string& str,
2623 const char separator,
2624 const bool fromBack)
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;
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 ) {
2640 else if ( fromBack ) {
2642 string::size_type pos = str.rfind( separator );
2643 if ( pos != str.npos ) {
2644 ptr2 = & str[ pos ];
2646 if ( separator != ' ')
2647 while ( *nextBeg && *nextBeg == separator ) ++nextBeg;
2650 else if ( separator == ' ' ) {
2651 while ( *ptr2 && !isspace(*ptr2)) ++ptr2;
2652 if ( *ptr2 ) nextBeg = ptr2 + 1;
2655 while ( *ptr2 && *ptr2 != separator ) ++ptr2;
2658 while ( *nextBeg && *nextBeg == separator ) ++nextBeg;
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 );
2672 //================================================================================
2674 * \brief divide a string into parts, return nb of parts
2676 //================================================================================
2678 int _ASCIIFileReader::split(const string& str,
2679 list<string> & parts,
2680 const char separator,
2681 const bool fromBack)
2687 const char* ptr1 = str.c_str();
2688 const char* back = ptr1 + str.size();
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;
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 ));
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;
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 ));
2727 ptr1 = sepPtr + int( sepPtr < back );
2733 //================================================================================
2735 * \brief check if a string contains only digits
2737 //================================================================================
2739 bool _ASCIIFileReader::isDigit(const string& str, const bool real)
2741 const char* s = str.c_str();
2745 if ( !isdigit(c) && c!='-' && c!='+' && c!='.' && c!=',' && c!='E' && c!='e')
2751 if ( !isdigit( *s++ ))
2758 //================================================================================
2760 * \brief Destructor of ASCIIFileReader closes its file
2762 //================================================================================
2764 _ASCIIFileReader::~_ASCIIFileReader()
2773 //================================================================================
2775 * \brief Constructor of _BinaryFileReader opens its file
2777 //================================================================================
2779 _BinaryFileReader::_BinaryFileReader(const string& fileName) throw (MEDEXCEPTION)
2780 : _exception(STRING("Unexpected EOF ") << fileName), _mySwapBytes(false)
2783 _file = ::_open (fileName.c_str(), _O_RDONLY|_O_BINARY);
2785 _file = ::open (fileName.c_str(), O_RDONLY);
2789 throw MEDEXCEPTION(STRING("Can't read from ") << fileName);
2791 _maxPos = ::lseek( _file, 0, SEEK_END);
2792 _pos = ::lseek( _file, 0, SEEK_SET);
2795 //================================================================================
2797 * \brief Destructor of _BinaryFileWriter closes its file
2799 //================================================================================
2801 _BinaryFileReader::~_BinaryFileReader()
2807 //================================================================================
2809 * \brief rewind the file backward
2811 //================================================================================
2813 void _BinaryFileReader::rewind()
2815 _pos = ::lseek( _file, 0, SEEK_SET);
2818 //================================================================================
2820 * \brief size of not read file portion in sizeof(int)
2822 //================================================================================
2824 int _BinaryFileReader::moreValuesAvailable() const
2826 return (_maxPos - _pos) / sizeof(int); // believe that sizeof(int) == sizeof(float)
2829 //================================================================================
2831 * \brief Return true if the whole file has been read
2833 //================================================================================
2835 bool _BinaryFileReader::eof()
2837 return _pos >= _maxPos;
2840 //================================================================================
2842 * \brief Skips given nb of bytes
2844 //================================================================================
2846 void _BinaryFileReader::skip(int size) throw (MEDEXCEPTION)
2848 if ( _pos + size > _maxPos )
2850 off_t newPos = ::lseek( _file, size, SEEK_CUR);
2851 if ( newPos < _pos + size )
2856 //================================================================================
2858 * \brief Read lines until TIME_STEP_BEG encounters
2860 //================================================================================
2862 void _BinaryFileReader::skipTimeStepBeginning() throw (MEDEXCEPTION)
2864 bool tsbReached = false;
2865 while ( !tsbReached ) {
2866 TStrOwner line( getLine() );
2867 tsbReached = ( strncmp( line, TIME_STEP_BEG, TIME_STEP_BEG_LEN ) == 0 );
2871 //================================================================================
2873 * \brief Constructor of _BinaryFileWriter opens its file
2875 //================================================================================
2877 _BinaryFileWriter::_BinaryFileWriter(const string& fileName) throw (MEDEXCEPTION)
2878 : _exception(STRING("Can't write into ") << fileName)
2881 _file = ::_open (fileName.c_str(), _O_WRONLY|_O_BINARY|_O_TRUNC);
2883 _file = ::open (fileName.c_str(), O_WRONLY|O_TRUNC); //length shall be truncated to 0
2890 //================================================================================
2892 * \brief Destructor of _BinaryFileWriter closes its file
2894 //================================================================================
2896 _BinaryFileWriter::~_BinaryFileWriter()
2902 //================================================================================
2904 * \brief Write string
2906 //================================================================================
2908 void _BinaryFileWriter::addString(const char* str) throw (MEDEXCEPTION)
2910 size_t len = strlen( str );
2911 if ((int) len > MAX_LINE_LENGTH )
2913 (LOCALIZED(STRING("_BinaryFileWriter::addString(), too long string (>80):\n") << str));
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
2921 add( buffer.c_str(), MAX_LINE_LENGTH );
2926 } // end namespace MEDMEM_ENSIGHT