]> SALOME platform Git repositories - modules/med.git/blob - medtool/src/MEDLoader/SauvReader.cxx
Salome HOME
Move medtool folder to MED base repository
[modules/med.git] / medtool / src / MEDLoader / SauvReader.cxx
1 // Copyright (C) 2007-2015  CEA/DEN, EDF R&D
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, or (at your option) any later version.
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 // File      : SauvReader.cxx
20 // Created   : Tue Aug 16 13:57:42 2011
21 // Author    : Edward AGAPOV (eap)
22 //
23
24 #include "SauvReader.hxx"
25
26 #include "SauvMedConvertor.hxx"
27 #include "MEDCouplingAutoRefCountObjectPtr.hxx"
28 #include "NormalizedUnstructuredMesh.hxx"
29 #include "MEDCouplingRefCountObject.hxx"
30
31 #include <cstring>
32 #include <sstream>
33 #include <iostream>
34
35 using namespace ParaMEDMEM;
36 using namespace SauvUtilities;
37 using namespace std;
38
39 #define GIBI_EQUAL(var_str, stat_str) (strncmp (var_str, stat_str, strlen(stat_str)) == 0)
40
41 //================================================================================
42 /*!
43  * \brief Creates a reader of a given sauve file
44  */
45 //================================================================================
46
47 SauvReader* SauvReader::New(const std::string& fileName)
48 {
49   if ( fileName.empty() ) THROW_IK_EXCEPTION("Invalid file name");
50
51   ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr< SauvUtilities::FileReader> parser;
52
53   // try to open as XRD
54   parser = new XDRReader( fileName.c_str() );
55   if ( parser->open() )
56     {
57       SauvReader* reader = new SauvReader;
58       reader->_fileReader = parser.retn();
59       return reader;
60     }
61
62   // try to open as ASCII
63   parser = new ASCIIReader( fileName.c_str() );
64   if ( parser->open() )
65     {
66       SauvReader* reader = new SauvReader;
67       reader->_fileReader = parser.retn();
68       return reader;
69     }
70
71   THROW_IK_EXCEPTION("Unable to open file |"<< fileName << "|");
72 }
73 //================================================================================
74 /*!
75  * \brief Destructor
76  */
77 //================================================================================
78
79 SauvReader::~SauvReader()
80 {
81   _fileReader->decrRef();
82 }
83
84 std::size_t SauvReader::getHeapMemorySizeWithoutChildren() const
85 {
86   return 0;
87 }
88
89 std::vector<const BigMemoryObject *> SauvReader::getDirectChildrenWithNull() const
90 {
91   return std::vector<const BigMemoryObject *>();
92 }
93
94 //================================================================================
95 /*!
96  * \brief Return current line of ASCII file to report an error
97  */
98 //================================================================================
99
100 std::string SauvReader::lineNb() const
101 {
102   if ( isASCII() )
103     return string(" (line #") + SauvUtilities::toString
104       ( static_cast<SauvUtilities::ASCIIReader*>( _fileReader )->lineNb() ) + ")";
105
106   return "";
107 }
108
109 //================================================================================
110 /*!
111  * \brief Reads contents of the sauve file and convert it to MEDFileData
112  */
113 //================================================================================
114
115 ParaMEDMEM::MEDFileData * SauvReader::loadInMEDFileDS()
116 {
117   SauvUtilities::IntermediateMED iMed; // intermadiate DS
118   _iMed = &iMed;
119
120   char* line; // a current line
121   const char* enregistrement_type=" ENREGISTREMENT DE TYPE";
122
123   while ( getNextLine(line, /*raiseOEF=*/false)) // external loop looking for "ENREGISTREMENT DE TYPE"
124     {
125       if ( isASCII() && !GIBI_EQUAL( line, enregistrement_type ))
126         continue; // "ENREGISTREMENT DE TYPE" not found -> read the next line
127
128       // read the number of a record
129       int recordNumber;
130       if ( isASCII() )
131         recordNumber = atoi( line + strlen(enregistrement_type) + 1 );
132       else
133         recordNumber = getInt();
134
135       // read the record
136       if ( recordNumber == 2 )
137         readRecord2();
138       else if (recordNumber == 4 )
139         readRecord4();
140       else if (recordNumber == 7 )
141         readRecord7();
142       else if (recordNumber == 5 )
143         break; // stop reading
144       else
145         if ( !isASCII() )
146           THROW_IK_EXCEPTION("XDR : ENREGISTREMENT DE TYPE " << recordNumber << " not implemented!!!");
147     }
148
149   ParaMEDMEM::MEDFileData* medFileData = iMed.convertInMEDFileDS();
150
151   return medFileData;
152 }
153
154 //================================================================================
155 /*!
156  * \brief Reads "ENREGISTREMENT DE TYPE 4"
157  */
158 //================================================================================
159
160 void SauvReader::readRecord4()
161 {
162   if ( !isASCII() )
163     {
164       getInt(); // skip NIVEAU
165       getInt(); // skip ERREUR
166       _iMed->_spaceDim = getInt();
167       getFloat(); // skip DENSITE
168     }
169   else
170     {
171       char* line;
172       getNextLine(line);
173       const char* s = " NIVEAU  15 NIVEAU ERREUR   0 DIMENSION";
174       _iMed->_spaceDim = atoi( line + strlen( s ) + 1 );
175       if ( !GIBI_EQUAL( line, " NIVEAU" ))
176         THROW_IK_EXCEPTION( "Could not read space dimension" << lineNb() );
177       }
178   if ( _iMed->_spaceDim < 1 )
179     THROW_IK_EXCEPTION( "Invalid space dimension:" << _iMed->_spaceDim );
180 }
181
182 //================================================================================
183 /*!
184  * \brief Reads "ENREGISTREMENT DE TYPE 7"
185  */
186 //================================================================================
187
188 void SauvReader::readRecord7()
189 {
190   if ( !isASCII() )
191     {
192       getInt(); // skip NOMBRE INFO CASTEM2000
193       getInt(); // skip IFOUR
194       getInt(); // skip NIFOUR
195       getInt(); // skip IFOMOD
196       getInt(); // skip IECHO
197       getInt(); // skip IIMPI
198       getInt(); // skip IOSPI
199       getInt(); // skip ISOTYP
200       getInt(); // skip NSDPGE
201     }
202   else
203     {
204       // skip 3 lines:
205       // NOMBRE INFO CASTEM2000   8
206       // IFOUR   2 NIFOUR   0 IFOMOD   2 IECHO   1 IIMPI   0 IOSPI   0 ISOTYP   1
207       // NSDPGE     0
208       char* line;
209       getNextLine(line);
210       getNextLine(line);
211       getNextLine(line);
212     }
213 }
214
215 //================================================================================
216 /*!
217  * \brief Reads the pile number, nb of objects and nb named of objects
218  */
219 //================================================================================
220
221 int SauvReader::readPileNumber(int& nbNamedObjects, int& nbObjects)
222 {
223   // FORMAT(' PILE NUMERO',I4,'NBRE ObjectS NOMMES',I8,'NBRE ObjectS',I8)
224   int pileNumber;
225   if ( !isASCII() )
226     {
227       initIntReading(3);
228       pileNumber     = getInt(); next();
229       nbNamedObjects = getInt(); next();
230       nbObjects      = getInt(); next();
231     }
232   else
233     {
234       char* line;
235       getNextLine(line);
236       const char *s1 = " PILE NUMERO", *s2 = "NBRE ObjectS NOMMES", *s3 = "NBRE ObjectS";
237       if ( ! GIBI_EQUAL( line, s1 ) )
238         THROW_IK_EXCEPTION("Could not read the pile number " << lineNb() );
239       line           = line + strlen(s1);
240       pileNumber     = atoi( line );
241       line           = line + 4 + strlen(s2);
242       nbNamedObjects = atoi( line );
243       line           = line + 8 + strlen(s3);
244       nbObjects      = atoi( line );
245     }
246   if ( nbNamedObjects<0 )
247     THROW_IK_EXCEPTION("Invalid nb of named objects: " << nbNamedObjects  << lineNb() );
248   if ( nbObjects<0)
249     THROW_IK_EXCEPTION("Invalid nb of objects: " << nbObjects  << lineNb() );
250   // It appears to be a valid case
251   // if ( nbObjects<nbNamedObjects)
252   //   THROW_IK_EXCEPTION("In PILE " << pileNumber <<
253   //                      " nb of objects is less than nb of named objects"  << lineNb() );
254   return pileNumber;
255 }
256
257 //================================================================================
258 /*!
259  * \brief Reads "ENREGISTREMENT DE TYPE 2"
260  */
261 //================================================================================
262
263 void SauvReader::readRecord2()
264 {
265   if ( _iMed->_spaceDim == 0 )
266     THROW_IK_EXCEPTION("Missing ENREGISTREMENT DE TYPE   4");
267
268   // read a pile number
269   int pileNumber, nbNamedObjects, nbObjects;
270   pileNumber = readPileNumber(nbNamedObjects, nbObjects);
271
272   if ( !_encounteredPiles.insert( pileNumber ).second && // piles may repeat
273        isASCII())
274     return;
275
276   // read object names and their indices
277   vector<string> objectNames(nbNamedObjects);
278   for ( initNameReading( nbNamedObjects ); more(); next() )
279     objectNames[ index() ] = getName();
280
281   vector<int> nameIndices(nbNamedObjects);
282   for ( initIntReading( nbNamedObjects ); more(); next() )
283     nameIndices[ index() ] = getInt();
284
285   switch ( pileNumber )
286     {
287     case PILE_SOUS_MAILLAGE:
288       read_PILE_SOUS_MAILLAGE(nbObjects, objectNames, nameIndices);
289       break;
290     case PILE_NODES_FIELD:
291       read_PILE_NODES_FIELD(nbObjects, objectNames, nameIndices);
292       break;
293     case PILE_TABLES:
294       read_PILE_TABLES(nbObjects, objectNames, nameIndices);
295       break;
296     case PILE_LREEL:
297       read_PILE_LREEL(nbObjects, objectNames, nameIndices);
298       break;
299     case PILE_LOGIQUES:
300       read_PILE_LOGIQUES(nbObjects, objectNames, nameIndices);
301       break;
302     case PILE_FLOATS:
303       read_PILE_FLOATS(nbObjects, objectNames, nameIndices);
304       break;
305     case PILE_INTEGERS:
306       read_PILE_INTEGERS(nbObjects, objectNames, nameIndices);
307       break;
308     case PILE_STRINGS:
309       read_PILE_STRINGS(nbObjects, objectNames, nameIndices);
310       break;
311     case PILE_LMOTS:
312       read_PILE_LMOTS(nbObjects, objectNames, nameIndices);
313       break;
314     case PILE_NOEUDS:
315       read_PILE_NOEUDS(nbObjects, objectNames, nameIndices);
316       break;
317     case PILE_COORDONNEES:
318       read_PILE_COORDONNEES(nbObjects, objectNames, nameIndices);
319       break;
320     case PILE_MODL:
321       read_PILE_MODL(nbObjects, objectNames, nameIndices);
322       break;
323     case PILE_FIELD:
324       read_PILE_FIELD(nbObjects, objectNames, nameIndices);
325       break;
326     default:
327       if ( !isASCII() )
328         THROW_IK_EXCEPTION("XDR : reading PILE " << pileNumber << " not implemented !!!");
329     }
330 }
331
332 //================================================================================
333 /*!
334  * \brief Reads "PILE NUMERO   1": gibi sub-meshes that are converted into med groups
335  */
336 //================================================================================
337
338 void SauvReader::read_PILE_SOUS_MAILLAGE(const int                 nbObjects,
339                                          std::vector<std::string>& objectNames,
340                                          std::vector<int>&         nameIndices)
341 {
342   _iMed->_groups.reserve(nbObjects*2); // fields may add some groups
343
344   char* line;
345   map<int,int> strangeGroupType;
346   int i;
347
348   for (int object=0; object!=nbObjects; ++object) // loop on sub-groups
349     {
350       initIntReading( 5 );
351       int castemCellType = getIntNext();
352       int nbSubGroups    = getIntNext();
353       int nbReferences   = getIntNext();
354       int nbNodesPerElem = getIntNext();
355       int nbElements     = getIntNext();
356
357       _iMed->_groups.push_back(Group());
358       SauvUtilities::Group & group = _iMed->_groups.back();
359
360       // Issue 0021311. Allocate places for names of referring groups
361       // that will be possibly filled after reading long names from
362       // PILE_TABLES and PILE_STRINGS
363       group._refNames.resize( nbReferences );
364
365       // castemCellType=0 corresponds to a sub-mesh composed of other sub-meshes
366       if (castemCellType==0 && nbSubGroups>0)
367         {
368           group._groups.resize( nbSubGroups );
369           for ( initIntReading( nbSubGroups ); more(); next() )
370             group._groups[ index() ] = & _iMed->_groups[ getInt() - 1 ];
371           //std::sort( group._groups.begin(), group._groups.end() ); // for _groups comparison in getFieldSupport()
372         }
373       // skip references
374       if ( isASCII() )
375         for ( i = 0; i < nbReferences; i += 10 ) // FORMAT(10I8)
376           getNextLine(line);
377       else
378         for (initIntReading(nbReferences); more(); next());
379
380       // skip colors
381       if ( isASCII() )
382         for ( i = 0; i < nbElements; i += 10 )
383           getNextLine(line);
384       else
385         for (initIntReading(nbElements); more(); next());
386
387       // not a composite group
388       if (castemCellType>0 && nbSubGroups==0)
389         {
390           group._cellType = SauvUtilities::gibi2medGeom(castemCellType);
391
392           initIntReading( nbElements * nbNodesPerElem );
393           if ( group._cellType == INTERP_KERNEL::NORM_ERROR ) // look for group end
394             {
395               for ( ; more();  next());
396               strangeGroupType.insert( make_pair( object, castemCellType ));
397             }
398           else
399             {
400               // if ( group._cellType == MED_POINT1 ) group._cellType = NORM_ERROR; // issue 21199
401
402               // read connectivity of elements of a group
403               SauvUtilities::Cell ma( nbNodesPerElem );
404               SauvUtilities::Node* pNode;
405               group._cells.resize( nbElements );
406               for ( i = 0; i < nbElements; ++i )
407                 {
408                   ma.init();
409                   for ( int n = 0; n < nbNodesPerElem; ++n )
410                     {
411                       int nodeID = getIntNext();
412                       pNode = _iMed->getNode( nodeID );
413                       ma._nodes[n] = pNode;
414                       _iMed->_nbNodes += ( !pNode->isUsed() );
415                       pNode->_number = nodeID;
416                     }
417                   ma._number = _iMed->getNbCellsOfType( group._cellType ) + 1;
418                   group._cells[i] = _iMed->insert( group._cellType, ma );
419                 }
420             }
421         }
422     } // loop on groups
423
424   // set group names
425   for (i=0; i!=(int)objectNames.size(); ++i)
426     {
427       int grpID = nameIndices[i];
428       SauvUtilities::Group & grp = _iMed->_groups[ grpID-1 ];
429       if ( !grp._name.empty() ) // a group has several names
430         { // create a group with subgroup grp and named grp.name
431           SauvUtilities::Group* newGroup = _iMed->addNewGroup();
432           newGroup->_groups.push_back( &_iMed->_groups[ grpID-1 ]);
433           newGroup->_name = grp._name;
434         }
435       grp._name=objectNames[i];
436 #ifdef _DEBUG
437       map<int,int>::iterator it = strangeGroupType.find( grpID - 1 );
438       if ( it != strangeGroupType.end() )
439         cout << "Skip " << grp._name << " of not supported CASTEM type: " << it->second << endl;
440 #endif
441     }
442 } // read_PILE_SOUS_MAILLAGE()
443
444 //================================================================================
445 /*!
446  * \brief Skip "PILE NUMERO  18" of XDR file
447  */
448 //================================================================================
449
450 void SauvReader::read_PILE_LREEL (const int nbObjects, std::vector<std::string>&, std::vector<int>&)
451 {
452   if ( isXRD() )
453     {
454       for (int object=0; object!=nbObjects; ++object) // pour chaque Group
455         {
456           initIntReading(1);
457           int nb_vals = getIntNext();
458           initDoubleReading(nb_vals);
459           for(int i=0; i<nb_vals; i++) next();
460         }
461     }
462 }
463
464 //================================================================================
465 /*!
466  * \brief Skip "PILE NUMERO  24" of XDR file
467  */
468 //================================================================================
469
470 void SauvReader::read_PILE_LOGIQUES (const int, std::vector<std::string>&, std::vector<int>&)
471 {
472   if ( isXRD() )
473     {
474       initIntReading(1);
475       int nb_vals = getIntNext();
476       initIntReading(nb_vals);
477       for(int i=0; i<nb_vals; i++) next();
478     }
479 }
480
481 //================================================================================
482 /*!
483  * \brief Skip "PILE NUMERO  25" of XDR file
484  */
485 //================================================================================
486
487 void SauvReader::read_PILE_FLOATS (const int, std::vector<std::string>&, std::vector<int>&)
488 {
489   if ( isXRD() )
490     {
491       initIntReading(1);
492       int nb_vals = getIntNext();
493       initDoubleReading(nb_vals);
494       for(int i=0; i<nb_vals; i++) next();
495     }
496 }
497
498 //================================================================================
499 /*!
500  * \brief Skip "PILE NUMERO  26" of XDR file
501  */
502 //================================================================================
503
504 void SauvReader::read_PILE_INTEGERS (const int, std::vector<std::string>&, std::vector<int>&)
505 {
506   if ( isXRD() )
507     {
508       initIntReading(1);
509       int nb_vals = getIntNext();
510       initIntReading(nb_vals);
511       for(int i=0; i<nb_vals; i++) next();
512     }
513 }
514
515 //================================================================================
516 /*!
517  * \brief Skip "PILE NUMERO  29" of XDR file
518  */
519 //================================================================================
520
521 void SauvReader::read_PILE_LMOTS (const int nbObjects, std::vector<std::string>&, std::vector<int>&)
522 {
523   if ( isXRD() )
524     {
525       for (int object=0; object!=nbObjects; ++object) // pour chaque Group
526         {
527           initIntReading(2);
528           int len = getIntNext();
529           int nb_vals = getIntNext();
530           int nb_char = len*nb_vals;
531           int nb_char_tmp = 0;
532           int fixed_length = 71;
533           while (nb_char_tmp < nb_char)
534             {
535               int remain_len = nb_char - nb_char_tmp;
536               int width;
537               if ( remain_len > fixed_length )
538                 {
539                   width = fixed_length;
540                 }
541               else
542                 {
543                   width = remain_len;
544                 }
545               initNameReading(1, width);
546               next();
547               nb_char_tmp += width;
548             }
549         }
550     }
551 }
552
553 //================================================================================
554 /*!
555  * \brief Skip "PILE NUMERO  38" of XDR file
556  */
557 //================================================================================
558
559 void SauvReader::read_PILE_MODL (const int nbObjects, std::vector<std::string>&, std::vector<int>&)
560 {
561   if ( isXRD() )
562     {
563       for (int object=0; object!=nbObjects; ++object) // pour chaque Group
564         {
565           // see wrmodl.eso
566           initIntReading(10);
567           int n1  = getIntNext();
568           int nm2 = getIntNext();
569           int nm3 = getIntNext();
570           int nm4 = getIntNext();
571           int nm5 = getIntNext();
572           int n45 = getIntNext();
573           /*int nm6 =*/ getIntNext();
574           /*int nm7 =*/ getIntNext();
575           next();
576           next();
577           int nm1 = n1 * n45;
578           int nm9 = n1 * 16;
579           for (initIntReading(nm1); more(); next());
580           for (initIntReading(nm9); more(); next());
581           for (initNameReading(nm5, 8); more(); next());
582           for (initNameReading(nm2, 8); more(); next());
583           for (initNameReading(nm3, 8); more(); next());
584           for (initIntReading(nm4); more(); next());
585         }
586     }
587 } // Fin case pile 38
588
589 //================================================================================
590 /*!
591  * \brief Read "PILE NUMERO  32": links to node coordinates
592  */
593 //================================================================================
594
595 void SauvReader::read_PILE_NOEUDS (const int nbObjects, std::vector<std::string>&, std::vector<int>&)
596 {
597   initIntReading(1);
598   int nb_indices = getIntNext();
599
600   if (nb_indices != nbObjects)
601     THROW_IK_EXCEPTION("Error of reading PILE NUMERO  " << PILE_NOEUDS << lineNb() );
602
603   for ( initIntReading( nbObjects ); more(); next() )
604     {
605       int coordID = getInt();
606       _iMed->getNode( index()+1 )->_coordID = coordID;
607     }
608 }
609
610 //================================================================================
611 /*!
612  * \brief Read "PILE NUMERO  33": node coordinates
613  */
614 //================================================================================
615
616 void SauvReader::read_PILE_COORDONNEES (const int nbObjects, std::vector<std::string>&, std::vector<int>&)
617 {
618   initIntReading(1);
619   int nbReals = getIntNext();
620
621   if ( nbReals < (int)(_iMed->_nbNodes*(_iMed->_spaceDim+1)) )
622     THROW_IK_EXCEPTION("Error of reading PILE NUMERO  " << PILE_COORDONNEES << lineNb() );
623
624   // there are coordinates + density for each node
625   _iMed->_coords.resize( nbReals - nbReals/(_iMed->_spaceDim+1));
626   double* coordPtr = &_iMed->_coords[0];
627
628   initDoubleReading( nbReals );
629   while ( more() )
630     {
631       for (unsigned j = 0; j < _iMed->_spaceDim; ++j, next())
632         *coordPtr++ = getDouble();
633       // skip density
634       getDouble();
635       next();
636     }
637 }
638
639 //================================================================================
640 /*!
641  * \brief Find or create a Group equal to a given field support
642  */
643 //================================================================================
644
645 void SauvReader::setFieldSupport(const vector<SauvUtilities::Group*>& supports,
646                                  SauvUtilities::DoubleField*          field)
647 {
648   SauvUtilities::Group* group = NULL;
649   set<SauvUtilities::Group*> sup_set( supports.begin(), supports.end() );
650   if ( sup_set.size() == 1 ) // one or equal supports
651     {
652       group = supports[0];
653     }
654   else
655     {
656       // check if sub-components are on cells of different types
657       map<int,int> nbGaussByCellType;
658       for ( size_t i = 0; i < supports.size(); ++i )
659         {
660           map<int,int>::iterator ct2ng = nbGaussByCellType.find( supports[i]->_cellType );
661           if ( ct2ng == nbGaussByCellType.end() )
662             nbGaussByCellType[ supports[i]->_cellType ] = field->_sub[i].nbGauss();
663           else if ( ct2ng->second != field->_sub[i].nbGauss() )
664             return;
665         }
666       bool isSameCellType = ( nbGaussByCellType.size() == 1 );
667       // try to find an existing composite group with the same sub-groups
668       if ( isSameCellType )
669         for ( size_t i = 0; i < _iMed->_groups.size() && !group; ++i )
670           {
671             Group & grp = _iMed->_groups[i];
672             if (sup_set.size() == grp._groups.size())
673               {
674                 bool sameOrder = true;
675                 for ( size_t j = 0; j < supports.size() && sameOrder; ++j )
676                   sameOrder = ( supports[j] == grp._groups[ j % grp._groups.size() ]);
677                 if ( sameOrder )
678                   group = & _iMed->_groups[i];
679               }
680           }
681       if ( !group ) // no such a group, add a new one
682         {
683           vector<SauvUtilities::Group*> newGroups( supports.begin(),
684                                                    supports.begin() + sup_set.size() );
685           // check if supports includes newGroups in the same order
686           bool sameOrder = true;
687           for ( size_t j = newGroups.size(); j < supports.size() && sameOrder; ++j )
688             sameOrder = ( supports[j] == newGroups[ j % newGroups.size() ]);
689           if ( sameOrder )
690             {
691               group = _iMed->addNewGroup( & newGroups );
692               group->_groups.swap( newGroups );
693             }
694         }
695       // sort field sub-components and supports by cell type
696       if ( group && !isSameCellType )
697         {
698           // sort groups
699           vector<SauvUtilities::Group*>& groups = group->_groups;
700           bool isModified = false, isSwapped = true;
701           while ( isSwapped )
702             {
703               isSwapped = false;
704               for ( size_t i = 1; i < groups.size(); ++i )
705                 {
706                   int nbN1 = groups[i-1]->empty() ? 0 : groups[i-1]->_cells[0]->_nodes.size();
707                   int nbN2 = groups[i  ]->empty() ? 0 : groups[i  ]->_cells[0]->_nodes.size();
708                   if ( nbN1 > nbN2 )
709                     {
710                       isSwapped = isModified = true;
711                       std::swap( groups[i], groups[i-1] );
712                     }
713                 }
714             }
715           // relocate sub-components according to a new order of groups
716           if ( isModified )
717             {
718               vector< DoubleField::_Sub_data > newSub   ( field->_sub.size() );
719               vector< vector< double > >       newValues( field->_comp_values.size() );
720               size_t iFromSub = 0, iNewSub = 0, iNewComp = 0;
721               for ( ; iFromSub < field->_sub.size(); iFromSub += groups.size() )
722                 {
723                   size_t iFromComp = iNewComp;
724                   for ( size_t iG = 0; iG < groups.size(); ++iG )
725                     {
726                       size_t iComp = iFromComp;
727                       for ( size_t iSub = iFromSub; iSub < field->_sub.size(); ++iSub )
728                         if ( field->_sub[ iSub ]._support == groups[ iG ] )
729                           {
730                             newSub[ iNewSub++ ] = field->_sub[ iSub ];
731                             int iC = 0, nbC = field->_sub[ iSub ].nbComponents();
732                             for ( ; iC < nbC; ++iC )
733                               newValues[ iNewComp++ ].swap( field->_comp_values[ iComp++ ]);
734                             break;
735                           }
736                         else
737                           {
738                             iComp += field->_sub[ iSub ].nbComponents();
739                           }
740                     }
741                 }
742               field->_sub.swap( newSub );
743               field->_comp_values.swap( newValues );
744             }
745         }
746     }
747   if ( group )
748     group->_isProfile = true;
749
750   field->_group = group;
751 }
752
753 //================================================================================
754 /*!
755  * \brief Set field names
756  */
757 //================================================================================
758
759 void SauvReader::setFieldNames(const vector<SauvUtilities::DoubleField* >& fields,
760                                const vector<string>&                       objets_nommes,
761                                const vector<int>&                          indices_objets_nommes)
762 {
763   unsigned i;
764   for ( i = 0; i < indices_objets_nommes.size(); ++i )
765     {
766       int fieldIndex = indices_objets_nommes[ i ];
767       if ( fields[ fieldIndex - 1 ] )
768         fields[ fieldIndex - 1 ]->_name = objets_nommes[ i ];
769     }
770 }
771
772 //================================================================================
773 /*!
774  * \brief Read "PILE NUMERO   2": NODE FIELDS
775  */
776 //================================================================================
777
778 void SauvReader::read_PILE_NODES_FIELD (const int                 nbObjects,
779                                         std::vector<std::string>& objectNames,
780                                         std::vector<int>&         nameIndices)
781 {
782   _iMed->_nodeFields.resize( nbObjects, (SauvUtilities::DoubleField*) 0 );
783   for (int object=0; object!=nbObjects; ++object) // loop on fields
784     {
785       // EXAMPLE ( with no values )
786
787       // (1)       4       7       2       1
788       // (2)     -88       0       3     -89       0       1     -90       0       2     -91
789       // (2)       0       1
790       // (3) FX   FY   FZ   FZ   FX   FY   FLX
791       // (4)       0       0       0       0       0       0       0
792       // (5)           cree  par  muc pri
793       // (6)
794       // (7)       2
795
796       // (1): nb subcomponents, nb components(total), IFOUR, nb attributes
797       int nb_sub, total_nb_comp, nb_attr;
798       int i_sub, i_comp;
799       initIntReading( 4 );
800       nb_sub        = getIntNext();
801       total_nb_comp = getIntNext();
802       next(); // ignore IFOUR
803       nb_attr       = getIntNext();
804       if ( nb_sub < 0 || total_nb_comp < 0 || nb_attr < 0 )
805         THROW_IK_EXCEPTION("Error of field reading " << lineNb());
806
807       // (2) loop on subcomponents of a field, for each read
808       // (a) support, (b) number of values and (c) number of components
809       vector<Group*> supports( nb_sub );
810       vector<int> nb_values  ( nb_sub );
811       vector<int> nb_comps   ( nb_sub );
812       int total_nb_values = 0;
813       initIntReading( nb_sub * 3 );
814       for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
815         {
816           int supId = -getIntNext(); // (a) reference to support
817           if ( supId < 1 || supId > (int)_iMed->_groups.size() )
818             THROW_IK_EXCEPTION("Wrong mesh reference: "<< supId << lineNb() );
819           supports[ i_sub ] = &_iMed->_groups[ supId-1 ]; // (a) reference to support
820
821           nb_values[ i_sub ] = getIntNext();    // (b) nb points
822           total_nb_values += nb_values[ i_sub ];
823           if ( nb_values[ i_sub ] < 0 )
824             THROW_IK_EXCEPTION(" Wrong nb of points: " << nb_values[ i_sub ]  << lineNb() );
825           nb_comps[ i_sub ] = getInt(); next();     // (c) nb of components in i_sub
826         }
827
828       // create a field if there are values
829       SauvUtilities::DoubleField* fdouble = 0;
830       if ( total_nb_values > 0 )
831         fdouble = new DoubleField( nb_sub, total_nb_comp );
832       _iMed->_nodeFields[ object ] = fdouble;
833
834       // (3) component names
835       initNameReading( total_nb_comp, 4 );
836       for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
837         {
838           // store support id and nb components of a sub
839           if ( fdouble )
840             fdouble->_sub[ i_sub ].setData( nb_comps[ i_sub ], supports[ i_sub ] );
841           for ( i_comp = 0; i_comp < nb_comps[ i_sub ]; ++i_comp, next() )
842             {
843               // store component name
844               string compName = getName();
845               if ( fdouble )
846                 fdouble->_sub[ i_sub ].compName( i_comp ) = compName;
847             }
848         }
849       // (4) nb harmonics ( ignored )
850       for ( initIntReading( total_nb_comp ); more(); next() );
851       // (5) TYPE ( ignored )
852       for (initNameReading(1, /*length=*/71); more(); next());
853       // (6) TITRE ( ignored )
854       for (initNameReading(1, /*length=*/71); more(); next());
855       // (7) attributes ( ignored )
856       for ( initIntReading( nb_attr ); more(); next() );
857
858       for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
859         {
860           // loop on components: read values
861           initDoubleReading( nb_values[ i_sub ] * nb_comps[ i_sub ] );
862           for ( i_comp = 0; i_comp < nb_comps[ i_sub ]; ++i_comp )
863             {
864               if ( fdouble )
865                 {
866                   vector<double>& vals = fdouble->addComponent( nb_values[ i_sub ] );
867                   for ( int i = 0; more() && i < nb_values[ i_sub ]; next(), ++i )
868                     vals[ i ] = getDouble();
869                 }
870               else
871                 {
872                   for ( int i = 0; i < nb_values[ i_sub ]; next(), ++i );
873                 }
874             }
875         } // loop on subcomponents of a field
876
877       // set a supporting group including all subs supports but only
878       // if all subs have the same components
879       if ( fdouble && fdouble->hasSameComponentsBySupport() )
880         setFieldSupport( supports, fdouble );
881       else
882         for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
883           fdouble->_sub[ i_sub ]._support->_isProfile = true;
884
885     } // end loop on field objects
886
887   // set field names
888   setFieldNames( _iMed->_nodeFields, objectNames, nameIndices );
889
890 }  // read_PILE_NODES_FIELD()
891
892 //================================================================================
893 /*!
894  * \brief Read "PILE NUMERO  39": FIELDS
895  */
896 //================================================================================
897
898 void SauvReader::read_PILE_FIELD (const int                 nbObjects,
899                                   std::vector<std::string>& objectNames,
900                                   std::vector<int>&         nameIndices)
901 {
902   // REAL EXAMPLE
903
904   // (1)        1       2       6      16
905   // (2)                                                         CARACTERISTIQUES
906   // (3)      -15  317773       4       0       0       0      -2       0       3
907   // (4)             317581
908   // (5)  0
909   // (6)   317767  317761  317755  317815
910   // (7)  YOUN     NU       H        SIGY
911   // (8)  REAL*8            REAL*8            REAL*8            REAL*8
912   // (9)        1       1       0       0
913   // (10)  2.00000000000000E+05
914   // (11)       1       1       0       0
915   // (12)  3.30000000000000E-01
916   // (13)       1       1       0       0
917   // (14)  1.00000000000000E+04
918   // (15)       6     706       0       0
919   // (16)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
920   // (17)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
921   // (18)  ...
922
923   _iMed->_cellFields.resize( nbObjects, (SauvUtilities::DoubleField*) 0 );
924   for (int object=0; object!=nbObjects; ++object) // pour chaque field
925     {
926       initIntReading( 4 );
927       int i_sub, nb_sub = getIntNext(); // (1) <nb_sub> 2 6 <title length>
928       next(); // skip "2"
929       next(); // skip "6"
930       int title_length = getIntNext(); // <title length>
931       if ( nb_sub < 1 )
932         THROW_IK_EXCEPTION("Error of field reading: wrong nb of subcomponents " << nb_sub << lineNb() );
933
934       string description;
935       if ( title_length )
936         {
937           if ( isXRD() )
938             {
939               initNameReading(1, title_length);
940               description = getName();
941               next();
942             }
943           else
944             {
945               char* line;
946               getNextLine( line ); // (2) title
947               const int len = 72; // line length
948               description = string(line + len - title_length, title_length); // title is right justified
949             }
950         }
951       // look for a line starting with '-' : <reference to support>
952       if ( isXRD() )
953         {
954           initIntReading( nb_sub * 9 );
955         }
956       else
957         {
958           do {
959             initIntReading( nb_sub * 9 );
960           } while ( getInt() >= 0 );
961         }
962       int total_nb_comp = 0;
963       vector<Group*> supports( nb_sub );
964       vector<int>     nb_comp( nb_sub );
965       for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
966         {                                    // (3)
967           int supportId     = -getIntNext(); // <reference to support>
968           next();                            // ignore <address>
969           nb_comp [ i_sub ] =  getIntNext(); // <nb of components in the sub>
970           for ( int i = 0; i < 6; ++i ) next();  // ignore 6 ints, in example "0 0 0 -2 0 3"
971
972           if ( supportId < 1 || supportId > (int)_iMed->_groups.size() )
973             THROW_IK_EXCEPTION("Error of field reading: wrong mesh reference "<< supportId << lineNb() );
974           if ( nb_comp[ i_sub ] < 0 )
975             THROW_IK_EXCEPTION("Error of field reading: wrong nb of components " <<nb_comp[ i_sub ] << lineNb() );
976
977           supports[ i_sub ] = &_iMed->_groups[ supportId-1 ];
978           total_nb_comp += nb_comp[ i_sub ];
979         }
980       // (4) dummy strings
981       for ( initNameReading( nb_sub, 17 ); more(); next() );
982       // (5) dummy strings
983       for ( initNameReading( nb_sub ); more(); next() );
984
985       // loop on subcomponents of a field, each of which refers to
986       // a certain support and has its own number of components;
987       // read component values
988       SauvUtilities::DoubleField* fdouble = 0;
989       for ( i_sub = 0; i_sub < nb_sub; ++ i_sub )
990         {
991           vector<string> comp_names( nb_comp[ i_sub ]), comp_type( nb_comp[ i_sub ]);
992           // (6) nb_comp addresses of MELVAL structure
993           for ( initIntReading( nb_comp[ i_sub ] ); more(); next() );
994           // (7) component names
995           for ( initNameReading( nb_comp[ i_sub ] ); more(); next() )
996             comp_names[ index() ] = getName();
997           // (8) component type
998           for ( initNameReading( nb_comp[ i_sub ], 17 ); more(); next() ) // 17 is name width
999             {
1000               comp_type[ index() ] = getName();
1001               // component types must be the same
1002               if ( index() > 0 && comp_type[ index() ] != comp_type[ index() - 1] )
1003                 THROW_IK_EXCEPTION( "Error of field reading: diff component types <"
1004                                     << comp_type[ index() ] << "> != <" << comp_type[ index() - 1 ]
1005                                     << ">" << lineNb() );
1006             }
1007           // now type is known, create a field, one for all subs
1008           bool isReal = (nb_comp[i_sub] > 0) ? (comp_type[0] == "REAL*8") : true;
1009           if ( !fdouble && total_nb_comp )
1010             {
1011               if ( !isReal )
1012                 cout << "Warning: read NOT REAL field, type <" << comp_type[0] << ">" << lineNb() << endl;
1013               _iMed->_cellFields[ object ] = fdouble = new SauvUtilities::DoubleField( nb_sub, total_nb_comp );
1014               fdouble->_description = description;
1015             }
1016           // store support id and nb components of a sub
1017           if ( fdouble )
1018             fdouble->_sub[ i_sub ].setData( nb_comp[ i_sub ], supports[ i_sub ]);
1019           // loop on components: read values
1020           for ( int i_comp = 0; i_comp < nb_comp[ i_sub ]; ++i_comp )
1021             {
1022               // (9) nb of values
1023               initIntReading( 4 );
1024               int nb_val_by_elem = getIntNext();
1025               int nb_values      = getIntNext();
1026               next();
1027               next();
1028               fdouble->_sub[ i_sub ]._nb_gauss[ i_comp ] = nb_val_by_elem;
1029
1030               // (10) values
1031               nb_values *= nb_val_by_elem;
1032               if ( fdouble )
1033                 {
1034                   vector<double> & vals = fdouble->addComponent( nb_values );
1035                   for ( isReal ? initDoubleReading( nb_values ) : initIntReading( nb_values ); more(); next())
1036                     vals[ index() ] = getDouble();
1037                   // store component name
1038                   fdouble->_sub[ i_sub ].compName( i_comp ) = comp_names[ i_comp ];
1039                 }
1040               else
1041                 {
1042                   for ( isReal ? initDoubleReading( nb_values ) : initIntReading( nb_values ); more(); next() ) ;
1043                 }
1044             }
1045         } // loop on subcomponents of a field
1046
1047       // set id of a group including all sub supports but only
1048       // if all subs have the same nb of components
1049       if ( fdouble && fdouble->hasSameComponentsBySupport() )
1050         setFieldSupport( supports, fdouble );
1051       else
1052         for ( i_sub = 0; i_sub < nb_sub; ++i_sub )
1053           fdouble->_sub[ i_sub ]._support->_isProfile = true;
1054
1055     } // end loop on field objects
1056
1057   // set field names
1058   setFieldNames( _iMed->_cellFields, objectNames, nameIndices );
1059
1060 } // read_PILE_FIELD()
1061
1062 //================================================================================
1063 /*!
1064  * \brief Read "PILE NUMERO  10": TABLES
1065  */
1066 //================================================================================
1067
1068 void SauvReader::read_PILE_TABLES (const int                 nbObjects,
1069                                    std::vector<std::string>& objectNames,
1070                                    std::vector<int>&         nameIndices)
1071 {
1072   // IMP 0020434: mapping GIBI names to MED names
1073
1074   string table_med_mail = "MED_MAIL";
1075   string table_med_cham = "MED_CHAM";
1076   string table_med_comp = "MED_COMP";
1077   int table_med_mail_id = -1;
1078   int table_med_cham_id = -1;
1079   int table_med_comp_id = -1;
1080   for (size_t iname = 0; iname < objectNames.size(); iname++)
1081     if      (objectNames[iname] == table_med_mail) table_med_mail_id = nameIndices[iname];
1082     else if (objectNames[iname] == table_med_cham) table_med_cham_id = nameIndices[iname];
1083     else if (objectNames[iname] == table_med_comp) table_med_comp_id = nameIndices[iname];
1084
1085   if ( isASCII() )
1086     if (table_med_mail_id < 0 && table_med_cham_id < 0 && table_med_comp_id < 0)
1087       return;
1088
1089   for (int itable = 1; itable <= nbObjects; itable++)
1090     {
1091       // read tables "MED_MAIL", "MED_CHAM" and "MED_COMP", that keeps correspondence
1092       // between GIBI names (8 symbols if any) and MED names (possibly longer)
1093       initIntReading(1);
1094       int nb_table_vals = getIntNext();
1095       if (nb_table_vals < 0)
1096         THROW_IK_EXCEPTION("Error of reading PILE NUMERO  10" << lineNb() );
1097
1098       int name_i_med_pile;
1099       initIntReading(nb_table_vals);
1100       for (int i = 0; i < nb_table_vals/4; i++)
1101         {
1102           if (itable == table_med_mail_id ||
1103               itable == table_med_cham_id ||
1104               itable == table_med_comp_id)
1105             {
1106               nameGIBItoMED name_i;
1107               name_i_med_pile  = getIntNext();
1108               name_i.med_id    = getIntNext();
1109               name_i.gibi_pile = getIntNext();
1110               name_i.gibi_id   = getIntNext();
1111
1112               if (name_i_med_pile != PILE_STRINGS)
1113                 {
1114                   // error: med(long) names are always kept in PILE_STRINGS
1115                 }
1116               if (itable == table_med_mail_id)
1117                 {
1118                   if (name_i.gibi_pile != PILE_SOUS_MAILLAGE) {
1119                     // error: must be PILE_SOUS_MAILLAGE
1120                   }
1121                   _iMed->_listGIBItoMED_mail.push_back(name_i);
1122                 }
1123               else if (itable == table_med_cham_id)
1124                 {
1125                   if (name_i.gibi_pile != PILE_FIELD &&
1126                       name_i.gibi_pile != PILE_NODES_FIELD)
1127                     {
1128                       // error: must be PILE_FIELD or PILE_NODES_FIELD
1129                     }
1130                   _iMed->_listGIBItoMED_cham.push_back(name_i);
1131                 }
1132               else if (itable == table_med_comp_id)
1133                 {
1134                   if (name_i.gibi_pile != PILE_STRINGS)
1135                     {
1136                       // error: gibi(short) names of components are kept in PILE_STRINGS
1137                     }
1138                   _iMed->_listGIBItoMED_comp.push_back(name_i);
1139                 }
1140             }
1141           else
1142             {
1143               // pass table
1144               for ( int ii = 0; ii < 4; ++ii ) next();
1145             }
1146         }
1147     } // for (int itable = 0; itable < nbObjects; itable++)
1148 }
1149
1150 //================================================================================
1151 /*!
1152  * \brief Read "PILE NUMERO  27"
1153  */
1154 //================================================================================
1155
1156 void SauvReader::read_PILE_STRINGS (const int                 nbObjects,
1157                                     std::vector<std::string>& objectNames,
1158                                     std::vector<int>&         nameIndices)
1159 {
1160   // IMP 0020434: mapping GIBI names to MED names
1161   initIntReading(2);
1162   int stringLen    = getIntNext();
1163   int nbSubStrings = getIntNext();
1164   if (nbSubStrings != nbObjects)
1165     THROW_IK_EXCEPTION("Error of reading PILE NUMERO  27" << lineNb() );
1166
1167   string aWholeString;
1168   if ( isXRD() )
1169     {
1170       const int fixedLength = 71;
1171       while ((int)aWholeString.length() < stringLen)
1172         {
1173           int remainLen = stringLen - aWholeString.length();
1174           int len;
1175           if ( remainLen > fixedLength )
1176             {
1177               len = fixedLength;
1178             }
1179           else
1180             {
1181               len = remainLen;
1182             }
1183           initNameReading(1, len);
1184           aWholeString += getName();
1185           next();
1186         }
1187     }
1188   else
1189     {
1190       char* line;
1191       const int fixedLength = 71;
1192       while ((int)aWholeString.length() < stringLen)
1193         {
1194           getNextLine( line );
1195           int remainLen = stringLen - aWholeString.length();
1196           if ( remainLen > fixedLength )
1197             {
1198               aWholeString += line + 1;
1199             }
1200           else
1201             {
1202               aWholeString += line + ( 72 - remainLen );
1203             }
1204         }
1205     }
1206   int prevOffset = 0;
1207   int currOffset = 0;
1208   initIntReading(nbSubStrings);
1209   for (int istr = 1; istr <= nbSubStrings; istr++, next())
1210     {
1211       currOffset = getInt();
1212       // fill mapStrings
1213       _iMed->_mapStrings[istr] = aWholeString.substr(prevOffset, currOffset - prevOffset);
1214       prevOffset = currOffset;
1215     }
1216 }