Salome HOME
Deal with eccentricity for structure element
[tools/medcoupling.git] / src / MEDLoader / MEDFileBlowStrEltUp.cxx
1 // Copyright (C) 2007-2017  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 // Author : Anthony Geay (EDF R&D)
20
21 #include "MEDFileBlowStrEltUp.hxx"
22 #include "MEDCouplingFieldDouble.hxx"
23 #include "MEDFileFieldVisitor.hxx"
24 #include "MEDCouplingPartDefinition.hxx"
25 #include "MCAuto.txx"
26
27 using namespace MEDCoupling;
28
29 const char MEDFileBlowStrEltUp::MED_BALL_STR[]="MED_BALL";
30
31 MEDFileBlowStrEltUp::MEDFileBlowStrEltUp(const MEDFileFields *fsOnlyOnSE, const MEDFileMeshes *ms, const MEDFileStructureElements *ses)
32 {
33   if(!fsOnlyOnSE || !ms || !ses)
34     throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp constructor : NULL input pointer !");
35   _ms.takeRef(ms); _ses.takeRef(ses);
36   std::vector< std::pair<std::string,std::string> > ps;
37   fsOnlyOnSE->getMeshSENames(ps);
38   std::size_t sz(ps.size());
39   _elts.resize(sz);
40   for(std::size_t i=0;i<sz;i++)
41     {
42       const std::pair<std::string,std::string>& p(ps[i]);
43       MCAuto<MEDFileFields> f(fsOnlyOnSE->partOfThisLyingOnSpecifiedMeshSEName(p.first,p.second));
44       _elts[i]=f;
45     }
46   for(std::size_t i=0;i<sz;i++)
47     {
48       const std::pair<std::string,std::string>& p(ps[i]);
49       MEDFileMesh *mesh(_ms->getMeshWithName(p.first));
50       if(!mesh)
51         throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp : NULL mesh !");
52       MEDFileUMesh *umesh(dynamic_cast<MEDFileUMesh *>(mesh));
53       if(!umesh)
54         throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp : Blow up of Stru Elt not managed yet for unstructured meshes !");
55     }
56 }
57
58 /*!
59  * \param [in] mesh - The mesh containing structure element called \a seName. After the call of this method the Structure elements parts will be removed.
60  * \param [out] mOut - the physical mesh of the structure element \a seName in mesh \a mesh
61  * \param [out] fsOut - the list of var attribute of structure element \a seName - \b WARNING no time steps here
62  */
63 MCAuto<MEDFileEltStruct4Mesh> MEDFileBlowStrEltUp::dealWithSEInMesh(const std::string& seName, MEDFileUMesh *mesh, MCAuto<MEDFileUMesh>& mOut, MCAuto<MEDFileFields>& fsOut) const
64 {
65   if(!mesh)
66     throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithSEInMesh : null pointer !");
67   if(seName==MED_BALL_STR)
68     {
69       MCAuto<MEDFileEltStruct4Mesh> ret(dealWithMEDBALLInMesh(mesh,mOut,fsOut));
70       mesh->killStructureElements();
71       return ret;
72     }
73   throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithSEInMesh : only MED_BALL is managed for the moment, but if you are interested please send spec to anthony.geay@edf.fr !");
74 }
75
76 MCAuto<MEDFileEltStruct4Mesh> MEDFileBlowStrEltUp::dealWithMEDBALLInMesh(const MEDFileUMesh *mesh, MCAuto<MEDFileUMesh>& mOut, MCAuto<MEDFileFields>& fsOut) const
77 {
78   mOut=MEDFileUMesh::New(); fsOut=MEDFileFields::New();
79   const std::vector< MCAuto<MEDFileEltStruct4Mesh> >& strs(mesh->getAccessOfUndergroundEltStrs());
80   MCAuto<MEDFileEltStruct4Mesh> zeStr;
81   for(std::vector< MCAuto<MEDFileEltStruct4Mesh> >::const_iterator it=strs.begin();it!=strs.end();it++)
82     {
83       if((*it)->getGeoTypeName()==MED_BALL_STR)
84         {
85           zeStr=*it;
86           break;
87         }
88     }
89   if(zeStr.isNull())
90     {
91       std::ostringstream oss; oss << "MEDFileBlowStrEltUp::dealWithMEDBALLInMesh : no geo type with name " <<  MED_BALL_STR << " in " << mesh->getName() << " !";
92       throw INTERP_KERNEL::Exception(oss.str());
93     }
94   const DataArrayDouble *coo(mesh->getCoords());
95   if(!coo)
96     throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLInMesh : null coords !");
97   MCAuto<DataArrayInt> conn(zeStr->getConn());
98   if(conn.isNull())
99     throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLInMesh : null connectivity !");
100   conn->checkAllocated();
101   if(conn->getNumberOfComponents()!=1)
102     throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLInMesh : excepted to be single compo !");
103   int nbCells(conn->getNumberOfTuples());
104   MCAuto<DataArrayDouble> connOut(coo->selectByTupleIdSafe(conn->begin(),conn->end()));
105   MCAuto<MEDCouplingUMesh> mcOut(MEDCouplingUMesh::Build0DMeshFromCoords(connOut));
106   mcOut->setName(BuildNewMeshName(mesh->getName(),MED_BALL_STR));
107   mOut->setMeshAtLevel(0,mcOut);
108   const DataArrayInt *ff1(mesh->getFamilyFieldAtLevel(1));
109   if(ff1)
110     {
111       MCAuto<DataArrayInt> ff1o(ff1->selectByTupleIdSafe(conn->begin(),conn->end()));
112       mOut->setFamilyFieldArr(1,ff1o);
113     }
114   const DataArrayInt *nf1(mesh->getNumberFieldAtLevel(1));
115   if(nf1)
116     {
117       MCAuto<DataArrayInt> nf1o(nf1->selectByTupleIdSafe(conn->begin(),conn->end()));
118       mOut->setRenumFieldArr(1,nf1o);
119     }
120   MCAuto<MEDFileUMeshPerTypeCommon> md(zeStr->getMeshDef());
121   const DataArrayInt *ff0(md->getFam());
122   if(ff0)
123     mOut->setFamilyFieldArr(0,const_cast<DataArrayInt *>(ff0));
124   const DataArrayInt *nf0(md->getNum());
125   if(nf0)
126     mOut->setRenumFieldArr(0,const_cast<DataArrayInt *>(nf0));
127   mOut->copyFamGrpMapsFrom(*mesh);
128   const std::vector< MCAuto<DataArray> >& vars(zeStr->getVars());
129   for(std::vector< MCAuto<DataArray> >::const_iterator it=vars.begin();it!=vars.end();it++)
130     {
131       const DataArray *elt(*it);
132       if(!elt)
133         continue;
134       {
135         const DataArrayDouble *eltC(dynamic_cast<const DataArrayDouble *>(elt));
136         if(eltC)
137           {
138             MCAuto<MEDFileFieldMultiTS> fmts(MEDFileFieldMultiTS::New());
139             MCAuto<MEDFileField1TS> f1ts(MEDFileField1TS::New());
140             MCAuto<MEDCouplingFieldDouble> f(MEDCouplingFieldDouble::New(ON_NODES));
141             f->setMesh(mcOut);
142             f->setArray(const_cast<DataArrayDouble *>(eltC));
143             f->setName(eltC->getName());
144             f1ts->setFieldNoProfileSBT(f);
145             fmts->pushBackTimeStep(f1ts);
146             fsOut->pushField(fmts);
147           }
148       }
149     }
150   return zeStr;
151 }
152
153 /*!
154  * \param [in] fs - fields lying all on same mesh and on same structure element
155  * \param [in] zeStr - ze structure of current structure element
156  * \param [in] varAtt - fields containing var att of current structure element. WARNING at this stage the number of iteration are equal to one for each field in \a varAtt
157  * \param [out] zeOutputs - ze fields that are the concatenation of fields in \a fs transformed and those in \a varAtt normalized in time space
158  */
159 void MEDFileBlowStrEltUp::dealWithSEInFields(const std::string& seName, const MEDFileFields *fs, const MEDFileEltStruct4Mesh *zeStr, const MEDFileFields *varAtt, MEDFileFields *zeOutputs) const
160 {
161   if(!fs)
162     throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithSEInFields : null pointer !");
163   if(seName==MED_BALL_STR)
164     {
165       dealWithMEDBALLSInFields(fs,zeStr,varAtt,zeOutputs);
166       return ;
167     }
168   throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithSEInFields : only MED_BALL is managed for the moment, but if you are interested please send spec to anthony.geay@edf.fr !");
169 }
170
171 void MEDFileBlowStrEltUp::dealWithMEDBALLSInFields(const MEDFileFields *fs, const MEDFileEltStruct4Mesh *zeStr, const MEDFileFields *varAtt, MEDFileFields *zeOutputs) const
172 {
173   int nbf(fs->getNumberOfFields());
174   std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > elts0;
175   std::vector< MEDFileAnyTypeFieldMultiTS * > elts1;
176   std::string zeMeshName;
177   for(int i=0;i<nbf;i++)
178     {
179       MCAuto<MEDFileAnyTypeFieldMultiTS> elt(fs->getFieldAtPos(i));
180       MCAuto<MEDFileAnyTypeFieldMultiTS> eltOut(elt->buildNewEmpty());
181       int nbTS(elt->getNumberOfTS());
182       for(int j=0;j<nbTS;j++)
183         {
184           MCAuto<MEDFileAnyTypeField1TS> eltt(elt->getTimeStepAtPos(j));
185           MCAuto<MEDFileAnyTypeField1TS> elttOut(eltt->deepCopy());
186           std::string meshName(eltt->getMeshName());
187           zeMeshName=BuildNewMeshName(meshName,MED_BALL_STR);
188           elttOut->setMeshName(zeMeshName);
189           elttOut->convertMedBallIntoClassic();
190           eltOut->pushBackTimeStep(elttOut);
191         }
192       elts0.push_back(eltOut); elts1.push_back(eltOut);
193     }
194   //
195   const MEDFileMesh *zeCurrentMesh(_ms->getMeshWithName(zeMeshName));
196   //
197   std::size_t ii(0);
198   std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > sp(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(elts1));
199   for(std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> >::const_iterator it0=sp.begin();it0!=sp.end();it0++,ii++)
200     {
201       std::vector< MCAuto<MEDFileFastCellSupportComparator> > fsc;
202       std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > sp2(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(*it0,zeCurrentMesh,fsc));
203       std::size_t jj(0);
204       for(std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> >::const_iterator it1=sp2.begin();it1!=sp2.end();it1++,jj++)
205         {
206           for(std::vector<MEDFileAnyTypeFieldMultiTS *>::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
207             zeOutputs->pushField(*it2);
208           // The most exciting part. Users that put profiles on struct elements part of fields. Reduce var att.
209           if((*it1).size()<1)
210             throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLSInFields : take a deep breath !");
211           MCAuto<MEDFileAnyTypeField1TS> zeGuideForPfl;// This var is the reference for pfl management.
212           {
213             if(!(*it1)[0])
214               throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLSInFields : take a deep breath 2 !");
215             int pdm((*it1)[0]->getNumberOfTS());
216             if(pdm<1)
217               throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLSInFields : take a deep breath 3 !");
218             zeGuideForPfl=(*it1)[0]->getTimeStepAtPos(0);
219           }
220           if(zeGuideForPfl.isNull())
221             throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLSInFields : take a deep breath 4 !");
222           std::vector<std::string> pfls(zeGuideForPfl->getPflsReallyUsed());
223           if(pfls.size()>=2)
224             throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLSInFields : drink less coffee");
225           MCAuto<DataArrayInt> pflMyLove;
226           if(pfls.size()==1)
227             pflMyLove.takeRef(zeGuideForPfl->getProfile(pfls[0]));
228           // Yeah we have pfls
229           std::vector<double> t2s;
230           std::vector< std::pair<int,int> > t1s((*it1)[0]->getTimeSteps(t2s));
231           std::size_t nbTS3(t2s.size());
232           int nbf2(varAtt->getNumberOfFields());
233           for(int i=0;i<nbf2;i++)
234             {
235               MCAuto<MEDFileAnyTypeFieldMultiTS> elt(varAtt->getFieldAtPos(i));
236               int nbTS2(elt->getNumberOfTS());
237               if(nbTS2!=1)
238                 throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::dealWithMEDBALLSInFields : internal error ! The dealWithMEDBALLInMesh is expected to return a single TS !");
239               MCAuto<MEDFileAnyTypeField1TS> elt2(elt->getTimeStepAtPos(0));
240               MCAuto<MEDFileAnyTypeFieldMultiTS> elt4(elt->buildNewEmpty());
241               for(std::size_t j=0;j<nbTS3;j++)
242                 {
243                   MCAuto<MEDFileAnyTypeField1TS> elt3(elt2->deepCopy());
244                   elt3->setTime(t1s[j].first,t1s[j].second,t2s[j]);
245                   elt3->setName(BuildVarAttName(ii,sp.size(),jj,sp2.size(),elt3->getName()));
246                   if(pflMyLove.isNotNull())
247                     elt3->makeReduction(INTERP_KERNEL::NORM_ERROR,ON_NODES,pflMyLove);
248                   elt4->pushBackTimeStep(elt3);
249                 }
250               zeOutputs->pushField(elt4);
251             }
252         }
253     }
254 }
255
256 void MEDFileBlowStrEltUp::generate(MEDFileMeshes *msOut, MEDFileFields *allZeOutFields)
257 {
258   for(std::vector< MCAuto<MEDFileFields> >::iterator elt=_elts.begin();elt!=_elts.end();elt++)
259     {
260       std::vector< std::pair<std::string,std::string> > ps;
261       (*elt)->getMeshSENames(ps);
262       if(ps.size()!=1)
263         throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::generateMeshes : internal error !");
264       MEDFileMesh *mesh(_ms->getMeshWithName(ps[0].first));
265       if(!mesh)
266         throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::generateMeshes : NULL mesh !");
267       MEDFileUMesh *umesh(dynamic_cast<MEDFileUMesh *>(mesh));
268       if(!umesh)
269         throw INTERP_KERNEL::Exception("MEDFileBlowStrEltUp::generateMeshes : Blow up of Stru Elt not managed yet for unstructured meshes !");
270       //
271       MCAuto<MEDFileFields> classicalSEFields(splitFieldsPerLoc(*elt,umesh,msOut,allZeOutFields));
272       if(classicalSEFields.isNotNull())
273         {
274           MCAuto<MEDFileUMesh> mOut;
275           MCAuto<MEDFileFields> fsOut1;
276           MCAuto<MEDFileEltStruct4Mesh> zeStr(dealWithSEInMesh(ps[0].second,umesh,mOut,fsOut1));
277           msOut->pushMesh(mOut);
278           dealWithSEInFields(ps[0].second,classicalSEFields,zeStr,fsOut1,allZeOutFields);
279         }
280     }
281 }
282
283 std::string MEDFileBlowStrEltUp::BuildNewMeshName(const std::string& meshName, const std::string& seName)
284 {
285   std::ostringstream mNameOut;
286   mNameOut << meshName << "_" << seName;
287   return mNameOut.str();
288 }
289
290 std::string MEDFileBlowStrEltUp::BuildVarAttName(std::size_t iPart, std::size_t totINbParts, std::size_t jPart, std::size_t totJNbParts, const std::string& name)
291 {
292   if(totINbParts==1 && totJNbParts==1)
293     return name;
294   std::ostringstream oss;
295   oss << name << "@" << iPart << "@" << jPart;
296   return oss.str();
297 }
298
299 void MEDFileBlowStrEltUp::DealWithSE(MEDFileFields *fs, MEDFileMeshes *ms, const MEDFileStructureElements *ses)
300 {
301   MCAuto<MEDFileFields> fsSEOnly(fs->partOfThisOnStructureElements());
302   fs->killStructureElements();
303   MEDFileBlowStrEltUp bu(fsSEOnly,ms,ses);
304   bu.generate(ms,fs);
305   fs->killStructureElementsInGlobs();
306 }
307
308 //
309
310 class FieldWalker2
311 {
312 public:
313   FieldWalker2(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd);
314   std::string getLoc() const { return _loc; }
315   std::string getPfl() const { return _pfl; }
316   bool isClassic() const { return _is_classic; }
317   bool operator!=(const FieldWalker2& other) const;
318   bool operator==(const FieldWalker2& other) const;
319   const SlicePartDefinition *getPartDef() const { return _pd; }
320 private:
321   std::string _loc;
322   std::string _pfl;
323   bool _is_classic;
324   MCAuto<SlicePartDefinition> _pd;
325 };
326
327 class LocInfo
328 {
329 public:
330   LocInfo() { }
331   LocInfo(const std::vector<FieldWalker2>& fw);
332   bool operator==(const LocInfo& other) const { return _locs==other._locs && _pfl==other._pfl; }
333   void push(const std::string& loc, const std::string& pfl) { checkUniqueLoc(loc); _locs.push_back(loc); _pfl.push_back(pfl); }
334   MCAuto<MEDFileUMesh> generateNonClassicalData(int zePos, const MEDFileUMesh *mesh, const MEDFileFieldGlobsReal *globs) const;
335   const PartDefinition *getPartDef() const { return _pd; }
336 private:
337   void checkUniqueLoc(const std::string& loc) const;
338   static MCAuto<DataArrayDouble> BuildMeshFromAngleVrille(INTERP_KERNEL::NormalizedCellType gt, const DataArrayDouble *angleDeVrille, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs);
339   static MCAuto<DataArrayDouble> BuildMeshFromEpaisseur(INTERP_KERNEL::NormalizedCellType gt, const DataArrayDouble *thikness, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs);
340   static MCAuto<DataArrayDouble> BuildMeshPipeSEG3(const DataArrayDouble *angle, const DataArrayDouble *scale, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs);
341   static MCAuto<MEDCouplingUMesh> BuildMeshCommon(INTERP_KERNEL::NormalizedCellType gt, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs, MCAuto<DataArrayDouble>& ptsForLoc);
342   static MCAuto<DataArrayDouble> BuildMeshFromStructure(INTERP_KERNEL::NormalizedCellType gt, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs);
343 public:
344   static const char ANGLE_DE_VRILLE[];
345   static const char ANGLE[];
346   static const char SCALE[];
347   static const char EPAISSEUR[];
348 private:
349   std::vector<std::string> _locs;
350   std::vector<std::string> _pfl;
351   MCAuto<PartDefinition> _pd;
352 };
353
354 const char LocInfo::ANGLE_DE_VRILLE[]="ANGLE DE VRILLE";
355
356 const char LocInfo::ANGLE[]="ANGLE";
357
358 const char LocInfo::SCALE[]="SCALE";
359
360 const char LocInfo::EPAISSEUR[]="EPAISSEUR";
361
362 LocInfo::LocInfo(const std::vector<FieldWalker2>& fw)
363 {
364   std::size_t sz(fw.size());
365   _locs.resize(sz); _pfl.resize(sz);
366   if(sz>0)
367     _pd=fw[0].getPartDef()->deepCopy();
368   for(std::size_t i=0;i<sz;i++)
369     {
370       _locs[i]=fw[i].getLoc();
371       _pfl[i]=fw[i].getPfl();
372       if(i>0)
373         _pd=(*_pd)+(*(fw[i].getPartDef()));
374     }
375 }
376
377 void LocInfo::checkUniqueLoc(const std::string& loc) const
378 {
379   if(std::find(_locs.begin(),_locs.end(),loc)!=_locs.end())
380     {
381       std::ostringstream oss; oss << "LocInfo::checkUniqueLoc : loc \"" << loc << "\" already exists !";
382       throw INTERP_KERNEL::Exception(oss.str());
383     }
384 }
385
386 MCAuto<MEDCouplingUMesh> LocInfo::BuildMeshCommon(INTERP_KERNEL::NormalizedCellType gt, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs, MCAuto<DataArrayDouble>& ptsForLoc)
387 {
388   MCAuto<DataArrayInt> conn(zeStr->getConn());
389   conn=conn->deepCopy(); conn->rearrange(1);
390   MCAuto<MEDCouplingUMesh> geoMesh;
391   {
392     MCAuto<MEDCoupling1SGTUMesh> umesh(MEDCoupling1SGTUMesh::New("",gt));
393     umesh->setCoords(mesh->getCoords());
394     umesh->setNodalConnectivity(conn);
395     geoMesh=umesh->buildUnstructured();
396   }
397   //
398   if(!pfl.empty())
399     {
400       const DataArrayInt *pflArr(globs->getProfile(pfl));
401       geoMesh=geoMesh->buildPartOfMySelf(pflArr->begin(),pflArr->end(),true);
402     }
403   //
404   MCAuto<MEDCouplingFieldDouble> fakeF(MEDCouplingFieldDouble::New(ON_GAUSS_PT));
405   fakeF->setMesh(geoMesh);
406   fakeF->setGaussLocalizationOnType(gt,loc.getRefCoords(),loc.getGaussCoords(),loc.getGaussWeights());
407   ptsForLoc=fakeF->getLocalizationOfDiscr();
408   //
409   return geoMesh;
410 }
411
412 MCAuto<DataArrayDouble> LocInfo::BuildMeshFromAngleVrille(INTERP_KERNEL::NormalizedCellType gt, const DataArrayDouble *angleDeVrille, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs)
413 {
414   MCAuto<DataArrayDouble> ptsForLoc;
415   MCAuto<MEDCouplingUMesh> geoMesh(BuildMeshCommon(gt,pfl,loc,zeStr,mesh,section,globs,ptsForLoc));
416   //
417   MCConstAuto<DataArrayDouble> angleVrille;
418   if(!pfl.empty())
419     {
420       const DataArrayInt *pflArr(globs->getProfile(pfl));
421       angleVrille=angleDeVrille->selectByTupleIdSafe(pflArr->begin(),pflArr->end());
422     }
423   else
424     angleVrille.takeRef(angleDeVrille);
425   //
426   MCAuto<MEDCouplingFieldDouble> dir(geoMesh->buildDirectionVectorField());
427   MCAuto<DataArrayDouble> rot(dir->getArray()->fromCartToSpher());
428   int nbCompo(ptsForLoc->getNumberOfComponents());
429   MCAuto<DataArrayDouble> secPts(section->getCoords()->changeNbOfComponents(nbCompo,0.));
430   int nbSecPts(secPts->getNumberOfTuples()),nbCells(geoMesh->getNumberOfCells()),nbg(loc.getGaussWeights().size());
431   {
432     const int TAB[3]={2,0,1};
433     std::vector<int> v(TAB,TAB+3);
434     secPts=secPts->keepSelectedComponents(v);
435   }
436   const double CENTER[3]={0.,0.,0.},AX0[3]={0.,0.,1.};
437   double AX1[3]; AX1[2]=0.;
438   std::vector< MCAuto<DataArrayDouble> > arrs(nbCells*nbg);
439   for(int j=0;j<nbCells;j++)
440     {
441       MCAuto<DataArrayDouble> p(secPts->deepCopy());
442       double ang0(rot->getIJ(j,2));
443       DataArrayDouble::Rotate3DAlg(CENTER,AX0,ang0,nbSecPts,p->begin(),p->getPointer());
444       AX1[0]=-sin(ang0); AX1[1]=cos(ang0);// rot Oy around OZ
445       double ang1(M_PI/2.-rot->getIJ(j,1));
446       DataArrayDouble::Rotate3DAlg(CENTER,AX1,-ang1,nbSecPts,p->begin(),p->getPointer());
447       DataArrayDouble::Rotate3DAlg(CENTER,dir->getArray()->begin()+j*3,angleVrille->getIJ(j,0),nbSecPts,p->begin(),p->getPointer());
448       for(int l=0;l<nbg;l++)
449         {
450           MCAuto<DataArrayDouble> p2(p->deepCopy());
451           for(int k=0;k<nbCompo;k++)
452             p2->applyLin(1.,ptsForLoc->getIJ(j*nbg+l,k),k);
453           arrs[j*nbg+l]=p2;
454         }
455     }
456   std::vector<const DataArrayDouble *> arrs2(VecAutoToVecOfCstPt(arrs));
457   MCAuto<DataArrayDouble> resu(DataArrayDouble::Aggregate(arrs2));
458   return resu;
459 }
460
461 MCAuto<DataArrayDouble> LocInfo::BuildMeshFromEpaisseur(INTERP_KERNEL::NormalizedCellType gt, const DataArrayDouble *thikness, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs)
462 {
463 #if __cplusplus >= 201103L
464   MCAuto<DataArrayDouble> ptsForLoc;
465   MCAuto<MEDCouplingUMesh> geoMesh(BuildMeshCommon(gt,pfl,loc,zeStr,mesh,section,globs,ptsForLoc));
466   int nbSecPts(section->getNumberOfNodes()),nbCells(geoMesh->getNumberOfCells()),nbg(loc.getGaussWeights().size());
467   MCConstAuto<DataArrayDouble> zeThikness;
468   if(!pfl.empty())
469     {
470       const DataArrayInt *pflArr(globs->getProfile(pfl));
471       zeThikness=thikness->selectByTupleIdSafe(pflArr->begin(),pflArr->end());
472     }
473   else
474     zeThikness.takeRef(thikness);
475   MCAuto<DataArrayDouble> orthoArr;
476   {
477     MCAuto<MEDCouplingFieldDouble> ortho(geoMesh->buildOrthogonalField());
478     orthoArr.takeRef(ortho->getArray());
479   }
480   int nbCompo(orthoArr->getNumberOfComponents());
481   MCAuto<DataArrayDouble> secPts(section->getCoords()->duplicateEachTupleNTimes(nbCompo));
482   secPts->rearrange(nbCompo);
483   std::vector< MCAuto<DataArrayDouble> > arrs(nbCells*nbg);
484   for(int j=0;j<nbCells;j++)
485     {
486       double thck(zeThikness->getIJ(j,0)),eccentricity(zeThikness->getIJ(j,1));
487       MCAuto<DataArrayDouble> fact(DataArrayDouble::New()),fact2(DataArrayDouble::New()); fact->alloc(1,nbCompo); fact2->alloc(1,nbCompo);
488       std::copy(orthoArr->begin()+j*nbCompo,orthoArr->begin()+(j+1)*nbCompo,fact->getPointer());
489       std::copy(orthoArr->begin()+j*nbCompo,orthoArr->begin()+(j+1)*nbCompo,fact2->getPointer());
490       std::for_each(fact2->rwBegin(),fact2->rwEnd(),[eccentricity](double& val){ val*=eccentricity; });
491       std::transform(fact->begin(),fact->end(),fact->getPointer(),[thck](const double& val){ return val*thck/2.; });
492       MCAuto<DataArrayDouble> p(DataArrayDouble::Multiply(secPts,fact));
493       p=DataArrayDouble::Add(p,fact2);
494       for(int l=0;l<nbg;l++)
495         {
496           MCAuto<DataArrayDouble> p2(p->deepCopy());
497           for(int k=0;k<nbCompo;k++)
498             p2->applyLin(1.,ptsForLoc->getIJ(j*nbg+l,k),k);
499           arrs[j*nbg+l]=p2;
500         }
501     }
502   std::vector<const DataArrayDouble *> arrs2(VecAutoToVecOfCstPt(arrs));
503   MCAuto<DataArrayDouble> resu(DataArrayDouble::Aggregate(arrs2));
504   return resu;
505 #else
506   throw INTERP_KERNEL::Exception("Broken news : 10% off for C++11 compiler :)");
507 #endif
508 }
509
510 MCAuto<DataArrayDouble> LocInfo::BuildMeshPipeSEG3(const DataArrayDouble *angle, const DataArrayDouble *scale, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs)
511 {
512   static const char MSG1[]="BuildMeshPipeSEG3 : not recognized pattern ! Send mail to anthony.geay@edf.fr with corresponding MED file !";
513   MCAuto<DataArrayDouble> ptsForLoc;
514   MCAuto<MEDCouplingUMesh> geoMesh(BuildMeshCommon(INTERP_KERNEL::NORM_SEG3,pfl,loc,zeStr,mesh,section,globs,ptsForLoc));
515   int nbSecPts(section->getNumberOfNodes()),nbCells(geoMesh->getNumberOfCells()),nbg(loc.getGaussWeights().size());
516   MCConstAuto<DataArrayDouble> zeAngle,zeScale;
517   if(!pfl.empty())
518     {
519       const DataArrayInt *pflArr(globs->getProfile(pfl));
520       zeAngle=angle->selectByTupleIdSafe(pflArr->begin(),pflArr->end());
521       zeScale=scale->selectByTupleIdSafe(pflArr->begin(),pflArr->end());
522     }
523   else
524     {
525       zeAngle.takeRef(angle);
526       zeScale.takeRef(scale);
527     }
528   if(zeAngle->getNumberOfComponents()!=3 || zeScale->getNumberOfComponents()!=2 || nbg!=3)
529     throw INTERP_KERNEL::Exception(MSG1);
530   MCAuto<MEDCouplingFieldDouble> dir;
531   {
532     MCAuto<MEDCouplingUMesh> geoMesh2(geoMesh->deepCopy());
533     geoMesh2->convertQuadraticCellsToLinear();
534     dir=geoMesh2->buildDirectionVectorField();
535   }
536   MCAuto<DataArrayDouble> rot(dir->getArray()->fromCartToSpher());
537   int nbCompo(ptsForLoc->getNumberOfComponents());
538   MCAuto<DataArrayDouble> secPts(section->getCoords()->changeNbOfComponents(nbCompo,0.));
539   {
540     const int TAB[3]={2,0,1};
541     std::vector<int> v(TAB,TAB+3);
542     secPts=secPts->keepSelectedComponents(v);
543   }
544   const double CENTER[3]={0.,0.,0.},AX0[3]={0.,0.,1.};
545   double AX1[3]; AX1[2]=0.;
546   std::vector< MCAuto<DataArrayDouble> > arrs(nbCells*nbg);
547   for(int j=0;j<nbCells;j++)
548     {
549       MCAuto<DataArrayDouble> p(secPts->deepCopy());
550       double ang0(rot->getIJ(j,2));
551       DataArrayDouble::Rotate3DAlg(CENTER,AX0,ang0,nbSecPts,p->begin(),p->getPointer());
552       AX1[0]=-sin(ang0); AX1[1]=cos(ang0);// rot Oy around OZ
553       double ang1(M_PI/2.-rot->getIJ(j,1));
554       DataArrayDouble::Rotate3DAlg(CENTER,AX1,-ang1,nbSecPts,p->begin(),p->getPointer());
555       for(int l=0;l<3;l++)
556         {
557           MCAuto<DataArrayDouble> p3(p->deepCopy());
558           DataArrayDouble::Rotate3DAlg(CENTER,dir->getArray()->begin()+j*3,zeAngle->getIJ(j,l),nbSecPts,p3->begin(),p3->getPointer());
559           MCAuto<DataArrayDouble> p2(p3->deepCopy());
560           for(int k=0;k<nbCompo;k++)
561             p2->applyLin(1.,ptsForLoc->getIJ(j*nbg+l,k),k);
562           arrs[j*nbg+l]=p2;
563         }
564     }
565   std::vector<const DataArrayDouble *> arrs2(VecAutoToVecOfCstPt(arrs));
566   MCAuto<DataArrayDouble> resu(DataArrayDouble::Aggregate(arrs2));
567   return resu;
568 }
569
570 MCAuto<DataArrayDouble> LocInfo::BuildMeshFromStructure(INTERP_KERNEL::NormalizedCellType gt, const std::string& pfl, const MEDFileFieldLoc& loc, const MEDFileEltStruct4Mesh *zeStr, const MEDFileUMesh *mesh, const MEDFileUMesh *section, const MEDFileFieldGlobsReal *globs)
571 {
572   static const char MSG1[]="BuildMeshFromStructure : not recognized pattern ! Send mail to anthony.geay@edf.fr with corresponding MED file !";
573   const std::vector< MCAuto<DataArray> >& vars(zeStr->getVars());
574   if(vars.size()==1)
575     {
576       MCAuto<DataArray> zeArr(vars[0]);
577       if(zeArr.isNull())
578         throw INTERP_KERNEL::Exception(MSG1);
579       MCAuto<DataArrayDouble> zeArr2(DynamicCast<DataArray,DataArrayDouble>(zeArr));
580       if(zeArr2.isNull())
581         throw INTERP_KERNEL::Exception(MSG1);
582       if(zeArr2->getName()==ANGLE_DE_VRILLE || zeArr2->getName()==ANGLE)
583         return BuildMeshFromAngleVrille(gt,zeArr2,pfl,loc,zeStr,mesh,section,globs);
584       if(zeArr2->getName()==EPAISSEUR || zeArr2->getName()==SCALE)
585         return BuildMeshFromEpaisseur(gt,zeArr2,pfl,loc,zeStr,mesh,section,globs);
586     }
587   if(vars.size()==2)
588     {
589       MCAuto<DataArray> zeArr0(vars[0]),zeArr1(vars[1]);
590       if(zeArr0.isNull() || zeArr1.isNull())
591         throw INTERP_KERNEL::Exception(MSG1);
592       MCAuto<DataArrayDouble> zeArr00(DynamicCastSafe<DataArray,DataArrayDouble>(zeArr0)),zeArr11(DynamicCastSafe<DataArray,DataArrayDouble>(zeArr1));
593       switch(gt)
594       {
595         case INTERP_KERNEL::NORM_SEG3:
596           {
597             MCAuto<DataArrayDouble> angle,scale;
598             if(zeArr00->getName()==ANGLE)
599               angle=zeArr00;
600             if(zeArr00->getName()==SCALE)
601               scale=zeArr00;
602             if(zeArr11->getName()==ANGLE)
603               angle=zeArr11;
604             if(zeArr11->getName()==SCALE)
605               scale=zeArr11;
606             if(angle.isNull() || scale.isNull())
607               throw INTERP_KERNEL::Exception(MSG1);
608             return BuildMeshPipeSEG3(angle,scale,pfl,loc,zeStr,mesh,section,globs);
609           }
610         default:
611           throw INTERP_KERNEL::Exception(MSG1);
612       }
613     }
614   throw INTERP_KERNEL::Exception(MSG1);
615 }
616
617 MCAuto<MEDFileUMesh> LocInfo::generateNonClassicalData(int zePos, const MEDFileUMesh *mesh, const MEDFileFieldGlobsReal *globs) const
618 {
619   static const char MSG1[]="LocInfo::generateNonClassicalData : no spec for GAUSS on StructureElement with more than one cell !";
620   std::size_t sz(_locs.size());
621   std::vector< MCAuto<DataArrayDouble> > arrs(sz);
622   for(std::size_t i=0;i<sz;i++)
623     {
624       const MEDFileFieldLoc& loc(globs->getLocalization(_locs[i]));
625       const MEDFileGTKeeper *gtk(loc.getUndergroundGTKeeper());
626       const MEDFileGTKeeperDyn *gtk2(dynamic_cast<const MEDFileGTKeeperDyn *>(gtk));
627       if(!gtk2)
628         throw INTERP_KERNEL::Exception("LocInfo::generateNonClassicalData : internal error !");
629       const MEDFileUMesh *meshLoc(gtk2->getMesh()),*section(gtk2->getSection());
630       const MEDFileStructureElement *se(gtk2->getSE());
631       INTERP_KERNEL::NormalizedCellType gt;
632       {
633         std::vector<int> nel(meshLoc->getNonEmptyLevels());
634         if(nel.size()!=1)
635           throw INTERP_KERNEL::Exception(MSG1);
636         if(nel[0]!=0)
637           throw INTERP_KERNEL::Exception(MSG1);
638         MCAuto<MEDCouplingUMesh> um(meshLoc->getMeshAtLevel(0));
639         if(um->getNumberOfCells()!=1)
640           throw INTERP_KERNEL::Exception(MSG1);
641         gt=um->getTypeOfCell(0);
642         std::vector<int> v;
643         um->getNodeIdsOfCell(0,v);
644         std::size_t sz2(v.size());
645         for(std::size_t j=0;j<sz2;j++)
646           if(v[j]!=j)
647             throw INTERP_KERNEL::Exception(MSG1);
648       }
649       const std::vector< MCAuto<MEDFileEltStruct4Mesh> >& strs(mesh->getAccessOfUndergroundEltStrs());
650       MCAuto<MEDFileEltStruct4Mesh> zeStr;
651       for(std::vector< MCAuto<MEDFileEltStruct4Mesh> >::const_iterator it=strs.begin();it!=strs.end();it++)
652         {
653           if((*it)->getGeoTypeName()==se->getName())
654             {
655               zeStr=*it;
656               break;
657             }
658         }
659       if(zeStr.isNull())
660         {
661           std::ostringstream oss; oss << "LocInfo::generateNonClassicalData :  : no geo type with name " <<  se->getName() << " in " << mesh->getName() << " !";
662           throw INTERP_KERNEL::Exception(oss.str());
663         }
664       arrs[i]=BuildMeshFromStructure(gt,_pfl[i],loc,zeStr,mesh,section,globs);
665     }
666   std::vector<const DataArrayDouble *> arrs2(VecAutoToVecOfCstPt(arrs));
667   MCAuto<DataArrayDouble> resu(DataArrayDouble::Aggregate(arrs2));
668   MCAuto<MEDFileUMesh> ret(MEDFileUMesh::New());
669   ret->setCoords(resu);
670   std::ostringstream meshName; meshName << mesh->getName() << "_on_" << sz << "_sections" << "_" << zePos;
671   ret->setName(meshName.str());
672   return ret;
673 }
674
675 FieldWalker2::FieldWalker2(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd)
676 {
677   _loc=pmptpd->getLocalization();
678   _pfl=pmptpd->getProfile();
679   _is_classic=pmptpd->getType()!=ON_GAUSS_PT;
680   _pd=SlicePartDefinition::New(pmptpd->getStart(),pmptpd->getEnd(),1);
681 }
682
683 bool FieldWalker2::operator!=(const FieldWalker2& other) const
684 {
685   return !((*this)==other);
686 }
687
688 bool FieldWalker2::operator==(const FieldWalker2& other) const
689 {
690   bool ret2(false);
691   {
692     std::string tmp;
693     ret2=_pd->isEqual(other._pd,tmp);
694   }
695   bool ret(_loc==other._loc && _pfl==other._pfl && _is_classic==other._is_classic && ret2);
696   return ret;
697 }
698
699 class FieldWalker1
700 {
701 public:
702   FieldWalker1(const MEDFileAnyTypeField1TSWithoutSDA *ts):_ts(ts),_pm_pt(0),_nb_mesh(0) { }
703   void newMeshEntry(const MEDFileFieldPerMesh *fpm);
704   void endMeshEntry(const MEDFileFieldPerMesh *fpm) { }
705   void newPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt);
706   void endPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt);
707   void newPerMeshPerTypePerDisc(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd);
708   void checkOK(const FieldWalker1& other) const;
709   bool isClassical() const;
710   std::vector<FieldWalker2> getNonClassicalData() const { return _fw; }
711 private:
712   const MEDFileAnyTypeField1TSWithoutSDA *_ts;
713   const MEDFileFieldPerMeshPerTypeDyn *_pm_pt;
714   std::vector<FieldWalker2> _fw;
715   int _nb_mesh;
716 };
717
718 void FieldWalker1::newMeshEntry(const MEDFileFieldPerMesh *fpm)
719 {
720   if(_nb_mesh++==1)
721     throw INTERP_KERNEL::Exception("FieldWalker1::newMeshEntry : multi mesh not supported !");
722 }
723
724 void FieldWalker1::newPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt)
725 {
726   if(_pm_pt)
727     throw INTERP_KERNEL::Exception("FieldWalker1::newPerMeshPerTypeEntry : multi SE loc not managed yet !");
728   const MEDFileFieldPerMeshPerTypeDyn *pmpt2(dynamic_cast<const MEDFileFieldPerMeshPerTypeDyn *>(pmpt));
729   if(!pmpt2)
730     throw INTERP_KERNEL::Exception("newPerMeshPerTypeEntry : internal error !");
731   _pm_pt=pmpt2;
732 }
733
734 void FieldWalker1::endPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *)
735 {
736   isClassical();
737 }
738
739 void FieldWalker1::newPerMeshPerTypePerDisc(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd)
740 {
741   _fw.push_back(FieldWalker2(pmptpd));
742 }
743
744 void FieldWalker1::checkOK(const FieldWalker1& other) const
745 {
746   std::size_t sz(_fw.size());
747   if(other._fw.size()!=sz)
748     throw INTERP_KERNEL::Exception("checkOK : not OK because size are not the same !");
749   for(std::size_t i=0;i<sz;i++)
750     if(_fw[i]!=other._fw[i])
751       throw INTERP_KERNEL::Exception("checkOK : not OK because an element mismatches !");
752 }
753
754 bool FieldWalker1::isClassical() const
755 {
756   if(_fw.empty())
757     throw INTERP_KERNEL::Exception("FieldWalker1::endPerMeshPerTypeEntry : internal error !");
758   std::size_t ic(0),inc(0);
759   for(std::vector<FieldWalker2>::const_iterator it=_fw.begin();it!=_fw.end();it++)
760     {
761       if((*it).isClassic())
762         ic++;
763       else
764         inc++;
765     }
766   if(ic!=0 && inc!=0)
767     throw INTERP_KERNEL::Exception("FieldWalker1::endPerMeshPerTypeEntry : mix is not allowed yet !");
768   return inc==0;
769 }
770
771 class FieldWalker
772 {
773 public:
774   FieldWalker(const MEDFileAnyTypeFieldMultiTSWithoutSDA *f):_f(f) { }
775   void newTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts);
776   void endTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts);
777   void newMeshEntry(const MEDFileFieldPerMesh *fpm);
778   void endMeshEntry(const MEDFileFieldPerMesh *fpm);
779   void newPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt);
780   void endPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt);
781   void newPerMeshPerTypePerDisc(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd);
782 public:
783   bool isEmpty() const;
784   bool isClassical() const;
785   const MEDFileAnyTypeFieldMultiTSWithoutSDA *field() const { return _f; }
786   std::vector<FieldWalker2> getNonClassicalData() const { return _fw_prev->getNonClassicalData(); }
787 private:
788   const MEDFileAnyTypeFieldMultiTSWithoutSDA *_f;
789   mutable INTERP_KERNEL::AutoCppPtr<FieldWalker1> _fw;
790   mutable INTERP_KERNEL::AutoCppPtr<FieldWalker1> _fw_prev;
791 };
792
793 bool FieldWalker::isEmpty() const
794 {
795   return _fw_prev.isNull();
796 }
797
798 bool FieldWalker::isClassical() const
799 {
800   return _fw_prev->isClassical();
801 }
802
803 void FieldWalker::newTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts)
804 {
805   _fw=new FieldWalker1(ts);
806 }
807
808 void FieldWalker::endTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts)
809 {
810   if(_fw_prev.isNull())
811     _fw_prev=new FieldWalker1(*_fw);
812   else
813     _fw_prev->checkOK(*_fw);
814   _fw=0;
815 }
816
817 void FieldWalker::newMeshEntry(const MEDFileFieldPerMesh *fpm)
818 {
819   _fw->newMeshEntry(fpm);
820 }
821
822 void FieldWalker::endMeshEntry(const MEDFileFieldPerMesh *fpm)
823 {
824   _fw->endMeshEntry(fpm);
825 }
826
827 void FieldWalker::newPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt)
828 {
829   _fw->newPerMeshPerTypeEntry(pmpt);
830 }
831
832 void FieldWalker::endPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt)
833 {
834   _fw->endPerMeshPerTypeEntry(pmpt);
835 }
836
837 void FieldWalker::newPerMeshPerTypePerDisc(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd)
838 {
839   _fw->newPerMeshPerTypePerDisc(pmptpd);
840 }
841
842 // this class splits fields into same
843 class LocSpliter : public MEDFileFieldVisitor
844 {
845 public:
846   LocSpliter(const MEDFileFieldGlobsReal *globs):_globs(globs),_fw(0) { }
847   MCAuto<MEDFileFields> getClassical() const { return _classical; }
848   void generateNonClassicalData(const MEDFileUMesh *mesh, std::vector< MCAuto<MEDFileFields> >& outFields, std::vector< MCAuto<MEDFileUMesh> >& outMeshes) const;
849 private:
850   void newFieldEntry(const MEDFileAnyTypeFieldMultiTSWithoutSDA *field);
851   void endFieldEntry(const MEDFileAnyTypeFieldMultiTSWithoutSDA *field);
852   //
853   void newTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts);
854   void endTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts);
855   //
856   void newMeshEntry(const MEDFileFieldPerMesh *fpm);
857   void endMeshEntry(const MEDFileFieldPerMesh *fpm);
858   //
859   void newPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt);
860   void endPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt);
861   //
862   void newPerMeshPerTypePerDisc(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd);
863 private:
864   const MEDFileFieldGlobsReal *_globs;
865   std::vector< LocInfo > _locs;
866   std::vector< MCAuto<MEDFileFields> > _fields_on_locs;//size of _locs== size of _fields_on_locs
867   MCAuto<MEDFileFields> _classical;
868 private:
869   mutable INTERP_KERNEL::AutoCppPtr<FieldWalker> _fw;
870 };
871
872 void LocSpliter::newFieldEntry(const MEDFileAnyTypeFieldMultiTSWithoutSDA *field)
873 {
874   _fw=new FieldWalker(field);
875 }
876
877 void LocSpliter::endFieldEntry(const MEDFileAnyTypeFieldMultiTSWithoutSDA *field)
878 {
879   if(_fw->isEmpty())
880     return ;
881   MCAuto<MEDFileAnyTypeFieldMultiTS> f(MEDFileAnyTypeFieldMultiTS::BuildNewInstanceFromContent(const_cast<MEDFileAnyTypeFieldMultiTSWithoutSDA *>(field)));
882   if(_fw->isClassical())
883     {
884       if(_classical.isNull())
885         {
886           _classical=MEDFileFields::New();
887           _classical->shallowCpyGlobs(*_globs);
888         }
889       _classical->pushField(f);
890     }
891   else
892     {
893       std::vector<FieldWalker2> fw2(_fw->getNonClassicalData());
894       LocInfo elt(fw2);
895       std::vector< LocInfo >::iterator it(std::find(_locs.begin(),_locs.end(),elt));
896       if(it==_locs.end())
897         {
898           _locs.push_back(elt);
899           MCAuto<MEDFileFields> zeF(MEDFileFields::New());
900           zeF->shallowCpyGlobs(*_globs);
901           zeF->pushField(f);
902           _fields_on_locs.push_back(zeF);
903         }
904       else
905         {
906           MCAuto<MEDFileFields> zeF(_fields_on_locs[std::distance(_locs.begin(),it)]);
907           zeF->pushField(f);
908         }
909     }
910 }
911
912 void LocSpliter::generateNonClassicalData(const MEDFileUMesh *mesh, std::vector< MCAuto<MEDFileFields> >& outFields, std::vector< MCAuto<MEDFileUMesh> >& outMeshes) const
913 {
914   int i(0);
915   for(std::vector<LocInfo>::const_iterator it=_locs.begin();it!=_locs.end();it++,i++)
916     {
917       MCAuto<MEDFileUMesh> m((*it).generateNonClassicalData(i,mesh,_globs));
918       outMeshes.push_back(m);
919       MCAuto<MEDCouplingUMesh> mcm(MEDCouplingUMesh::Build0DMeshFromCoords(m->getCoords()));
920       mcm->setName(m->getName());
921       MCAuto<MEDFileFields> fs(_fields_on_locs[i]);
922       MCAuto<MEDFileFields> outFs(MEDFileFields::New());
923       for(int j=0;j<fs->getNumberOfFields();j++)
924         {
925           MCAuto<MEDFileAnyTypeFieldMultiTS> fmtsNC(fs->getFieldAtPos(j));
926           MCAuto<MEDFileFieldMultiTS> fmts(DynamicCastSafe<MEDFileAnyTypeFieldMultiTS,MEDFileFieldMultiTS>(fmtsNC));
927           MCAuto<MEDFileFieldMultiTS> outFmts(MEDFileFieldMultiTS::New());
928           for(int k=0;k<fmts->getNumberOfTS();k++)
929             {
930               MCAuto<MEDFileField1TS> outF1t(MEDFileField1TS::New());
931               MCAuto<MEDFileField1TS> f1ts(fmts->getTimeStepAtPos(k));
932               int t2,t3;
933               double t1(f1ts->getTime(t2,t3));
934               MCAuto<MEDCouplingFieldDouble> mcf(MEDCouplingFieldDouble::New(ON_NODES));
935               MCAuto<DataArrayDouble> arr,arr2;
936               arr.takeRef(f1ts->getUndergroundDataArray());
937               arr2=arr->selectPartDef((*it).getPartDef());
938               mcf->setArray(arr2);
939               mcf->setTime(t1,t2,t3);
940               mcf->setName(f1ts->getName());
941               mcf->setMesh(mcm);
942               outF1t->setFieldNoProfileSBT(mcf);
943               outFmts->pushBackTimeStep(outF1t);
944             }
945           outFs->pushField(outFmts);
946         }
947       outFields.push_back(outFs);
948     }
949 }
950
951 void LocSpliter::newTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts)
952 {
953   _fw->newTimeStepEntry(ts);
954 }
955
956 void LocSpliter::endTimeStepEntry(const MEDFileAnyTypeField1TSWithoutSDA *ts)
957 {
958   _fw->endTimeStepEntry(ts);
959 }
960
961 void LocSpliter::newMeshEntry(const MEDFileFieldPerMesh *fpm)
962 {
963   _fw->newMeshEntry(fpm);
964 }
965
966 void LocSpliter::endMeshEntry(const MEDFileFieldPerMesh *fpm)
967 {
968   _fw->endMeshEntry(fpm);
969 }
970
971 void LocSpliter::newPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt)
972 {
973   _fw->newPerMeshPerTypeEntry(pmpt);
974 }
975
976 void LocSpliter::endPerMeshPerTypeEntry(const MEDFileFieldPerMeshPerTypeCommon *pmpt)
977 {
978   _fw->endPerMeshPerTypeEntry(pmpt);
979 }
980
981 void LocSpliter::newPerMeshPerTypePerDisc(const MEDFileFieldPerMeshPerTypePerDisc *pmptpd)
982 {
983   _fw->newPerMeshPerTypePerDisc(pmptpd);
984 }
985
986 void MEDFileBlowStrEltUp::DealWithConflictNames(MEDFileAnyTypeFieldMultiTS *fmtsToAdd, const MEDFileFields *fs)
987 {
988   std::vector<std::string> fnames(fs->getFieldsNames());
989   for(int i=0;i<1000;i++)
990     {
991       std::ostringstream oss; oss << fmtsToAdd->getName();
992       if(i>=1)
993         oss << "_" << i-1;
994       if(std::find(fnames.begin(),fnames.end(),oss.str())==fnames.end())
995         {
996           fmtsToAdd->setName(oss.str());
997           return ;
998         }
999     }
1000   throw INTERP_KERNEL::Exception("DealWithConflictNames : Eh eh interesting !");
1001 }
1002
1003 MCAuto<MEDFileFields> MEDFileBlowStrEltUp::splitFieldsPerLoc(const MEDFileFields *fields, const MEDFileUMesh *mesh, MEDFileMeshes *msOut, MEDFileFields *allZeOutFields)
1004 {
1005   LocSpliter ls(fields);
1006   fields->accept(ls);
1007   std::vector< MCAuto<MEDFileFields> > outFields;
1008   std::vector< MCAuto<MEDFileUMesh> > outMeshes;
1009   ls.generateNonClassicalData(mesh,outFields,outMeshes);
1010   for(std::vector< MCAuto<MEDFileFields> >::iterator it=outFields.begin();it!=outFields.end();it++)
1011     {
1012       for(int j=0;j<(*it)->getNumberOfFields();j++)
1013         {
1014           MCAuto<MEDFileAnyTypeFieldMultiTS> fmts((*it)->getFieldAtPos(j));
1015           //DealWithConflictNames(fmts,allZeOutFields);// uncomment to have a writable data structure
1016           allZeOutFields->pushField(fmts);
1017         }
1018     }
1019   for(std::vector< MCAuto<MEDFileUMesh> >::iterator it=outMeshes.begin();it!=outMeshes.end();it++)
1020     msOut->pushMesh(*it);
1021   return ls.getClassical();
1022 }