Salome HOME
Merge branch V7_3_1_BR
[plugins/ghs3dprlplugin.git] / src / GHS3DPRLPlugin / GHS3DPRLPlugin_GHS3DPRL.cxx
1 // Copyright (C) 2007-2014  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
20 // ---
21 // File   : GHS3DPRLPlugin_GHS3DPRL.cxx
22 // Author : Christian VAN WAMBEKE (CEA) (from Hexotic plugin Lioka RAZAFINDRAZAKA)
23 // ---
24 //
25 #include "GHS3DPRLPlugin_GHS3DPRL.hxx"
26 #include "GHS3DPRLPlugin_Hypothesis.hxx"
27
28 #include <SMDS_MeshElement.hxx>
29 #include <SMDS_MeshNode.hxx>
30 #include <SMESH_subMesh.hxx>
31
32 #include <TopExp_Explorer.hxx>
33 #include <OSD_File.hxx>
34
35 #include "utilities.h"
36
37 #ifndef WIN32
38 #include <sys/sysinfo.h>
39 #endif
40
41 #ifdef _DEBUG_
42 #define DUMP(txt) \
43 //  cout << txt
44 #else
45 #define DUMP(txt)
46 #endif
47
48 #include <SMESH_Gen.hxx>
49 #include <SMESHDS_Mesh.hxx>
50 #include <SMESH_ControlsDef.hxx>
51
52 #include <list>
53 #include <Standard_ProgramError.hxx>
54 #include <TCollection_AsciiString.hxx>
55 #include <TopTools_MapOfShape.hxx>
56 #include <TopoDS.hxx>
57 #include <BRepGProp.hxx>
58 #include <GProp_GProps.hxx>
59
60 /* 
61 extern "C" {
62 #include <med.h>
63 //#include <med_config.h>
64 #include <med_utils.h>
65 //#include <med_misc.h>
66 #include <stdlib.h>
67 }
68 using namespace med_2_2;*/
69
70 static void removeFile( const TCollection_AsciiString& fileName )
71 {
72   try {
73     OSD_File( fileName ).Remove();
74   }
75   catch ( Standard_ProgramError ) {
76     MESSAGE("Can't remove file: " << fileName.ToCString() << " ; file does not exist or permission denied");
77   }
78 }
79
80 //=============================================================================
81 GHS3DPRLPlugin_GHS3DPRL::GHS3DPRLPlugin_GHS3DPRL(int hypId, int studyId, SMESH_Gen* gen)
82   : SMESH_3D_Algo(hypId, studyId, gen)
83 {
84   MESSAGE("GHS3DPRLPlugin_GHS3DPRL::GHS3DPRLPlugin_GHS3DPRL");
85   _name = "GHS3DPRL_3D";
86   _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID);// 1 bit /shape type
87   _countSubMesh=0;
88   _nodeRefNumber=0;
89   _compatibleHypothesis.push_back("GHS3DPRL_Parameters");
90 }
91
92 //=============================================================================
93 GHS3DPRLPlugin_GHS3DPRL::~GHS3DPRLPlugin_GHS3DPRL()
94 {
95   MESSAGE("GHS3DPRLPlugin_GHS3DPRL::~GHS3DPRLPlugin_GHS3DPRL");
96 }
97
98 //=============================================================================
99 bool GHS3DPRLPlugin_GHS3DPRL::CheckHypothesis
100                          (SMESH_Mesh& aMesh,
101                           const TopoDS_Shape& aShape,
102                           SMESH_Hypothesis::Hypothesis_Status& aStatus)
103 {
104   //MESSAGE("GHS3DPRLPlugin_GHS3DPRL::CheckHypothesis");
105   _hypothesis = NULL;
106
107   list<const SMESHDS_Hypothesis*>::const_iterator itl;
108   const SMESHDS_Hypothesis* theHyp;
109
110   const list<const SMESHDS_Hypothesis*>& hyps = GetUsedHypothesis(aMesh, aShape);
111   int nbHyp = hyps.size();
112   if (!nbHyp)
113   {
114     aStatus = SMESH_Hypothesis::HYP_OK;
115     return true;  // can work with no hypothesis
116   }
117
118   itl = hyps.begin();
119   theHyp = (*itl); // use only the first hypothesis
120
121   string hypName = theHyp->GetName();
122   if (hypName == "GHS3DPRL_Parameters")
123   {
124     _hypothesis = static_cast<const GHS3DPRLPlugin_Hypothesis*> (theHyp);
125     ASSERT(_hypothesis);
126     aStatus = SMESH_Hypothesis::HYP_OK;
127   }
128   else
129     aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
130
131   return aStatus == SMESH_Hypothesis::HYP_OK;
132 }
133
134 //=======================================================================
135 static bool writeGHS3DPRLFiles (const TCollection_AsciiString &  GHS3DPRL_In,
136                                 SMESHDS_Mesh *                   theMesh,
137                                 map <int,int> &                  theSmdsToGHS3DPRLIdMap,
138                                 map <int,const SMDS_MeshNode*> & theGHS3DPRLIdToNodeMap)
139 {
140    bool Ok;
141    int ifam=0;
142    TCollection_AsciiString namefile(GHS3DPRL_In);
143    namefile+=".points";
144    removeFile(namefile);
145    ofstream theFile;
146    theFile.open(namefile.ToCString(),ios::out);
147 #ifdef WIN32
148    Ok = theFile.is_open();
149 #else
150    Ok=theFile.rdbuf()->is_open();
151 #endif
152    if (!Ok)
153    {
154       INFOS("Can't write into "<<namefile.ToCString());
155       return false;
156    }
157    cout<<endl<<"writeGHS3DPRLFiles version 2.1 "<<endl;
158    cout<<endl<<"Creating GHS3DPRL processed mesh file : "<<namefile<<endl;
159
160    int nbVertices=theMesh->NbNodes();
161    int nbFaces=theMesh->NbFaces();        //triangles or quadrangles
162    const char* space="  ";
163    const int dummyint=1;                  //nrs,nsd,refnum=1 (for wrap)
164
165    // Writing SMESH points into GHS3DPRL File.points
166    theFile<<nbVertices<<endl;
167
168    int aSmdsNodeID=1;
169    const SMDS_MeshNode* node_2;
170    SMDS_NodeIteratorPtr itOnNode=theMesh->nodesIterator();
171    //int ifam=100;//test famille
172    theFile.precision(15); theFile.setf(ios::scientific,ios::floatfield);
173    //cout<<"set precision 15 on float\n";
174    while (itOnNode->more())
175    {
176       node_2 = itOnNode->next();
177       theSmdsToGHS3DPRLIdMap.insert(map <int,int>::value_type(node_2->GetID(),aSmdsNodeID));
178       theGHS3DPRLIdToNodeMap.insert(map <int,const SMDS_MeshNode*>::value_type(aSmdsNodeID,node_2));
179       theFile<<node_2->X()<<space<<node_2->Y()<<space<<node_2->Z()<<space<<ifam<<endl;
180       aSmdsNodeID++;
181       //if (aSmdsNodeID==11) ifam++;
182    }
183    //no specified points;
184    theFile.close();
185
186    namefile=GHS3DPRL_In+".faces";
187    removeFile(namefile);
188    theFile.open(namefile.ToCString(),ios::out);
189 #ifdef WIN32
190    Ok=theFile.is_open();
191 #else
192    Ok=theFile.rdbuf()->is_open();
193 #endif
194    if (!Ok)
195    {
196       INFOS("Can't write into "<<namefile.ToCString());
197       return false;
198    }
199    cout<<endl<<"Creating GHS3DPRL processed mesh file : "<<namefile<<endl;
200
201    // Writing SMESH faces into GHS3DPRL File.faces
202    theFile<<nbFaces<<" 0"<<endl;   //NB_ELEMS DUMMY_INT
203                                    //" 0" is a reserved parameter
204
205    const SMDS_MeshElement* aFace;
206    map<int,int>::const_iterator itOnSmdsNode;
207    SMDS_ElemIteratorPtr itOnFaceNode;
208    SMDS_FaceIteratorPtr itOnSmdsFace = theMesh->facesIterator();
209    long nbNoTriangles=0;
210    int ifaces=0;
211    //ifam=300;
212    while (itOnSmdsFace->more())
213    {
214       aFace=itOnSmdsFace->next();
215       itOnFaceNode=aFace->nodesIterator();
216       const int nbNodes=aFace->NbNodes();
217       if (nbNodes!=3) nbNoTriangles++;
218       ifaces++;
219       theFile<<nbNodes<<space;        // NB_NODES
220       while (itOnFaceNode->more())
221       {
222           aSmdsNodeID=itOnFaceNode->next()->GetID();
223           itOnSmdsNode=theSmdsToGHS3DPRLIdMap.find(aSmdsNodeID);
224           ASSERT(itOnSmdsNode!=theSmdsToGHS3DPRLIdMap.end());
225           theFile<<space<<(*itOnSmdsNode).second; //NODE_1 NODE_2 ...
226       }
227       //(NB_NODES+1) times: DUMMY_INT
228       //if (ifaces==11) ifam++;
229       theFile<<space<<ifam;
230       for ( int i=1; i<=nbNodes; i++) theFile<<space<<200+i;
231       theFile<<endl;
232    }
233    theFile.close();
234
235    cout<<"Processed mesh files created, they contains :\n";
236    cout<<"    "<<nbVertices<<" vertices\n";
237    if (nbNoTriangles==0)
238       cout<<"    "<<nbFaces<<" faces\n\n";
239    else
240       cout<<"    "<<nbFaces<<" faces with "<<nbNoTriangles<<"faces no triangles\n\n";
241    return true;
242 }
243
244 //=======================================================================
245 static bool getInt( int & theValue, char * & theLine )
246 {
247   char *ptr;
248   theValue = strtol( theLine, &ptr, 10 );
249   if ( ptr == theLine ||
250       // there must not be neither '.' nor ',' nor 'E' ...
251       (*ptr != ' ' && *ptr != '\n' && *ptr != '\0'))
252     return false;
253
254   DUMP( "  " << theValue );
255   theLine = ptr;
256   return true;
257 }
258
259 //=======================================================================
260 static bool getDouble( double & theValue, char * & theLine )
261 {
262   char *ptr;
263   theValue = strtod( theLine, &ptr );
264   if ( ptr == theLine )
265     return false;
266
267   DUMP( "   " << theValue );
268   theLine = ptr;
269   return true;
270 }
271
272 //=======================================================================
273
274 #define GHS3DPRLPlugin_BUFLENGTH 256
275 #define GHS3DPRLPlugin_ReadLine(aPtr,aBuf,aFile,aLineNb) \
276 {  aPtr = fgets( aBuf, GHS3DPRLPlugin_BUFLENGTH - 2, aFile ); aLineNb++; DUMP(endl); }
277
278 //=======================================================================
279 static bool readResult(FILE *                           theFile,
280                        SMESHDS_Mesh *                   theMesh,
281                        const TopoDS_Shape &             theShape,
282                        map <int,const SMDS_MeshNode*> & theGHS3DPRLIdToNodeMap,
283                        const TCollection_AsciiString &  GHS3DPRL_Out,
284                        int &                            nodeRefNumber)
285 {
286   // ---------------------------------
287   // Read generated elements and nodes
288   // ---------------------------------
289
290   cout << "Reading GHS3DPRL output file : " << GHS3DPRL_Out << endl;
291   cout << endl;
292
293   char aBuffer[ GHS3DPRLPlugin_BUFLENGTH ];
294   char * aPtr;
295   int aLineNb = 0;
296   int shapeID = theMesh->ShapeToIndex( theShape );
297
298   int line = 1, EndOfFile = 0, nbElem = 0, nField = 10, nbRef = 0, aGHS3DPRLNodeID = 0;
299   const char * theField;
300
301   vector<const char*> tabField = vector<const char*>(nField);
302   vector<int> tabRef = vector<int>(nField);
303
304   tabField[0] = "MeshVersionFormatted";    tabRef[0] = 0;
305   tabField[1] = "Dimension";               tabRef[1] = 0;
306   tabField[2] = "Vertices";                tabRef[2] = 3;
307   tabField[3] = "Edges";                   tabRef[3] = 2;
308   tabField[4] = "Triangles";               tabRef[4] = 3;
309   tabField[5] = "Quadrilaterals";          tabRef[5] = 4;
310   tabField[6] = "Hexahedra";               tabRef[6] = 8;
311   tabField[7] = "Corners";                 tabRef[7] = 1;
312   tabField[8] = "Ridges";                  tabRef[0] = 1;
313   tabField[9] = "End";                     tabRef[0] = 0;
314
315   nodeRefNumber += theMesh->NbNodes();
316
317   SMDS_NodeIteratorPtr itOnGHS3DPRLInputNode = theMesh->nodesIterator();
318   while ( itOnGHS3DPRLInputNode->more() )
319       theMesh->RemoveNode( itOnGHS3DPRLInputNode->next() );
320
321   while ( EndOfFile == 0  ) {
322       GHS3DPRLPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
323       for ( int iField = 0; iField < nField; iField++ ) {
324           stringstream theMessage;
325           theField = tabField[iField];
326           if ( strncmp(aPtr, theField, strlen(theField)) == 0 ) {
327               if ( strcmp(theField, "End") == 0 ) {
328                   EndOfFile = 1;
329                   theMessage << "End of GHS3DPRL output file has been reached";
330               }
331               else {
332                   GHS3DPRLPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
333                   line++;
334                   getInt( nbElem, aPtr );
335
336                   if ( strcmp(theField, "MeshVersionFormatted") == 0 )
337                       theMessage << "GHS3DPRL mesh descriptor : " << theField << " " << nbElem;
338                   else if ( strcmp(theField, "Dimension") == 0 )
339                       theMessage << "GHS3DPRL mesh of " << nbElem << "D dimension";
340                   else if ( strcmp(theField, "Vertices")       == 0 ||
341                             strcmp(theField, "Edges")          == 0 ||
342                             strcmp(theField, "Quadrilaterals") == 0 ||
343                             strcmp(theField, "Hexahedra")      == 0 ) {
344                       nbRef = tabRef[iField];
345                       GHS3DPRLPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );           // read blank line
346
347                       if ( strcmp(theField, "Vertices") == 0 ) {
348                           int aGHS3DPRLID;
349                           vector<double> coord = vector<double>(nbRef);
350                           SMDS_MeshNode * aGHS3DPRLNode;
351
352                           for ( int iElem = 0; iElem < nbElem; iElem++ ) {
353                               aGHS3DPRLID = iElem + 1 + nodeRefNumber;
354                               GHS3DPRLPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );   // read file lines
355                               for ( int iCoord = 0; iCoord < 3; iCoord++ )
356                                   getDouble ( coord[ iCoord ], aPtr );
357                               aGHS3DPRLNode = theMesh->AddNode(coord[0], coord[1], coord[2]);
358                               theMesh->SetNodeInVolume( aGHS3DPRLNode, shapeID );
359                               theGHS3DPRLIdToNodeMap[ aGHS3DPRLID ] = aGHS3DPRLNode;
360                           }
361                       }
362                       else {
363                           vector<const SMDS_MeshNode*> node = vector<const SMDS_MeshNode*>(nbRef);
364                           SMDS_MeshElement* aGHS3DPRLElement;
365                           map <int,const SMDS_MeshNode*>::iterator itOnGHS3DPRLNode;
366
367                           for ( int iElem = 0; iElem < nbElem; iElem++ ) {
368                               GHS3DPRLPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );   // read file lines
369                               for ( int iRef = 0; iRef < nbRef; iRef++ ) {
370                                   getInt ( aGHS3DPRLNodeID, aPtr );                         // read nbRef aGHS3DPRLNodeID
371                                   aGHS3DPRLNodeID += nodeRefNumber;
372                                   itOnGHS3DPRLNode = theGHS3DPRLIdToNodeMap.find( aGHS3DPRLNodeID );
373                                   node[ iRef ] = itOnGHS3DPRLNode->second;
374                               }
375
376                               if ( strcmp(theField, "Edges") == 0 )                        // create an element
377                                   aGHS3DPRLElement = theMesh->AddEdge( node[0], node[1] );
378                               else if ( strcmp(theField, "Quadrilaterals") == 0 )
379                                   aGHS3DPRLElement = theMesh->AddFace( node[0], node[1], node[2], node[3] );
380                               else if ( strcmp(theField, "Hexahedra") == 0 )
381                                   aGHS3DPRLElement = theMesh->AddVolume( node[0], node[1], node[2], node[3], node[4], node[5], node[6], node[7] );
382
383                               theMesh->SetMeshElementOnShape( aGHS3DPRLElement, shapeID );
384                           }
385                       }
386                       theMessage << nbElem << " " << theField << " created";
387                   }
388               }
389               if ( theMessage.str().size() != 0 ) {
390                   cout << theMessage.str() << endl;
391                   break;
392               }
393           }
394       }
395   }
396   cout << endl;
397   return true;
398 }
399
400 //=============================================================================
401 // Pass parameters to GHS3DPRL
402 void GHS3DPRLPlugin_GHS3DPRL::SetParameters(const GHS3DPRLPlugin_Hypothesis* hyp)
403 {
404   if (hyp) {
405     MESSAGE("GHS3DPRLPlugin_GHS3DPRL::SetParameters");
406     _MEDName = hyp->GetMEDName();  //"DOMAIN\0"
407     _NbPart = hyp->GetNbPart();
408     _KeepFiles = hyp->GetKeepFiles();
409     _Background = hyp->GetBackground();
410     _ToMeshHoles = hyp->GetToMeshHoles();
411   }
412 }
413
414 //=======================================================================
415 //before launching salome
416 //SALOME_TMP_DIR (for keep tepal intermediates files) could be set in user's directories
417 static TCollection_AsciiString getTmpDir()
418 {
419   TCollection_AsciiString aTmpDir;
420   char *Tmp_dir = getenv("SALOME_TMP_DIR");
421   if (Tmp_dir == NULL) Tmp_dir = getenv("TMP");
422   if(Tmp_dir != NULL)
423   {
424     aTmpDir = Tmp_dir;
425 #ifdef WIN32
426       if(aTmpDir.Value(aTmpDir.Length()) != '\\') aTmpDir+='\\';
427 #else
428       if(aTmpDir.Value(aTmpDir.Length()) != '/') aTmpDir+='/';
429 #endif
430   }
431   else
432   {
433 #ifdef WIN32
434       aTmpDir = TCollection_AsciiString("C:\\");
435 #else
436       aTmpDir = TCollection_AsciiString("/tmp/");
437 #endif
438   }
439   return aTmpDir;
440 }
441
442 //=============================================================================
443 // Here we are going to use the GHS3DPRL mesher
444 bool GHS3DPRLPlugin_GHS3DPRL::Compute(SMESH_Mesh& theMesh,
445                                      const TopoDS_Shape& theShape)
446 {
447    bool Ok;
448    TCollection_AsciiString pluginerror("ghs3dprl: ");
449    SMESHDS_Mesh* meshDS = theMesh.GetMeshDS();
450    //cout<<"GetMeshDS done\n";
451    if (_countSubMesh==0){
452       MESSAGE("GHS3DPRLPlugin_GHS3DPRL::Compute");
453       _countTotal=0;
454       TopExp_Explorer expf(meshDS->ShapeToMesh(), TopAbs_SOLID);
455       for ( ; expf.More(); expf.Next() ) _countTotal++;
456    }
457    _countSubMesh++;
458    //cout<<"Compute _countSubMesh "<<_countSubMesh<<endl;
459    //no stuff if multiples submesh, multiple call compute
460    //mesh all in one pass tepal (the last)
461    if (_countSubMesh != _countTotal ) return true;
462
463    //stuff on last call 
464    if (_hypothesis==NULL){
465       pluginerror += "No existing parameters/hypothesis for GHS3DPRL";
466       cout <<"\n"<<pluginerror<<"\n\n";
467       error(COMPERR_ALGO_FAILED, pluginerror.ToCString() );
468       return false;
469    }
470    SetParameters(_hypothesis);
471    cout << "\n" << _name << " parameters :\n" << endl;
472    cout << "     generic path/name of MED Files = " << _MEDName << endl;
473    cout << "     number of partitions = " << _NbPart << endl;
474    cout << "     mesh holes = " << _ToMeshHoles << endl;
475    cout << "     keep intermediates files (from tepal) = " << _KeepFiles << endl;
476    cout << "     background (from tepal) = " << _Background << "\n\n";
477
478       //string tmpDir=getTmpDir_new();
479       TCollection_AsciiString
480          tmpDir=getTmpDir(),
481          GHS3DPRL_In,GHS3DPRL_Out,GHS3DPRL_Outxml,
482          run_GHS3DPRL("tepal2med "),rm("rm "),run_nokeep_files,
483          NbPart,fileskinmed(""),path,casenamemed;  //_MEDName.c_str());
484
485       casenamemed += (char *)_MEDName.c_str();
486       int n=casenamemed.SearchFromEnd('/');
487       if (n>0) {
488          path=casenamemed.SubString(1,n);
489          casenamemed=casenamemed.SubString(n+1,casenamemed.Length()); 
490       }
491       else
492          path=tmpDir;
493
494       if (casenamemed.Length()>20){
495          casenamemed=casenamemed.SubString(1,20);
496          cerr<<"MEDName truncated (no more 20 characters) = "<<casenamemed<<endl;
497       }
498       cout<<"path="<<path<<endl;
499       cout<<"casenamemed="<<casenamemed<<endl;
500
501       map <int,int> aSmdsToGHS3DPRLIdMap;
502       map <int,const SMDS_MeshNode*> aGHS3DPRLIdToNodeMap;
503       GHS3DPRL_In = path + "GHS3DPRL";
504       GHS3DPRL_Out = path + casenamemed;
505       GHS3DPRL_Outxml = path + casenamemed + ".xml"; //master file
506       NbPart=_NbPart;
507       //tepal2med --casename=/home/whoami/tmp/GHS3DPRL --number=5 --medname=DOMAIN 
508       //          --limitswap=1000 --verbose=0 --test=no --menu=no --launchtepal=yes;
509       
510       run_GHS3DPRL = run_GHS3DPRL +
511                      " --casename=" + GHS3DPRL_In + 
512                      " --number=" + NbPart + 
513                      " --medname=" + GHS3DPRL_Out +
514                      " --launchtepal=yes";
515       if (_ToMeshHoles) run_GHS3DPRL += " --meshholes=yes"; else run_GHS3DPRL += " --meshholes=no";
516       if (_Background) run_GHS3DPRL += " --background=yes"; else run_GHS3DPRL += " --background=no";
517       run_nokeep_files = rm +GHS3DPRL_In + "* " + path + "tepal.log";
518       system( run_nokeep_files.ToCString() ); //clean files
519       run_nokeep_files = rm + GHS3DPRL_In + "* ";
520
521       cout<<"GHS3DPRL command : "<<run_GHS3DPRL.ToCString()<<endl;
522       cout<<"Write files .faces .point ...";
523       GHS3DPRL_Out = path + casenamemed;
524       removeFile( GHS3DPRL_Outxml ); //only the master xml file
525       Ok=writeGHS3DPRLFiles(GHS3DPRL_In, meshDS, aSmdsToGHS3DPRLIdMap, aGHS3DPRLIdToNodeMap);
526       if (Ok) {cout<<" ...done\n";}
527       else {
528          cout<<" ...NOT done\n";
529          pluginerror = pluginerror + "problem writing input tepal files " + GHS3DPRL_In + "[.faces|.points]";
530       }
531
532       //Ecriture dans un fichier MED ?v2.? meme si not Ok
533       //create empty file -> avoid warning message
534       //med_idt fid=MEDouvrir((const char *)fileskinmed.ToCString(),MED_CREATION);
535       //med_err ret=MEDfermer(fid);
536       //fileskinmed=fileskinmed + "cp /home/wambeke/empty.med "+ path + "GHS3DPRL_skin.med";
537       //system( fileskinmed.ToCString() );
538       fileskinmed=path + "GHS3DPRL_skin.med";
539       cout<<"Write file "<<fileskinmed<<"...";
540       theMesh.ExportMED(fileskinmed.ToCString(),"SKIN_INITIAL",true,1);
541       cout<<" ...done\n";
542
543       if (!Ok) {
544          error(COMPERR_ALGO_FAILED, pluginerror.ToCString());
545          return false; //pb abandonne
546       }
547
548       //sometimes it is better to wait flushing files on slow filesystem...
549       system( "sleep 3" );
550       system( run_GHS3DPRL.ToCString() );
551       system( "sleep 3" );
552
553       //cout<<"!!!reprise plugin!!!"<<endl;
554
555       if (_Background) {
556          pluginerror = pluginerror + "backgrounding... plugin is not waiting for output files "+ casenamemed + "_*.med";
557          cout<<pluginerror<<endl;
558          error(COMPERR_ALGO_FAILED, pluginerror.ToCString());
559          return false; //but it is not a problem but if true my message is overwritten
560          //return true; //but it is not a problem, 
561       }
562
563       // read a result, GHS3DPRL_Out is the name of master file (previous xml format)
564       FILE * aResultFile = fopen( GHS3DPRL_Outxml.ToCString(), "r" );
565       if (aResultFile){
566           //Ok = readResult( aResultFile, meshDS, theShape, aGHS3DPRLIdToNodeMap, GHS3DPRL_Out, _nodeRefNumber );
567           Ok = true;
568           Ok = false; //but it is not a problem but if true my message is overwritten
569           fclose(aResultFile);
570           cout<<"GHS3DPRL OK output master file "<<casenamemed<<".xml exist !\n\n";
571           pluginerror = pluginerror + "new tetaedra not in memory, but stored in files "+ casenamemed + "_*.med";
572           cout<<pluginerror<<endl;
573           error(COMPERR_ALGO_FAILED, pluginerror.ToCString());
574           if (!_KeepFiles) system( run_nokeep_files.ToCString() );
575       }
576       else{
577           Ok = false; //it is a problem AND my message is NOT overwritten
578           pluginerror = pluginerror + "output master file " + casenamemed + ".xml do not exist";
579           cout<<pluginerror<<endl;
580           error(COMPERR_ALGO_FAILED, pluginerror.ToCString() );
581           cout<<"GHS3DPRL KO output files "<<GHS3DPRL_Out<<" do not exist ! see intermediates files keeped:\n";
582           TCollection_AsciiString run_list_files("ls -alt ");
583           run_list_files +=  GHS3DPRL_Out + "* " + GHS3DPRL_In + "* " + path + "tepal.log";
584           system( run_list_files.ToCString() );
585           cout<<endl;
586       }
587       _countSubMesh=0;
588
589
590     return Ok;
591    /*pid_t pid = fork();
592    if (pid > 0) {
593       //Processus père
594       cout<<"le pere est la\n";
595       system("echo le pere > lepere.tmp");
596       system("sleep 10");
597       system("ps -edf |grep wambeke > pslepere.tmp");
598       cout<<"le pere return 0\n";
599       return 0; //ok
600       //exit(0);
601    } else if (pid == 0) {
602       //Processus fils
603       cout<<"le fils est la\n";
604       //On rend le fils indépendant de tout terminal
605       setsid();
606       system("sleep 20");
607       system("echo le fils > lefils.tmp");
608       system("sleep 20");
609       system("ps -edf |grep wambeke > pslefils.tmp");
610       cout<<"le fils return 0\n";
611       return 0; //ok
612    } else {
613       //Traitement d'erreur
614       cout<<"ya probleme sur fork()\n";
615       return 1; //ko
616    }*/
617 }
618
619 //=============================================================================
620 ostream & GHS3DPRLPlugin_GHS3DPRL::SaveTo(ostream & save)
621 {
622   return save;
623 }
624
625 //=============================================================================
626 istream & GHS3DPRLPlugin_GHS3DPRL::LoadFrom(istream & load)
627 {
628   return load;
629 }
630
631 //=============================================================================
632 ostream & operator << (ostream & save, GHS3DPRLPlugin_GHS3DPRL & hyp)
633 {
634   return hyp.SaveTo( save );
635 }
636
637 //=============================================================================
638 istream & operator >> (istream & load, GHS3DPRLPlugin_GHS3DPRL & hyp)
639 {
640   return hyp.LoadFrom( load );
641 }
642
643
644 //=============================================================================
645 /*!
646  *  
647  */
648 //=============================================================================
649 bool GHS3DPRLPlugin_GHS3DPRL::Evaluate(SMESH_Mesh& aMesh,
650                                        const TopoDS_Shape& aShape,
651                                        MapShapeNbElems& aResMap)
652 {
653   int nbtri = 0, nbqua = 0;
654   double fullArea = 0.0;
655   for (TopExp_Explorer exp(aShape, TopAbs_FACE); exp.More(); exp.Next()) {
656     TopoDS_Face F = TopoDS::Face( exp.Current() );
657     SMESH_subMesh *sm = aMesh.GetSubMesh(F);
658     MapShapeNbElemsItr anIt = aResMap.find(sm);
659     if( anIt==aResMap.end() ) {
660       SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
661       smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,
662                                             "Submesh can not be evaluated",this));
663       return false;
664     }
665     std::vector<int> aVec = (*anIt).second;
666     nbtri += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
667     nbqua += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
668     GProp_GProps G;
669     BRepGProp::SurfaceProperties(F,G);
670     double anArea = G.Mass();
671     fullArea += anArea;
672   }
673
674   // collect info from edges
675   int nb0d_e = 0, nb1d_e = 0;
676   bool IsQuadratic = false;
677   bool IsFirst = true;
678   TopTools_MapOfShape tmpMap;
679   for (TopExp_Explorer exp(aShape, TopAbs_EDGE); exp.More(); exp.Next()) {
680     TopoDS_Edge E = TopoDS::Edge(exp.Current());
681     if( tmpMap.Contains(E) )
682       continue;
683     tmpMap.Add(E);
684     SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
685     MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
686     std::vector<int> aVec = (*anIt).second;
687     nb0d_e += aVec[SMDSEntity_Node];
688     nb1d_e += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
689     if(IsFirst) {
690       IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
691       IsFirst = false;
692     }
693   }
694   tmpMap.Clear();
695
696   double ELen = sqrt(2.* ( fullArea/(nbtri+nbqua*2) ) / sqrt(3.0) );
697
698   GProp_GProps G;
699   BRepGProp::VolumeProperties(aShape,G);
700   double aVolume = G.Mass();
701   double tetrVol = 0.1179*ELen*ELen*ELen;
702   double CoeffQuality = 0.9;
703   int nbVols = (int)aVolume/tetrVol/CoeffQuality;
704   int nb1d_f = (nbtri*3 + nbqua*4 - nb1d_e) / 2;
705   int nb1d_in = (int) ( nbVols*6 - nb1d_e - nb1d_f ) / 5;
706   std::vector<int> aVec(SMDSEntity_Last);
707   for(int i=0; i<SMDSEntity_Last; i++) aVec[i]=0;
708   if( IsQuadratic ) {
709     aVec[SMDSEntity_Node] = nb1d_in/6 + 1 + nb1d_in;
710     aVec[SMDSEntity_Quad_Tetra] = nbVols - nbqua*2;
711     aVec[SMDSEntity_Quad_Pyramid] = nbqua;
712   }
713   else {
714     aVec[SMDSEntity_Node] = nb1d_in/6 + 1;
715     aVec[SMDSEntity_Tetra] = nbVols - nbqua*2;
716     aVec[SMDSEntity_Pyramid] = nbqua;
717   }
718   SMESH_subMesh *sm = aMesh.GetSubMesh(aShape);
719   aResMap.insert(std::make_pair(sm,aVec));
720
721   return true;
722 }