ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur_auto/src/remailleur.cpp
Revision: 283
Committed: Tue Sep 13 21:11:20 2011 UTC (13 years, 8 months ago) by francois
Original Path: magic/lib/mailleur/src/remailleur.cpp
File size: 52958 byte(s)
Log Message:
structure de l'écriture

File Contents

# User Rev Content
1 francois 283 //---------------------------------------------------------------------------
2    
3     #include "gestionversion.h"
4     #include <math.h>
5     #include "remailleur.h"
6     #include "mg_gestionnaire.h"
7     #include "mailleur3d.h"
8     #include "mailleur0d.h"
9     #include "mailleur1d.h"
10     #include "mailleur2d.h"
11     #include "ot_mathematique.h"
12     #include "mg_geometrie_outils.h"
13     //---------------------------------------------------------------------------
14     #pragma package(smart_init)
15    
16    
17    
18    
19    
20     REMAILLEUR::REMAILLEUR(class MG_GESTIONNAIRE* g1,MG_GESTIONNAIRE* g2,MG_GEOMETRIE* geo1,MG_GEOMETRIE* geo2,class FCT_TAILLE* fct_taille1,FCT_TAILLE* fct_taille2,MG_MAILLAGE* mori,MG_MAILLAGE* mmodi,VCT_COMPARAISON_RESULTAT& cmp):gestorigine(g1),gestmodifie(g2),carteori(fct_taille1),cartemod(fct_taille2),rescmp(cmp),maiorigine(mori),maimodifie(mmodi),geoorigine(geo1),geomodifie(geo2)
21     {
22     int nbvolume=geomodifie->get_nb_mg_volume();
23     lsttrifront=new TPL_LISTE_ENTITE<MG_TRIANGLE*>[nbvolume];
24     }
25    
26     REMAILLEUR::~REMAILLEUR()
27     {
28     //TPL_MAP_ENTITE<class MG_SEGMENT_FRONTIERE*>::ITERATEUR it;
29     //for (MG_SEGMENT_FRONTIERE* seg=lstsegfront.get_premier(it);seg!=NULL;seg=lstsegfront.get_suivant(it) ) delete seg;
30     delete octree;
31     delete [] lsttrifront;
32     }
33    
34    
35     void REMAILLEUR::maille(MG_GROUPE_TOPOLOGIQUE* mggt)
36     {
37     maille(20);
38     }
39    
40     void REMAILLEUR::maille(int etape)
41     {
42     //etape 0 initialisation
43     initialise_octree();
44     //etape 1 Destruction des entites autour des entites disparues
45     if (etape<1) return;
46     affiche((char*)" Destruction autour des disparitions");
47     TPL_MAP_ENTITE<MG_ELEMENT_MAILLAGE*> lstdetruire;
48     int nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_DISPARUE);
49     for (int i=0;i<nb;i++)
50     {
51     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_DISPARUE,i);
52     int dim=ele->get_dimension();
53     if (dim==0)
54     {
55     MG_SOMMET* som=(MG_SOMMET*)ele;
56     int nblien=som->get_lien_maillage()->get_nb();
57     for (int j=0;j<nblien;j++)
58     {
59     MG_NOEUD* no=(MG_NOEUD*)som->get_lien_maillage()->get(j);
60     if (maiorigine->get_mg_noeudid(no->get_id())!=NULL) lstdetruire.ajouter(no);
61     }
62     }
63     if (dim==1)
64     {
65     MG_ARETE* are=(MG_ARETE*)ele;
66     int nblien=are->get_lien_maillage()->get_nb();
67     for (int j=0;j<nblien;j++)
68     {
69     MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(j);
70     if (maiorigine->get_mg_segmentid(seg->get_id())!=NULL)
71     {
72     lstdetruire.ajouter(seg->get_noeud1());
73     lstdetruire.ajouter(seg->get_noeud2());
74     }
75     }
76     }
77     if (dim==2)
78     {
79     MG_FACE* face=(MG_FACE*)ele;
80     int nblien=face->get_lien_maillage()->get_nb();
81     for (int j=0;j<nblien;j++)
82     {
83     MG_TRIANGLE* tri=(MG_TRIANGLE*)face->get_lien_maillage()->get(j);
84     if (maiorigine->get_mg_triangleid(tri->get_id())!=NULL)
85     {
86     lstdetruire.ajouter(tri->get_noeud1());
87     lstdetruire.ajouter(tri->get_noeud2());
88     lstdetruire.ajouter(tri->get_noeud3());
89     }
90     }
91     }
92     }
93     TPL_MAP_ENTITE<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
94     for (MG_NOEUD* no=(MG_NOEUD*)lstdetruire.get_premier(it);no!=NULL;no=(MG_NOEUD*)lstdetruire.get_suivant(it))
95     {
96     double ereelle=0.;
97     int nbseg=no->get_lien_segment()->get_nb();
98     for (int i=0;i<nbseg;i++)
99     ereelle=ereelle+no->get_lien_segment()->get(i)->get_longueur();
100     ereelle=ereelle/nbseg;
101     detruit_noeud(no,2*ereelle);
102     }
103     //etape 2 : transfert des noeuds de sommets;
104     if (etape<2) return;
105     affiche((char*)" Transfert sommet");
106     nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE);
107     for (int i=0;i<nb;i++)
108     {
109     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE,i);
110     int dim=ele->get_dimension();
111     if (dim==0)
112     {
113     MG_SOMMET* som=(MG_SOMMET*)ele;
114     int nblien=som->get_lien_maillage()->get_nb();
115     for (int j=0;j<nblien;j++)
116     {
117     MG_NOEUD* no=(MG_NOEUD*)som->get_lien_maillage()->get(j);
118     if ((maiorigine->get_mg_noeudid(no->get_id())!=NULL) && (no->get_nouveau_numero()==CONSERVE))
119     {
120     MG_SOMMET* nvsom=geomodifie->get_mg_sommetid(rescmp.get_liste_correspondance_modifie(som->get_id()));
121     transfert_noeud(no,nvsom);
122     }
123     }
124    
125     }
126     }
127     // etape 3 maillage 0d
128     if (etape<3) return;
129     affiche((char*)" MAILLAGE 0D");
130     MAILLEUR0D m0d(maimodifie,geomodifie);
131     m0d.adapte();
132     // etape 4 destruction autour des nouveaux sommets
133     if (etape<4) return;
134     affiche((char*)" Destruction autour des nouveaux sommets");
135     lstdetruire.vide();
136     nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::MODIFIE_APPARUE);
137     for (int i=0;i<nb;i++)
138     {
139     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::MODIFIE_APPARUE,i);
140     int dim=ele->get_dimension();
141     if (dim==0)
142     {
143     MG_SOMMET* som=(MG_SOMMET*)ele;
144     int nblien=som->get_lien_maillage()->get_nb();
145     for (int j=0;j<nblien;j++)
146     {
147     MG_NOEUD* no=(MG_NOEUD*)som->get_lien_maillage()->get(j);
148     if (maimodifie->get_mg_noeudid(no->get_id())!=NULL) lstdetruire.ajouter(no);
149     }
150     }
151     }
152     for (MG_ELEMENT_MAILLAGE* ele=lstdetruire.get_premier(it);ele!=NULL;ele=lstdetruire.get_suivant(it))
153     {
154     detruit_noeud(ele,-1,1);
155     }
156     // etape 5 transferts des segments d arete
157     if (etape<5) return;
158     /*affiche(" Transfert arete");
159     nb=rescmp.get_nb_liste_topologie(ORIGINE_CONSERVEE);
160     for (int i=0;i<nb;i++)
161     {
162     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(ORIGINE_CONSERVEE,i);
163     int dim=ele->get_dimension();
164     if (dim==1)
165     {
166     MG_ARETE* are=(MG_ARETE*)ele;
167     int nblien=are->get_lien_maillage()->get_nb();
168     for (int j=0;j<nblien;j++)
169     {
170     MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(j);
171     if ((maiorigine->get_mg_segmentid(seg->get_id())!=NULL))
172     {
173     MG_NOEUD* n1=seg->get_noeud1();
174     MG_NOEUD* n2=seg->get_noeud2();
175     if (n1->get_nouveau_numero()==CONSERVE)
176     if (n2->get_nouveau_numero()==CONSERVE)
177     {
178     MG_ARETE* nvare=geomodifie->get_mg_areteid(rescmp.get_liste_correspondance(are->get_id()));
179     MG_NOEUD* nv1=get_liste_correspondance(n1);
180     MG_NOEUD* nv2=get_liste_correspondance(n2);
181     if (nv1==NULL) nv1=transfert_noeud(n1,nvare);
182     if (nv2==NULL) nv2=transfert_noeud(n2,nvare);
183     MG_SEGMENT *seg=new MG_SEGMENT(nvare,nv1,nv2,IMPOSE);
184     maimodifie->ajouter_mg_segment(seg);
185     }
186     }
187     }
188    
189     }
190     }
191     */
192     affiche((char*)" Creation des blocs de maille 1D");
193     int nbarete=geomodifie->get_nb_mg_arete();
194     nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE);
195     for (int i=0;i<nb;i++)
196     {
197     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE,i);
198     int dim=ele->get_dimension();
199     if (dim==1)
200     {
201     multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > lstseg;
202     MG_ARETE* arete=(MG_ARETE*)ele;
203     TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR ittet;
204     for (MG_ELEMENT_MAILLAGE* ele=arete->get_lien_maillage()->get_premier(ittet);ele!=NULL;ele=arete->get_lien_maillage()->get_suivant(ittet))
205     {
206     MG_SEGMENT* seg=(MG_SEGMENT*)ele;
207     if ((maiorigine->get_mg_segmentid(seg->get_id())!=NULL))
208     {
209     MG_NOEUD* n1=seg->get_noeud1();
210     MG_NOEUD* n2=seg->get_noeud2();
211     if (n1->get_nouveau_numero()==CONSERVE)
212     if (n2->get_nouveau_numero()==CONSERVE)
213     {
214     seg->change_nouveau_numero(NONTRAITE);
215     pair<unsigned long,MG_SEGMENT*> tmp(n1->get_id(),seg);
216     lstseg.insert(tmp);
217     pair<unsigned long,MG_SEGMENT*> tmp2(n2->get_id(),seg);
218     lstseg.insert(tmp2);
219     }
220     }
221     }
222    
223     cree_bloc_maille(lstseg,arete);
224     }
225     }
226     // etape 6 Maillage des aretes
227     if (etape<6) return;
228     affiche((char*)" MAILLAGE 1D");
229     int nbare=geomodifie->get_nb_mg_arete();
230     for (int i=0;i<nbare;i++)
231     {
232     MG_ARETE* arete=geomodifie->get_mg_arete(i);
233     TPL_LISTE_ENTITE<double> param;
234     int typeface=arete->get_courbe()->get_type_geometrique(param);
235     multimap<unsigned long,BLOC_MAILLE_1D*,less <unsigned long> >::iterator it=lstb1d.begin();
236     while (it!=lstb1d.end())
237     {
238     BLOC_MAILLE_1D *b1d=(*it).second;
239     if ((b1d->etat==NONATTACHE) && (b1d->type==typeface))
240     {
241     MG_NOEUD* no1=b1d->lst[0]->get_noeud1();
242     MG_NOEUD* no2=b1d->lst[0]->get_noeud2();
243     double xyz1[3],xyz2[3],xyz[3];
244     transfert_coord(no1->get_coord(),xyz1);
245     transfert_coord(no2->get_coord(),xyz2);
246     double t1,t2;
247     arete->inverser(t1,xyz1);
248     arete->inverser(t2,xyz2);
249     int sensidentique=1;
250     OT_VECTEUR_3D vec(xyz1,xyz2);
251     vec.norme();
252     double dxyz1[3];
253     arete->deriver(t1,dxyz1);
254     OT_VECTEUR_3D dvec(dxyz1);
255     dvec.norme();
256     double ps=dvec*vec;
257     if (ps <0) sensidentique=0;
258     if (sensidentique==0)
259     {
260     double t=t1;
261     t1=t2;
262     t2=t;
263     double xyz[3];
264     xyz[0]=xyz1[0];
265     xyz[1]=xyz1[1];
266     xyz[2]=xyz1[2];
267     xyz1[0]=xyz2[0];
268     xyz1[1]=xyz2[1];
269     xyz1[2]=xyz2[2];
270     xyz2[0]=xyz[0];
271     xyz2[1]=xyz[1];
272     xyz2[2]=xyz[2];
273     }
274     if (arete->get_courbe()->est_periodique()==1)
275     {
276     //if (t1<arete->get_tmin()) t1=t1+arete->get_courbe()->get_periode();
277     if (t2<t1+1e-6) t2=t2+arete->get_courbe()->get_periode();
278     }
279     double xyz1tmp[3],xyz2tmp[3];
280     arete->evaluer(t1,xyz1tmp);
281     arete->evaluer(t2,xyz2tmp);
282     OT_VECTEUR_3D vec1(xyz1,xyz1tmp);
283     OT_VECTEUR_3D vec2(xyz2,xyz2tmp);
284     double metrique[9];
285     cartemod->evaluer(xyz1,metrique);
286     double eps=(1./sqrt(metrique[0])*1e-6);
287     double t=0.5*(t1+t2);
288     if ((vec1.get_longueur()<eps)&&(vec2.get_longueur()<eps))
289     {
290     if (arete->get_courbe()->est_periodique()==1)
291     if (t<arete->get_tmin()) t=t+arete->get_courbe()->get_periode();
292     if ((t>=arete->get_tmin()) && (t<=arete->get_tmax()))
293     {
294     int nbseg=b1d->lst.size();
295     b1d->etat=ATTACHE;
296     for (int i=0;i<nbseg;i++)
297     {
298     MG_SEGMENT *seg=b1d->lst[i];
299     MG_NOEUD* no1=seg->get_noeud1();
300     MG_NOEUD* no2=seg->get_noeud2();
301     double xyz1[3],xy2[3];
302     transfert_coord(no1->get_coord(),xyz1);
303     transfert_coord(no2->get_coord(),xyz2);
304     MG_NOEUD* nv1=get_liste_correspondance(no1);
305     MG_NOEUD* nv2=get_liste_correspondance(no2);
306     if (nv1==NULL) nv1=transfert_noeud(no1,arete);
307     if (nv2==NULL) nv2=transfert_noeud(no2,arete);
308     MG_SEGMENT *nvseg;
309     if (sensidentique==1)
310     nvseg=new MG_SEGMENT(arete,nv1,nv2,IMPOSE);
311     else
312     nvseg=new MG_SEGMENT(arete,nv2,nv1,IMPOSE);
313     maimodifie->ajouter_mg_segment(nvseg);
314     }
315     }
316     }
317     }
318     it++;
319     }
320     MAILLEUR1D m1d(maimodifie,geomodifie,arete,cartemod);
321     m1d.adapte();
322     }
323     //cree_quadtree(maimodifie,quad,&lstsegfront);
324     // etape 7 Destruction autour des nouveaux segments
325     if (etape<7) return;
326     affiche((char*)" Destruction autour des nouvelles aretes");
327     lstdetruire.vide();
328     /*nb=rescmp.get_nb_liste_topologie(MODIFIE_APPARUE);
329     for (int i=0;i<nb;i++)
330     {
331     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(MODIFIE_APPARUE,i);
332     int dim=ele->get_dimension();
333     if (dim==1)
334     {
335     MG_ARETE* are=(MG_ARETE*)ele;
336     int nblien=are->get_lien_maillage()->get_nb();
337     for (int j=0;j<nblien;j++)
338     {
339     MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(j);
340     if (maimodifie->get_mg_segmentid(seg->get_id())!=NULL)
341     {
342     lstnoeuddetruire.ajouter(seg->get_noeud1());
343     lstnoeuddetruire.ajouter(seg->get_noeud2());
344     }
345     }
346     }
347     }*/
348     LISTE_MG_SEGMENT::iterator itnvseg;
349     for (MG_SEGMENT* seg=maimodifie->get_premier_segment(itnvseg);seg!=NULL;seg=maimodifie->get_suivant_segment(itnvseg))
350     {
351     if (seg->get_origine()==MAILLEUR_AUTO)
352     {
353     lstdetruire.ajouter(seg);
354     }
355     }
356     for (MG_ELEMENT_MAILLAGE* ele=lstdetruire.get_premier(it);ele!=NULL;ele=lstdetruire.get_suivant(it))
357     {
358     detruit_noeud(ele,-1,2);
359     }
360    
361     // etape 8 Creation des blocs de mailles 2D
362     if (etape<8) return;
363     affiche((char*)" Creation des blocs de maille 2D");
364     nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE);
365     for (int i=0;i<nb;i++)
366     {
367     MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE,i);
368     int dim=ele->get_dimension();
369     if (dim==2)
370     {
371     multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > lsttri;
372     MG_FACE* face=(MG_FACE*)ele;
373     TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR ittet;
374     int compteur=0,compteur2=0;
375     for (MG_ELEMENT_MAILLAGE* ele=face->get_lien_maillage()->get_premier(ittet);ele!=NULL;ele=face->get_lien_maillage()->get_suivant(ittet))
376     {
377     MG_TRIANGLE* tri=(MG_TRIANGLE*)ele;
378     if ((maiorigine->get_mg_triangleid(tri->get_id())!=NULL))
379     {
380     MG_NOEUD* n1=tri->get_noeud1();
381     MG_NOEUD* n2=tri->get_noeud2();
382     MG_NOEUD* n3=tri->get_noeud3();
383     compteur2++;
384     if (n1->get_nouveau_numero()==CONSERVE)
385     if (n2->get_nouveau_numero()==CONSERVE)
386     if (n3->get_nouveau_numero()==CONSERVE)
387     {
388     tri->change_nouveau_numero(NONTRAITE);
389     pair<unsigned long,MG_TRIANGLE*> tmp(n1->get_id(),tri);
390     lsttri.insert(tmp);
391     pair<unsigned long,MG_TRIANGLE*> tmp2(n2->get_id(),tri);
392     lsttri.insert(tmp2);
393     pair<unsigned long,MG_TRIANGLE*> tmp3(n3->get_id(),tri);
394     lsttri.insert(tmp3);
395     compteur++;
396     }
397     }
398     }
399    
400     cree_bloc_maille(lsttri,face);
401     }
402     }
403     // etape 9 MAILLAGE 2D
404     if (etape<9) return;
405     affiche((char*)" MAILLAGE 2D");
406     nbfacemod=geomodifie->get_nb_mg_face();
407     for (int i=0;i<nbfacemod;i++)
408     {
409     MG_GEOMETRIE_OUTILS ot;
410     TPL_LISTE_ENTITE<MG_TRIANGLE*> lsttricontraint;
411     MG_FACE* face=geomodifie->get_mg_face(i);
412     TPL_LISTE_ENTITE<double> param;
413     int typeface=face->get_surface()->get_type_geometrique(param);
414     multimap<unsigned long,BLOC_MAILLE_2D*,less <unsigned long> >::iterator it=lstb2d.begin();
415     while (it!=lstb2d.end())
416     {
417     BLOC_MAILLE_2D *b2d=(*it).second;
418     if ((b2d->etat==NONATTACHE) && (b2d->type==typeface))
419     {
420     MG_NOEUD* no1=b2d->lst[0]->get_noeud1();
421     MG_NOEUD* no2=b2d->lst[0]->get_noeud2();
422     MG_NOEUD* no3=b2d->lst[0]->get_noeud3();
423     double uv[2],uv1[2],uv2[2],uv3[2];
424     double xyz1[3],xyz2[3],xyz3[3];
425     transfert_coord(no1->get_coord(),xyz1);
426     transfert_coord(no2->get_coord(),xyz2);
427     transfert_coord(no3->get_coord(),xyz3);
428     face->inverser(uv1,xyz1);
429     face->inverser(uv2,xyz2);
430     face->inverser(uv3,xyz3);
431     double xyz1tmp[3];
432     face->evaluer(uv1,xyz1tmp);
433     double xyz2tmp[3];
434     face->evaluer(uv2,xyz2tmp);
435     double xyz3tmp[3];
436     face->evaluer(uv3,xyz3tmp);
437     OT_VECTEUR_3D vec1(xyz1,xyz1tmp);
438     OT_VECTEUR_3D vec2(xyz2,xyz2tmp);
439     OT_VECTEUR_3D vec3(xyz3,xyz3tmp);
440     double metrique[9];
441     cartemod->evaluer(xyz1,metrique);
442     double eps=(1./sqrt(metrique[0])*1e-6);
443     if (vec1.get_longueur()<eps)
444     if (vec2.get_longueur()<eps)
445     if (vec3.get_longueur()<eps)
446     {
447     OT_DECALAGE_PARAMETRE deca(face->get_surface()->get_periode_u(),face->get_surface()->get_periode_v());
448     double du=deca.calcul_decalage_parametre_u(uv1[0]);
449     double dv=deca.calcul_decalage_parametre_v(uv1[1]);
450     uv[0]=0.333333333333333333*(deca.decalage_parametre_u(uv1[0],du)+deca.decalage_parametre_u(uv2[0],du)+deca.decalage_parametre_u(uv3[0],du));
451     uv[1]=0.333333333333333333*(deca.decalage_parametre_v(uv1[1],dv)+deca.decalage_parametre_v(uv2[1],dv)+deca.decalage_parametre_v(uv3[1],dv));
452     uv[0]=deca.decalage_parametre_u(uv[0],-du);
453     uv[1]=deca.decalage_parametre_v(uv[1],-dv);
454     double dist=ot.calcule_distance_contour_face_uv(uv,face);
455     if (dist>0.)
456     {
457     int sensidentique=1;
458     OT_VECTEUR_3D vec12(xyz1,xyz2);
459     OT_VECTEUR_3D vec13(xyz1,xyz3);
460     double normal[3];
461     face->calcul_normale_unitaire(uv1,normal);
462     OT_VECTEUR_3D nor=vec13 & vec12;
463     nor.norme();
464     if ((nor*normal)<0.) sensidentique=0;
465     int nbtri=b2d->lst.size();
466     b2d->etat=ATTACHE;
467     //cout << " attachment du bloc " << b2d->face->get_id();
468     for (int i=0;i<nbtri;i++)
469     {
470     MG_TRIANGLE *tri=b2d->lst[i];
471     MG_NOEUD* no1=tri->get_noeud1();
472     MG_NOEUD* no2=tri->get_noeud2();
473     MG_NOEUD* no3=tri->get_noeud3();
474     double xyz1[3],xy2[3],xyz3[3];
475     transfert_coord(no1->get_coord(),xyz1);
476     transfert_coord(no2->get_coord(),xyz2);
477     transfert_coord(no3->get_coord(),xyz3);
478     MG_NOEUD* nvno1=new MG_NOEUD(NULL,xyz1[0],xyz1[1],xyz1[2],IMPOSE);
479     MG_NOEUD* nvno2=new MG_NOEUD(NULL,xyz2[0],xyz2[1],xyz2[2],IMPOSE);
480     MG_NOEUD* nvno3=new MG_NOEUD(NULL,xyz3[0],xyz3[1],xyz3[2],IMPOSE);
481     MG_TRIANGLE* nvtri;
482     if (sensidentique==1)
483     nvtri=new MG_TRIANGLE(NULL,nvno1,nvno2,nvno3,NULL,NULL,NULL,IMPOSE);
484     else
485     nvtri=new MG_TRIANGLE(NULL,nvno2,nvno1,nvno3,NULL,NULL,NULL,IMPOSE);
486     lsttricontraint.ajouter(nvtri);
487     nvno1->change_nouveau_numero(no1->get_id());
488     nvno2->change_nouveau_numero(no2->get_id());
489     nvno3->change_nouveau_numero(no3->get_id());
490     }
491     }
492     }
493     }
494     it++;
495     }
496     MAILLEUR2D m2d(maimodifie,geomodifie,face,cartemod);
497     char mess[50];
498     sprintf(mess," Face %i identificateur %lu",i,face->get_id());
499     affiche(mess);
500     m2d.maille(face,NULL,&(lsttricontraint));
501     int nbtricon=lsttricontraint.get_nb();
502     for (int i=0;i<nbtricon;i++)
503     {
504     MG_TRIANGLE* tri=lsttricontraint.get(0);
505     MG_NOEUD* noeud1=tri->get_noeud1();
506     MG_NOEUD* noeud2=tri->get_noeud2();
507     MG_NOEUD* noeud3=tri->get_noeud3();
508     lsttricontraint.supprimer(tri);
509     delete tri;
510     delete noeud1;
511     delete noeud2;
512     delete noeud3;
513     }
514     }
515     multimap<unsigned long,BLOC_MAILLE_2D*,less <unsigned long> >::iterator itfin2=lstb2d.begin();
516     while (itfin2!=lstb2d.end())
517     {
518     delete (*itfin2).second;
519     itfin2++;
520     }
521     cree_liste_frontiere(maimodifie);
522     // etape 10 Destruction autour des nouveaux triangles
523     if (etape<10) return;
524     affiche((char*)" Destruction autour des nouvelles faces");
525     lstdetruire.vide();
526     LISTE_MG_TRIANGLE::iterator itnvtri;
527     for (MG_TRIANGLE* tri=maimodifie->get_premier_triangle(itnvtri);tri!=NULL;tri=maimodifie->get_suivant_triangle(itnvtri))
528     {
529     if (tri->get_origine()==MAILLEUR_AUTO)
530     {
531     lstdetruire.ajouter(tri);
532     }
533     }
534     for (MG_ELEMENT_MAILLAGE* ele=lstdetruire.get_premier(it);ele!=NULL;ele=lstdetruire.get_suivant(it))
535     {
536     detruit_noeud(ele,-1,3);
537     }
538     // etape 11 Creation des blocs de mailles 3D
539     if (etape<11) return;
540     affiche((char*)" Creation des blocs de maille 3D");
541     int nbvolume=geoorigine->get_nb_mg_volume();
542     for (int i=0;i<nbvolume;i++)
543     {
544     multimap<unsigned long,MG_TETRA*,less <unsigned long> > lsttet;
545     MG_VOLUME* vol=geoorigine->get_mg_volume(i);
546     int nblien=vol->get_lien_maillage()->get_nb();
547     TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR ittet;
548     for (MG_ELEMENT_MAILLAGE* ele=vol->get_lien_maillage()->get_premier(ittet);ele!=NULL;ele=vol->get_lien_maillage()->get_suivant(ittet))
549     {
550     MG_TETRA* tet=(MG_TETRA*)ele;
551     if ((maiorigine->get_mg_tetraid(tet->get_id())!=NULL))
552     {
553     MG_NOEUD* n1=tet->get_noeud1();
554     MG_NOEUD* n2=tet->get_noeud2();
555     MG_NOEUD* n3=tet->get_noeud3();
556     MG_NOEUD* n4=tet->get_noeud4();
557     if (n1->get_nouveau_numero()==CONSERVE)
558     if (n2->get_nouveau_numero()==CONSERVE)
559     if (n3->get_nouveau_numero()==CONSERVE)
560     if (n4->get_nouveau_numero()==CONSERVE)
561     {
562     tet->change_nouveau_numero(NONTRAITE);
563     pair<unsigned long,MG_TETRA*> tmp(n1->get_id(),tet);
564     lsttet.insert(tmp);
565     pair<unsigned long,MG_TETRA*> tmp2(n2->get_id(),tet);
566     lsttet.insert(tmp2);
567     pair<unsigned long,MG_TETRA*> tmp3(n3->get_id(),tet);
568     lsttet.insert(tmp3);
569     pair<unsigned long,MG_TETRA*> tmp4(n4->get_id(),tet);
570     lsttet.insert(tmp4);
571     }
572     }
573    
574     }
575     cree_bloc_maille(lsttet,vol);
576     }
577     // etape 12 Maillage 3D
578     if (etape<12) return;
579     affiche((char*)" MAILLEUR 3D");
580     nbvolume=geomodifie->get_nb_mg_volume();
581     for (int i=0;i<nbvolume;i++)
582     {
583     TPL_LISTE_ENTITE<MG_TETRA*> lsttetcontraint;
584     MG_VOLUME* vol=geomodifie->get_mg_volume(i);
585     multimap<unsigned long,BLOC_MAILLE_3D*,less <unsigned long> >::iterator it=lstb3d.begin();
586     while (it!=lstb3d.end())
587     {
588     BLOC_MAILLE_3D *b3d=(*it).second;
589     if (b3d->etat==NONATTACHE)
590     {
591     MG_NOEUD* no1=b3d->lst[0]->get_noeud1();
592     MG_NOEUD* no2=b3d->lst[0]->get_noeud2();
593     MG_NOEUD* no3=b3d->lst[0]->get_noeud3();
594     MG_NOEUD* no4=b3d->lst[0]->get_noeud4();
595     double xyz1[3],xyz2[3],xyz3[3],xyz4[3],xyz[3];
596     transfert_coord(no1->get_coord(),xyz1);
597     transfert_coord(no2->get_coord(),xyz2);
598     transfert_coord(no3->get_coord(),xyz3);
599     transfert_coord(no4->get_coord(),xyz4);
600     xyz[0]=0.25*(xyz1[0]+xyz2[0]+xyz3[0]+xyz4[0]);
601     xyz[1]=0.25*(xyz1[1]+xyz2[1]+xyz3[1]+xyz4[1]);
602     xyz[2]=0.25*(xyz1[2]+xyz2[2]+xyz3[2]+xyz4[2]);
603     int res=point_appartient_volume(xyz,i);
604     if (res==1)
605     {
606     int nbtet=b3d->lst.size();
607     b3d->etat=ATTACHE;
608     //cout << "bloc recupere " << b3d->get_id() << endl;
609     for (int j=0;j<nbtet;j++)
610     {
611     MG_TETRA *tet=b3d->lst[j];
612     MG_NOEUD* no1=tet->get_noeud1();
613     MG_NOEUD* no2=tet->get_noeud2();
614     MG_NOEUD* no3=tet->get_noeud3();
615     MG_NOEUD* no4=tet->get_noeud4();
616     double xyz1[3],xy2[3],xyz3[3],xyz4[3];
617     transfert_coord(no1->get_coord(),xyz1);
618     transfert_coord(no2->get_coord(),xyz2);
619     transfert_coord(no3->get_coord(),xyz3);
620     transfert_coord(no4->get_coord(),xyz4);
621     MG_NOEUD* nvno1=new MG_NOEUD(NULL,xyz1[0],xyz1[1],xyz1[2],IMPOSE);
622     MG_NOEUD* nvno2=new MG_NOEUD(NULL,xyz2[0],xyz2[1],xyz2[2],IMPOSE);
623     MG_NOEUD* nvno3=new MG_NOEUD(NULL,xyz3[0],xyz3[1],xyz3[2],IMPOSE);
624     MG_NOEUD* nvno4=new MG_NOEUD(NULL,xyz4[0],xyz4[1],xyz4[2],IMPOSE);
625     MG_TETRA* nvtet=new MG_TETRA(NULL,nvno1,nvno2,nvno3,nvno4,NULL,NULL,NULL,NULL,IMPOSE);
626     lsttetcontraint.ajouter(nvtet);
627     nvno1->change_nouveau_numero(no1->get_id());
628     nvno2->change_nouveau_numero(no2->get_id());
629     nvno3->change_nouveau_numero(no3->get_id());
630     nvno4->change_nouveau_numero(no4->get_id());
631     nvtet->change_nouveau_numero(tet->get_id());
632     }
633     }
634     }
635     it++;
636     }
637     MAILLEUR3D m3d(maimodifie,geomodifie,vol,cartemod);
638     char mess[50];
639     sprintf(mess," Volume %i identificateur %lu",i,vol->get_id());
640     affiche(mess);
641     m3d.active_affichage(affiche);
642     m3d.maille(vol,NULL,&(lsttetcontraint));
643     int nbtetcon=lsttetcontraint.get_nb();
644     for (int i=0;i<nbtetcon;i++)
645     {
646     MG_TETRA* tet=lsttetcontraint.get(0);
647     MG_NOEUD* noeud1=tet->get_noeud1();
648     MG_NOEUD* noeud2=tet->get_noeud2();
649     MG_NOEUD* noeud3=tet->get_noeud3();
650     MG_NOEUD* noeud4=tet->get_noeud4();
651     lsttetcontraint.supprimer(tet);
652     delete tet;
653     delete noeud1;
654     delete noeud2;
655     delete noeud3;
656     delete noeud4;
657     }
658     }
659     multimap<unsigned long,BLOC_MAILLE_3D*,less <unsigned long> >::iterator itfin3=lstb3d.begin();
660     while (itfin3!=lstb3d.end())
661     {
662     delete (*itfin3).second;
663     itfin3++;
664     }
665     // etape 13 Recuperation des correspondances
666     if (etape<13) return;
667     affiche((char*)" Sauvegarde des correspondances");
668     LISTE_MG_NOEUD::iterator itnoeud;
669     for (MG_NOEUD* noeud=maimodifie->get_premier_noeud(itnoeud);noeud!=NULL;noeud=maimodifie->get_suivant_noeud(itnoeud))
670     {
671     if (noeud->get_origine()==IMPOSE)
672     rescmp.ajouter_correspondance(noeud->get_nouveau_numero(),noeud->get_id());
673     }
674     LISTE_MG_TETRA::iterator ittet;
675     for (MG_TETRA* tet=maimodifie->get_premier_tetra(ittet);tet!=NULL;tet=maimodifie->get_suivant_tetra(ittet))
676     {
677     if (tet->get_origine()==IMPOSE)
678     rescmp.ajouter_correspondance(tet->get_nouveau_numero(),tet->get_id());
679     }
680    
681     }
682    
683    
684     void REMAILLEUR::initialise_octree(void)
685     {
686     FCT_GENERATEUR_3D<4> *carte2=(FCT_GENERATEUR_3D<4> *)carteori;
687     int nb_cellule=carte2->get_nb_cellule();
688     double param[32];
689     carte2->get_information(0,0,param);
690     double xmin=param[0]+0.005*(param[1]-param[0]);
691     double ymin=param[8]+0.005*(param[10]-param[8]);
692     double zmin=param[16]+0.005*(param[20]-param[16]);
693     carte2->get_information(nb_cellule-1,0,param);
694     double xmax=param[1]-0.005*(param[1]-param[0]);
695     double ymax=param[10]-0.005*(param[10]-param[8]);
696     double zmax=param[20]-0.005*(param[20]-param[16]);
697     octree=new TPL_OCTREE_FCT<MG_NOEUD*,FCT_GENERATEUR_3D<4> >;
698     octree->initialiser(*carte2,xmin,ymin,zmin,xmax,ymax,zmax);
699     LISTE_MG_NOEUD::iterator it;
700     for (MG_NOEUD* noeud=maiorigine->get_premier_noeud(it);noeud!=NULL;noeud=maiorigine->get_suivant_noeud(it))
701     {
702     noeud->change_nouveau_numero(CONSERVE);
703     octree->inserer(noeud);
704     }
705    
706     }
707    
708     void REMAILLEUR::detruit_noeud(MG_ELEMENT_MAILLAGE* elebase,double distance,int type)
709     {
710     if (elebase->get_origine()!=MAILLEUR_AUTO) return;
711     double xyz[3];
712     double distmin;
713     OT_VECTEUR_3D direction;
714     if (type==0)
715     {
716     MG_NOEUD* no=(MG_NOEUD*)elebase;
717     OT_VECTEUR_3D vecnoeud(no->get_coord());
718     //OT_VECTEUR_3D vecnoeud2=rescmp.change_coord_mod_ori(vecnoeud);
719     OT_VECTEUR_3D vecnoeud2=vecnoeud;
720     xyz[0]=vecnoeud2.get_x();
721     xyz[1]=vecnoeud2.get_y();
722     xyz[2]=vecnoeud2.get_z();
723     }
724     if (type==1)
725     {
726     MG_NOEUD* no=(MG_NOEUD*)elebase;
727     OT_VECTEUR_3D vecnoeud(no->get_coord());
728     OT_VECTEUR_3D vecnoeud2=rescmp.change_coord_mod_ori(vecnoeud);
729     xyz[0]=vecnoeud2.get_x();
730     xyz[1]=vecnoeud2.get_y();
731     xyz[2]=vecnoeud2.get_z();
732     }
733     if (type==2)
734     {
735     MG_SEGMENT* seg=(MG_SEGMENT*)elebase;
736     OT_VECTEUR_3D vecnoeud1(seg->get_noeud1()->get_coord());
737     OT_VECTEUR_3D vecnoeud2(seg->get_noeud2()->get_coord());
738     OT_VECTEUR_3D vecnoeud=0.5*(vecnoeud1+vecnoeud2);
739     OT_VECTEUR_3D vecnoeud3=rescmp.change_coord_mod_ori(vecnoeud);
740     xyz[0]=vecnoeud3.get_x();
741     xyz[1]=vecnoeud3.get_y();
742     xyz[2]=vecnoeud3.get_z();
743     direction=vecnoeud2-vecnoeud1;
744     distmin=direction.get_longueur()/2.;
745     direction.norme();
746     }
747     if (type==3)
748     {
749     MG_TRIANGLE* tri=(MG_TRIANGLE*)elebase;
750     OT_VECTEUR_3D vecnoeud1(tri->get_noeud1()->get_coord());
751     OT_VECTEUR_3D vecnoeud2(tri->get_noeud2()->get_coord());
752     OT_VECTEUR_3D vecnoeud3(tri->get_noeud3()->get_coord());
753     OT_VECTEUR_3D vecnoeud=0.3333333333333333*(vecnoeud1+vecnoeud2+vecnoeud3);
754     OT_VECTEUR_3D vecnoeud4=rescmp.change_coord_mod_ori(vecnoeud);
755     xyz[0]=vecnoeud4.get_x();
756     xyz[1]=vecnoeud4.get_y();
757     xyz[2]=vecnoeud4.get_z();
758     distmin=std::max((vecnoeud-vecnoeud1).get_longueur(),(vecnoeud-vecnoeud2).get_longueur());
759     distmin=std::max((vecnoeud-vecnoeud3).get_longueur(),distmin);
760     direction=(vecnoeud3-vecnoeud1)&(vecnoeud2-vecnoeud1);
761     direction.norme();
762     }
763     TPL_MAP_ENTITE<MG_NOEUD*> lst;
764     double Ea,Eb;
765     if (type!=0)
766     {
767     double m[9];
768     carteori->evaluer(xyz,m);
769     Ea=1./sqrt(m[0]);
770     }
771     octree->rechercher(xyz[0],xyz[1],xyz[2],2*Ea,lst);
772     double alpha=1.05;
773     distmin=distmin*alpha;
774     LISTE_MG_NOEUD::iterator it;
775     for (MG_NOEUD* noeud=lst.get_premier(it);noeud!=NULL;noeud=lst.get_suivant(it))
776     {
777     double *xyztmp=noeud->get_coord();
778     OT_VECTEUR_3D vec(xyz,xyztmp);
779     if (type!=0)
780     {
781     double m[9];
782     carteori->evaluer(xyztmp,m);
783     Eb=1./sqrt(m[0]);
784     distance=std::max(1.5*Ea,1.5*Eb);
785     distance=std::max(distance,2.6*fabs(Ea-Eb));
786     }
787     if (type==2)
788     {
789     double cosinus=vec*direction/vec.get_longueur();
790     double sinus=sqrt(1-cosinus*cosinus);
791     double facteur;
792     if (OPERATEUR::egal(cosinus,0.,0.000001)) facteur=distance;
793     else
794     {
795     double delta=sinus*sinus+4.*distance*distance/distmin/distmin*cosinus*cosinus;
796     facteur=(-sinus+sqrt(delta))/2./distance/cosinus/cosinus*distmin*distmin;
797     distance=facteur;
798     }
799     }
800     if (type==3)
801     {
802     double sinus=vec*direction/vec.get_longueur();
803     if (sinus>-0.0000001)
804     {
805     double cosinus=sqrt(1-sinus*sinus);
806     double facteur;
807     if (OPERATEUR::egal(cosinus,0.,0.000001)) facteur=distance;
808     else
809     {
810     double delta=sinus*sinus+4*distance*distance/distmin/distmin*cosinus*cosinus;
811     facteur=(-sinus+sqrt(delta))/2./distance/cosinus/cosinus*distmin*distmin;
812     distance=facteur;
813     }
814     }
815    
816     }
817     if (vec.get_longueur()<distance)
818     {
819     noeud->change_nouveau_numero(DETRUIT);
820     }
821    
822     }
823    
824     }
825    
826    
827     /*void REMAILLEUR::cree_quadtree(MG_MAILLAGE* mai,TPL_QUADTREE<MG_SEGMENT_FRONTIERE*,MG_NOEUD*> **quad,TPL_MAP_ENTITE<MG_SEGMENT_FRONTIERE*> *lstsegfront)
828     {
829     int nbface=mai->get_mg_geometrie()->get_nb_mg_face();
830     for (int i=0;i<nbface;i++)
831     {
832     ((quad[i]))=new TPL_QUADTREE<MG_SEGMENT_FRONTIERE*,MG_NOEUD*>;
833     TPL_LISTE_ENTITE<MG_NOEUD*> lstn;
834     double umax=-1e300;
835     double vmax=-1e300;
836     double umin=1e300;
837     double vmin=1e300;
838     MG_FACE* face=mai->get_mg_geometrie()->get_mg_face(i);
839     int nbboucle=face->get_nb_mg_boucle();
840     int nbsegfravant=lstsegfront->get_nb();
841     for (int j=0;j<nbboucle;j++)
842     {
843     MG_BOUCLE* bou=face->get_mg_boucle(j);
844     int nbarete=bou->get_nb_mg_coarete();
845     for (int k=0;k<nbarete;k++)
846     {
847     MG_ARETE* are=bou->get_mg_coarete(k)->get_arete();
848     int nbseg=are->get_lien_maillage()->get_nb();
849     for (int l=0;l<nbseg;l++)
850     {
851     MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(l);
852     if (mai->get_mg_segmentid(seg->get_id())==NULL) continue;
853     double *xyz=seg->get_noeud1()->get_coord();
854     double uv[2];
855     face->inverser(uv,xyz);
856     seg->get_noeud1()->change_u(uv[0]);
857     seg->get_noeud1()->change_v(uv[1]);
858     xyz=seg->get_noeud2()->get_coord();
859     face->inverser(uv,xyz);
860     seg->get_noeud2()->change_u(uv[0]);
861     seg->get_noeud2()->change_v(uv[1]);
862     lstn.ajouter(seg->get_noeud2());
863     if (uv[0]<umin) umin=uv[0];
864     if (uv[0]>umax) umax=uv[0];
865     if (uv[1]<vmin) vmin=uv[1];
866     if (uv[1]>vmax) vmax=uv[1];
867     MG_SEGMENT_FRONTIERE* segfr=new MG_SEGMENT_FRONTIERE(seg);
868     lstsegfront->ajouter(segfr);
869     }
870     }
871     }
872     double periodeenu=face->get_surface()->get_periode_u();
873     double periodeenv=face->get_surface()->get_periode_v();
874     if (periodeenu!=0.0)
875     {
876     umin=0.;
877     umax=periodeenu;
878     }
879     else
880     {
881     double diff=umax-umin;
882     umin=umin-0.125*diff;
883     umax=umax+0.125*diff;
884     }
885     if (periodeenv!=0.0)
886     {
887     vmin=0.;
888     vmax=periodeenv;
889     }
890     else
891     {
892     double diff=vmax-vmin;
893     vmin=vmin-0.125*diff;
894     vmax=vmax+0.125*diff;
895     }
896     (quad[i])->initialiser(&lstn,1,umin,vmin,umax,vmax,face->get_surface()->get_periode_u(),face->get_surface()->get_periode_v());
897     int nbsegfr=lstsegfront->get_nb();
898     for (int j=nbsegfravant;j<nbsegfr;j++)
899     {
900     double *uv=lstsegfront->get(j)->get_uv1();
901     uv[0]=lstsegfront->get(j)->get_segment()->get_noeud1()->get_u();
902     uv[1]=lstsegfront->get(j)->get_segment()->get_noeud1()->get_v();
903     uv=lstsegfront->get(j)->get_uv2();
904     uv[0]=lstsegfront->get(j)->get_segment()->get_noeud2()->get_u();
905     uv[1]=lstsegfront->get(j)->get_segment()->get_noeud2()->get_v();
906     (quad[i])->inserer( lstsegfront->get(j));
907     }
908     }
909    
910     }
911     */
912     void REMAILLEUR::cree_liste_frontiere(MG_MAILLAGE* mai)
913     {
914     int nbvolume=mai->get_mg_geometrie()->get_nb_mg_volume();
915     for (int i=0;i<nbvolume;i++)
916     {
917     MG_VOLUME* vol=mai->get_mg_geometrie()->get_mg_volume(i);
918     int nbco=vol->get_nb_mg_coquille();
919     for (int j=0;j<nbco;j++)
920     {
921     MG_COQUILLE* coq=vol->get_mg_coquille(j);
922     int nbface=coq->get_nb_mg_coface();
923     for (int k=0;k<nbface;k++)
924     {
925     MG_FACE* face=coq->get_mg_coface(k)->get_face();
926     int nbseg=face->get_lien_maillage()->get_nb();
927     for (int l=0;l<nbseg;l++)
928     {
929     MG_TRIANGLE* tri=(MG_TRIANGLE*)face->get_lien_maillage()->get(l);
930     if (mai->get_mg_triangleid(tri->get_id())==NULL) continue;
931     lsttrifront[i].ajouter(tri);
932     }
933     }
934     }
935    
936    
937     }
938     }
939    
940    
941    
942     void REMAILLEUR::ajouter_liste_correspondance_noeud(MG_NOEUD* n1,MG_NOEUD* n2)
943     {
944     CORRESPONDANCENOEUD corr(n1,n2);
945     pair<unsigned long,CORRESPONDANCENOEUD> tmp(n1->get_id(),corr);
946     lstcorrnoeud.insert(tmp);
947     }
948    
949     MG_NOEUD* REMAILLEUR::get_liste_correspondance(MG_NOEUD* n1)
950     {
951     MG_NOEUD *n2=NULL;
952     std::map<unsigned long,CORRESPONDANCENOEUD,std::less<unsigned long> >::iterator j=lstcorrnoeud.find(n1->get_id());
953     if (j!=lstcorrnoeud.end())
954     n2=(*j).second.noeudmod;
955     return n2;
956     }
957    
958    
959     MG_NOEUD* REMAILLEUR::transfert_noeud(MG_NOEUD* no,MG_ELEMENT_TOPOLOGIQUE* ele)
960     {
961     double *xyz=no->get_coord();
962     OT_VECTEUR_3D pt(xyz);
963     OT_VECTEUR_3D nvpt=rescmp.change_coord_ori_mod(pt);
964     MG_NOEUD* noeud=new MG_NOEUD(ele,nvpt.get_x(),nvpt.get_y(),nvpt.get_z(),IMPOSE);
965     maimodifie->ajouter_mg_noeud(noeud);
966     ajouter_liste_correspondance_noeud(no,noeud);
967     noeud->change_nouveau_numero(no->get_id());
968     return noeud;
969     }
970    
971     void REMAILLEUR::transfert_coord(double *xyz,double *xyz2)
972     {
973     OT_VECTEUR_3D pt(xyz);
974     OT_VECTEUR_3D nvpt=rescmp.change_coord_ori_mod(pt);
975     xyz2[0]=nvpt.get_x();
976     xyz2[1]=nvpt.get_y();
977     xyz2[2]=nvpt.get_z();
978     return;
979     }
980    
981     void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > &lst,MG_ARETE* arete)
982     {
983     while (lst.size()>0)
984     {
985     multimap<unsigned long,MG_SEGMENT*,less <unsigned long> >::iterator it=lst.begin();
986     MG_SEGMENT* seg=(*it).second;
987     lst.erase(it);
988     TPL_LISTE_ENTITE<double> param;
989     BLOC_MAILLE_1D *b1d=new BLOC_MAILLE_1D(arete,arete->get_courbe()->get_type_geometrique(param));
990     pair<unsigned long,BLOC_MAILLE_1D*> tmp(b1d->get_id(),b1d);
991     lstb1d.insert(tmp);
992     b1d->lst.push_back(seg);
993     seg->change_nouveau_numero(TRAITE);
994     char mess[255];
995     sprintf(mess," Arete %lu , Bloc %d",arete->get_id(),lstb1d.size());
996     //affiche(mess);
997     int num=0;
998     while (num!=b1d->lst.size())
999     {
1000     MG_SEGMENT* itseg=b1d->lst[num++];
1001     MG_SEGMENT* seg=trouve_segment(lst,itseg->get_noeud1());
1002     while (seg!=NULL)
1003     {
1004     if (seg->get_nouveau_numero()==NONTRAITE)
1005     {
1006     b1d->lst.push_back(seg);
1007     seg->change_nouveau_numero(TRAITE);
1008     }
1009     seg=trouve_segment(lst,itseg->get_noeud1());
1010     }
1011     seg=trouve_segment(lst,itseg->get_noeud2());
1012     while (seg!=NULL)
1013     {
1014     if (seg->get_nouveau_numero()==NONTRAITE)
1015     {
1016     b1d->lst.push_back(seg);
1017     seg->change_nouveau_numero(TRAITE);
1018     }
1019     seg=trouve_segment(lst,itseg->get_noeud2());
1020     }
1021     }
1022    
1023     sprintf(mess," Segment %d ",b1d->lst.size());
1024     //affiche(mess);
1025     }
1026     }
1027    
1028     void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > &lst,MG_FACE* face)
1029     {
1030     while (lst.size()>0)
1031     {
1032     multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> >::iterator it=lst.begin();
1033     MG_TRIANGLE* tri=(*it).second;
1034     lst.erase(it);
1035     TPL_LISTE_ENTITE<double> param;
1036     BLOC_MAILLE_2D *b2d=new BLOC_MAILLE_2D(geoorigine->get_mg_faceid(face->get_id()),face->get_surface()->get_type_geometrique(param));
1037     pair<unsigned long,BLOC_MAILLE_2D*> tmp(b2d->get_id(),b2d);
1038     lstb2d.insert(tmp);
1039     b2d->lst.push_back(tri);
1040     tri->change_nouveau_numero(TRAITE);
1041     char mess[255];
1042     sprintf(mess," Face %lu , Bloc %d",face->get_id(),lstb2d.size());
1043     //affiche(mess);
1044     int num=0;
1045     while (num!=b2d->lst.size())
1046     {
1047     MG_TRIANGLE* ittri=b2d->lst[num++];
1048     MG_TRIANGLE* tri=trouve_triangle(lst,ittri->get_noeud1());
1049     while (tri!=NULL)
1050     {
1051     if (tri->get_nouveau_numero()==NONTRAITE)
1052     {
1053     b2d->lst.push_back(tri);
1054     tri->change_nouveau_numero(TRAITE);
1055     }
1056     tri=trouve_triangle(lst,ittri->get_noeud1());
1057     }
1058     tri=trouve_triangle(lst,ittri->get_noeud2());
1059     while (tri!=NULL)
1060     {
1061     if (tri->get_nouveau_numero()==NONTRAITE)
1062     {
1063     b2d->lst.push_back(tri);
1064     tri->change_nouveau_numero(TRAITE);
1065     }
1066     tri=trouve_triangle(lst,ittri->get_noeud2());
1067     }
1068     tri=trouve_triangle(lst,ittri->get_noeud3());
1069     while (tri!=NULL)
1070     {
1071     if (tri->get_nouveau_numero()==NONTRAITE)
1072     {
1073     b2d->lst.push_back(tri);
1074     tri->change_nouveau_numero(TRAITE);
1075     }
1076     tri=trouve_triangle(lst,ittri->get_noeud3());
1077     }
1078     }
1079    
1080     sprintf(mess," Triangle %d ",b2d->lst.size());
1081     //affiche(mess);
1082     }
1083     }
1084    
1085     void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_TETRA*,less <unsigned long> > &lst,MG_VOLUME* vol)
1086     {
1087     while (lst.size()>0)
1088     {
1089     multimap<unsigned long,MG_TETRA*,less <unsigned long> >::iterator it=lst.begin();
1090     MG_TETRA* tet=(*it).second;
1091     lst.erase(it);
1092     BLOC_MAILLE_3D *b3d=new BLOC_MAILLE_3D(geoorigine->get_mg_volumeid(vol->get_id()));
1093     pair<unsigned long,BLOC_MAILLE_3D*> tmp(b3d->get_id(),b3d);
1094     lstb3d.insert(tmp);
1095     b3d->lst.push_back(tet);
1096     tet->change_nouveau_numero(TRAITE);
1097     char mess[255];
1098     sprintf(mess," Volume %lu , Bloc %d",vol->get_id(),lstb3d.size());
1099     //affiche(mess);
1100     int num=0;
1101     while (num!=b3d->lst.size())
1102     {
1103     MG_TETRA* ittet=b3d->lst[num++];
1104     MG_TETRA* tet=trouve_tetra(lst,ittet->get_noeud1());
1105     while (tet!=NULL)
1106     {
1107     if (tet->get_nouveau_numero()==NONTRAITE)
1108     {
1109     b3d->lst.push_back(tet);
1110     tet->change_nouveau_numero(TRAITE);
1111     }
1112     tet=trouve_tetra(lst,ittet->get_noeud1());
1113     }
1114     tet=trouve_tetra(lst,ittet->get_noeud2());
1115     while (tet!=NULL)
1116     {
1117     if (tet->get_nouveau_numero()==NONTRAITE)
1118     {
1119     b3d->lst.push_back(tet);
1120     tet->change_nouveau_numero(TRAITE);
1121     }
1122     tet=trouve_tetra(lst,ittet->get_noeud2());
1123     }
1124     tet=trouve_tetra(lst,ittet->get_noeud3());
1125     while (tet!=NULL)
1126     {
1127     if (tet->get_nouveau_numero()==NONTRAITE)
1128     {
1129     b3d->lst.push_back(tet);
1130     tet->change_nouveau_numero(TRAITE);
1131     }
1132     tet=trouve_tetra(lst,ittet->get_noeud3());
1133     }
1134     tet=trouve_tetra(lst,ittet->get_noeud4());
1135     while (tet!=NULL)
1136     {
1137     if (tet->get_nouveau_numero()==NONTRAITE)
1138     {
1139     b3d->lst.push_back(tet);
1140     tet->change_nouveau_numero(TRAITE);
1141     }
1142     tet=trouve_tetra(lst,ittet->get_noeud4());
1143     }
1144     }
1145    
1146     sprintf(mess," Tetra %d ",b3d->lst.size());
1147     //affiche(mess);
1148     }
1149     }
1150    
1151    
1152     MG_SEGMENT* REMAILLEUR::trouve_segment(multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > &lst,MG_NOEUD* no)
1153     {
1154     MG_SEGMENT* seg=NULL;
1155     unsigned long id=no->get_id();
1156     multimap<unsigned long,MG_SEGMENT*,less <unsigned long> >::iterator it=lst.find(id);
1157     if (it!=lst.end())
1158     {
1159     seg=(*it).second;
1160     lst.erase(it);
1161     }
1162     return seg;
1163     }
1164    
1165     MG_TRIANGLE* REMAILLEUR::trouve_triangle(multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > &lst,MG_NOEUD* no)
1166     {
1167     MG_TRIANGLE* tri=NULL;
1168     unsigned long id=no->get_id();
1169     multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> >::iterator it=lst.find(id);
1170     if (it!=lst.end())
1171     {
1172     tri=(*it).second;
1173     lst.erase(it);
1174     }
1175     return tri;
1176     }
1177    
1178     MG_TETRA* REMAILLEUR::trouve_tetra(multimap<unsigned long,MG_TETRA*,less <unsigned long> > &lst,MG_NOEUD* no)
1179     {
1180     MG_TETRA* tet=NULL;
1181     unsigned long id=no->get_id();
1182     multimap<unsigned long,MG_TETRA*,less <unsigned long> >::iterator it=lst.find(id);
1183     if (it!=lst.end())
1184     {
1185     tet=(*it).second;
1186     lst.erase(it);
1187     }
1188     return tet;
1189     }
1190    
1191     int REMAILLEUR::point_appartient_volume(double *xyz,int numvol)
1192     {
1193     int ok=0;
1194     int i=0;
1195     multimap<double,double,less<double> > intersection;
1196     OT_VECTEUR_3D vecteur_dir;
1197     while (ok==0)
1198    
1199     {
1200     int ok2=0;
1201     while (ok2==0)
1202     {
1203     MG_TRIANGLE* tri=lsttrifront[numvol].get(i);
1204     MG_NOEUD* no1=tri->get_noeud1();
1205     MG_NOEUD* no2=tri->get_noeud2();
1206     MG_NOEUD* no3=tri->get_noeud3();
1207     double xb=(no1->get_x()+no2->get_x()+no3->get_x())/3.;
1208     double yb=(no1->get_y()+no2->get_y()+no3->get_y())/3.;
1209     double zb=(no1->get_z()+no2->get_z()+no3->get_z())/3.;
1210     OT_VECTEUR_3D directeur(xb-xyz[0],yb-xyz[1],zb-xyz[2]);
1211     directeur.norme();
1212     vecteur_dir=directeur;
1213     OT_VECTEUR_3D n1n3(no1->get_coord(),no3->get_coord());
1214     OT_VECTEUR_3D n1n2(no1->get_coord(),no2->get_coord());
1215     OT_VECTEUR_3D n=n1n3&n1n2;
1216     n.norme();
1217     double ps1=n*directeur;
1218     if (fabs(ps1)<1e-6) i++;
1219     else ok2=1;
1220     }
1221     ok2=0;
1222     intersection.clear();
1223     for (int j=0;j<lsttrifront[numvol].get_nb();j++)
1224     {
1225     MG_TRIANGLE* tri=lsttrifront[numvol].get(j);
1226    
1227     double t;
1228     int res=inter_droite_triangle(xyz,vecteur_dir,tri,&t);
1229     if (res==2) {
1230     ok2=1;
1231     break;
1232     }
1233     if (res==1)
1234     {
1235     pair<double,double> tmp(t,t);
1236     intersection.insert(tmp);
1237     }
1238     }
1239     if (ok2==0) ok=1;
1240     i++;
1241    
1242     }
1243     multimap<double,double,less <double> >::iterator it=intersection.begin();
1244     int nbinterneg=0;
1245     double tavant=1e308;
1246     while (it!=intersection.end())
1247     {
1248     if ((*it).second<0.)
1249     {
1250     double t=(*it).second;
1251     if (!(OPERATEUR::egal(t,tavant,1e-6)))
1252     nbinterneg++;
1253     tavant=t;
1254     }
1255     it++;
1256     }
1257     if (nbinterneg%2==1) return 1;
1258     return 0;
1259     }
1260    
1261    
1262     int REMAILLEUR::inter_droite_triangle(double *xyz,double *dir,MG_TRIANGLE* tri,double *t)
1263     {
1264     int inter=0;
1265     MG_NOEUD *noeud1=tri->get_noeud1();
1266     MG_NOEUD *noeud2=tri->get_noeud2();
1267     MG_NOEUD *noeud3=tri->get_noeud3();
1268     OT_VECTEUR_3D g1(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
1269     OT_VECTEUR_3D g2(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
1270     OT_VECTEUR_3D g3(dir);
1271     OT_MATRICE_3D systeme(g1,g2,g3);
1272     double det=systeme.get_determinant();
1273     double eps=0.333333333333*(g1.diff()+g2.diff()+g3.diff());
1274     eps=eps*eps*0.0018;
1275     if (OPERATEUR::egal(det,0.0,eps)==true)
1276     {
1277     OT_VECTEUR_3D g3b(xyz[0]-noeud1->get_x(),xyz[1]-noeud1->get_y(),xyz[2]-noeud1->get_z());
1278     OT_MATRICE_3D systeme2(g1,g2,g3b);
1279     double det2=systeme2.get_determinant();
1280     double eps2=0.333333333333*(g1.diff()+g2.diff()+g3b.diff());
1281     eps2=18.*eps2*eps2*eps2*1e-6;
1282     if (OPERATEUR::egal(det2,0.0,eps)==true) // cas 2D
1283     {
1284     return 2;
1285     }
1286     else return 0;
1287     }
1288     else
1289     {
1290     double x,y,z;
1291     x=1.0/det*(g2.get_y()*g3.get_z()-g2.get_z()*g3.get_y());
1292     y=1.0/det*(g3.get_x()*g2.get_z()-g2.get_x()*g3.get_z());
1293     z=1.0/det*(g2.get_x()*g3.get_y()-g2.get_y()*g3.get_x());
1294     OT_VECTEUR_3D g1b(x,y,z);
1295     x=1.0/det*(g3.get_y()*g1.get_z()-g1.get_y()*g3.get_z());
1296     y=1.0/det*(g1.get_x()*g3.get_z()-g3.get_x()*g1.get_z());
1297     z=1.0/det*(g3.get_x()*g1.get_y()-g1.get_x()*g3.get_y());
1298     OT_VECTEUR_3D g2b(x,y,z);
1299     x=1.0/det*(g1.get_y()*g2.get_z()-g1.get_z()*g2.get_y());
1300     y=1.0/det*(g2.get_x()*g1.get_z()-g1.get_x()*g2.get_z());
1301     z=1.0/det*(g1.get_x()*g2.get_y()-g1.get_y()*g2.get_x());
1302     OT_VECTEUR_3D g3b(x,y,z);
1303     OT_VECTEUR_3D n1n4(xyz[0]-noeud1->get_x(),xyz[1]-noeud1->get_y(),xyz[2]-noeud1->get_z());
1304     double alpha1=n1n4*g1b;
1305     double alpha2=n1n4*g2b;
1306     double alpha3=-(n1n4*g3b);
1307     double alpha4=1-alpha1-alpha2;
1308     double eps=0.000001;
1309     if ((alpha1>-eps) && (alpha1<1.+eps))
1310     if ((alpha2>-eps) && (alpha2<1.+eps))
1311     if ((alpha4>-eps) && (alpha4<1.+eps)) {
1312     *t=alpha3;
1313     return 1;
1314     }
1315     return 0;
1316     }
1317    
1318     }