ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/aster/src/mgopt_posttraitement.cpp
Revision: 603
Committed: Fri Nov 21 18:52:21 2014 UTC (10 years, 5 months ago) by francois
File size: 66954 byte(s)
Log Message: 
ajout du calcul en elasticite 2D contrainte plane et deformation plane a partir d'une coque CAO plane en XY + OT_PARAMETRE gestion des parametres par default

File Contents

# User Rev Content
1 francois 222 #include "gestionversion.h"
2 francois 460 #include "mgopt_posttraitement.h"
3 francois 222 #include "fem_solution.h"
4     #include "fem_maillage.h"
5     #include "mg_maillage.h"
6     #include "mg_gestionnaire.h"
7     #include "mg_triangle_peau.h"
8 francois 335 #include "mg_geometrie_outils.h"
9 francois 412 #include "fem_maillage_outils.h"
10 francois 388 #include "tpl_fonctions_generiques.h"
11 cuillier 438 #include "ot_algorithme_geometrique.h"
12 francois 343 #include "mg_export.h"
13 francois 222 #include <stdio.h>
14 picher 230 #include <stdlib.h>
15     #include <time.h>
16 francois 222 #include <string.h>
17 francois 224 #include <math.h>
18 francois 222 //---------------------------------------------------------------------------
19    
20    
21 francois 224
22 francois 460 MGOPT_POSTTRAITEMENT::MGOPT_POSTTRAITEMENT():affichageactif(0)
23 francois 222 {
24 francois 224 for (int i=0;i<9;i++) etat[i]=0;
25 francois 457 params.ajouter("decoupage",0.,OT_PARAMETRES::DOUBLE,"0. découpage en conservant les mailles entieres 1. découpage par isodensite");
26     params.ajouter("seuil",0.8,OT_PARAMETRES::DOUBLE,"Valeur du seuil de densité pour les 2 découpages");
27     params.ajouter("nomdensextra","densiteextra.sol",OT_PARAMETRES::STRING,"Nom du fichier comportant la solution de la densité extrapolé au noeud");
28     params.ajouter("consimpose",1.,OT_PARAMETRES::DOUBLE,"0. Maillage non design non inclus dans le decoupage par mailles entieres 1. Maillage design inclus");
29     params.ajouter("consoptimise",1.,OT_PARAMETRES::DOUBLE,"0. Maillage optimise non inclus dans le decoupage par mailles entieres 1. Maillage optimise inclus");
30     params.ajouter("consnooptimise",0.,OT_PARAMETRES::DOUBLE,"0. Maillage non optimise non inclus dans le decoupage par mailles entieres 1. Maillage non optimise inclus");
31     params.ajouter("reactivation",1.,OT_PARAMETRES::DOUBLE,"1. Reactivation de mailles isolées sur le long des aretes 0. non reactivation");
32     params.ajouter("chen2005_debut",0.,OT_PARAMETRES::DOUBLE,"Numero d'itération du debut d'application deu lissage de chen2005 0. méthode non utilisee");
33     params.ajouter("chen2005_itermax",10.,OT_PARAMETRES::DOUBLE,"Nombre d'itérations maxiamles du lissage de chen2005");
34     params.ajouter("chen2005_epsilon",0.01,OT_PARAMETRES::DOUBLE,"Critere de convergence du lissage de chen2005");
35 francois 458 params.ajouter("chen2005_sigma",0.1,OT_PARAMETRES::DOUBLE,"Écart type du filtre de lissage de chen2005");
36     params.ajouter("chen2005_filtre",1.,OT_PARAMETRES::DOUBLE,"Choix du filtre de lissage de chen2005 1. Gaussian 2.Laplacien 3.elfallahford");
37 francois 457 params.ajouter("chen2008_debut",0.,OT_PARAMETRES::DOUBLE,"Numero d'itération du debut d'application deu lissage de chen2008 0. méthode non utilisee");
38     params.ajouter("chen2008_itermax",10.,OT_PARAMETRES::DOUBLE,"Nombre d'itérations maxiamles du lissage de chen2008");
39     params.ajouter("chen2008_epsilon",0.01,OT_PARAMETRES::DOUBLE,"Critere de convergence du lissage de chen2008");
40 francois 458 params.ajouter("chen2008_sigma",0.1,OT_PARAMETRES::DOUBLE,"Écart type du filtre de lissage de chen2008");
41 francois 457 params.ajouter("chen2008_gamma",1.,OT_PARAMETRES::DOUBLE,"Vitesse d'avancement de déplacement du lissage de chen2008");
42 francois 458 params.ajouter("chen2008_filtre",1.,OT_PARAMETRES::DOUBLE,"Choix du filtre de lissage de chen2008 1. Gaussian 2.Laplacien 3.elfallahford");
43 francois 457 params.ajouter("jiao2012_debut",0.,OT_PARAMETRES::DOUBLE,"Numero d'itération du debut d'application deu lissage de jiao2012 0. méthode non utilisee");
44     params.ajouter("jiao2012_itermax",10.,OT_PARAMETRES::DOUBLE,"Nombre d'itérations maxiamles du lissage de jioa2012");
45 francois 222 }
46    
47    
48 francois 460 MGOPT_POSTTRAITEMENT::~MGOPT_POSTTRAITEMENT()
49 francois 222 {
50 francois 224 for (int i=0;i<lst_peau.size();i++)
51     delete lst_peau[i];
52 francois 222 }
53 picher 233
54 francois 460 void MGOPT_POSTTRAITEMENT::active_affichage(void (*fonc)(char*))
55 francois 232 {
56     affiche=fonc;
57     affichageactif=1;
58     }
59 francois 222
60 picher 233
61 francois 457
62 francois 461 void MGOPT_POSTTRAITEMENT::posttraite(char *fichier)
63 francois 457 {
64     affiche((char*)"");
65     affiche((char*)"***********************************************************");
66     affiche((char*)"Post-traitement d'un calcul d'optimisation de topologie");
67     affiche((char*)"***********************************************************");
68     affiche((char*)"");
69     affiche((char*)"");
70     affiche((char*)"");
71     affiche((char*)"Écriture d'un fichier de parametres par défaut");
72     params.enregistrer(fichier);
73     }
74    
75    
76    
77 francois 460 void MGOPT_POSTTRAITEMENT::conserve(int origine)
78 francois 224 {
79     etat[(origine-1000)/10]=1;
80     }
81    
82 picher 233
83 francois 460 void MGOPT_POSTTRAITEMENT::copieorigine(FEM_MAILLAGE* mai)
84 picher 231 {
85 francois 309 LISTE_FEM_ELEMENT3::iterator it_tetra;
86     for (FEM_ELEMENT3* tet=mai->get_premier_element3(it_tetra);tet!=NULL;tet=mai->get_suivant_element3(it_tetra))
87 picher 233 {
88     MG_TETRA* mgtet=(MG_TETRA*)tet->get_mg_element_maillage();
89     mgtet->change_origine(mgtet->get_origine());
90     }
91     }
92 francois 457
93    
94 francois 460 void MGOPT_POSTTRAITEMENT::gain_poids(MG_MAILLAGE* mg_mai,MG_GESTIONNAIRE& gest2)
95 picher 233 {
96     LISTE_MG_TETRA::iterator it_tetra;
97     double volume_initial = 0.;
98     double volume_final = 0.;
99     for (MG_TETRA* mgtet=mg_mai->get_premier_tetra(it_tetra);mgtet!=NULL;mgtet=mg_mai->get_suivant_tetra(it_tetra))
100     {
101     int originetet=mgtet->get_origine();
102     MG_NOEUD* noeud1 = mgtet->get_noeud1();
103     double x1 = noeud1->get_x();
104     double y1 = noeud1->get_y();
105     double z1 = noeud1->get_z();
106     MG_NOEUD* noeud2 = mgtet->get_noeud2();
107     double x2 = noeud2->get_x();
108     double y2 = noeud2->get_y();
109     double z2 = noeud2->get_z();
110     MG_NOEUD* noeud3 = mgtet->get_noeud3();
111     double x3 = noeud3->get_x();
112     double y3 = noeud3->get_y();
113     double z3 = noeud3->get_z();
114     MG_NOEUD* noeud4 = mgtet->get_noeud4();
115     double x4 = noeud4->get_x();
116     double y4 = noeud4->get_y();
117     double z4 = noeud4->get_z();
118     double volume_tet = (x2-x1)*((y3-y1)*(z4-z1)-(z3-z1)*(y4-y1))-(y2-y1)*((x3-x1)*(z4-z1)-(z3-z1)*(x4-x1)) + (z2-z1)*((x3-x1)*(y4-y1)-(y3-y1)*(x4-x1));
119     volume_tet = volume_tet/6.;
120     volume_initial = volume_initial + volume_tet;
121     if (originetet != MAILLEUR_AUTO)
122     {
123     volume_final = volume_final + volume_tet;
124     }
125     }
126     }
127 francois 457
128 francois 460 void MGOPT_POSTTRAITEMENT::reactivation(MG_MAILLAGE* mg_mai,MG_GESTIONNAIRE& gest2)
129 picher 233 {
130 francois 457 cout << " Reactivation d'elements manquants sur aretes vives et faces planes" << endl;
131 picher 233 LISTE_MG_TETRA::iterator it_tetra;
132     for (MG_TETRA* mgtet=mg_mai->get_premier_tetra(it_tetra);mgtet!=NULL;mgtet=mg_mai->get_suivant_tetra(it_tetra))
133     {
134     int originetet=mgtet->get_origine();
135     int compteur = 0;
136     if (originetet == MAILLEUR_AUTO)
137     {
138     MG_TRIANGLE* tri1 = mgtet->get_triangle1();
139     MG_TRIANGLE* tri2 = mgtet->get_triangle2();
140     MG_TRIANGLE* tri3 = mgtet->get_triangle3();
141     MG_TRIANGLE* tri4 = mgtet->get_triangle4();
142     int nb_tet1 = tri1->get_lien_tetra()->get_nb();
143     for (int k = 0;k<nb_tet1;k++)
144     {
145     MG_TETRA* tet1 = tri1->get_lien_tetra()->get(k);
146     int ori1 = tet1->get_origine();
147     if ((tet1 != mgtet) && (ori1 == MAILLEUR_AUTO)) compteur++;
148     }
149     int nb_tet2 = tri2->get_lien_tetra()->get_nb();
150     for (int k = 0;k<nb_tet2;k++)
151     {
152     MG_TETRA* tet2 = tri2->get_lien_tetra()->get(k);
153     int ori2 = tet2->get_origine();
154     if ((tet2 != mgtet) && (ori2 == MAILLEUR_AUTO)) compteur++;
155     }
156     int nb_tet3 = tri3->get_lien_tetra()->get_nb();
157     for (int k = 0;k<nb_tet3;k++)
158     {
159     MG_TETRA* tet3 = tri3->get_lien_tetra()->get(k);
160     int ori3 = tet3->get_origine();
161     if ((tet3 != mgtet) && (ori3 == MAILLEUR_AUTO)) compteur++;
162     }
163     int nb_tet4 = tri4->get_lien_tetra()->get_nb();
164     for (int k = 0;k<nb_tet4;k++)
165     {
166     MG_TETRA* tet4 = tri4->get_lien_tetra()->get(k);
167     int ori4 = tet4->get_origine();
168     if ((tet4 != mgtet) && (ori4 == MAILLEUR_AUTO)) compteur++;
169     }
170 picher 248 if (compteur == 0) mgtet->change_origine(OPTIMISE);
171 picher 233 }
172     }
173     }
174    
175    
176 francois 343
177    
178 francois 460 MG_MAILLAGE* MGOPT_POSTTRAITEMENT::extract_skin_par_decoupage(FEM_SOLUTION* sol,double limit,MG_GESTIONNAIRE& gest2,std::string nom)
179 francois 343 {
180     affiche((char*)"Extraction de l'enveloppe par découpage");
181     sol->active_solution(0);
182     FEM_MAILLAGE* mai=sol->get_maillage();
183     affiche((char*)" Extrapolation de la densité aux noeuds");
184     LISTE_FEM_NOEUD::iterator it;
185     for (FEM_NOEUD* no=mai->get_premier_noeud(it);no!=NULL;no=mai->get_suivant_noeud(it))
186     {
187     double nume=0.;
188     double deno=0.;
189     int nb=no->get_lien_element3()->get_nb();
190     int passe=0;
191     for (int i=0;i<nb;i++)
192     {
193     FEM_ELEMENT3* tet=no->get_lien_element3()->get(i);
194     if (tet->get_mg_element_maillage()->get_origine()!=IMPOSE)
195     {
196     double jac[9];
197     int col,ligne;
198     double volume=tet->get_jacobien(jac,jac,ligne,col,1.);
199     passe=1;
200     nume=nume+tet->get_solution()*volume;
201     deno=deno+volume;
202     }
203     }
204     if (passe==1) no->change_solution(nume/deno); else no->change_solution(1.);
205     }
206     if (nom!="")
207     {
208     affiche((char*)" Enregistrement de la densité aux noeuds");
209     MG_GESTIONNAIRE *gest=mai->get_mg_maillage()->get_gestionnaire();
210     std::string chemin=nom+".sol";
211 francois 383 FEM_SOLUTION* femdens=new FEM_SOLUTION(mai,1,(char*)chemin.c_str(),mai->get_nb_fem_element3(),"Extrapolation aux noeuds",MAGIC::ENTITE_SOLUTION::ENTITE_ELEMENT3_NOEUD,MAGIC::TYPE_SOLUTION::SCALAIRE);
212 francois 343 gest->ajouter_fem_solution(femdens);
213     femdens->change_legende(0,"Densite");
214     LISTE_FEM_ELEMENT3::iterator it3;
215     int i=0;
216     for (FEM_ELEMENT3* tet=mai->get_premier_element3(it3);tet!=NULL;tet=mai->get_suivant_element3(it3))
217     {
218     MG_TETRA* tet2=(MG_TETRA*)tet->get_mg_element_maillage();
219     if (tet2->get_origine()==IMPOSE)
220     {
221 francois 375 femdens->ecrire(1.,i,0,0);
222     femdens->ecrire(1.,i,0,1);
223     femdens->ecrire(1.,i,0,2);
224     femdens->ecrire(1.,i,0,3);
225 francois 343 }
226     else
227     {
228 francois 375 femdens->ecrire(tet->get_fem_noeud(0)->get_solution(),i,0,0);
229     femdens->ecrire(tet->get_fem_noeud(1)->get_solution(),i,0,1);
230     femdens->ecrire(tet->get_fem_noeud(2)->get_solution(),i,0,2);
231     femdens->ecrire(tet->get_fem_noeud(3)->get_solution(),i,0,3);
232 francois 343 }
233     i++;
234     }
235     MG_EXPORT exp;
236     exp.gmsh(mai,nom);
237     gest->supprimer_fem_solution(gest->get_nb_fem_solution()-1);
238     }
239     MG_MAILLAGE* mgmai=new MG_MAILLAGE(NULL);
240     gest2.ajouter_mg_maillage(mgmai);
241    
242     MG_MAILLAGE* mgfem=mai->get_mg_maillage();
243     LISTE_MG_TRIANGLE::iterator itmg;
244     for (MG_TRIANGLE* tri=mgfem->get_premier_triangle(itmg);tri!=NULL;tri=mgfem->get_suivant_triangle(itmg))
245     tri->change_nouveau_numero(0);
246     affiche((char*)" Decoupage des tetra optimises");
247     LISTE_FEM_ELEMENT3::iterator it3;
248     for (FEM_ELEMENT3* tet=mai->get_premier_element3(it3);tet!=NULL;tet=mai->get_suivant_element3(it3))
249     {
250     MG_TETRA* tet2=(MG_TETRA*)tet->get_mg_element_maillage();
251     if (tet2->get_origine()==IMPOSE)
252     {
253     tet2->get_triangle1()->change_origine(IMPOSE);
254     tet2->get_triangle2()->change_origine(IMPOSE);
255     tet2->get_triangle3()->change_origine(IMPOSE);
256     tet2->get_triangle4()->change_origine(IMPOSE);
257     tet2->get_noeud1()->change_solution(tet->get_fem_noeud(0)->get_solution());
258     tet2->get_noeud2()->change_solution(tet->get_fem_noeud(1)->get_solution());
259     tet2->get_noeud3()->change_solution(tet->get_fem_noeud(2)->get_solution());
260     tet2->get_noeud4()->change_solution(tet->get_fem_noeud(3)->get_solution());
261     tet2->get_noeud1()->change_nouveau_numero(tet->get_fem_noeud(0)->get_id());
262     tet2->get_noeud2()->change_nouveau_numero(tet->get_fem_noeud(1)->get_id());
263     tet2->get_noeud3()->change_nouveau_numero(tet->get_fem_noeud(2)->get_id());
264     tet2->get_noeud4()->change_nouveau_numero(tet->get_fem_noeud(3)->get_id());
265     if (tet2->get_triangle1()->get_nouveau_numero()==0)
266     tet2->get_triangle1()->change_nouveau_numero(tet2->get_noeud4()->get_id());
267     else
268     tet2->get_triangle1()->change_nouveau_numero(-1);
269     if (tet2->get_triangle2()->get_nouveau_numero()==0)
270     tet2->get_triangle2()->change_nouveau_numero(tet2->get_noeud3()->get_id());
271     else
272     tet2->get_triangle2()->change_nouveau_numero(-1);
273     if (tet2->get_triangle3()->get_nouveau_numero()==0)
274     tet2->get_triangle3()->change_nouveau_numero(tet2->get_noeud1()->get_id());
275     else
276     tet2->get_triangle3()->change_nouveau_numero(-1);
277     if (tet2->get_triangle4()->get_nouveau_numero()==0)
278     tet2->get_triangle4()->change_nouveau_numero(tet2->get_noeud2()->get_id());
279     else
280     tet2->get_triangle4()->change_nouveau_numero(-1);
281     }
282     }
283    
284     TPL_LISTE_ENTITE<MG_TRIANGLE*> tri_impose_interne;
285     for (MG_TRIANGLE* tri=mgfem->get_premier_triangle(itmg);tri!=NULL;tri=mgfem->get_suivant_triangle(itmg))
286     {
287     if (tri->get_origine()==IMPOSE)
288     if (tri->get_lien_topologie()->get_dimension()==3)
289     if (tri->get_nouveau_numero()>0)
290     tri_impose_interne.ajouter(tri);
291     }
292     for (FEM_ELEMENT3* tet=mai->get_premier_element3(it3);tet!=NULL;tet=mai->get_suivant_element3(it3))
293     {
294     if (tet->get_mg_element_maillage()->get_origine()!=IMPOSE)
295     {
296     FEM_NOEUD* no1=tet->get_fem_noeud(0);
297     FEM_NOEUD* no2=tet->get_fem_noeud(1);
298     FEM_NOEUD* no3=tet->get_fem_noeud(2);
299     FEM_NOEUD* no4=tet->get_fem_noeud(3);
300     std::vector<MG_NOEUD*> tab;
301     interpole_segment(no1,no2,&tab,limit,mgmai);
302     interpole_segment(no1,no3,&tab,limit,mgmai);
303     interpole_segment(no1,no4,&tab,limit,mgmai);
304     interpole_segment(no2,no3,&tab,limit,mgmai);
305     interpole_segment(no2,no4,&tab,limit,mgmai);
306     interpole_segment(no3,no4,&tab,limit,mgmai);
307     int nb=tab.size();
308     FEM_NOEUD* noext;
309     if (nb>0)
310     {
311     if (no1->get_solution()<limit) noext=no1;
312     if (no2->get_solution()<limit) noext=no2;
313     if (no3->get_solution()<limit) noext=no3;
314     if (no4->get_solution()<limit) noext=no4;
315     }
316     if (nb==3)
317     {
318     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[1],tab[2],mgmai,TRIANGULATION);
319     oriente_tri(tri,noext->get_coord());
320     }
321     if (nb==4)
322     {
323     if (test_du_point_milieu(tab[0],tab[1],tet)==1)
324     {
325     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[1],tab[2],mgmai,TRIANGULATION);
326     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[0],tab[1],tab[3],mgmai,TRIANGULATION);
327     oriente_tri(tri,noext->get_coord());
328     oriente_tri(tri2,noext->get_coord());
329     }
330     else if (test_du_point_milieu(tab[0],tab[2],tet)==1)
331     {
332     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[2],tab[1],mgmai,TRIANGULATION);
333     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[0],tab[2],tab[3],mgmai,TRIANGULATION);
334     oriente_tri(tri,noext->get_coord());
335     oriente_tri(tri2,noext->get_coord());}
336     else if (test_du_point_milieu(tab[0],tab[3],tet)==1)
337     {
338     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[3],tab[1],mgmai,TRIANGULATION);
339     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[0],tab[3],tab[2],mgmai,TRIANGULATION);
340     oriente_tri(tri,noext->get_coord());
341     oriente_tri(tri2,noext->get_coord());}
342     else if (test_du_point_milieu(tab[1],tab[2],tet)==1)
343     {
344     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[1],tab[2],tab[0],mgmai,TRIANGULATION);
345     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[1],tab[2],tab[3],mgmai,TRIANGULATION);
346     oriente_tri(tri,noext->get_coord());
347     oriente_tri(tri2,noext->get_coord());
348     }
349     else if (test_du_point_milieu(tab[1],tab[3],tet)==1)
350     {
351     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[1],tab[3],tab[0],mgmai,TRIANGULATION);
352     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[1],tab[3],tab[2],mgmai,TRIANGULATION);
353     oriente_tri(tri,noext->get_coord());
354     oriente_tri(tri2,noext->get_coord());}
355     else
356     {
357     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[2],tab[3],tab[0],mgmai,TRIANGULATION);
358     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[2],tab[3],tab[1],mgmai,TRIANGULATION);
359     oriente_tri(tri,noext->get_coord());
360     oriente_tri(tri2,noext->get_coord());
361     }
362    
363     }
364     }
365     }
366     affiche((char*)" Decoupage des triangles non design interne au volume");
367     for (int i=0;i<tri_impose_interne.get_nb();i++)
368     {
369     MG_TRIANGLE* tri=tri_impose_interne.get(i);
370     MG_NOEUD* no1=tri->get_noeud1();
371     MG_NOEUD* no2=tri->get_noeud2();
372     MG_NOEUD* no3=tri->get_noeud3();
373     std::vector<MG_NOEUD*> tab;
374     int num1=-1;
375     int num2=-1;
376     int num3=-1;
377     if (no1->get_solution()<limit)
378     {
379     MG_NOEUD* no=get_noeud_peau(mai->get_fem_noeudid(no1->get_nouveau_numero()),mgmai);
380     tab.push_back(no);
381     num1=1;
382     }
383     if (no2->get_solution()<limit)
384     {
385     MG_NOEUD* no=get_noeud_peau(mai->get_fem_noeudid(no2->get_nouveau_numero()),mgmai);
386     tab.push_back(no);
387     num2=1;
388     }
389     if (no3->get_solution()<limit)
390     {
391     MG_NOEUD* no=get_noeud_peau(mai->get_fem_noeudid(no3->get_nouveau_numero()),mgmai);
392     tab.push_back(no);
393     num3=1;
394     }
395     if (num1*num2<0) interpole_segment(mai->get_fem_noeudid(no1->get_nouveau_numero()),mai->get_fem_noeudid(no2->get_nouveau_numero()),&tab,limit,mgmai);
396     if (num1*num3<0) interpole_segment(mai->get_fem_noeudid(no1->get_nouveau_numero()),mai->get_fem_noeudid(no3->get_nouveau_numero()),&tab,limit,mgmai);
397     if (num2*num3<0) interpole_segment(mai->get_fem_noeudid(no2->get_nouveau_numero()),mai->get_fem_noeudid(no3->get_nouveau_numero()),&tab,limit,mgmai);
398     MG_NOEUD* noint=mgfem->get_mg_noeudid(tri->get_nouveau_numero());
399     int nb=tab.size();
400     if (nb==3)
401     {
402     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[1],tab[2],mgmai,IMPOSE);
403     oriente_tri(tri,noint->get_coord());
404     tri->inverse_sens();
405     }
406     if (nb==4)
407     {
408     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[1],tab[2],mgmai,IMPOSE);
409     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[1],tab[2],tab[3],mgmai,IMPOSE);
410     oriente_tri(tri,noint->get_coord());
411     tri->inverse_sens();
412     oriente_tri(tri2,noint->get_coord());
413     tri2->inverse_sens();
414     }
415     }
416 francois 345 affiche((char*)" Decoupage des triangles non design sur la frontiere du volume");
417 francois 343 LISTE_FEM_ELEMENT2::iterator it2;
418     for (FEM_ELEMENT2* tri=mai->get_premier_element2(it2);tri!=NULL;tri=mai->get_suivant_element2(it2))
419     {
420     std::vector<MG_NOEUD*> tab;
421     FEM_NOEUD *no1=tri->get_fem_noeud(0);
422     FEM_NOEUD *no2=tri->get_fem_noeud(1);
423     FEM_NOEUD *no3=tri->get_fem_noeud(2);
424     int ori=((MG_FACE*)tri->get_lien_topologie())->get_mg_coface(0)->get_orientation();
425     OT_VECTEUR_3D n1n2(no1->get_coord(),no2->get_coord());
426     OT_VECTEUR_3D n1n3(no1->get_coord(),no3->get_coord());
427     OT_VECTEUR_3D nor=n1n3&n1n2;
428     double xyzext[3];
429     xyzext[0]=no1->get_x()+nor.get_x()*ori;
430     xyzext[1]=no1->get_y()+nor.get_y()*ori;
431     xyzext[2]=no1->get_z()+nor.get_z()*ori;
432     if (tri->get_mg_element_maillage()->get_origine()==IMPOSE)
433     {
434     MG_NOEUD *mgno1=get_noeud_peau(no1,mgmai);
435     MG_NOEUD *mgno2=get_noeud_peau(no2,mgmai);
436     MG_NOEUD *mgno3=get_noeud_peau(no3,mgmai);
437     int num1=interpole_segment(no1,no2,&tab,limit,mgmai,0);
438     int num2=interpole_segment(no1,no3,&tab,limit,mgmai,0);
439     int num3=interpole_segment(no2,no3,&tab,limit,mgmai,0);
440     int nb=tab.size();
441     if (nb==0)
442     {
443     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno2,mgno3,mgmai,IMPOSE);
444     oriente_tri(tri,xyzext);
445     }
446     if (nb==1)
447     {
448     if (num1==1)
449     {
450     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno3,tab[0],mgmai,IMPOSE);
451     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,mgno2,mgno3,tab[0],mgmai,IMPOSE);
452     oriente_tri(tri,xyzext);
453     oriente_tri(tri2,xyzext);
454     }
455     if (num2==1)
456     {
457     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno2,tab[0],mgmai,IMPOSE);
458     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,mgno2,mgno3,tab[0],mgmai,IMPOSE);
459     oriente_tri(tri,xyzext);
460     oriente_tri(tri2,xyzext);
461     }
462     if (num3==1)
463     {
464     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno2,tab[0],mgmai,IMPOSE);
465     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,mgno1,mgno3,tab[0],mgmai,IMPOSE);
466     oriente_tri(tri,xyzext);
467     oriente_tri(tri2,xyzext);
468     }
469     }
470     if (nb==2)
471     {
472     if (num1==0)
473     {
474     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno2,tab[0],mgmai,IMPOSE);
475     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,mgno2,tab[0],tab[1],mgmai,IMPOSE);
476     MG_TRIANGLE_PEAU* tri3=insere_triangle(NULL,mgno3,tab[0],tab[1],mgmai,IMPOSE);
477     oriente_tri(tri,xyzext);
478     oriente_tri(tri2,xyzext);
479     oriente_tri(tri3,xyzext);
480     }
481     if (num2==0)
482     {
483     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno3,tab[0],mgmai,IMPOSE);
484     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,mgno3,tab[0],tab[1],mgmai,IMPOSE);
485     MG_TRIANGLE_PEAU* tri3=insere_triangle(NULL,mgno2,tab[0],tab[1],mgmai,IMPOSE);
486     oriente_tri(tri,xyzext);
487     oriente_tri(tri2,xyzext);
488     oriente_tri(tri3,xyzext);
489     }
490     if (num2==0)
491     {
492     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno2,mgno3,tab[0],mgmai,IMPOSE);
493     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,mgno3,tab[0],tab[1],mgmai,IMPOSE);
494     MG_TRIANGLE_PEAU* tri3=insere_triangle(NULL,mgno1,tab[0],tab[1],mgmai,IMPOSE);
495     oriente_tri(tri,xyzext);
496     oriente_tri(tri2,xyzext);
497     oriente_tri(tri3,xyzext);
498     }
499     }
500     }
501     else if ((no1->get_solution()>=limit) && (no2->get_solution()>=limit) && (no3->get_solution()>=limit))
502     {
503     MG_NOEUD *mgno1=get_noeud_peau(no1,mgmai);
504     MG_NOEUD *mgno2=get_noeud_peau(no2,mgmai);
505     MG_NOEUD *mgno3=get_noeud_peau(no3,mgmai);
506     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,mgno1,mgno2,mgno3,mgmai,TRIANGULATION);
507     oriente_tri(tri,xyzext);
508     }
509     else if (!((no1->get_solution()<limit) && (no2->get_solution()<limit) && (no3->get_solution()<limit)))
510     {
511     std::vector<MG_NOEUD*> tab;
512     int num1=-1;
513     int num2=-1;
514     int num3=-1;
515     if (no1->get_solution()>=limit)
516     {
517     MG_NOEUD* no=get_noeud_peau(no1,mgmai);
518     tab.push_back(no);
519     num1=1;
520     }
521     if (no2->get_solution()>=limit)
522     {
523     MG_NOEUD* no=get_noeud_peau(no2,mgmai);
524     tab.push_back(no);
525     num2=1;
526     }
527     if (no3->get_solution()>=limit)
528     {
529     MG_NOEUD* no=get_noeud_peau(no3,mgmai);
530     tab.push_back(no);
531     num3=1;
532     }
533     if (num1*num2<0) interpole_segment(no1,no2,&tab,limit,mgmai);
534     if (num1*num3<0) interpole_segment(no1,no3,&tab,limit,mgmai);
535     if (num2*num3<0) interpole_segment(no2,no3,&tab,limit,mgmai);
536     int nb=tab.size();
537     if (nb==3)
538     {
539     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[1],tab[2],mgmai,TRIANGULATION);
540     oriente_tri(tri,xyzext);
541     }
542     if (nb==4)
543     {
544     MG_TRIANGLE_PEAU* tri=insere_triangle(NULL,tab[0],tab[1],tab[2],mgmai,TRIANGULATION);
545     MG_TRIANGLE_PEAU* tri2=insere_triangle(NULL,tab[1],tab[2],tab[3],mgmai,TRIANGULATION);
546     oriente_tri(tri,xyzext);
547     oriente_tri(tri2,xyzext);
548     }
549     }
550     }
551     LISTE_MG_TRIANGLE::iterator it_tri;
552     for (MG_TRIANGLE* mgtri=mgmai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mgmai->get_suivant_triangle(it_tri))
553     {
554     MG_TRIANGLE_PEAU* tripeau=(MG_TRIANGLE_PEAU*)mgtri;
555     MG_TRIANGLE_PEAU* voisin1=recherche_voisin(tripeau->get_noeud1(),tripeau->get_noeud2(),tripeau);
556     MG_TRIANGLE_PEAU* voisin2=recherche_voisin(tripeau->get_noeud2(),tripeau->get_noeud3(),tripeau);
557     MG_TRIANGLE_PEAU* voisin3=recherche_voisin(tripeau->get_noeud3(),tripeau->get_noeud1(),tripeau);
558     tripeau->change_voisin1(voisin1);
559     tripeau->change_voisin2(voisin2);
560     tripeau->change_voisin3(voisin3);
561     tripeau->change_nouveau_numero(0);
562     }
563     int fin;
564     do
565     {
566     fin=1;
567     for (MG_TRIANGLE* mgtri=mgmai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mgmai->get_suivant_triangle(it_tri))
568     {
569     MG_TRIANGLE_PEAU *tripeau=(MG_TRIANGLE_PEAU*)mgtri;
570     if (tripeau->get_nouveau_numero()==0)
571     {
572     fin=0;
573     std::vector<MG_TRIANGLE_PEAU*> *peau=new std::vector<MG_TRIANGLE_PEAU*>;
574     lst_peau.push_back(peau);
575     tripeau->change_nouveau_numero(1);
576     peau->push_back(tripeau);
577     determine_peau(peau);
578     }
579     }
580     }
581     while (fin==0);
582 francois 457 return mgmai;
583     }
584 francois 460 void MGOPT_POSTTRAITEMENT::suppression_peaux_isoles(MG_MAILLAGE* mgmai)
585 francois 457 {
586     affiche((char*)"Suppression des peaux isolées");
587 francois 343
588     char message[500];
589     for (int cas=0;cas<2;cas++)
590     {
591 francois 457 if (cas==0) affiche((char*)" Analyse initiale des peaux");
592     if (cas==1) affiche((char*)" Analyse finale des peaux");
593     int nbisole=0;
594 francois 343 for (int i=0;i<lst_peau.size();i++)
595     {
596     int isole=1;
597     for (int j=0;j<lst_peau[i]->size();j++)
598     if ((*lst_peau[i])[j]->get_origine()==IMPOSE) {isole=0;break;}
599     if (isole==1)
600     {
601 francois 457 nbisole++;
602 francois 343 for (int j=0;j<lst_peau[i]->size();j++)
603     {
604     mgmai->supprimer_mg_triangleid((*lst_peau[i])[j]->get_id());
605     }
606     lst_peau[i]->clear();
607     }
608     }
609 francois 457 sprintf(message," %d peaux, %d non isoles, %d isoles",(int)lst_peau.size(),(int)lst_peau.size()-nbisole,nbisole);
610     affiche(message);
611 francois 343 for (std::vector<std::vector<MG_TRIANGLE_PEAU*> *>::iterator j=lst_peau.end()-1;j>lst_peau.begin();j--)
612     if ((*j)->size()==0)
613     lst_peau.erase(j);
614     }
615 francois 457 LISTE_MG_TRIANGLE::iterator itmg;
616 francois 343 for (MG_TRIANGLE* tri=mgmai->get_premier_triangle(itmg);tri!=NULL;tri=mgmai->get_suivant_triangle(itmg))
617     {
618     if (tri->get_origine()==IMPOSE)
619     {
620     tri->get_noeud1()->change_origine(IMPOSE);
621     tri->get_noeud2()->change_origine(IMPOSE);
622     tri->get_noeud3()->change_origine(IMPOSE);
623     }
624     }
625 cuillier 438
626 francois 343 }
627    
628 francois 460 void MGOPT_POSTTRAITEMENT::oriente_tri(MG_TRIANGLE_PEAU* tri,double *xyz)
629 francois 343 {
630     MG_NOEUD* no1=tri->get_noeud1();
631     MG_NOEUD* no2=tri->get_noeud2();
632     MG_NOEUD* no3=tri->get_noeud3();
633     OT_VECTEUR_3D n1n2(no1->get_coord(),no2->get_coord());
634     OT_VECTEUR_3D n1n3(no1->get_coord(),no3->get_coord());
635     OT_VECTEUR_3D normal=n1n2&n1n3;
636     OT_VECTEUR_3D dir(no1->get_coord(),xyz);
637     dir.norme();
638     normal.norme();
639     double ps=normal*dir;
640     if (ps<0) tri->inverse_sens();
641     }
642    
643    
644    
645 francois 460 MG_NOEUD* MGOPT_POSTTRAITEMENT::get_noeud_peau(FEM_NOEUD* no,MG_MAILLAGE* mai)
646 francois 343 {
647     static std::map<std::string,MG_NOEUD*> correspond;
648     unsigned long id=no->get_id();
649     char message[255];
650     sprintf(message,"_%lu_",id);
651     std::string key=message;
652     MG_NOEUD* mgno=correspond[key];
653     if (mgno==NULL)
654     {
655     mgno=new MG_NOEUD(NULL,no->get_x(),no->get_y(),no->get_z(),TRIANGULATION);
656     mai->ajouter_mg_noeud(mgno);
657     correspond[key]=mgno;
658     }
659     return mgno;
660     }
661    
662    
663 francois 460 int MGOPT_POSTTRAITEMENT::test_du_point_milieu(MG_NOEUD* no1,MG_NOEUD* no2,FEM_ELEMENT3* tet)
664 francois 343 {
665     double *xyz1=no1->get_coord();
666     double *xyz2=no2->get_coord();
667     double xyz[3];
668     xyz[0]=0.5*(xyz1[0]+xyz2[0]);
669     xyz[1]=0.5*(xyz1[1]+xyz2[1]);
670     xyz[2]=0.5*(xyz1[2]+xyz2[2]);
671 francois 412 FEM_MAILLAGE_OUTILS ot;
672 francois 388 if (((ot.estdansletetra(tet,xyz[0],xyz[1],xyz[2])>>1)&1)==1) return 1;
673     //if (((MG_MAILLAGE_OUTILS::estdansletetra(tet,xyz[0],xyz[1],xyz[2])>>1)&1)==1) return 1;
674 francois 343 return 0;
675     }
676    
677    
678    
679 francois 460 int MGOPT_POSTTRAITEMENT::interpole_segment(FEM_NOEUD* no1,FEM_NOEUD* no2,std::vector<MG_NOEUD*> *tab,double limit,MG_MAILLAGE* mai,int creation)
680 francois 343 {
681     static std::map<std::string,MG_NOEUD*> correspond;
682     unsigned long id1=no1->get_id();
683     unsigned long id2=no2->get_id();
684     char message[255];
685     sprintf(message,"_%lu_%lu_",id1,id2);
686     std::string key1=message;
687     sprintf(message,"_%lu_%lu_",id2,id1);
688     std::string key2=message;
689     double sol1=no1->get_solution();
690     double sol2=no2->get_solution();
691 francois 345 if (fabs(sol1-limit)<0.0000001)
692     {
693     sprintf(message,"%lu",id1);
694     std::string key=message;
695     MG_NOEUD* no=correspond[key];
696     if (no==NULL)
697     {
698     double *xyz1=no1->get_coord();
699     no=new MG_NOEUD(NULL,xyz1[0],xyz1[1],xyz1[2],TRIANGULATION);
700     mai->ajouter_mg_noeud(no);
701     correspond[key]=no;
702     }
703     int present=0;
704     for (int i=0;i<tab->size();i++)
705     if ((*tab)[i]==no) present=1;
706     if (present==0)
707     {
708     tab->push_back(no);
709     return 1;
710     }
711     else return 0;
712     }
713     if (fabs(sol2-limit)<0.0000001)
714     {
715     sprintf(message,"%lu",id2);
716     std::string key=message;
717     MG_NOEUD* no=correspond[key];
718     if (no==NULL)
719     {
720     double *xyz2=no2->get_coord();
721     no=new MG_NOEUD(NULL,xyz2[0],xyz2[1],xyz2[2],TRIANGULATION);
722     mai->ajouter_mg_noeud(no);
723     correspond[key]=no;
724     }
725     int present=0;
726     for (int i=0;i<tab->size();i++)
727     if ((*tab)[i]==no) present=1;
728     if (present==0)
729     {
730     tab->push_back(no);
731     return 1;
732     }
733     else return 0;
734    
735     }
736    
737     if (((sol1<limit) && (sol2>limit))||((sol1>limit) && (sol2<limit)))
738 francois 343 {
739     MG_NOEUD* no=correspond[key1];
740     if (no==NULL) no=correspond[key2];
741     if ((no==NULL) && (creation==1))
742     {
743     double t=(limit-sol1)/(sol2-sol1);
744     double xyz[3];
745     double *xyz1=no1->get_coord();
746     double *xyz2=no2->get_coord();
747     xyz[0]=xyz1[0]+t*(xyz2[0]-xyz1[0]);
748     xyz[1]=xyz1[1]+t*(xyz2[1]-xyz1[1]);
749     xyz[2]=xyz1[2]+t*(xyz2[2]-xyz1[2]);
750     no=new MG_NOEUD(NULL,xyz[0],xyz[1],xyz[2],TRIANGULATION);
751     mai->ajouter_mg_noeud(no);
752     correspond[key1]=no;
753     }
754     if (no!=NULL)
755     {
756     tab->push_back(no);
757     return 1;
758     }
759     }
760     return 0;
761     }
762    
763 francois 460 void MGOPT_POSTTRAITEMENT::lire_params(char *fichier)
764 francois 457 {
765     params.lire(fichier);
766     }
767 francois 343
768 francois 461 void MGOPT_POSTTRAITEMENT::posttraite(FEM_SOLUTION* sol, MG_GESTIONNAIRE& gest2,char *nomparam)
769 picher 233 {
770 francois 457 if (nomparam!=NULL) lire_params(nomparam);
771 francois 272 affiche((char*)"Extraction du maillage de surface");
772 francois 457 FEM_MAILLAGE* fem=sol->get_maillage();
773     int decoup=params.get_valeur("decoupage");
774     MG_MAILLAGE* mgmai;
775     if (decoup==0)
776     {
777     affiche((char*)" Découpage par mailles 3D entieres");
778     adapte_seuil(fem,sol);
779     mgmai=extract_skin_maille_entiere(fem,gest2);
780     }
781     if (decoup==1)
782     {
783     affiche((char*)" Découpage par isodensité");
784     double seuil=params.get_valeur("seuil");
785     std::string nomfich=params.get_nom("nomdensextra");
786     mgmai=extract_skin_par_decoupage(sol,seuil,gest2,nomfich.c_str());
787     }
788     suppression_peaux_isoles(mgmai);
789     affiche((char*)"Procedure de lissage");
790     int chen2005_debut=params.get_valeur("chen2005_debut");
791     int chen2008_debut=params.get_valeur("chen2008_debut");
792     int jiao2012_debut=params.get_valeur("jiao2012_debut");
793     int chen2005fait=0;
794     if (chen2005_debut==0) chen2005fait=1;
795     int chen2008fait=0;
796     if (chen2008_debut==0) chen2008fait=1;
797     int jiao2012fait=0;
798     if (jiao2012_debut==0) jiao2012fait=1;
799     int fin=0;
800     int iteration=1;
801     do
802     {
803     if ((iteration>=chen2005_debut) && (chen2005fait==0))
804     {
805     char message[300];
806     sprintf(message, " Iteration %d : lissage de Chen 2005",iteration);
807     affiche(message);
808     double epsilon=params.get_valeur("chen2005_epsilon");
809     double sigma=params.get_valeur("chen2005_sigma");
810     int iter_max=params.get_valeur("chen2005_itermax");
811     int nbiter=lissage_chen2005(mgmai,gest2,epsilon,sigma,iter_max);
812     if (nbiter==iter_max) sprintf(message," Arrêt de la procédure de lissage de Chen2005 après %d itérations",nbiter);
813     else sprintf(message," Convergence de la procédure de lissage de Chen2005 après %d itérations",nbiter);
814     affiche(message);
815     chen2005fait=1;
816     iteration=iteration+nbiter;
817     }
818     else if ((iteration>=chen2008_debut) && (chen2008fait==0))
819     {
820     char message[300];
821     sprintf(message, " Iteration %d : lissage de Chen 2008",iteration);
822     affiche(message);
823     double epsilon=params.get_valeur("chen2008_epsilon");
824     double sigma=params.get_valeur("chen2008_sigma");
825     double gamma=params.get_valeur("chen2008_gamma");
826     int iter_max=params.get_valeur("chen2008_itermax");
827     int nbiter=lissage_chen2008(mgmai,gest2,sigma,gamma,epsilon,iter_max);
828     if (nbiter==iter_max) sprintf(message," Arrêt de la procédure de lissage de Chen2008 après %d itérations",nbiter);
829     else sprintf(message," Convergence de la procédure de lissage de Chen2008 après %d itérations",nbiter);
830     affiche(message);
831     chen2008fait=1;
832     iteration=iteration+nbiter;
833     }
834     else if ((iteration>=jiao2012_debut) && (jiao2012fait==0))
835     {
836     char message[300];
837     sprintf(message, " Iteration %d : lissage de Jiao 2012",iteration);
838     affiche(message);
839     int iter_max=params.get_valeur("jiao2012_itermax");
840     int nbiter=lissage_jiao2012(mgmai,gest2,iter_max);
841     if (nbiter==iter_max) sprintf(message," Arrêt de la procédure de lissage de Jiao2012 après %d itérations",nbiter);
842     else sprintf(message," Convergence de la procédure de lissage de Jiao2012 après %d itérations",nbiter);
843     affiche(message);
844     jiao2012fait=1;
845     iteration=iteration+nbiter;
846     }
847     else iteration++;
848     if (chen2005fait==1)
849     if (chen2008fait==1)
850     if (jiao2012fait==1)
851     fin=1;
852     }
853     while (fin==0);
854     }
855    
856    
857 francois 460 MG_MAILLAGE* MGOPT_POSTTRAITEMENT::extract_skin_maille_entiere(FEM_MAILLAGE* mai,MG_GESTIONNAIRE &gest2)
858 francois 457 {
859     //Menu de la methode de lissage
860     affiche((char*)"Extraction du maillage de surface");
861 picher 231 int coderesu = 0;
862 picher 233 int mai2_id;
863 francois 457 int imp=params.get_valeur("consimpose");
864     int opti=params.get_valeur("consoptimise");
865     int m_auto=params.get_valeur("consnooptimise");
866 francois 232 if (opti==1) conserve(OPTIMISE);
867     if (imp==1) conserve(IMPOSE);
868 picher 233 if (m_auto==1) conserve(MAILLEUR_AUTO);
869 francois 457 copieorigine(mai);
870 picher 233 MG_MAILLAGE* mg_mai = (MG_MAILLAGE*)mai->get_mg_maillage();
871 francois 457 //gain_poids(mg_mai,gest2);
872     int reactiv=params.get_valeur("reactivation");
873 picher 233 if (reactiv == 1)
874     {
875     reactivation(mg_mai,gest2);
876     }
877 francois 457 affiche((char*)" Analyse des cas non manifold");
878 picher 231 do
879     {
880 francois 457 int nbmaniare,nbmanino,nbpeau;
881     coderesu = extract_skin(mg_mai,gest2,nbpeau,nbmaniare,nbmanino,&mai2_id);
882     char message[500];
883     sprintf(message," %d peaux. %d manifold par arete. %d manifold par noeud",nbpeau,nbmaniare,nbmanino);
884     affiche(message);
885     // gain_poids(mg_mai,gest2);
886 picher 231 }
887     while (coderesu == 0);
888 picher 233
889     MG_MAILLAGE* mg_mai2=gest2.get_mg_maillageid(mai2_id);
890 francois 457
891     return mg_mai2;
892    
893 picher 231 }
894    
895 francois 460 void MGOPT_POSTTRAITEMENT::adapte_seuil(class FEM_MAILLAGE* fem,FEM_SOLUTION* solution)
896 picher 230 {
897 francois 457 double seuil=params.get_valeur("seuil");
898     solution->active_solution(0);
899     LISTE_FEM_ELEMENT3::iterator it;
900     for (FEM_ELEMENT3 *tet=fem->get_premier_element3(it);tet!=NULL;tet=fem->get_suivant_element3(it))
901 picher 230 {
902 francois 457 if (((MG_TETRA*)tet->get_mg_element_maillage())->get_origine()!=IMPOSE)
903     if (tet->get_solution()>seuil)
904     ((MG_TETRA*)tet->get_mg_element_maillage())->change_origine(OPTIMISE);
905     else
906     ((MG_TETRA*)tet->get_mg_element_maillage())->change_origine(MAILLEUR_AUTO);
907    
908 picher 230 }
909     }
910 francois 224
911 francois 457
912    
913    
914     //---------------------------- Lissage Chen 2005 Gilles Philippe -------------------------------------------------------------------------
915 francois 460 int MGOPT_POSTTRAITEMENT::lissage_chen2005(MG_MAILLAGE* mg_mai,MG_GESTIONNAIRE& gest2, double epsilon, double sigma, int iter_max)
916 picher 230 {
917 francois 457 //Methode de lissage inspiree de Chen(2005) par Gilles Philipe Picher Martel
918 picher 230 double un_sur_pi = 1./M_PI;
919     int compteur = 0;
920     int fin = 0;
921 francois 224
922 francois 457
923 picher 230 do
924     {
925     TPL_LISTE_ENTITE<OT_VECTEUR_3D> liste_normales;
926 picher 233 TPL_LISTE_ENTITE<OT_VECTEUR_3D> liste_normales2;
927 picher 230 TPL_LISTE_ENTITE<double> liste_wij;
928     LISTE_MG_TRIANGLE::iterator it_tri;
929     int k = 0; //pour identifier les triangles pour liste_normales et liste_wij
930     for (MG_TRIANGLE* mgtri=mg_mai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mg_mai->get_suivant_triangle(it_tri))
931     {
932     MG_TRIANGLE_PEAU* mgtri_i = (MG_TRIANGLE_PEAU*)mgtri;
933     OT_VECTEUR_3D normal_f_i = mgtri_i->calcul_normal();
934 picher 233 normal_f_i.norme();
935     liste_normales2.ajouter(normal_f_i);
936 picher 230 //Remplissage de la liste des voisins du triangle i
937     TPL_MAP_ENTITE<MG_TRIANGLE_PEAU*> liste_voisins;
938     MG_NOEUD* noeud1 = mgtri_i->get_noeud1();
939     double nb_voisins1 = noeud1->get_lien_triangle()->get_nb();
940     for (int j = 0;j<nb_voisins1;j++)
941     {
942     MG_TRIANGLE_PEAU* mgtri_1 = (MG_TRIANGLE_PEAU*) noeud1->get_lien_triangle()->get(j);
943     liste_voisins.ajouter(mgtri_1);
944     }
945     MG_NOEUD* noeud2 = mgtri_i->get_noeud2();
946     double nb_voisins2 = noeud2->get_lien_triangle()->get_nb();
947     for (int j = 0;j<nb_voisins2;j++)
948     {
949     MG_TRIANGLE_PEAU* mgtri_2 = (MG_TRIANGLE_PEAU*) noeud2->get_lien_triangle()->get(j);
950     liste_voisins.ajouter(mgtri_2);
951     }
952     MG_NOEUD* noeud3 = mgtri_i->get_noeud3();
953     double nb_voisins3 = noeud3->get_lien_triangle()->get_nb();
954     for (int j = 0;j<nb_voisins3;j++)
955     {
956     MG_TRIANGLE_PEAU* mgtri_3 = (MG_TRIANGLE_PEAU*) noeud3->get_lien_triangle()->get(j);
957     liste_voisins.ajouter(mgtri_3);
958     }
959     liste_voisins.supprimer(mgtri_i);
960     int nb_voisins = liste_voisins.get_nb();
961     double w_ij = 1./nb_voisins;
962     double phi_i_min = 10.;
963     double s_i = 0.0;
964     double phi_im = 0.0;
965 picher 231 double *phi_ij = new double[nb_voisins];
966 picher 230 OT_VECTEUR_3D normal_f_i_mean(0.,0.,0.);
967     OT_VECTEUR_3D eta_i(0.,0.,0.);
968     TPL_MAP_ENTITE<MG_TRIANGLE_PEAU*>::ITERATEUR it;
969     int j = 0;
970     for(MG_TRIANGLE_PEAU* mgtri_j=liste_voisins.get_premier(it);mgtri_j!=NULL;mgtri_j=liste_voisins.get_suivant(it))
971     {
972     OT_VECTEUR_3D normal_f_j = mgtri_j->calcul_normal();
973     //1-Calculer la normale moyenne pour chaque triangle
974     normal_f_i_mean = normal_f_i_mean + w_ij*normal_f_j;
975 francois 457 //2.1-On calcul l'angle entre normal_f_i et normal_f_j pour j allant de 1 a Nb_voisins
976 picher 230 double prod_scalaire = normal_f_i*normal_f_j;
977 picher 233 if (prod_scalaire > 1.)
978     {
979     prod_scalaire = 1.;
980     }
981     if (prod_scalaire < -1.)
982     {
983     prod_scalaire = -1.;
984     }
985 picher 230 phi_ij[j] = acos(prod_scalaire)*un_sur_pi;
986     //2.2-On trouve le plus petit des angles et la normale heta_i correspondante
987     if (phi_ij[j] < phi_i_min)
988     {
989     phi_i_min = phi_ij[j];
990     eta_i = normal_f_j;
991     }
992     //3.1-On calcul l'angle moyen phi_im
993     phi_im = phi_im + w_ij*phi_ij[j];
994     j++;
995     }
996     normal_f_i_mean.norme();
997     j = 0;
998     for(MG_TRIANGLE_PEAU* mgtri_j=liste_voisins.get_premier(it);mgtri_j!=NULL;mgtri_j=liste_voisins.get_suivant(it))
999     {
1000     //3.2-Calcul de s_i selon la variance
1001     s_i = s_i + w_ij*pow((phi_ij[j] - phi_im),2);
1002     j++;
1003     }
1004 picher 231 delete[] phi_ij;
1005 picher 230 //4-On calcule une nouvelle normale pour chaque triangle
1006 francois 458 double pond;
1007     int num=params.get_valeur("chen2005_filtre");
1008     if (num==1) pond=ponderation_gaussian(s_i,sigma);
1009     else if (num==2) pond=ponderation_laplacian(s_i,sigma);
1010     else if (num==3) pond=ponderation_elfallahford(s_i,sigma);
1011     OT_VECTEUR_3D normal_f_i_new = pond*normal_f_i_mean + (1. - pond)*eta_i;
1012 picher 230 normal_f_i_new.norme();
1013     liste_normales.ajouter(normal_f_i_new);
1014     liste_wij.ajouter(w_ij);
1015     mgtri->change_nouveau_numero(k);
1016     k++;
1017     }
1018    
1019     LISTE_MG_NOEUD::iterator it_no;
1020     for (MG_NOEUD* noeud_i=mg_mai->get_premier_noeud(it_no);noeud_i!=NULL;noeud_i=mg_mai->get_suivant_noeud(it_no))
1021     {
1022     int nb_voisins_j = noeud_i->get_lien_triangle()->get_nb();
1023     double w_ij_prime = 0.0;
1024     OT_VECTEUR_3D v_temp(0.,0.,0.);
1025     OT_VECTEUR_3D v_i(noeud_i->get_x(),noeud_i->get_y(),noeud_i->get_z());
1026     for(int j=0;j<nb_voisins_j;j++)
1027     {
1028     MG_TRIANGLE_PEAU* mgtri_j = (MG_TRIANGLE_PEAU*) noeud_i->get_lien_triangle()->get(j);
1029 francois 457 //On calcule le centroide cj du triangle mgtri_j
1030 picher 230 MG_NOEUD* n1 = mgtri_j->get_noeud1();
1031     MG_NOEUD* n2 = mgtri_j->get_noeud2();
1032     MG_NOEUD* n3 = mgtri_j->get_noeud3();
1033     double cj_x = 0.333333333333333*(n1->get_x() + n2->get_x() + n3->get_x());
1034     double cj_y = 0.333333333333333*(n1->get_y() + n2->get_y() + n3->get_y());
1035     double cj_z = 0.333333333333333*(n1->get_z() + n2->get_z() + n3->get_z());
1036     //On forme le vecteur vi_cj
1037     OT_VECTEUR_3D vi_cj(cj_x - noeud_i->get_x(),cj_y - noeud_i->get_y(),cj_z - noeud_i->get_z());
1038     OT_VECTEUR_3D normal_f_i_new = liste_normales.get(mgtri_j->get_nouveau_numero());
1039 francois 457 // w_ij_prime correspond a la somme des aires des triangles voisins du noeuds
1040 picher 233 OT_VECTEUR_3D AB(n2->get_x() - n1->get_x(),n2->get_y() - n1->get_y(),n2->get_z() - n1->get_z());
1041     OT_VECTEUR_3D AC(n3->get_x() - n1->get_x(),n3->get_y() - n1->get_y(),n3->get_z() - n1->get_z());
1042     OT_VECTEUR_3D prodvect = AB&AC;
1043     double w_ij = 0.5*prodvect.get_longueur();
1044 picher 230 w_ij_prime = w_ij_prime + w_ij;
1045     v_temp = v_temp + w_ij*(vi_cj*normal_f_i_new)*normal_f_i_new;
1046     }
1047 francois 457 //5-On met a jour la position des noeuds
1048 picher 230 v_i = v_i + v_temp/w_ij_prime;
1049 picher 233 int origine = noeud_i->get_origine();
1050     if (origine != IMPOSE)
1051     {
1052     noeud_i->change_x(v_i.get_x());
1053     noeud_i->change_y(v_i.get_y());
1054     noeud_i->change_z(v_i.get_z());
1055     }
1056 picher 230 }
1057 francois 457 //Critere d'arret de l'algorithme
1058 picher 230 int l=0;
1059     int nb_tri = mg_mai->get_nb_mg_triangle();
1060     for (MG_TRIANGLE* mgtri=mg_mai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mg_mai->get_suivant_triangle(it_tri))
1061     {
1062     MG_TRIANGLE_PEAU* mgtri_i = (MG_TRIANGLE_PEAU*)mgtri;
1063 picher 233 OT_VECTEUR_3D normal_f_i = liste_normales2.get(mgtri_i->get_nouveau_numero());
1064 picher 230 OT_VECTEUR_3D normal_f_i_new = liste_normales.get(mgtri_i->get_nouveau_numero());
1065 picher 248 double critere = 1. - normal_f_i*normal_f_i_new;
1066 picher 230 if (critere <= epsilon) l++;
1067     }
1068 picher 248 double tolerance = 0.01*nb_tri;
1069     if (nb_tri - l <= 0) fin = 1;
1070 picher 230 compteur++;
1071     }
1072     while ((fin == 0) && (compteur < iter_max));
1073    
1074 francois 457 return compteur;
1075     }
1076 picher 230
1077 francois 457 //---------------------------- Lissage Jiao 2012 -------------------------------------------------------------------------
1078 francois 460 int MGOPT_POSTTRAITEMENT::lissage_jiao2012(MG_MAILLAGE* mg_mai,MG_GESTIONNAIRE& gest2, int iter_max)
1079 picher 248 {
1080 francois 457 // Lissage global avec methode de lissage inspiree de JIAO IMR2012 par Jean Christophe sept 2013
1081 picher 248 double un_sur_pi = 1./M_PI;
1082     int compteur = 0;
1083    
1084 cuillier 438 do
1085     {
1086     vector<double> nouv_position_x;
1087     vector<double> nouv_position_y;
1088     vector<double> nouv_position_z;
1089     LISTE_MG_NOEUD::iterator it_no;
1090     for (MG_NOEUD* noeud_i=mg_mai->get_premier_noeud(it_no);noeud_i!=NULL;noeud_i=mg_mai->get_suivant_noeud(it_no))
1091     {
1092     int origine = noeud_i->get_origine();
1093     if (origine != IMPOSE)
1094     {
1095     int nb_voisins_j = noeud_i->get_lien_triangle()->get_nb();
1096     double wij = 0.0;
1097     double wij_sommej = 0.0;
1098     OT_VECTEUR_3D v_temp(0.,0.,0.);
1099     OT_VECTEUR_3D v_i(noeud_i->get_x(),noeud_i->get_y(),noeud_i->get_z());
1100     for(int j=0;j<nb_voisins_j;j++)
1101     {
1102     MG_TRIANGLE_PEAU* mgtri_j = (MG_TRIANGLE_PEAU*) noeud_i->get_lien_triangle()->get(j);
1103 francois 457 //On calcule le centroide cj du triangle mgtri_j
1104 cuillier 438 MG_NOEUD* n1 = mgtri_j->get_noeud1();
1105     MG_NOEUD* n2 = mgtri_j->get_noeud2();
1106     MG_NOEUD* n3 = mgtri_j->get_noeud3();
1107     double cj_x = 0.333333333333333*(n1->get_x() + n2->get_x() + n3->get_x());
1108     double cj_y = 0.333333333333333*(n1->get_y() + n2->get_y() + n3->get_y());
1109     double cj_z = 0.333333333333333*(n1->get_z() + n2->get_z() + n3->get_z());
1110     //On forme le vecteur vi_cj et le vecteur cj
1111     OT_VECTEUR_3D cj(cj_x ,cj_y,cj_z);
1112     OT_VECTEUR_3D vi_cj(cj_x - noeud_i->get_x(),cj_y - noeud_i->get_y(),cj_z - noeud_i->get_z());
1113     wij=vi_cj.get_longueur();
1114     wij_sommej=wij_sommej+wij;
1115     v_temp = v_temp + (wij*cj);
1116     }
1117     //5-On calcule la nouvelle position des noeuds et on affecte au tableau temporaire
1118     v_i =v_temp/wij_sommej;
1119     nouv_position_x.push_back(v_i.get_x());
1120     nouv_position_y.push_back(v_i.get_y());
1121     nouv_position_z.push_back(v_i.get_z());
1122     }
1123     }
1124     //On actualise la position des noeuds
1125     int ind_noeud=0;
1126     for (MG_NOEUD* noeud_i=mg_mai->get_premier_noeud(it_no);noeud_i!=NULL;noeud_i=mg_mai->get_suivant_noeud(it_no))
1127     {
1128     int origine = noeud_i->get_origine();
1129     if (origine != IMPOSE)
1130     {
1131     noeud_i->change_x(nouv_position_x[ind_noeud]);
1132     noeud_i->change_y(nouv_position_y[ind_noeud]);
1133     noeud_i->change_z(nouv_position_z[ind_noeud]);
1134     ind_noeud++;
1135     }
1136     }
1137     //Critere d'arret de l'algorithme
1138     compteur++;
1139     }
1140     while (compteur < iter_max);
1141    
1142 francois 457 return compteur;
1143 cuillier 438 }
1144    
1145 francois 457 //---------------------------- Lissage Chen 2008 -------------------------------------------------------------------------
1146 francois 460 int MGOPT_POSTTRAITEMENT::lissage_chen2008(MG_MAILLAGE* mg_mai,MG_GESTIONNAIRE& gest2, double sigma, double gamma, double epsilon, int iter_max)
1147 cuillier 438 {
1148 francois 457 //Lissage global inspire de Chen(2008) corrige par Jean Christophe sept 2013
1149 cuillier 438 double un_sur_pi = 1./M_PI;
1150     int compteur = 0;
1151     int fin = 0;
1152 francois 457 compteur = 0;
1153 cuillier 438
1154     do
1155     {
1156     vector<double> nouv_position_x;
1157     vector<double> nouv_position_y;
1158     vector<double> nouv_position_z;
1159     TPL_LISTE_ENTITE<OT_VECTEUR_3D> liste_normales;
1160     TPL_LISTE_ENTITE<OT_VECTEUR_3D> liste_normales2;
1161     TPL_LISTE_ENTITE<double> liste_wij;
1162     LISTE_MG_TRIANGLE::iterator it_tri;
1163     int k = 0; //pour identifier les triangles pour liste_normales et liste_wij
1164     for (MG_TRIANGLE* mgtri=mg_mai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mg_mai->get_suivant_triangle(it_tri))
1165     {
1166     MG_TRIANGLE_PEAU* mgtri_i = (MG_TRIANGLE_PEAU*)mgtri;
1167     OT_VECTEUR_3D normal_f_i = mgtri_i->calcul_normal();
1168     normal_f_i.norme();
1169     liste_normales2.ajouter(normal_f_i);
1170     //Remplissage de la liste des voisins du triangle i
1171     TPL_MAP_ENTITE<MG_TRIANGLE_PEAU*> liste_voisins;
1172     MG_NOEUD* noeud1 = mgtri_i->get_noeud1();
1173     double nb_voisins1 = noeud1->get_lien_triangle()->get_nb();
1174     for (int j = 0;j<nb_voisins1;j++)
1175     {
1176     MG_TRIANGLE_PEAU* mgtri_1 = (MG_TRIANGLE_PEAU*) noeud1->get_lien_triangle()->get(j);
1177     liste_voisins.ajouter(mgtri_1);
1178     }
1179     MG_NOEUD* noeud2 = mgtri_i->get_noeud2();
1180     double nb_voisins2 = noeud2->get_lien_triangle()->get_nb();
1181     for (int j = 0;j<nb_voisins2;j++)
1182     {
1183     MG_TRIANGLE_PEAU* mgtri_2 = (MG_TRIANGLE_PEAU*) noeud2->get_lien_triangle()->get(j);
1184     liste_voisins.ajouter(mgtri_2);
1185     }
1186     MG_NOEUD* noeud3 = mgtri_i->get_noeud3();
1187     double nb_voisins3 = noeud3->get_lien_triangle()->get_nb();
1188     for (int j = 0;j<nb_voisins3;j++)
1189     {
1190     MG_TRIANGLE_PEAU* mgtri_3 = (MG_TRIANGLE_PEAU*) noeud3->get_lien_triangle()->get(j);
1191     liste_voisins.ajouter(mgtri_3);
1192     }
1193     liste_voisins.supprimer(mgtri_i);
1194     int nb_voisins = liste_voisins.get_nb();
1195     double w_ij = 1./nb_voisins;
1196     double phi_i_min = 10.;
1197     double s_i = 0.0;
1198     double phi_im = 0.0;
1199     double *phi_ij = new double[nb_voisins];
1200     OT_VECTEUR_3D normal_f_i_mean(0.,0.,0.);
1201     normal_f_i_mean = normal_f_i;
1202     OT_VECTEUR_3D eta_i(0.,0.,0.);
1203     TPL_MAP_ENTITE<MG_TRIANGLE_PEAU*>::ITERATEUR it;
1204     int j = 0;
1205     for(MG_TRIANGLE_PEAU* mgtri_j=liste_voisins.get_premier(it);mgtri_j!=NULL;mgtri_j=liste_voisins.get_suivant(it))
1206     {
1207     OT_VECTEUR_3D normal_f_j = mgtri_j->calcul_normal();
1208     //1-Calculer la normale moyenne pour chaque triangle
1209     normal_f_i_mean = normal_f_i_mean + normal_f_j;
1210 francois 457 //2.1-On calcule l'angle entre normal_f_i et normal_f_j pour j allant de 1 a Nb_voisins
1211 cuillier 438 double prod_scalaire = normal_f_i*normal_f_j;
1212     if (prod_scalaire > 1.)
1213     {
1214     prod_scalaire = 1.;
1215     }
1216     if (prod_scalaire < -1.)
1217     {
1218     prod_scalaire = -1.;
1219     }
1220     phi_ij[j] = acos(prod_scalaire)*un_sur_pi;
1221     //2.2-On trouve le plus petit des angles et la normale heta_i correspondante
1222     if (phi_ij[j] < phi_i_min)
1223     {
1224     phi_i_min = phi_ij[j];
1225     eta_i = normal_f_j;
1226     }
1227 francois 457 //3.1-On calcule l'angle moyen phi_im
1228 cuillier 438 phi_im = phi_im + w_ij*phi_ij[j];
1229     j++;
1230     }
1231     normal_f_i_mean.norme();
1232     j = 0;
1233     for(MG_TRIANGLE_PEAU* mgtri_j=liste_voisins.get_premier(it);mgtri_j!=NULL;mgtri_j=liste_voisins.get_suivant(it))
1234     {
1235     //3.2-Calcul de s_i selon la variance
1236     s_i = s_i + w_ij*pow((phi_ij[j] - phi_im),2);
1237     j++;
1238     }
1239     delete[] phi_ij;
1240     //4-On calcule une nouvelle normale pour chaque triangle
1241 francois 458 double pond;
1242     int num=params.get_valeur("chen2008_filtre");
1243     if (num==1) pond=ponderation_gaussian(s_i,sigma);
1244     else if (num==2) pond=ponderation_laplacian(s_i,sigma);
1245     else if (num==3) pond=ponderation_elfallahford(s_i,sigma);
1246     OT_VECTEUR_3D normal_f_i_new = pond*normal_f_i_mean + (1. - pond)*eta_i;
1247 cuillier 438 normal_f_i_new.norme();
1248     liste_normales.ajouter(normal_f_i_new);
1249     liste_wij.ajouter(w_ij);
1250     mgtri->change_nouveau_numero(k);
1251     k++;
1252     }
1253    
1254     LISTE_MG_NOEUD::iterator it_no;
1255     for (MG_NOEUD* noeud_i=mg_mai->get_premier_noeud(it_no);noeud_i!=NULL;noeud_i=mg_mai->get_suivant_noeud(it_no))
1256     {
1257     int origine = noeud_i->get_origine();
1258     if (origine != IMPOSE)
1259     {
1260     int nb_voisins_j = noeud_i->get_lien_triangle()->get_nb();
1261     double w_ij_prime = 0.0;
1262     OT_VECTEUR_3D v_temp(0.,0.,0.);
1263     OT_VECTEUR_3D v_i(noeud_i->get_x(),noeud_i->get_y(),noeud_i->get_z());
1264     for(int j=0;j<nb_voisins_j;j++)
1265     {
1266     MG_TRIANGLE_PEAU* mgtri_j = (MG_TRIANGLE_PEAU*) noeud_i->get_lien_triangle()->get(j);
1267 francois 457 //On calcule le centroide cj du triangle mgtri_j
1268 cuillier 438 MG_NOEUD* n1 = mgtri_j->get_noeud1();
1269     MG_NOEUD* n2 = mgtri_j->get_noeud2();
1270     MG_NOEUD* n3 = mgtri_j->get_noeud3();
1271     double cj_x = 0.333333333333333*(n1->get_x() + n2->get_x() + n3->get_x());
1272     double cj_y = 0.333333333333333*(n1->get_y() + n2->get_y() + n3->get_y());
1273     double cj_z = 0.333333333333333*(n1->get_z() + n2->get_z() + n3->get_z());
1274     //On forme le vecteur vi_cj
1275     OT_VECTEUR_3D vi_cj(cj_x - noeud_i->get_x(),cj_y - noeud_i->get_y(),cj_z - noeud_i->get_z());
1276     OT_VECTEUR_3D normal_f_i_new = liste_normales.get(mgtri_j->get_nouveau_numero());
1277     v_temp = v_temp + (vi_cj*normal_f_i_new)*normal_f_i_new;
1278     }
1279 francois 457 //5-On met a jour la position des noeuds
1280 cuillier 438 v_i = v_i + (gamma/(2*nb_voisins_j))*v_temp;
1281     nouv_position_x.push_back(v_i.get_x());
1282     nouv_position_y.push_back(v_i.get_y());
1283     nouv_position_z.push_back(v_i.get_z());
1284     }
1285     }
1286    
1287     //On actualise la position des noeuds
1288     int ind_noeud=0;
1289     for (MG_NOEUD* noeud_i=mg_mai->get_premier_noeud(it_no);noeud_i!=NULL;noeud_i=mg_mai->get_suivant_noeud(it_no))
1290     {
1291     int origine = noeud_i->get_origine();
1292     if (origine != IMPOSE)
1293     {
1294     noeud_i->change_x(nouv_position_x[ind_noeud]);
1295     noeud_i->change_y(nouv_position_y[ind_noeud]);
1296     noeud_i->change_z(nouv_position_z[ind_noeud]);
1297     ind_noeud++;
1298     }
1299     }
1300    
1301     //Critere d'arret de l'algorithme
1302     int l=0;
1303     int nb_tri = mg_mai->get_nb_mg_triangle();
1304     for (MG_TRIANGLE* mgtri=mg_mai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mg_mai->get_suivant_triangle(it_tri))
1305     {
1306     MG_TRIANGLE_PEAU* mgtri_i = (MG_TRIANGLE_PEAU*)mgtri;
1307     OT_VECTEUR_3D normal_f_i = liste_normales2.get(mgtri_i->get_nouveau_numero());
1308     OT_VECTEUR_3D normal_f_i_new = liste_normales.get(mgtri_i->get_nouveau_numero());
1309     double critere = 1. - normal_f_i*normal_f_i_new;
1310 picher 248 if (critere <= epsilon) l++;
1311     }
1312     double tolerance = 0.01*nb_tri;
1313     if (nb_tri - l <= 0) fin = 1;
1314     compteur++;
1315     }
1316     while ((fin == 0) && (compteur < iter_max));
1317    
1318 francois 457 return compteur;
1319    
1320 picher 248 }
1321    
1322 cuillier 438
1323 francois 457 //---------------------------- Extract skin -------------------------------------------------------------------------
1324 francois 460 int MGOPT_POSTTRAITEMENT::extract_skin(MG_MAILLAGE* mg_mai,MG_GESTIONNAIRE& gest2,int &nbpeau,int &nbmaniare,int &nbmanino, int *mai2_id)
1325 cuillier 438 {
1326    
1327 francois 457 //etape 0 - On commence par mettre a zero tous les nouveau_numero des triangles et des noeuds du maillage
1328 francois 222 LISTE_MG_TRIANGLE::iterator it_tri;
1329     for (MG_TRIANGLE* mgtri=mg_mai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mg_mai->get_suivant_triangle(it_tri))
1330     {
1331     mgtri->change_nouveau_numero(0);
1332 francois 224 mgtri->change_origine(MAILLEUR_AUTO);
1333 francois 222 }
1334     LISTE_MG_NOEUD::iterator it_noeud;
1335     for (MG_NOEUD* mgnoeud=mg_mai->get_premier_noeud(it_noeud);mgnoeud!=NULL;mgnoeud=mg_mai->get_suivant_noeud(it_noeud))
1336     {
1337     mgnoeud->change_nouveau_numero(0);
1338     }
1339    
1340 francois 457 //etape 1 - On boucle ensuite tous les tetraedres pour ajouter un compteur du nombre de fois qu'un triangle appartient a 1 tetra
1341 picher 233 LISTE_MG_TETRA::iterator it_tetra;
1342     for (MG_TETRA* mgtet=mg_mai->get_premier_tetra(it_tetra);mgtet!=NULL;mgtet=mg_mai->get_suivant_tetra(it_tetra))
1343 francois 222 {
1344     int origine = mgtet->get_origine();
1345 francois 224 if (origine==IMPOSE)
1346 francois 222 {
1347 francois 224 mgtet->get_triangle1()->change_origine(IMPOSE);
1348     mgtet->get_triangle2()->change_origine(IMPOSE);
1349     mgtet->get_triangle3()->change_origine(IMPOSE);
1350     mgtet->get_triangle4()->change_origine(IMPOSE);
1351     }
1352     if (origine==OPTIMISE)
1353     {
1354     if (mgtet->get_triangle1()->get_origine()!=IMPOSE) mgtet->get_triangle1()->change_origine(OPTIMISE);
1355     if (mgtet->get_triangle2()->get_origine()!=IMPOSE) mgtet->get_triangle2()->change_origine(OPTIMISE);
1356     if (mgtet->get_triangle3()->get_origine()!=IMPOSE) mgtet->get_triangle3()->change_origine(OPTIMISE);
1357     if (mgtet->get_triangle4()->get_origine()!=IMPOSE) mgtet->get_triangle4()->change_origine(OPTIMISE);
1358     }
1359    
1360     if (((origine == OPTIMISE) && (etat[(OPTIMISE-1000)/10]==1) ) ||
1361     ((origine == IMPOSE) && (etat[(IMPOSE-1000)/10]==1) ) ||
1362     ((origine == MAILLEUR_AUTO) && (etat[(MAILLEUR_AUTO-1000)/10]==1) ) ||
1363     ((origine == TRIANGULATION) && (etat[(TRIANGULATION-1000)/10]==1) ) ||
1364     ((origine == MODIFICATION) && (etat[(MODIFICATION-1000)/10]==1) ) ||
1365     ((origine == DUPLIQUER) && (etat[(DUPLIQUER-1000)/10]==1) ) )
1366    
1367     {
1368 francois 222 int num1 = mgtet->get_triangle1()->get_nouveau_numero();
1369     int num2 = mgtet->get_triangle2()->get_nouveau_numero();
1370     int num3 = mgtet->get_triangle3()->get_nouveau_numero();
1371     int num4 = mgtet->get_triangle4()->get_nouveau_numero();
1372     num1++;
1373     num2++;
1374     num3++;
1375     num4++;
1376     mgtet->get_triangle1()->change_nouveau_numero(num1);
1377     mgtet->get_triangle2()->change_nouveau_numero(num2);
1378     mgtet->get_triangle3()->change_nouveau_numero(num3);
1379     mgtet->get_triangle4()->change_nouveau_numero(num4);
1380     }
1381     }
1382    
1383 francois 457 //etape 2 - On boucle l'ensemble des triangles identifies a l'etape 1 pour identifier les noeuds leur appartenant
1384 francois 222 for (MG_TRIANGLE* mgtri=mg_mai->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mg_mai->get_suivant_triangle(it_tri))
1385     {
1386     int num = mgtri->get_nouveau_numero();
1387     if (num == 1)
1388     {
1389     MG_NOEUD* noeud1 = mgtri->get_noeud1();
1390     MG_NOEUD* noeud2 = mgtri->get_noeud2();
1391     MG_NOEUD* noeud3 = mgtri->get_noeud3();
1392     noeud1->change_nouveau_numero(1);
1393     noeud2->change_nouveau_numero(1);
1394     noeud3->change_nouveau_numero(1);
1395 picher 248 if (mgtri->get_origine()==IMPOSE)
1396     {
1397     noeud1->change_origine(IMPOSE);
1398     noeud2->change_origine(IMPOSE);
1399     noeud3->change_origine(IMPOSE);
1400     }
1401 francois 222 }
1402     }
1403    
1404 francois 457 //etape 3 - On cree un nouveau maillage pour la peau
1405 francois 222 MG_MAILLAGE* mai2 = new MG_MAILLAGE(NULL);
1406     gest2.ajouter_mg_maillage(mai2);
1407    
1408 francois 457 //etape 4 - On boucle l'ensemble des noeuds identifies a l'etape 2 pour les recreer dans le second maillage
1409 francois 222 for (MG_NOEUD* mgnoeud=mg_mai->get_premier_noeud(it_noeud);mgnoeud!=NULL;mgnoeud=mg_mai->get_suivant_noeud(it_noeud))
1410     {
1411     int num = mgnoeud->get_nouveau_numero();
1412     if (num == 1)
1413     {
1414     double x = mgnoeud->get_x();
1415     double y = mgnoeud->get_y();
1416     double z = mgnoeud->get_z();
1417 francois 224 int origine=TRIANGULATION;
1418     if (mgnoeud->get_origine()==IMPOSE) origine=IMPOSE;
1419     MG_NOEUD* noeud1 = new MG_NOEUD(NULL,x,y,z,origine);
1420 francois 222 mai2->ajouter_mg_noeud(noeud1);
1421     mgnoeud->change_nouveau_numero(noeud1->get_id());
1422 francois 234 noeud1->change_nouveau_numero(mgnoeud->get_id());
1423 francois 222 }
1424     }
1425    
1426 francois 457 // etape 5 - On boucle l'ensemble des triangles identifies a l'etape 1 pour les recreer dans le maillage 2
1427 picher 233 for (MG_TETRA* mgtet=mg_mai->get_premier_tetra(it_tetra);mgtet!=NULL;mgtet=mg_mai->get_suivant_tetra(it_tetra))
1428 francois 222 {
1429 francois 224 int originetet=mgtet->get_origine();
1430     if (((originetet == OPTIMISE) && (etat[(OPTIMISE-1000)/10]==1) ) ||
1431     ((originetet == IMPOSE) && (etat[(IMPOSE-1000)/10]==1) ) ||
1432     ((originetet == MAILLEUR_AUTO) && (etat[(MAILLEUR_AUTO-1000)/10]==1) ) ||
1433     ((originetet == TRIANGULATION) && (etat[(TRIANGULATION-1000)/10]==1) ) ||
1434     ((originetet == MODIFICATION) && (etat[(MODIFICATION-1000)/10]==1) ) ||
1435     ((originetet == DUPLIQUER) && (etat[(DUPLIQUER-1000)/10]==1) ) )
1436 francois 222 {
1437 francois 224 MG_NOEUD* noeud1 = mgtet->get_noeud1();
1438     MG_NOEUD* noeud2 = mgtet->get_noeud2();
1439     MG_NOEUD* noeud3 = mgtet->get_noeud3();
1440     MG_NOEUD* noeud4 = mgtet->get_noeud4();
1441 francois 222 MG_NOEUD* node1 = mai2->get_mg_noeudid(noeud1->get_nouveau_numero());
1442     MG_NOEUD* node2 = mai2->get_mg_noeudid(noeud2->get_nouveau_numero());
1443     MG_NOEUD* node3 = mai2->get_mg_noeudid(noeud3->get_nouveau_numero());
1444 francois 224 MG_NOEUD* node4 = mai2->get_mg_noeudid(noeud4->get_nouveau_numero());
1445     MG_TRIANGLE* tri1=mgtet->get_triangle1();
1446     MG_TRIANGLE* tri2=mgtet->get_triangle2();
1447     MG_TRIANGLE* tri3=mgtet->get_triangle3();
1448     MG_TRIANGLE* tri4=mgtet->get_triangle4();
1449     if (tri1->get_nouveau_numero()==1)
1450     {
1451     int origine=TRIANGULATION;
1452     if (tri1->get_origine()==IMPOSE) origine=IMPOSE;
1453     MG_TRIANGLE_PEAU* tripeau = insere_triangle(NULL,node1,node2,node3,mai2,origine);
1454     tripeau->change_nouveau_numero(0);
1455     tri1->change_nouveau_numero(0);
1456     }
1457     if (tri2->get_nouveau_numero()==1)
1458     {
1459     int origine=TRIANGULATION;
1460     if (tri2->get_origine()==IMPOSE) origine=IMPOSE;
1461     MG_TRIANGLE_PEAU* tripeau = insere_triangle(NULL,node1,node4,node2,mai2,origine);
1462     tripeau->change_nouveau_numero(0);
1463     tri2->change_nouveau_numero(0);
1464     }
1465     if (tri3->get_nouveau_numero()==1)
1466     {
1467     int origine=TRIANGULATION;
1468     if (tri3->get_origine()==IMPOSE) origine=IMPOSE;
1469     MG_TRIANGLE_PEAU* tripeau = insere_triangle(NULL,node2,node4,node3,mai2,origine);
1470     tripeau->change_nouveau_numero(0);
1471     tri3->change_nouveau_numero(0);
1472     }
1473     if (tri4->get_nouveau_numero()==1)
1474     {
1475     int origine=TRIANGULATION;
1476     if (tri4->get_origine()==IMPOSE) origine=IMPOSE;
1477     MG_TRIANGLE_PEAU* tripeau = insere_triangle(NULL,node1,node3,node4,mai2,origine);
1478     tripeau->change_nouveau_numero(0);
1479     tri4->change_nouveau_numero(0);
1480     }
1481 francois 222 }
1482 francois 224 }
1483     // etape 6 recherche des vosins
1484     for (MG_TRIANGLE* mgtri=mai2->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mai2->get_suivant_triangle(it_tri))
1485     {
1486     MG_TRIANGLE_PEAU* tripeau=(MG_TRIANGLE_PEAU*)mgtri;
1487     MG_TRIANGLE_PEAU* voisin1=recherche_voisin(tripeau->get_noeud1(),tripeau->get_noeud2(),tripeau);
1488     MG_TRIANGLE_PEAU* voisin2=recherche_voisin(tripeau->get_noeud2(),tripeau->get_noeud3(),tripeau);
1489     MG_TRIANGLE_PEAU* voisin3=recherche_voisin(tripeau->get_noeud3(),tripeau->get_noeud1(),tripeau);
1490     tripeau->change_voisin1(voisin1);
1491     tripeau->change_voisin2(voisin2);
1492     tripeau->change_voisin3(voisin3);
1493     }
1494 francois 457 // etape 7 - On recherche les peaux
1495 francois 224 int fin;
1496     do
1497     {
1498     fin=1;
1499     for (MG_TRIANGLE* mgtri=mai2->get_premier_triangle(it_tri);mgtri!=NULL;mgtri=mai2->get_suivant_triangle(it_tri))
1500     {
1501     MG_TRIANGLE_PEAU *tripeau=(MG_TRIANGLE_PEAU*)mgtri;
1502     if (tripeau->get_nouveau_numero()==0)
1503     {
1504     fin=0;
1505     std::vector<MG_TRIANGLE_PEAU*> *peau=new std::vector<MG_TRIANGLE_PEAU*>;
1506     lst_peau.push_back(peau);
1507     tripeau->change_nouveau_numero(1);
1508     peau->push_back(tripeau);
1509     determine_peau(peau);
1510     }
1511     }
1512     }
1513     while (fin==0);
1514 francois 457 nbpeau=lst_peau.size();
1515 francois 224 //test de manifold
1516     LISTE_MG_SEGMENT::iterator itseg;
1517     for (MG_SEGMENT* seg=mai2->get_premier_segment(itseg);seg!=NULL;seg=mai2->get_suivant_segment(itseg))
1518     seg->change_nouveau_numero(0);
1519     TPL_MAP_ENTITE<MG_SEGMENT*> nonmanifoldarete;
1520     TPL_MAP_ENTITE<MG_NOEUD*> nonmanifoldnoeud;
1521 francois 234 TPL_MAP_ENTITE<MG_NOEUD*> nonmanifoldnoeuddepuisarete;
1522 francois 457 nbpeau=lst_peau.size();
1523 francois 224 for (int i=0;i<nbpeau;i++)
1524     {
1525     int nbtri=lst_peau[i]->size();
1526     for (int j=0;j<nbtri;j++)
1527     {
1528     MG_TRIANGLE_PEAU* tripeau=(*lst_peau[i])[j];
1529     tripeau->get_segment1()->change_nouveau_numero(tripeau->get_segment1()->get_nouveau_numero()+1);
1530     tripeau->get_segment2()->change_nouveau_numero(tripeau->get_segment2()->get_nouveau_numero()+1);
1531     tripeau->get_segment3()->change_nouveau_numero(tripeau->get_segment3()->get_nouveau_numero()+1);
1532     if (tripeau->get_segment1()->get_nouveau_numero()>2)
1533     nonmanifoldarete.ajouter(tripeau->get_segment1());
1534     if (tripeau->get_segment2()->get_nouveau_numero()>2)
1535     nonmanifoldarete.ajouter(tripeau->get_segment2());
1536     if (tripeau->get_segment3()->get_nouveau_numero()>2)
1537     nonmanifoldarete.ajouter(tripeau->get_segment3());
1538     }
1539     }
1540 picher 233
1541 francois 224 int nbnonmanifoldarete=nonmanifoldarete.get_nb();
1542     for (int i=0;i<nbnonmanifoldarete;i++)
1543     {
1544     MG_SEGMENT* seg=nonmanifoldarete.get(i);
1545     MG_NOEUD* n1=seg->get_noeud1();
1546     MG_NOEUD* n2=seg->get_noeud2();
1547 francois 234 nonmanifoldnoeuddepuisarete.ajouter(n1);
1548     nonmanifoldnoeuddepuisarete.ajouter(n2);
1549 francois 224 }
1550 francois 234 LISTE_MG_NOEUD::iterator itnoeud;
1551     for (MG_NOEUD* no=mai2->get_premier_noeud(itnoeud);no!=NULL;no=mai2->get_suivant_noeud(itnoeud))
1552     {
1553     if (nonmanifoldnoeuddepuisarete.existe(no)) continue;
1554     if (est_non_manifold(no)) nonmanifoldnoeud.ajouter(no);
1555     }
1556 francois 224
1557 francois 457 nbmaniare = nonmanifoldarete.get_nb();
1558     nbmanino = nonmanifoldnoeud.get_nb();
1559 picher 230
1560 francois 457 //etape 8 - Traitement des cas de non manifold
1561     //non manifold par arete
1562 picher 233 for (int i=0;i<nbnonmanifoldarete;i++)
1563     {
1564     MG_SEGMENT* segment=nonmanifoldarete.get(i);
1565 francois 234 MG_NOEUD* noeud1 = mg_mai->get_mg_noeudid(segment->get_noeud1()->get_nouveau_numero());
1566     MG_NOEUD* noeud2 = mg_mai->get_mg_noeudid(segment->get_noeud2()->get_nouveau_numero());
1567     int nb_tetra = noeud1->get_lien_tetra()->get_nb();
1568     for (int j=0;j<nb_tetra;j++)
1569 picher 233 {
1570 francois 234 MG_TETRA* mgtet =noeud1->get_lien_tetra()->get(j);
1571     int originetet=mgtet->get_origine();
1572     if (originetet == MAILLEUR_AUTO)
1573 picher 233 {
1574 francois 457 //On reactive le tetra si l'autre noeud du segment lui appartient aussi
1575 francois 234 MG_NOEUD* no1 = mgtet->get_noeud1();
1576     MG_NOEUD* no2 = mgtet->get_noeud2();
1577     MG_NOEUD* no3 = mgtet->get_noeud3();
1578     MG_NOEUD* no4 = mgtet->get_noeud4();
1579     if((no1 == noeud2) ||(no2 == noeud2) || (no3 == noeud2) || (no4 == noeud2))
1580     mgtet->change_origine(OPTIMISE);
1581 picher 233 }
1582     }
1583     }
1584    
1585     //non manifold par noeud
1586 picher 231 int nbnonmanifoldnoeud=nonmanifoldnoeud.get_nb();
1587     for (int i=0;i<nbnonmanifoldnoeud;i++)
1588     {
1589 francois 234 MG_NOEUD* noeud=mg_mai->get_mg_noeudid(nonmanifoldnoeud.get(i)->get_nouveau_numero());
1590     int nb_tetra = noeud->get_lien_tetra()->get_nb();
1591     for (int j = 0;j<nb_tetra;j++)
1592 francois 345 {
1593 francois 234 MG_TETRA* mgtet =noeud->get_lien_tetra()->get(j);
1594     int originetet=mgtet->get_origine();
1595     if (originetet == MAILLEUR_AUTO)
1596 picher 231 {
1597 francois 234 mgtet->change_origine(OPTIMISE);
1598 picher 231 }
1599     }
1600     }
1601 francois 234 *mai2_id = mai2->get_id();
1602     if ((nbnonmanifoldarete != 0) || (nbnonmanifoldnoeud != 0))
1603     {
1604     for (int i=0;i<lst_peau.size();i++)
1605     {
1606     delete lst_peau[i];
1607     }
1608     lst_peau.clear();
1609     gest2.supprimer_mg_maillageid(*mai2_id);
1610     return 0;
1611     }
1612 picher 230
1613 francois 234 return 1;
1614 francois 222 }
1615    
1616 francois 457 //---------------------------- Determine peau -------------------------------------------------------------------------
1617 francois 460 void MGOPT_POSTTRAITEMENT::determine_peau(std::vector<MG_TRIANGLE_PEAU*> * peau)
1618 francois 222 {
1619 francois 224 int num=0;
1620     while (num!=peau->size())
1621     {
1622     MG_TRIANGLE_PEAU* p=(*peau)[num];
1623     if (p->get_voisin1()->get_nouveau_numero()==0)
1624     {
1625     peau->push_back(p->get_voisin1());
1626     p->get_voisin1()->change_nouveau_numero(1);
1627     }
1628     if (p->get_voisin2()->get_nouveau_numero()==0)
1629     {
1630     peau->push_back(p->get_voisin2());
1631     p->get_voisin2()->change_nouveau_numero(1);
1632     }
1633     if (p->get_voisin3()->get_nouveau_numero()==0)
1634     {
1635     peau->push_back(p->get_voisin3());
1636     p->get_voisin3()->change_nouveau_numero(1);
1637     }
1638     num++;
1639     }
1640     }
1641    
1642 francois 457 //---------------------------- Insere triangle -------------------------------------------------------------------------
1643 francois 460 MG_TRIANGLE_PEAU* MGOPT_POSTTRAITEMENT::insere_triangle(MG_ELEMENT_TOPOLOGIQUE* topo,class MG_NOEUD *mgnoeud1,class MG_NOEUD *mgnoeud2,class MG_NOEUD *mgnoeud3,MG_MAILLAGE* mg_maillage,int origine)
1644 francois 224 {
1645 francois 222 MG_TRIANGLE_PEAU* voisin1=NULL,*voisin2=NULL,*voisin3=NULL;
1646     MG_SEGMENT* mgsegment1=mg_maillage->get_mg_segment(mgnoeud1->get_id(),mgnoeud2->get_id());
1647     MG_SEGMENT* mgsegment2=mg_maillage->get_mg_segment(mgnoeud2->get_id(),mgnoeud3->get_id());
1648     MG_SEGMENT* mgsegment3=mg_maillage->get_mg_segment(mgnoeud3->get_id(),mgnoeud1->get_id());
1649     if (mgsegment1==NULL)
1650 francois 224 mgsegment1=mg_maillage->ajouter_mg_segment(topo,mgnoeud1,mgnoeud2,origine);
1651 francois 222 if (mgsegment2==NULL)
1652 francois 224 mgsegment2=mg_maillage->ajouter_mg_segment(topo,mgnoeud2,mgnoeud3,origine);
1653 francois 222 if (mgsegment3==NULL)
1654 francois 224 mgsegment3=mg_maillage->ajouter_mg_segment(topo,mgnoeud3,mgnoeud1,origine);
1655     MG_TRIANGLE_PEAU* mtriangle=new MG_TRIANGLE_PEAU(topo,mgnoeud1,mgnoeud2,mgnoeud3,mgsegment1,mgsegment2,mgsegment3,origine);
1656 francois 222 mg_maillage->ajouter_mg_triangle(mtriangle);
1657     return mtriangle;
1658     }
1659    
1660 francois 457 //---------------------------- Recherche voisin -------------------------------------------------------------------------
1661 francois 460 MG_TRIANGLE_PEAU* MGOPT_POSTTRAITEMENT::recherche_voisin(MG_NOEUD* mg_noeud1,MG_NOEUD* mg_noeud2,MG_TRIANGLE_PEAU* triref)
1662 francois 222 {
1663 francois 224 MG_TRIANGLE_PEAU* trisol=NULL;
1664     double angleref=4.*M_PI;
1665 francois 222 int nb1=mg_noeud1->get_lien_triangle()->get_nb();
1666     int nb2=mg_noeud2->get_lien_triangle()->get_nb();
1667     for (int i=0;i<nb1;i++)
1668     for (int j=0;j<nb2;j++)
1669 francois 224 if (mg_noeud1->get_lien_triangle()->get(i)==mg_noeud2->get_lien_triangle()->get(j))
1670     if ((MG_TRIANGLE_PEAU*)mg_noeud1->get_lien_triangle()->get(i)!=triref)
1671     {
1672     double angle=calcul_angle(triref,(MG_TRIANGLE_PEAU*)mg_noeud1->get_lien_triangle()->get(i));
1673     if (angle<angleref)
1674     {
1675     angleref=angle;
1676     trisol=(MG_TRIANGLE_PEAU*)mg_noeud1->get_lien_triangle()->get(i);
1677     }
1678     }
1679     return trisol;
1680 francois 222 }
1681 francois 223
1682 francois 460 int MGOPT_POSTTRAITEMENT::est_non_manifold(MG_NOEUD* no)
1683 francois 234 {
1684     static int compteur=0;
1685     compteur++;
1686     int nb_tri=no->get_lien_triangle()->get_nb();
1687     TPL_MAP_ENTITE<MG_TRIANGLE*> lst;
1688     for (int i=0;i<nb_tri;i++)
1689     lst.ajouter(no->get_lien_triangle()->get(i));
1690     MG_TRIANGLE_PEAU* premier_triangle=(MG_TRIANGLE_PEAU*)lst.get(0);
1691     lst.supprimer(premier_triangle);
1692     MG_TRIANGLE_PEAU* tricour=(MG_TRIANGLE_PEAU*)premier_triangle;
1693     int i=0;
1694     do
1695     {
1696     if (lst.existe(tricour->get_voisin1()) || ((tricour->get_voisin1()==premier_triangle)&& (i>1)))
1697     {
1698     tricour=tricour->get_voisin1();
1699     lst.supprimer(tricour);
1700     }
1701     else if (lst.existe(tricour->get_voisin2()) || ((tricour->get_voisin2()==premier_triangle)&& (i>1)))
1702     {
1703     tricour=tricour->get_voisin2();
1704     lst.supprimer(tricour);
1705     }
1706     else if (lst.existe(tricour->get_voisin3()) || ((tricour->get_voisin3()==premier_triangle)&& (i>1)))
1707     {
1708     tricour=tricour->get_voisin3();
1709     lst.supprimer(tricour);
1710     }
1711     i++;
1712     }
1713     while (tricour!=premier_triangle);
1714     if (lst.get_nb()>0)
1715     return 1;
1716     return 0;
1717     }
1718 francois 224
1719 francois 457 //---------------------------- Calcule angle -------------------------------------------------------------------------
1720 francois 460 double MGOPT_POSTTRAITEMENT::calcul_angle(MG_TRIANGLE_PEAU* ft1,MG_TRIANGLE_PEAU* ft2)
1721 francois 223 {
1722 francois 388 MG_NOEUD* noeud1=ft1->get_noeud1();
1723     MG_NOEUD* noeud2=ft1->get_noeud2();
1724     MG_NOEUD* noeud3=ft1->get_noeud3();
1725     MG_NOEUD* noeud4=ft2->get_noeud1();
1726     MG_NOEUD* noeud5=ft2->get_noeud2();
1727     MG_NOEUD* noeud6=ft2->get_noeud3();
1728     double angle=get_angle_par_noeud<MG_NOEUD*>(noeud1,noeud2,noeud3,noeud4,noeud5,noeud6);
1729     return angle;
1730 francois 223 }
1731 picher 230
1732 francois 457 //---------------------------- Ponderation gaussian -------------------------------------------------------------------------
1733 francois 460 double MGOPT_POSTTRAITEMENT::ponderation_gaussian(double s,double sigma)
1734 picher 230 {
1735     double w_s;
1736     w_s = exp(-(s*s)/(2.*sigma*sigma));
1737     return w_s;
1738     }
1739    
1740 francois 457 //---------------------------- Ponderation Laplacian -------------------------------------------------------------------------
1741 francois 460 double MGOPT_POSTTRAITEMENT::ponderation_laplacian(double s,double sigma)
1742 picher 230 {
1743     double w_s;
1744     w_s = exp(-(sqrt(2.)*fabs(s))/sigma);
1745     return w_s;
1746     }
1747    
1748 francois 457 //---------------------------- Ponderation elfallahford -------------------------------------------------------------------------
1749 francois 460 double MGOPT_POSTTRAITEMENT::ponderation_elfallahford(double s,double sigma)
1750 picher 230 {
1751     double w_s;
1752     w_s = 1./sqrt(1+pow((s/sigma),2));
1753     return w_s;
1754     }