ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/CAD4FE/src/CAD4FE_MakeLoops.cpp
Revision: 1158
Committed: Thu Jun 13 22:18:49 2024 UTC (11 months ago) by francois
File size: 14976 byte(s)
Log Message:
compatibilité Ubuntu 22.04
Suppression des refeences à Windows
Ajout d'une banière

File Contents

# User Rev Content
1 francois 1158 //####//------------------------------------------------------------
2     //####//------------------------------------------------------------
3     //####// MAGiC
4     //####// Jean Christophe Cuilliere et Vincent FRANCOIS
5     //####// Departement de Genie Mecanique - UQTR
6     //####//------------------------------------------------------------
7     //####// MAGIC est un projet de recherche de l equipe ERICCA
8     //####// du departement de genie mecanique de l Universite du Quebec a Trois Rivieres
9     //####// http://www.uqtr.ca/ericca
10     //####// http://www.uqtr.ca/
11     //####//------------------------------------------------------------
12     //####//------------------------------------------------------------
13     //####//
14     //####// CAD4FE_MakeLoops.cpp
15     //####//
16     //####//------------------------------------------------------------
17     //####//------------------------------------------------------------
18     //####// COPYRIGHT 2000-2024
19     //####// jeu 13 jun 2024 11:58:56 EDT
20     //####//------------------------------------------------------------
21     //####//------------------------------------------------------------
22 foucault 27
23     #include <sstream>
24     #include <string>
25    
26     #pragma hdrstop
27    
28     #include "gestionversion.h"
29    
30     #include "ot_mathematique.h"
31     #include "ot_algorithme_geometrique.h"
32     #include "mg_coarete.h"
33     #include "mg_face.h"
34     #include "mg_arete.h"
35     #include "mg_boucle.h"
36     #include "mg_maillage.h"
37    
38 foucault 569 #include "CAD4FE_Geometric_Tools.h"
39 foucault 27
40     #include "CAD4FE_MakeLoops.h"
41     #include "CAD4FE_MCVertex.h"
42     #include "CAD4FE_MCEdge.h"
43     #include "CAD4FE_MCFace.h"
44    
45    
46     #pragma package(smart_init)
47    
48     using namespace CAD4FE;
49    
50    
51     MakeLoops::CoEdge::CoEdge(MG_FACE * __face, MG_ARETE * __e, int __sense):f(__face),sense(__sense),e(__e),prev(0),next(0){}
52     MakeLoops::CoEdge::~CoEdge()
53     {
54    
55     }
56     int MakeLoops::CoEdge::Sense(){return sense;}
57     MG_ARETE * MakeLoops::CoEdge::Edge(){return e;}
58     MG_FACE * MakeLoops::CoEdge::Face(){return f;}
59     bool MakeLoops::CoEdge::IsInverse(CoEdge & __mcEdgeSense)
60     {
61     return (__mcEdgeSense.e == e && __mcEdgeSense.Sense() == -Sense() ) ;
62     }
63     MG_SOMMET * MakeLoops::CoEdge::StartVertex()
64     {
65     if (Sense()==1)
66     return e->get_cosommet1()->get_sommet();
67     else
68     return e->get_cosommet2()->get_sommet();
69     }
70     MG_SOMMET * MakeLoops::CoEdge::EndVertex()
71     {
72     if (Sense()==1)
73     return e->get_cosommet2()->get_sommet();
74     else
75     return e->get_cosommet1()->get_sommet();
76     }
77     MG_COSOMMET * MakeLoops::CoEdge::StartCoVertex()
78     {
79     if (Sense()==1)
80     return e->get_cosommet1();
81     else
82     return e->get_cosommet2();
83     }
84     MG_COSOMMET * MakeLoops::CoEdge::EndCoVertex()
85     {
86     if (Sense()==1)
87     return e->get_cosommet2();
88     else
89     return e->get_cosommet1();
90     }
91     bool MakeLoops::CoEdge::IsAfter(CoEdge & __mcEdgeSense)
92     {
93     return this->StartVertex() == __mcEdgeSense.EndVertex();
94     }
95     bool MakeLoops::CoEdge::IsBefore(CoEdge & __mcEdgeSense)
96     {
97     return this->EndVertex() == __mcEdgeSense.StartVertex();
98     }
99     OT_VECTEUR_3D MakeLoops::CoEdge::StartDir()
100     {
101     OT_VECTEUR_3D dir;
102     double t=StartCoVertex()->get_t();
103     e->deriver(t,dir);
104     dir=Sense()*dir;
105     dir.norme();
106     return dir;
107     }
108    
109     OT_VECTEUR_3D MakeLoops::CoEdge::EndDir()
110     {
111     OT_VECTEUR_3D dir;
112     double t=EndCoVertex()->get_t();
113     e->deriver(t,dir);
114     dir=Sense()*dir;
115     dir.norme();
116     return dir;
117     }
118    
119     OT_VECTEUR_3D MakeLoops::CoEdge::StartMeshDir(MG_MAILLAGE * __mesh)
120     {
121     OT_VECTEUR_3D result;
122     MG_SOMMET * v1=StartVertex(), *v2=EndVertex();
123     int meshSense = 0; // if meshSense = sense Then use segment which direction is equal to sense
124     // if meshSense = -sense Then use segment which direction is opposite to sense
125     // if meshSense = 0 Then use any segment to compute mesh derivative
126     if (v1==v2)
127     meshSense = Sense();
128     else
129     meshSense = 0;
130    
131     result=MeshDir(__mesh, v1, meshSense);
132     if (result.get_longueur2()>1E100)
133     {
134     result = StartDir();
135     }
136    
137     return result;
138     }
139     OT_VECTEUR_3D MakeLoops::CoEdge::EndMeshDir(MG_MAILLAGE * __mesh)
140     {
141     OT_VECTEUR_3D result;
142     MG_SOMMET * v1=StartVertex(), *v2=EndVertex();
143     int meshSense = 0;
144     if (v1==v2) // if meshSense = sense Then use segment which direction is equal to sense
145     // if meshSense = -sense Then use segment which direction is opposite to sense
146     // if meshSense = 0 Then use any segment to compute mesh derivative
147     meshSense = -Sense();
148     else
149     meshSense = 0;
150    
151     result=MeshDir(__mesh, v2, meshSense);
152     if (result.get_longueur2()>1E100)
153     {
154     result = EndDir();
155     }
156    
157     return result;
158     }
159     OT_VECTEUR_3D MakeLoops::CoEdge::MeshDir(MG_MAILLAGE * __mesh, MG_SOMMET * __v, int __meshSense)
160     {
161     OT_VECTEUR_3D dir;
162     TPL_SET < MG_ELEMENT_MAILLAGE *> * mesh = __v->get_lien_maillage();
163     TPL_SET < MG_ELEMENT_MAILLAGE *>::ITERATEUR it;
164     MG_NOEUD * n;
165     MG_SEGMENT * vertexSeg = NULL;
166     double t1,t2,*xyz1,*xyz2, t1o, o, t2o;
167     double period=e->get_courbe()->get_periode();
168     for (n = (MG_NOEUD*)mesh->get_premier(it); n; n = (MG_NOEUD*)mesh->get_suivant(it))
169     if (!__mesh || __mesh->contient(n))
170     {
171     int nSeg = n->get_lien_segment()->get_nb();
172     TPL_LISTE_ENTITE< MG_SEGMENT * > * lstSeg = n->get_lien_segment();
173     for (int i=0; i<nSeg; i++)
174     {
175    
176     MG_SEGMENT * seg = lstSeg->get(i);
177     if (seg->get_lien_topologie() == e)
178     {
179     bool valid_sense;
180    
181     if (__meshSense == +1 || __meshSense == -1)
182     {
183     if (seg->get_noeud1()->get_lien_topologie()==__v)
184     {
185     xyz1=seg->get_noeud1()->get_coord();
186     xyz2=seg->get_noeud2()->get_coord();
187     }
188     else
189     {
190     xyz1=seg->get_noeud2()->get_coord();
191     xyz2=seg->get_noeud1()->get_coord();
192     }
193    
194     e->inverser(t1, xyz1);
195     e->inverser(t2, xyz2);
196    
197     t1o=.5*period;
198     o=t1o-t1;
199     t2o=t2+o;
200     if (t2o>period)
201     t2o -= period;
202     if (t2o<0)
203     t2o += period;
204    
205     if (__meshSense == +1)
206     valid_sense = (t1o<t2o);
207     else if (__meshSense == -1)
208     valid_sense = (t1o>t2o);
209     }
210     else
211     valid_sense = true;
212    
213     if (valid_sense)
214     {
215     if (vertexSeg == NULL || seg->get_longueur() > vertexSeg->get_longueur())
216     vertexSeg = seg;
217     }
218     }
219     }
220     }
221    
222     MG_NOEUD * vertexSegN[2];
223    
224     if (vertexSeg == NULL)
225     {
226     return OT_VECTEUR_3D(1E308,1E308,1E308);
227     }
228    
229     xyz1=vertexSeg->get_noeud1()->get_coord();
230     xyz2=vertexSeg->get_noeud2()->get_coord();
231     e->inverser(t1, xyz1);
232     e->inverser(t2, xyz2);
233    
234     t1o=.5*e->get_courbe()->get_periode();
235     o=t1o-t1;
236     t2o=t2+o;
237     if (t2o>period)
238     t2o -= period;
239     if (t2o<0)
240     t2o += period;
241    
242     dir = OT_VECTEUR_3D(xyz2);
243     dir -= OT_VECTEUR_3D(xyz1);
244     dir /= Sense()*(t2o-t1o);
245     return dir;
246     }
247    
248     double MakeLoops::CoEdge::MeshAngleInPlane(CoEdge & __coEdge, OT_VECTEUR_3D & __normal, MG_MAILLAGE * __mesh)
249     {
250     OT_MATRICE_3D repereNormal, transform3DToRepereNormal;
251     OT_VECTEUR_3D dir1_3D, dir2_3D, dir1_2D, dir2_2D, normal_2D;
252     if (IsInverse(__coEdge))
253     return -M_PI;
254     if (IsAfter(__coEdge))
255     {
256     dir2_3D=StartMeshDir(__mesh);
257     dir1_3D=__coEdge.EndMeshDir(__mesh);
258     }
259     else if (IsBefore(__coEdge))
260     {
261     dir1_3D=EndMeshDir(__mesh);
262     dir2_3D=__coEdge.StartMeshDir(__mesh);
263     }
264     else return 1E308;
265     double angle = AngleInPlane(dir1_3D, dir2_3D, __normal);
266 francois 1075 printf("Angle Between %lu and %lu = %f\n", e->get_id(), __coEdge.e->get_id(), angle);
267 foucault 27 return angle;
268     }
269    
270     double MakeLoops::CoEdge::AngleInPlane(CoEdge & __coEdge, OT_VECTEUR_3D & __normal)
271     {
272     OT_VECTEUR_3D dir1_3D, dir2_3D;
273     if (IsInverse(__coEdge))
274     return -M_PI;
275     if (IsAfter(__coEdge))
276     {
277     dir2_3D=StartDir();
278     dir1_3D=__coEdge.EndDir();
279     }
280     else if (IsBefore(__coEdge))
281     {
282     dir1_3D=EndDir();
283     dir2_3D=__coEdge.StartDir();
284     }
285     else return 1E308;
286     double angle = AngleInPlane(dir1_3D, dir2_3D, __normal);
287 francois 1075 printf("Angle Between %lu and %lu = %f\n", e->get_id(), __coEdge.e->get_id(), angle);
288 foucault 27 return angle;
289     }
290    
291     double MakeLoops::CoEdge::AngleInPlane(OT_VECTEUR_3D & dir1_3D, OT_VECTEUR_3D & dir2_3D, OT_VECTEUR_3D & __normal)
292     {
293     OT_VECTEUR_3D dir1_2D, dir2_2D;
294     OT_MATRICE_3D repereNormal, transform3DToRepereNormal;
295     repereNormal = GeometricTools::GetPlaneFrame(__normal);
296     repereNormal.transpose(transform3DToRepereNormal);
297     dir1_2D=transform3DToRepereNormal*dir1_3D;
298     dir2_2D=transform3DToRepereNormal*dir2_3D;
299     dir2_2D[2]=dir1_2D[2]=0;
300     dir1_2D = dir1_2D;
301     double x=1/(dir1_2D.get_longueur()*dir2_2D.get_longueur());
302     double cs = x*(dir1_2D*dir2_2D);
303     OT_VECTEUR_3D Y=OT_VECTEUR_3D(0,0,1)&dir1_2D;
304     double sn = dir2_2D*Y;
305     if (cs>1)cs=1;else if (cs<-1)cs=-1;
306     double angle = acos (cs);
307     if (sn<0) angle = -angle;
308     return angle;
309     }
310    
311     MakeLoops::MakeLoops(std::vector <MG_FACE*> __faceList, std::vector <MG_ARETE*> __edgeList, std::vector <int> __senseList)
312     {
313     for (unsigned i=0; i<__edgeList.size(); i++)
314     {
315     MG_FACE * face = __faceList[i];
316     MG_ARETE * edge = __edgeList[i];
317     int sense = __senseList[i];
318     if (_mapFaceCoEdge.find(face) == _mapFaceCoEdge.end())
319     {
320     std::set < CoEdge * > faceCoEdges;
321     faceCoEdges.insert(new CoEdge(face, edge, sense));
322     _mapFaceCoEdge[face] = faceCoEdges;
323     }
324     else
325     {
326     _mapFaceCoEdge[face].insert(new CoEdge(face, edge, sense));
327     }
328     }
329     }
330     MakeLoops::~MakeLoops()
331     {
332     for (std::map < MG_FACE * , std::set < CoEdge * > >::iterator itLst = _mapFaceCoEdge.begin();
333     itLst != _mapFaceCoEdge.end();
334     itLst++)
335     {
336     std::set < CoEdge * > & lst = itLst->second;
337     for ( std::set < CoEdge * >::iterator itCoEdge = lst.begin(); itCoEdge != lst.end(); itCoEdge++)
338     {
339     CoEdge * coedge = *itCoEdge;
340     delete coedge;
341     }
342     }
343     }
344     std::string MakeLoops::PrintFaceNormalAtVertices(MG_FACE * __face)
345     {
346     std::stringstream out;
347     unsigned N=0;
348     MCFace * mcFace = (MCFace*) __face;
349     out << "BaseColor { rgb 0.5 0 0 }\n";
350     out << "\n Coordinate3 {\n point [ \n";
351     std::set < CoEdge * > & lst = _mapFaceCoEdge[__face];
352     for (std::set < CoEdge * >::iterator itCoEdge = lst.begin();
353     itCoEdge != lst.end();
354     itCoEdge++)
355     {
356     CoEdge * current = *itCoEdge;
357     MCVertex * vertex = (MCVertex *)current->StartVertex();
358     OT_VECTEUR_3D normal(0,0,0);
359     int nbRefFaces;
360     mcFace->calcul_normale_unitaire(vertex,normal,&nbRefFaces);
361     double xyz[3],xyz2[3];
362     vertex->get_point()->evaluer(xyz);
363     for (int i=0; i<3; i++) xyz2[i] = xyz[i]+normal[i]*.005;
364     out << xyz[0] <<" "<< xyz[1] <<" "<< xyz[2]<<",\n";
365     out << xyz2[0] <<" "<< xyz2[1] <<" "<< xyz2[2]<<",\n";
366     N+=2;
367     }
368     out << "\n]\n}\n";
369     out << "\nIndexedLineSet {";
370     out << "\ncoordIndex\n [ \n";
371     for (unsigned int j=0; j+1<N; j+=2)
372     {
373     out << j << ", ";
374     out << j+1 << ", ";
375     out << "-1,\n";
376     }
377     out << "] \n}\n";
378    
379     return out.str();
380     }
381     void MakeLoops::GetFaceLoops(MG_FACE* __face, std::vector < std::vector < CoEdge *> > & __loops )
382     {
383     std::set < CoEdge * > & lst = _mapFaceCoEdge[__face];
384     std::set < CoEdge * > unvisited = lst;
385    
386 foucault 64 if (lst.size()==0) return;
387    
388 foucault 27 bool isMCT;
389     {
390     std::string MCSTR("MC");
391     std::string idorig = __face->get_idoriginal();
392     string::size_type loc = idorig.find( MCSTR, 0 );
393     isMCT = ( loc == 0 );
394     }
395    
396     std::vector<CoEdge*> currentLoop;
397     CoEdge * current = 0;
398    
399     while (unvisited.size())
400     {
401     if (current == 0 || current->next != 0 )
402     {
403     if (currentLoop.size())
404     {
405     __loops.push_back(currentLoop);
406     currentLoop.clear();
407     }
408 francois 1075 printf("Loop %lu of Face %lu\n", __loops.size(), __face->get_id());
409 foucault 27 current = *(unvisited.begin());
410     }
411     currentLoop.push_back(current);
412     unvisited.erase(current);
413 francois 1075 printf("Edge %lu\n", current->e->get_id());
414 foucault 27
415     double score_max = -(+M_PI+.001);
416     OT_VECTEUR_3D normal;
417     if (isMCT)
418     {
419     int nbRefFaceNormal;
420     MCVertex * mcVertex = (MCVertex*) current->EndVertex();
421     MCFace * mcFace = (MCFace*)__face;
422     mcFace->calcul_normale_unitaire(mcVertex, normal, &nbRefFaceNormal);
423     }
424     else
425     {
426     MG_SOMMET * vertex = current->EndVertex();
427     double xyzVertex[3]; vertex->get_point()->evaluer(xyzVertex);
428     double uvVertex[2]; __face->inverser(uvVertex, xyzVertex);
429     __face->calcul_normale_unitaire(uvVertex, normal);
430     }
431    
432    
433     for (std::set<CoEdge*>::const_iterator it = unvisited.begin(); it != unvisited.end(); it++)
434     {
435     CoEdge * candidate = *it;
436     if (candidate->IsAfter(*current))
437     {
438     double score = current->MeshAngleInPlane(*candidate, normal,0);
439     if (score > score_max)
440     {
441     score_max = score;
442     current->next = candidate;
443     }
444     }
445     }
446    
447     current = current->next;
448     }
449    
450     __loops.push_back(currentLoop);
451     int N=0;
452     for (unsigned k=0; k<__loops.size(); k++)
453     N+=__loops[k].size();
454 francois 1075 printf("N=%d lst=%lu\n", N, lst.size());
455 foucault 27 }